On 06/18/2013 12:09 PM, Joseph Rushton Wakeling wrote:
> ... we get different sequences. However, the .save as implemented needs
> rewriting (I haven't touched it yet) as currently it simply returns a
> reference
> copy. I'll get onto that. :-)
Tweaked version attached. This adds a new constructor which copies a PRNG of
the same type:
this(typeof(this) source)
{
this.mt[] = source.mt[];
this.mti = source.mti;
this._y = source._y;
}
... and .save then makes use of this to return a duplicate:
@property typeof(this) save()
{
return new typeof(this)(this);
}
With this in hand, we can do,
auto gen = new MtClass19937(unpredictableSeed);
auto t1 = gen.take(5);
auto t2 = gen.take(5);
writeln(t1);
writeln(t2);
... and get two independent sequences, whereas,
auto t3 = gen.take(5);
auto t4 = gen.save.take(5);
writeln(t3);
writeln(t4);
... produces two identical sequences.
This still requires some intelligence/virtue on behalf of the programmer,
though: if we do,
auto t5 = gen.save.take(5);
writeln(t5);
... and don't also take 5 from gen itself, then subsequent uses of gen will be
correlated with t5, probably unintentionally. However, I think that's a fair
tradeoff for the ability to .save at all.
It's also now trivial to rewrite things like my SimpleRandomRange such that,
auto gen = new MtClass19937(unpredictableSeed);
auto r = simpleRandomRange(0.0L, 1.0L, gen);
auto t1 = r5.take(5);
auto t2 = r5.take(5);
writeln(t1);
writeln(t2);
... produces different sequences, whereas,
auto t3 = r5.take(5);
auto t4 = r5.save.take(5);
writeln(t3);
writeln(t4);
... produce the same. Ditto for "real" random ranges, such as
std.random.RandomSample.
So we have a solution to the first set of problems identified, that is much
nicer than the one I had in mind. Very good!
However, at least one more niggle remains -- to be followed up on.
module mtclass;
import std.range, std.traits;
final class MersenneTwisterClass(UIntType, size_t w, size_t n, size_t m, size_t r,
UIntType a, size_t u, size_t s,
UIntType b, size_t t,
UIntType c, size_t l)
if(isUnsigned!UIntType)
{
///Mark this as a Rng
enum bool isUniformRandom = true;
/**
Parameter for the generator.
*/
enum size_t wordSize = w;
enum size_t stateSize = n;
enum size_t shiftSize = m;
enum size_t maskBits = r;
enum UIntType xorMask = a;
enum UIntType temperingU = u;
enum size_t temperingS = s;
enum UIntType temperingB = b;
enum size_t temperingT = t;
enum UIntType temperingC = c;
enum size_t temperingL = l;
/// Smallest generated value (0).
enum UIntType min = 0;
/// Largest generated value.
enum UIntType max =
w == UIntType.sizeof * 8 ? UIntType.max : (1u << w) - 1;
/// The default seed value.
enum UIntType defaultSeed = 5489u;
static assert(1 <= m && m <= n);
static assert(0 <= r && 0 <= u && 0 <= s && 0 <= t && 0 <= l);
static assert(r <= w && u <= w && s <= w && t <= w && l <= w);
static assert(0 <= a && 0 <= b && 0 <= c);
static assert(a <= max && b <= max && c <= max);
/**
Constructs a MersenneTwisterEngine object.
*/
this(UIntType value)
{
seed(value);
}
this(typeof(this) source)
{
this.mt[] = source.mt[];
this.mti = source.mti;
this._y = source._y;
}
/**
Seeds a MersenneTwisterEngine object.
Note:
This seed function gives 2^32 starting points. To allow the RNG to be started in any one of its
internal states use the seed overload taking an InputRange.
*/
void seed()(UIntType value = defaultSeed)
{
static if (w == UIntType.sizeof * 8)
{
mt[0] = value;
}
else
{
static assert(max + 1 > 0);
mt[0] = value % (max + 1);
}
for (mti = 1; mti < n; ++mti)
{
mt[mti] =
cast(UIntType)
(1812433253UL * (mt[mti-1] ^ (mt[mti-1] >> (w - 2))) + mti);
/* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
/* In the previous versions, MSBs of the seed affect */
/* only MSBs of the array mt[]. */
/* 2002/01/09 modified by Makoto Matsumoto */
//mt[mti] &= ResultType.max;
/* for >32 bit machines */
}
popFront();
}
/**
Seeds a MersenneTwisterEngine object using an InputRange.
Throws:
$(D Exception) if the InputRange didn't provide enough elements to seed the generator.
The number of elements required is the 'n' template parameter of the MersenneTwisterEngine struct.
Examples:
----------------
Mt19937 gen;
gen.seed(map!((a) => unpredictableSeed)(repeat(0)));
----------------
*/
void seed(T)(T range) if(isInputRange!T && is(Unqual!(ElementType!T) == UIntType))
{
size_t j;
for(j = 0; j < n && !range.empty; ++j, range.popFront())
{
mt[j] = range.front;
}
mti = n;
if(range.empty && j < n)
{
throw new Exception(format("MersenneTwisterEngine.seed: Input range didn't provide enough"
" elements: Need %s elemnets.", n));
}
popFront();
}
/**
Advances the generator.
*/
void popFront()
{
if (mti == size_t.max) seed();
enum UIntType
upperMask = ~((cast(UIntType) 1u <<
(UIntType.sizeof * 8 - (w - r))) - 1),
lowerMask = (cast(UIntType) 1u << r) - 1;
static immutable UIntType mag01[2] = [0x0UL, a];
ulong y = void;
if (mti >= n)
{
/* generate N words at one time */
int kk = 0;
const limit1 = n - m;
for (; kk < limit1; ++kk)
{
y = (mt[kk] & upperMask)|(mt[kk + 1] & lowerMask);
mt[kk] = cast(UIntType) (mt[kk + m] ^ (y >> 1)
^ mag01[cast(UIntType) y & 0x1U]);
}
const limit2 = n - 1;
for (; kk < limit2; ++kk)
{
y = (mt[kk] & upperMask)|(mt[kk + 1] & lowerMask);
mt[kk] = cast(UIntType) (mt[kk + (m -n)] ^ (y >> 1)
^ mag01[cast(UIntType) y & 0x1U]);
}
y = (mt[n -1] & upperMask)|(mt[0] & lowerMask);
mt[n - 1] = cast(UIntType) (mt[m - 1] ^ (y >> 1)
^ mag01[cast(UIntType) y & 0x1U]);
mti = 0;
}
y = mt[mti++];
/* Tempering */
y ^= (y >> temperingU);
y ^= (y << temperingS) & temperingB;
y ^= (y << temperingT) & temperingC;
y ^= (y >> temperingL);
_y = cast(UIntType) y;
}
/**
Returns the current random value.
*/
@property UIntType front()
{
if (mti == size_t.max) seed();
return _y;
}
///
@property typeof(this) save()
{
return new typeof(this)(this);
}
/**
Always $(D false).
*/
enum bool empty = false;
private UIntType mt[n];
private size_t mti = size_t.max; /* means mt is not initialized */
UIntType _y = UIntType.max;
}
alias MersenneTwisterClass!(uint, 32, 624, 397, 31, 0x9908b0df, 11, 7,
0x9d2c5680, 15, 0xefc60000, 18)
MtClass19937;
