1000/1000

Hot
Most Recent

Submitted Successfully!

To reward your contribution, here is a gift for you: A free trial for our video production service.

Thank you for your contribution! You can also upload a video entry or images related to this topic.

Do you have a full video?

Are you sure to Delete?

Cite

If you have any further questions, please contact Encyclopedia Editorial Office.

HandWiki. Multiply-with-carry. Encyclopedia. Available online: https://encyclopedia.pub/entry/27990 (accessed on 14 August 2024).

HandWiki. Multiply-with-carry. Encyclopedia. Available at: https://encyclopedia.pub/entry/27990. Accessed August 14, 2024.

HandWiki. "Multiply-with-carry" *Encyclopedia*, https://encyclopedia.pub/entry/27990 (accessed August 14, 2024).

HandWiki. (2022, September 29). Multiply-with-carry. In *Encyclopedia*. https://encyclopedia.pub/entry/27990

HandWiki. "Multiply-with-carry." *Encyclopedia*. Web. 29 September, 2022.

Copy Citation

In computer science, multiply-with-carry (MWC) is a method invented by George Marsaglia for generating sequences of random integers based on an initial set from two to many thousands of randomly chosen seed values. The main advantages of the MWC method are that it invokes simple computer integer arithmetic and leads to very fast generation of sequences of random numbers with immense periods, ranging from around 260 to 22000000. As with all pseudorandom number generators, the resulting sequences are functions of the supplied seed values.

multiply-with-carry
pseudorandom number
random numbers

An MWC generator is a special form of Lehmer random number generator *x _{n}* =

Normal Lehmer generator implementations choose a modulus close to the machine word size. An MWC generator instead maintains its state in base *b*, so multiplying by *b* is done implicitly by shifting one word. The MWC multiplier *a* and lag *r* determine the modulus *p* = *ab ^{r}*−1.

The base *b* is typically chosen to equal the computer's word size, as this makes arithmetic modulo *b* trivial. This may vary from *b* = 2^{8} for a microcontroller to *b* = 2^{64}. (This article uses *b* = 2^{32} for examples.)

Sometimes an odd base is preferred, in which case *b* = 2^{k} − 1 can be used, which is almost as simple to implement.

The period of a lag-*r* MWC generator is the order of *b* in the multiplicative group of numbers modulo *p* = *ab ^{r}* − 1. While it is theoretically possible to choose a non-prime modulus, a prime modulus eliminates the possibility of the initial seed sharing a common divisor with the modulus, which would reduce the generator's period.

Because 2 is a quadratic residue of numbers of the form 8*k*±1, *b* = 2^{k} cannot be a primitive root of *p = ab ^{r}* − 1. Therefore MWC generators with base 2

The basic form of an MWC generator has parameters *a*, *b* and *r*, and *r*+1 words of state. The state consists of *r* residues modulo *b*

- 0 ≤
*x*_{0},*x*_{1},*x*_{2},...,*x*_{r−1}< b,

and a carry *c*_{r−1} < *a*.

The initial state ("seed") values are arbitrary, except that they must not be all zero, nor all at the maximum permitted values (*x _{i}* =

The lag-*r* MWC sequence is then a sequence of pairs *x _{n}*,

- [math]\displaystyle{ x_n=(ax_{n-r}+c_{n-1})\,\bmod\,b,\ c_n=\left\lfloor\frac{ax_{n-r}+c_{n-1}}{b}\right\rfloor }[/math]

and the MWC generator output is the sequence of *x*'s,

*x*,_{r}*x*_{r+1},*x*_{r+2}, ...

Although the theory of MWC generators permits *a* > *b*, *a* is almost always chosen smaller for convenience of implementation.

The state transformation function of an MWC generator is one step of Montgomery reduction modulo *p*. The state is a large integer with most significant word *c*_{n−1} and least significant word *x*_{n−r}. Each step, *x*_{n−r}·(*ab ^{r}*−1) is added to this integer. This is done in two parts: −1·

So far, this has simply added a multiple of *p* to the state, resulting in a different representative of the same residue class modulo *p*. But finally, the state is shifted down one word, dividing by *b*. This discards the least significant word of zero (which, in practice, is never computed at all) and effectively multiplies the state by *b*^{−1} (mod *p*).

Thus, a multiply-with-carry generator *is* a Lehmer generator with modulus *p* and multiplier *b*^{−1} (mod *p*). This is the same as a generator with multiplier *b*, but producing output in reverse order, which does not affect the quality of the resultant pseudorandom numbers.

This in turn means that the mathematical techniques developed for such generators (such as the spectral test) can be applied to multiply-with-carry generators.

A linear congruential generator with base *b* = 2^{32} is implemented as

- [math]\displaystyle{ x_{n+1}=(ax_n+c)\ \bmod\,2^{32}, }[/math]

where *c* is a constant. If *a* ≡ 1 (mod 4) and *c* is odd, the resulting base-2^{32} congruential sequence will have period 2^{32}.^{[1]}

This can be computed using only the low 32 bits of the product of *a* and the current *x*. However, many microprocessors can compute a full 64-bit product in almost the same time as the low 32 bits. Indeed, many compute the 64-bit product and then ignore the high half.

A lag-1 multiply-with-carry generator allows us to make the period nearly 2^{63} by using almost the same computer operations, except that the top half of the 64-bit product is not ignored after the product is computed. Instead, a 64-bit sum is computed, and the top half is used as a *new* carry value *c* rather than the fixed additive constant of the standard congruential sequence: Compute *ax*+*c* in 64 bits, then use the top half as the new *c*, and bottom half as the new *x*.

With multiplier *a* specified, each pair of input values *x*, *c* is converted to a new pair,

- [math]\displaystyle{ x\leftarrow (ax+c)\,\bmod\,2^{32},\ \ c\leftarrow \left\lfloor\frac{ax+c}{2^{32}}\right\rfloor. }[/math]

If *x* and *c* are not both zero, then the period of the resulting multiply-with-carry sequence will be the order of *b* = 2^{32} in the multiplicative group of residues modulo *ab ^{r}* − 1, that is, the smallest

If *p* = *ab ^{r}* − 1 is prime, then Fermat's little theorem guarantees that the order of any element must divide

More specifically, in such a case, the order of *any* element divides *p* − 1, and there are only four possible divisors: 1, 2, *ab ^{r}*/2 − 1, or

Following are some maximal values of *a* for computer applications which satisfy the above safe prime condition, for lag-1 generators:

Bits in a |
b |
Maximum a such that ab−1 is a safe prime |
Period |
---|---|---|---|

15 | 2^{16} |
2^{15}−50 = 32,718 |
1,072,103,423 |

16 | 2^{16} |
2^{16}−352 = 65,184 |
2,135,949,311 |

31 | 2^{32} |
2^{31}−563 = 2,147,483,085 |
4,611,684,809,394,094,079 |

32 | 2^{32} |
2^{32}−178 = 4,294,967,118 |
9,223,371,654,602,686,463 |

64 | 2^{64} |
2^{64}−742 = 18,446,744,073,709,550,874 |
2^{63}(2^{64}−742)−1 ≈ 1.701×10^{38} |

128 | 2^{128} |
2^{128}−10,408 |
2^{127}(2^{128}−10,408)−1 ≈ 5.790×10^{76} |

256 | 2^{256} |
2^{256}−9166 |
2^{255}(2^{256}−9166)−1 ≈ 6.704×10^{153} |

512 | 2^{512} |
2^{512}−150,736 |
2^{511}(2^{512}−150,736)−1 ≈ 8.988×10^{307} |

While a safe prime ensures that almost *any* element of the group has a large order, the period of the generator is specifically the order of *b*. For small moduli, more computationally expensive methods can be used to find multipliers *a* where the period is *ab*/2 − 1. The following are again maximum values of *a* of various sizes.

Bits in a |
b^{r} |
Maximum a such thatb has order ab/2−1^{r} |
Period |
---|---|---|---|

8 | 2^{8} |
2^{8}−7 = 249 |
31,871 |

8 | 2^{82} = 2^{16} |
2^{8}−32 = 224 |
7,340,031 |

15 | 2^{16} |
2^{15}−29 = 32,739 |
1,072,791,551 |

16 | 2^{16} |
2^{16}−22 = 65,514 |
2,146,762,751 |

8 | 2^{84} = 2^{32} |
2^{8}−64 = 192 |
412,316,860,415 |

15 | 2^{162} = 2^{32} |
2^{15}−26 = 32,742 |
70,312,909,602,815 |

16 | 2^{162} = 2^{32} |
2^{16}−2 = 65,534 |
140,733,193,388,031 |

31 | 2^{32} |
2^{31}−68 = 2,147,483,580 |
4,611,685,872,398,499,839 |

32 | 2^{32} |
2^{32}−76 = 4,294,967,220 |
9,223,371,873,646,018,559 |

8 | 2^{88} = 2^{64} |
2^{8}−41 = 215 |
2^{63}(2^{8}−41)−1 ≈ 1.983×10^{21} |

15 | 2^{164} = 2^{64} |
2^{15}−50 = 32,718 |
2^{63}(2^{15}−50)−1 ≈ 3.018×10^{23} |

16 | 2^{164} = 2^{64} |
2^{16}−56 = 65,480 |
2^{63}(2^{16}−56)−1 ≈ 6.039×10^{23} |

31 | 2^{322} = 2^{64} |
2^{31}−38 = 2,147,483,610 |
2^{63}(2^{31}−38)−1 ≈ 1.981×10^{28} |

32 | 2^{322} = 2^{64} |
2^{32}−43 = 4,294,967,253 |
2^{63}(2^{32}−43)−1 ≈ 3.961×10^{28} |

63 | 2^{64} |
2^{63}−140 = 9,223,372,036,854,775,668 |
2^{63}(2^{63}−140)−1 ≈ 8.507×10^{37} |

64 | 2^{64} |
2^{64}−116 = 18,446,744,073,709,551,500 |
2^{63}(2^{64}−116)−1 ≈ 1.701×10^{38} |

The output of a multiply-with-carry generator is equivalent to the radix-*b* expansion of a fraction with denominator *p* = *ab ^{r}* − 1. Here is an example for the simple case of

Starting with [math]\displaystyle{ c_0=1,x_0=0 }[/math], the MWC sequence

- [math]\displaystyle{ x_n=(7x_{n-1}+c_{n-1})\,\bmod\,10,\ c_n=\left\lfloor\frac{7x_{n-1}+c_{n-1}}{10}\right\rfloor, }[/math]

produces this sequence of states:

- 10,01,07,49,67,55,40,04,28,58,61,13,22,16,43,25,37,52,19,64,34,31, 10,01,07,...

with period 22. Consider just the sequence of *x _{i}:*

- 0,1,7,9,7,5,0,4,8,8,1,3,2,6,3,5,7,2,9,4,4,1, 0,1,7,9,7,5,0,...

Notice that if those repeated segments of *x* values are put *in reverse order*:

- [math]\displaystyle{ 1449275\cdots 9710\,1449275\cdots 9710\,144\cdots }[/math]

we get the expansion *j*/(*ab*−1) with *a*=7, *b*=10, *j=10*:

- [math]\displaystyle{ \frac{10}{69}=0.1449275362318840579710\,14492753623\ldots }[/math]

This is true in general: The sequence of *x*s produced by a lag-*r* MWC generator:

- [math]\displaystyle{ x_n=(ax_{n-r}+c_{n-1})\bmod\,b\,,\ \ c_n=\left\lfloor\frac{ax_{n-r}+c_{n-1}}{b}\right\rfloor, }[/math]

when put in reverse order, will be the radix-*b* expansion of a fraction *j*/(*ab ^{r}* − 1) for some 0 <

Continuing the above example, if we start with [math]\displaystyle{ x_0=34 }[/math], and generate the ordinary congruential sequence

- [math]\displaystyle{ x_n=7x_{n-1}\,\bmod\,69 }[/math],

we get the period 22 sequence

- 31,10,1,7,49,67,55,40,4,28,58,61,13,22,16,43,25,37,52,19,64,34, 31,10,1,7,...

and that sequence, reduced mod 10, is

- 1,0,1,7,9,7,5,0,4,8,8,1,3,2,6,3,5,7,2,9,4,4, 1,0,1,7,9,7,5,0,...

the same sequence of *x'*s resulting from the MWC sequence.

This is true in general, (but apparently only for lag-1 MWC sequences):

Given initial values [math]\displaystyle{ x_0,c_0 }[/math], the sequence [math]\displaystyle{ x_1,x_2,\ldots }[/math] resulting from the lag-1 MWC sequence

- [math]\displaystyle{ x_n=(ax_{n-1}+c_{n-1})\,\bmod b\,,\ \ c_n=\left\lfloor\frac{ax_{n-1}+c_{n-1}}{b}\right\rfloor }[/math]

is exactly the Lehmer random number generator output sequence *y _{n}* =

Choosing a different initial value *y*_{0} merely rotates the cycle of *x'*s.

Establishing the period of a lag-*r* MWC generator usually entails choosing multiplier *a* so that *p*=*ab ^{r}* − 1 is prime. Then

However, the procedure also permits a *negative* multiplier. This leads to a slight modification of the MWC procedure, and produces a modulus *p* = |−*ab ^{r}* − 1| =

The modified procedure is called complementary-multiply-with-carry (CMWC), and the setup is the same as that for lag-*r* MWC: multiplier *a*, base *b*, and *r* + 1 seeds,

*x*_{0},*x*_{1},*x*_{2}, ...,*x*_{r−1}, and*c*_{r−1}.

The modification is to the generation of a new pair (*x*, *c*). Rearranging the computation to avoid negative numbers, the new *x* value is complemented by subtracting it from *b* − 1:

- [math]\displaystyle{ x_n=(b-1)-(ax_{n-r}+c_{n-1})\,\bmod\,b,\ c_n=\left\lfloor\frac{ax_{n-r}+c_{n-1}}{b}\right\rfloor. }[/math]

The resulting sequence of *x*'s produced by the CMWC RNG will have period the order of *b* in the multiplicative group of residues modulo *ab ^{r}*+1, and the output

Use of lag-*r* CMWC makes it much easier to find periods for *r'*s as large as 512, 1024, 2048, etc. (Making *r* a power of 2 makes it slightly easier (and faster) to access elements in the array containing the *r* most recent *x'*s.)

Another advantage of this modified procedure is that the period is a multiple of *b*, so the output is exactly equidistributed mod *b*.^{[2]} (The ordinary MWC has a slight bias over its full period, which is negligible if the period is long enough.) produces one fewer zero over its full period.)

One *disadvantage* of the CMWC construction is that, with a power-of-two base, the maximum achievable period is less than for a similar-sized MWC generator; you lose several bits. Thus, an MWC generator is usually preferable for small lags. This can remedied by using *b* = 2^{k}−1, or choosing a lag one word longer to compensate.

Some examples: With *b*=2^{32}, and *a*=109111 or 108798 or 108517, the period of the lag-1024 CMWC

- [math]\displaystyle{ x_n=(b-1)-(ax_{n-1024}+c_{n-1})\,\bmod\,b,\ c_n=\left\lfloor\frac{ax_{n-1024}+c_{n-1}}{b}\right\rfloor. }[/math]

will be *a*·2^{32762} = *ab*^{1024}/64, about 10^{9867}.

With *b* = 2^{32} and *a* = 3636507990, *p* = *ab*^{1359} − 1 is a safe prime, so the MWC sequence based on that *a* has period 3636507990·2^{43487} [math]\displaystyle{ \approx }[/math] 10^{13101}.

With *b* = 2^{32}, a CMWC RNG with near record period may be based on the prime *p* = 15455296*b*^{42658} + 1. The order of *b* for that prime is 241489·2^{1365056}, about 10^{410928}.

The MWC modulus of *ab ^{r}*−1 is chosen to make computation particularly simple, but brings with it some disadvantages, notably that the period is at most half the modulus. There are several ways to generalize this, at the cost of more multiplications per iteration.

First, it is possible to add additional terms to the product, producing a modulus of the form *a _{r}b^{r}*+

However, this does not fix the period issue, which depends on the low bits of the modulus. Fortunately, the Montgomery reduction algorithm permits other moduli, as long as they are relatively prime to the base *b*, and this can be applied to permit a modulus of the form *a _{r}b^{r}*−

- [math]\displaystyle{ t = c_{n-1} + \sum_1^r a_i x_{n-i}, \ c_n = \left\lfloor\frac{t}{b}\right\rfloor, \ x_n = a_0^{-1}t \bmod b. }[/math]

In the common case that *b* = 2^{k}, *a*_{0} must be odd for the inverse to exist.

The following is an implementation of the CMWC algorithm in the C programming language. Also, included in the program is a sample initialization function. In this implementation the base is 2^{32}−1 and lag *r*=4096. The period of the resulting generator is about [math]\displaystyle{ 2^{131104} }[/math].

// C99 Complementary Multiply With Carry generator #include <stdint.h> #include <stdio.h> #include <stdlib.h> #include <time.h> // CMWC working parts #define CMWC_CYCLE 4096 // as Marsaglia recommends #define CMWC_C_MAX 809430660 // as Marsaglia recommends struct cmwc_state { uint32_t Q[CMWC_CYCLE]; uint32_t c; // must be limited with CMWC_C_MAX unsigned i; }; // Make 32 bit random number (some systems use 16 bit RAND_MAX [Visual C 2012 uses 15 bits!]) uint32_t rand32(void) { uint32_t result = rand(); return result << 16 | rand(); } // Init the state with seed void initCMWC(struct cmwc_state *state, unsigned int seed) { srand(seed); for (int i = 0; i < CMWC_CYCLE; i++) state->Q[i] = rand32(); do state->c = rand32(); while (state->c >= CMWC_C_MAX); state->i = CMWC_CYCLE - 1; } // CMWC engine uint32_t randCMWC(struct cmwc_state *state) //EDITED parameter *state was missing { uint64_t const a = 18782; // as Marsaglia recommends uint32_t const m = 0xfffffffe; // as Marsaglia recommends uint64_t t; uint32_t x; state->i = (state->i + 1) & (CMWC_CYCLE - 1); t = a * state->Q[state->i] + state->c; /* Let c = t / 0xfffffff, x = t mod 0xffffffff */ state->c = t >> 32; x = t + state->c; if (x < state->c) { x++; state->c++; } return state->Q[state->i] = m - x; } int main() { struct cmwc_state cmwc; unsigned int seed = time(NULL); initCMWC(&cmwc, seed); printf("Random CMWC: %u\n", randCMWC(&cmwc)); }

Because of simplicity, speed, quality (it passes statistical tests very well) and astonishing period, CMWC is known to be used in game development, particularly in modern roguelike games. It is informally known as the Mother of All PRNGs. In libtcod, CMWC4096 replaced MT19937 as the default PRNG.^{[4]}

- Hull, T. E.; Dobell, A. R. (July 1962). "Random Number Generators". SIAM Review 4 (3): 230–254. doi:10.1137/1004061. http://chagall.med.cornell.edu/BioinfoCourse/PDFs/Lecture4/random_number_generator.pdf. Retrieved 2016-06-26.
- Couture, Raymond; L'Ecuyer, Pierre (April 1997). "Distribution properties of Multiply-with-carry random number generators". Mathematics of Computation 66 (218): 591–607. doi:10.1090/S0025-5718-97-00827-2. http://www.ams.org/mcom/1997-66-218/S0025-5718-97-00827-2/S0025-5718-97-00827-2.pdf. "We shall see that, for the complementary MWC, each bit of the output value is fair, that is, the two binary digits will appear equally often in a full period, a property not shared by MWC generators.".
- Goresky, Mark; Klapper, Andrew (October 2003). "Efficient multiply-with-carry random number generators with maximal period". ACM Transactions on Modeling and Computer Simulation 13 (4): 310–321. doi:10.1145/945511.945514. https://www.math.ias.edu/~goresky/pdf/p1-goresky.pdf.
- "Random number generator - RogueBasin". http://www.roguebasin.com/index.php?title=Random_number_generator#Complementary_multiply-with-carry.

More

Information

Subjects:
Others

Contributor
MDPI registered users' name will be linked to their SciProfiles pages. To register with us, please refer to https://encyclopedia.pub/register
:

View Times:
1.4K

Entry Collection:
HandWiki

Update Date:
29 Sep 2022