Skip to content

Instantly share code, notes, and snippets.

@andanteyk
Last active October 4, 2022 03:40
Show Gist options
  • Save andanteyk/971f8fb5d2aa55f2d685413368d02078 to your computer and use it in GitHub Desktop.
Save andanteyk/971f8fb5d2aa55f2d685413368d02078 to your computer and use it in GitHub Desktop.
Seiran128 sample implementation in C# (.NET 5)
using System;
using System.Buffers;
using System.Collections;
using System.Collections.Generic;
using System.Linq;
using System.Runtime.InteropServices;
using System.Security.Cryptography;
// This code is licensed under the terms of the MIT license.
// For details, see https://gist.github.com/andanteyk/d08ab296665b3fc68df58beff3ea39cb .
namespace AndanteSoft
{
/// <summary>
/// Seiran128 pseudorandom number generator.
/// </summary>
/// <seealso cref="https://github.com/andanteyk/prng-seiran"/>
public sealed class Seiran
{
private ulong State0;
private ulong State1;
[ThreadStatic]
private static Seiran StaticInstance = null;
/// <summary>
/// Gets a static instance of <see cref="Seiran"/>.
/// </summary>
/// <remarks>
/// If you don't need reproducibility, we recommend using this.
/// This instance is thread static.
/// </remarks>
public static Seiran I => StaticInstance ??= new Seiran();
/// <summary>
/// Initializes a new instance of the <see cref="Seiran"/> class using a cryptographically secure pseudo-random number generator.
/// </summary>
/// <remarks>
/// Basically you should use this constructor.
/// </remarks>
public Seiran()
{
Span<ulong> state = stackalloc ulong[2];
var stateBytes = MemoryMarshal.Cast<ulong, byte>(state);
using (var csprng = new RNGCryptoServiceProvider())
{
do
{
csprng.GetBytes(stateBytes);
} while (state[0] == 0 && state[1] == 0);
}
State0 = state[0];
State1 = state[1];
}
/// <summary>
/// Initializes a new instance of the <see cref="Seiran"/> class using the specified <paramref name="seed"/> value.
/// </summary>
/// <param name="seed">
/// The seed value used to initialize the <see cref="Seiran"/> class.
/// If the same value is specified, the same sequence of random numbers will be generated.
/// </param>
/// <remarks>
/// This constructor is provided for compatibility with <see cref="System.Random"/>.
/// You should use <see cref="Seiran()"/> if you want a completely random instance,
/// and <see cref="Seiran(ReadOnlySpan{ulong})"/> if you need to reproduce a sequence of random numbers.
/// </remarks>
public Seiran(long seed)
{
ulong x = (ulong)seed;
static (ulong value, ulong state) splitMix(ulong state)
{
ulong z = state += 0x9e3779b97f4a7c15;
z = (z ^ z >> 30) * 0xbf58476d1ce4e5b9;
z = (z ^ z >> 27) * 0x94d049bb133111eb;
return (z ^ z >> 31, state);
}
(State0, x) = splitMix(x);
(State1, _) = splitMix(x);
}
/// <summary>
/// Initializes a new instance of the <see cref="Seiran"/> class using the specified <paramref name="state"/>.
/// </summary>
/// <param name="state">
/// Internal state obtained by <see cref="GetState"/>.
/// </param>
/// <exception cref="ArgumentOutOfRangeException">
/// The length of <paramref name="state"/> is less than 2, or all elements of <paramref name="state"/> are 0.
/// </exception>
public Seiran(ReadOnlySpan<ulong> state)
{
if (state.Length < 2)
throw new ArgumentOutOfRangeException(nameof(state));
if (state[0] == 0 && state[1] == 0)
throw new ArgumentOutOfRangeException(nameof(state));
State0 = state[0];
State1 = state[1];
}
/// <summary>
/// Initialize a new instance of the <see cref="Seiran"/> class that will generate the same sequence of random numbers as the <paramref name="original"/> instance.
/// </summary>
/// <param name="original">Copy source instance.</param>
/// <exception cref="ArgumentNullException"><paramref name="original"/> is null.</exception>
public Seiran(Seiran original)
{
if (original == null)
throw new ArgumentNullException(nameof(original));
State0 = original.State0;
State1 = original.State1;
}
/// <summary>
/// Gets the internal state of this <see cref="Seiran"/> instance.
/// </summary>
/// <returns>Internal state. The element length is 2.</returns>
/// <remarks>
/// You can use <see cref="Seiran(ReadOnlySpan{ulong})"/> to restore the current state of this instance.
/// </remarks>
public ulong[] GetState() => new ulong[] { State0, State1 };
/// <summary>
/// Generates a 64-bit random number.
/// </summary>
/// <returns>[0, <see cref="ulong.MaxValue"/>]</returns>
public ulong Next()
{
static ulong rotl(ulong x, int k) => x << k | x >> -k;
ulong s0 = State0, s1 = State1;
ulong result = rotl((s0 + s1) * 9, 29) + s0;
State0 = s0 ^ rotl(s1, 29);
State1 = s0 ^ s1 << 9;
return result;
}
public uint NextUintStrict(uint maxExclusive)
{
var hi = Math.BigMul(Next(), maxExclusive, out var lo);
if (lo < maxExclusive)
{
ulong lowerBound = (0ul - maxExclusive) % maxExclusive;
while (lo < lowerBound)
{
hi = Math.BigMul(Next(), maxExclusive, out lo);
}
}
return (uint)hi;
}
/// <summary>
/// Generates a non-negative random number.
/// </summary>
/// <returns>[0, <see cref="int.MaxValue"/>]</returns>
public int NextInt() => (int)(Next() >> 33);
/// <summary>
/// Generates a random number smaller than the specified maximum.
/// </summary>
/// <param name="maxExclusive">
/// An exclusive upper limit for the generated random numbers.
/// <paramref name="maxExclusive"/> must be greater than 0.
/// </param>
/// <returns>[0, <paramref name="maxExclusive"/>)</returns>
/// <exception cref="ArgumentOutOfRangeException">
/// <paramref name="maxExclusive"/> is less than or equal to 0.
/// </exception>
public int NextInt(int maxExclusive)
{
if (maxExclusive <= 0)
throw new ArgumentOutOfRangeException(nameof(maxExclusive));
return (int)NextUintStrict((uint)maxExclusive);
}
/// <summary>
/// Generates a random number in the specified range.
/// </summary>
/// <param name="minInclusive"> An inclusive lower limit for the generated random numbers. </param>
/// <param name="maxExclusive"> An exclusive upper limit for the generated random numbers. <paramref name="maxExclusive"/> must be greater than <paramref name="minInclusive"/>. </param>
/// <returns>[<paramref name="minInclusive"/>, <paramref name="maxExclusive"/>)</returns>
/// <exception cref="ArgumentOutOfRangeException">
/// <paramref name="minInclusive"/> is greater than or equal to <paramref name="maxExclusive"/>.
/// </exception>
public int NextInt(int minInclusive, int maxExclusive)
{
if (minInclusive >= maxExclusive)
throw new ArgumentOutOfRangeException(nameof(minInclusive));
return (int)NextUintStrict((uint)(maxExclusive - minInclusive)) + minInclusive;
}
/// <summary>
/// Generates a floating-point random number greater than or equal to 0.0 and less than 1.0.
/// </summary>
/// <returns>[0.0, 1.0)</returns>
/// <remarks>
/// The actual minimum value greater than 0 is 1.1102230246251565E-16 and the actual maximum value is 0.99999999999999989.
/// </remarks>
public double NextDouble() => (Next() >> 11) * (1.0 / (1ul << 53));
/// <summary>
/// Generates a floating-point random number in the specified range.
/// </summary>
/// <param name="minInclusive"> An inclusive lower limit for the generated random numbers. </param>
/// <param name="maxExclusive"> An exclusive upper limit for the generated random numbers. <paramref name="maxExclusive"/> must be greater than <paramref name="minInclusive"/>. </param>
/// <returns>[<paramref name="minInclusive"/>, <paramref name="maxExclusive"/>)</returns>
/// <exception cref="ArgumentOutOfRangeException">
/// <paramref name="minInclusive"/> is greater than or equal to <paramref name="maxExclusive"/>.
/// </exception>
public double NextDouble(double minInclusive, double maxExclusive)
{
if (minInclusive >= maxExclusive)
throw new ArgumentOutOfRangeException(nameof(minInclusive));
double t = NextDouble();
double ret = minInclusive * (1 - t) + maxExclusive * t;
if (ret == maxExclusive && double.IsFinite(ret))
ret = Math.BitDecrement(ret);
return ret;
}
/// <summary>
/// Fills the elements of the specified array with random numbers.
/// </summary>
/// <param name="bytes">An array to be filled with random numbers.</param>
/// <exception cref="ArgumentNullException">
/// <paramref name="bytes"/> is null.
/// </exception>
/// <remarks>
/// The range of random numbers is [0x00, 0xff].
/// </remarks>
public void NextBytes(byte[] bytes)
{
if (bytes == null)
throw new ArgumentNullException(nameof(bytes));
NextBytes(bytes.AsSpan());
}
/// <summary>
/// Fills the elements of the specified span with random numbers.
/// </summary>
/// <param name="buffer">A span to be filled with random numbers.</param>
/// <remarks>
/// The range of random numbers is [0x00, 0xff].
/// If you want to apply to another type, use <see cref="MemoryMarshal.Cast{TFrom, TTo}(Span{TFrom})"/>.
/// </remarks>
public void NextBytes(Span<byte> buffer)
{
var ulongSpan = MemoryMarshal.Cast<byte, ulong>(buffer);
for (int i = 0; i < ulongSpan.Length; i++)
ulongSpan[i] = Next();
if ((buffer.Length & 7) != 0)
{
ulong r = Next();
for (int i = buffer.Length & ~7; i < buffer.Length; i++)
{
buffer[i] = (byte)r;
r >>= 8;
}
}
}
/// <summary>
/// Generates independent pairs of random numbers that follow a normal distribution.
/// </summary>
/// <returns>Independent pairs of random numbers that follow a normal distribution with mean 0.0 and standard deviation 1.0.</returns>
public (double x, double y) NextGaussian()
{
double u, v, s;
do
{
u = NextDouble(-1, 1);
v = NextDouble(-1, 1);
s = u * u + v * v;
} while (s >= 1 || s == 0);
s = Math.Sqrt(-2 * Math.Log(s) / s);
return (u * s, v * s);
}
/// <summary>
/// Shuffles the order of the elements.
/// </summary>
/// <typeparam name="T">The type of the elements of <paramref name="list"/>.</typeparam>
/// <param name="list">List to be shuffled.</param>
/// <exception cref="ArgumentNullException"> <paramref name="list"/> is null. </exception>
public void ShuffleList<T>(IList<T> list)
{
if (list == null)
throw new ArgumentNullException(nameof(list));
for (int i = list.Count - 1; i > 0; i--)
{
int r = NextInt(i + 1);
(list[i], list[r]) = (list[r], list[i]);
}
}
/// <summary>
/// Shuffles the order of the elements.
/// </summary>
/// <typeparam name="T">The type of the elements of <paramref name="span"/>.</typeparam>
/// <param name="span">List to be shuffled.</param>
public void ShuffleSpan<T>(Span<T> span)
{
for (int i = span.Length - 1; i > 0; i--)
{
int r = NextInt(i + 1);
(span[i], span[r]) = (span[r], span[i]);
}
}
/// <summary>
/// Shuffles the order of the elements.
/// </summary>
/// <typeparam name="T">The type of the elements of <paramref name="source"/>.</typeparam>
/// <param name="source">Sequence to be shuffled.</param>
/// <returns>A sequence in which the order of the elements is randomly changed.</returns>
/// <exception cref="ArgumentNullException"> <paramref name="source"/> is null. </exception>
/// <remarks>
/// This method is implemented using deferred execution.
/// The <paramref name="source"/> is not evaluated and shuffled until GetEnumerator() is called.
/// </remarks>
public IEnumerable<T> Shuffle<T>(IEnumerable<T> source)
{
if (source == null)
throw new ArgumentNullException(nameof(source));
return new ShuffleIterator<T>(this, source);
}
public sealed class ShuffleIterator<T> : IEnumerable<T>, IEnumerator<T>
{
private readonly IEnumerable<T> Source;
private readonly Seiran RandomSource;
private int CurrentIndex;
private T[] Elements;
internal ShuffleIterator(Seiran randomSource, IEnumerable<T> source)
{
RandomSource = randomSource ?? throw new ArgumentNullException(nameof(randomSource));
Source = source ?? throw new ArgumentNullException(nameof(source));
CurrentIndex = -1;
}
public T Current
{
get
{
if ((uint)CurrentIndex < Elements?.Length)
return Elements[CurrentIndex];
else
return default;
}
}
object IEnumerator.Current => Current;
public void Dispose() { }
public bool MoveNext()
{
if (CurrentIndex >= Elements.Length)
return false;
return ++CurrentIndex < Elements.Length;
}
public void Reset() => throw new NotSupportedException();
public IEnumerator<T> GetEnumerator()
{
if (Elements == null)
{
Elements = Source.ToArray();
RandomSource.ShuffleSpan(Elements.AsSpan());
}
if (CurrentIndex == -1)
return this;
var clone = MemberwiseClone() as ShuffleIterator<T>;
clone.CurrentIndex = -1;
return clone;
}
IEnumerator IEnumerable.GetEnumerator() => GetEnumerator();
}
private void Jump(ReadOnlySpan<ulong> jumpPolynomial)
{
ulong t0 = 0, t1 = 0;
for (int i = 0; i < 2; i++)
{
for (int b = 0; b < 64; b++)
{
if (((jumpPolynomial[i] >> b) & 1) != 0)
{
t0 ^= State0;
t1 ^= State1;
}
Next();
}
}
State0 = t0;
State1 = t1;
}
/// <summary>
/// Advance the internal state by 2^32 steps.
/// </summary>
/// <remarks>
/// This method is equivalent to 2^32 <see cref="Next"/> calls.
/// It can be executed in the same amount of time as 128 <see cref="Next"/> calls.
/// </remarks>
public void Jump32() => Jump(stackalloc ulong[] { 0x40165CBAE9CA6DEB, 0x688E6BFC19485AB1 });
/// <summary>
/// Advance the internal state by 2^64 steps.
/// </summary>
/// <remarks>
/// This method is equivalent to 2^64 <see cref="Next"/> calls.
/// It can be executed in the same amount of time as 128 <see cref="Next"/> calls.
/// </remarks>
public void Jump64() => Jump(stackalloc ulong[] { 0xF4DF34E424CA5C56, 0x2FE2DE5C2E12F601 });
/// <summary>
/// Advance the internal state by 2^96 steps.
/// </summary>
/// <remarks>
/// This method is equivalent to 2^96 <see cref="Next"/> calls.
/// It can be executed in the same amount of time as 128 <see cref="Next"/> calls.
/// </remarks>
public void Jump96() => Jump(stackalloc ulong[] { 0x185F4DF8B7634607, 0x95A98C7025F908B2 });
/// <summary>
/// Rewinds the internal state to the previous state.
/// </summary>
public void Prev()
{
static ulong rotl(ulong x, int k) => x << k | x >> -k;
ulong t1 = rotl(State0 ^ State1, 64 - 29);
t1 ^= t1 << 44 ^ (t1 & ~0xffffful) << 24;
t1 ^= (t1 & (0x7ffful << 40)) << 4;
t1 ^= (t1 & (0x7fful << 40)) << 8;
t1 ^= (t1 & (0x7ul << 40)) << 16;
t1 ^= (t1 & (0xffffful << 35)) >> 20;
t1 ^= (t1 & (0x7ffful << 20)) >> 20;
State0 ^= rotl(t1, 29);
State1 = t1;
}
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment