Search All of the Math Forum:

Views expressed in these public forums are not endorsed by NCTM or The Math Forum.

Notice: We are no longer accepting new posts, but the forums will continue to be readable.

Topic: The CMWC4827 RNG, an improvement on MWC4691.
Replies: 16   Last Post: Oct 30, 2010 7:15 PM

 Messages: [ Previous | Next ]
 geo Posts: 38 Registered: 6/19/08
Re: The CMWC4827 RNG, an improvement on MWC4691.
Posted: Oct 26, 2010 6:43 PM

On Oct 25, 7:21 pm, "christian.bau" <christian....@cbau.wanadoo.co.uk>
wrote:
> On Oct 22, 12:18 pm, geo <gmarsag...@gmail.com> wrote:
>

> > CMWC4827 RNG snipped
> > George Marsaglia

>
> I tried to figure out the maths behind this, and I ended up with a
> just slightly different algorithm:
>
> The important number is p = 4095 * (2^32)^4827 + 1, which is a prime.
> Let M = (2^32)^4827, then p = 4095 * M + 1. We can start with any
> number 1 <= x < p, then replace x with (4095 x) modulo p, and again
> and again. Each value x produces 4827 random integers of 32 bits. The
> result x will never be 0.
>
> A slight change: Instead of storing x, we store x-1. So instead of
> replacing x with (4095 x) modulo p, we add 1 to x, calculate (4095 x)
> modulo p, and subtract 1. So we calculate ((4095 (x+1)) modulo p) - 1
> = (4095 x + 4094) modulo p. Really the same thing, but we now have 0
> <= x < 4095 * M.
>
> Now assume x = X + c * M, where 0 <= X < M, and 0 <= c < 4095. Then
>
>     (4095 x + 4094) modulo p =
>     (4095 * (X + c * M) + 4094) modulo p =
>     (4095 * X + 4094 + c * (4095 M)) modulo p =
>     (4095 * X + 4094 + c * (p - 1)) modulo p =
>     (4095 * X + 4094 - c) modulo p =
>     4095 * X + (4094 - c).
>
> We store X in an array of 4827 32-bit integers and keep c separate.
> Initialiase X to random integers and c to a random value from 0 to
> 4094. The random number calculation then is
>
> unsigned long   CMWC4827(void)
> {
>     unsigned long   t,
>                     x;
>     static int      j = 4827;
>     if (j < 4826) ++j; else { j = 0; c = 4094 - c; }
>     x = Q[j];
>     t = (x << 12) + carry;
>     carry = (x >> 20) - (t < x);
>     return (Q[j] = (t - x));
>
> }
>
> So the change is that we replace c with 4094-c when we restart at the
> beginning of the array, and don't do the "not" operation in the
> calculation of the next component of x. You could change this to use
> 64 bit arithmetic and write
>
> uint32_t CMWC4827 (void)
> {
>     static long j = 4827;
>     static uint64_t carry = 1271;
>     if (++j >= 4827) { j = 0; carry = 4094 - carry; }
>     uint64_t x = Q [j] * (uint64_t) 4095 + carry;
>     carry = x >> 32;
>     return (Q [j] = x);
>
>
>
> }- Hide quoted text -
>
> - Show quoted text -

With b=2^32 and p=4095*b^4827+1, your suggestion
will generate a sequence of 4827 "digits" forming,
in reverse order, the base-b expansion of j1/(p-2)
the expansion of j2/(p-2) for some other j2, and so on.
Such "digits" may well serve as satisfactory 32-bit
random integers, but we cannot be sure of the period,
or that we will avoid partial overlapping with string

If you change your suggested method so as to keep
the old carry c---rather than reset c before
refilling the 4827 array---and go ahead with the
MWC rather than CMWC method, you will generate,
in reverse order, the full base-b MWC expansion of
j/m, with m=4095*b^4827-1 and period the order of b
for the composite m.

The CMWC method produces, in reverse order,
the full base-b=2^32 digits in the expansion of
some j/p, with p=4095*b^4827+1 and j determined by
the 4827 seed values and initial c in 0<=c<4095.

The order of b for the MWC method applied to m
probably differs from the CMWC order, 4095*2^154458,
by at most a few dozen powers of 10, and, with 10^46500,
we have 46,500 in hand.

But finding that IGO, (Industrial-Grade-Order), which
amounts to finding all, or all the up-to-30-digit
prime factors of m=4095*b^4827-1, say p1,p2,...,pn,
then
IGO = lcm(order(b,p1),order(b,p2)...,order(b,pn))
is difficult.

I have considered that approach, but using
m=a*b^4096-1, or even m=a*b^8192-1 with a=2^i-1,
so that one can more easily effect incrementing
the index of the array element.
With a couple of my PCs working on finding
up-to-30-digit factors of such composites,
I hope to report some interesting IGOs soon.

Can any of you, with time and factoring
programs available, find all---or all the
up-to-30-digit---prime factors of a*b^4096-1
or a*b^8192-1 with b=2^32 and a=2^i-1, for MWC?

Or the same for CMWC primes or composites,
a*b^4096+1 or a*b^8192+1?

George Marsaglia

Date Subject Author
10/22/10 geo
10/22/10 Dann Corbit
10/22/10 Jens Thoms Toerring
10/22/10 Ian
10/22/10 Jens Thoms Toerring
10/22/10 Ian
10/23/10 gnasher729
10/23/10 Ilmari Karonen
10/23/10 gnasher729
10/23/10 Jens Thoms Toerring
10/23/10 Jens Thoms Toerring
10/22/10 Ian
10/25/10 gnasher729
10/26/10 geo
10/30/10 gnasher729
10/27/10 John Smith
10/30/10 rossum