/* "CMR" stands for "Constant Multiply Rotate", "63" stands for the period: > 2^63 (almost 2^64). It's based on a publication by Mark A. Overton(2011): http://www. Fast, High-Quality, Parallel Random Number Generators CMR63: private static uint cmr() { z0 *= m0; z1 *= m1; z0 = z0 << s0 | z0 >> r0; z1 = z1 << s1 | z1 >> r1; return z0 ^ z1; } The constants can be found in: http://www.OvertonTables.zip (cmr-32.txt) A small part of the cmr-32 table: cmr-32: x = Constant*x; x = rotl(x,RL); Seed value can be 1, or SeedStart..(SeedStart+SeedCnt-1). SB = number of bits in seed (2^SB <= SeedCnt). SeedCnt SeedStart SB Period Constant RL Prime factors of period -------- ---------- -- ---------- ---------- -- ----------------------- 62973467 1377002680 25 4294966876 3563976171 16 2^2 1073741719 33525602 3202323436 24 4294965919 1422968075 16 307 757 18481 32022541 1180658772 24 4294966152 1977089609 19 2^3 3 41 1051 4153 31280659 554048116 24 4294966449 433149435 17 3 1259 1137137 2749407 2613580375 21 4294950337 272690735 19 4294950337 (prime) 1095230 3119024045 20 4294928147 64333559 18 4294928147 (prime) 998895 3255944955 19 4294915769 3152644205 13 4294915769 (prime) 664885 3993266363 19 4294881427 4031235431 15 4294881427 (prime) The full table contains 5586 data rows, 272 periods are prime. The smallest period is 4284340603. Thus the smallest period of two combined generators (with periods having no common factors) is at least: 4284340603^2 4284340603^2 > 2^63.99 ==> 2^63.99 < period < 2^64 I use seed values to get unique pairs of constants (multipliers), assuming random sequences generated by those unique pairs are not correlated. Constants belonging to the 272 prime periods can be combined with other constants (belonging to 272 of the 5586 - 272 non prime periods). A more refined scheme might get the number of combinations closer to 5586/2. On my pc (Athlon II X4 640, 3GHz, 2GB, XP) 10^9 iterations of CMR63 take 6.27 seconds, so the full period takes about 3600 years! 2^63.99 / 10^9 * 6.27 / (60 * 60 * 24 * 365) ~ 3600 The scaled version: "cmr_uint(u)" returns a value x ( 0 <= x <= u ). The time to generate one scaled random number is 21 ns, mean 14 ns, more accurate, mean maximum time: 21 ns, mean mean is mean: 14 ns. Lots of random bits are spoiled, it does not matter: 3600 years! CMR16: Uses rotation (shift) value 16. Periods of the sub generators ( z2 = m2 * ... , z3 = m3 * ... ) are prime. Result: 22 Divergent High Quality Random Number Sequences. */ class cmr63 { private static int s0 = 16, r0 = 16, s1 = 15, r1 = 17; // init_cmr(0); private static uint m0 = 3563976171, z0 = 4125873261, m1 = 4031235431, z1 = 3803445283; private static void init_cmr(uint seed) { byte[] s = { 16, 16, 19, 17, 19, 18, 13, 15 }; uint[] m = { 3563976171, 1422968075, 1977089609, 433149435, 272690735, 64333559, 3152644205, 4031235431}; seed &= 3; m0 = m[seed]; s0 = s[seed]; r0 = 32 - s0; z0 = 1; seed ^= 7; m1 = m[seed]; s1 = s[seed]; r1 = 32 - s1; z1 = 1; cmr(); } private static uint cmr() { z0 *= m0; z1 *= m1; z0 = z0 << s0 | z0 >> r0; z1 = z1 << s1 | z1 >> r1; return z0 ^ z1; } private static uint g_u = 0; private static int g_s = 0; private static uint cmr_uint(uint u) { uint x; if (u == 0) return 0; if (u != g_u) { g_u = u; g_s = 32 - nob(u); } do { z0 *= m0; z1 *= m1; z0 = z0 << s0 | z0 >> r0; z1 = z1 << s1 | z1 >> r1; x = (z0 ^ z1) >> g_s; } while (x > u); return x; } private static int nob(uint u) // number of bits, position of the highest set bit, u > 0 { return u < 1u << 16 ? u < 1u << 08 ? u < 1u << 04 ? u < 1u << 02 ? u < 1u << 01 ? 01 : 02 : u < 1u << 03 ? 03 : 04 : u < 1u << 06 ? u < 1u << 05 ? 05 : 06 : u < 1u << 07 ? 07 : 08 : u < 1u << 12 ? u < 1u << 10 ? u < 1u << 09 ? 09 : 10 : u < 1u << 11 ? 11 : 12 : u < 1u << 14 ? u < 1u << 13 ? 13 : 14 : u < 1u << 15 ? 15 : 16 : u < 1u << 24 ? u < 1u << 20 ? u < 1u << 18 ? u < 1u << 17 ? 17 : 18 : u < 1u << 19 ? 19 : 20 : u < 1u << 22 ? u < 1u << 21 ? 21 : 22 : u < 1u << 23 ? 23 : 24 : u < 1u << 28 ? u < 1u << 26 ? u < 1u << 25 ? 25 : 26 : u < 1u << 27 ? 27 : 28 : u < 1u << 30 ? u < 1u << 29 ? 29 : 30 : u < 1u << 31 ? 31 : 32; } static void Main() { // http://www.TestU01.BigCrush ? } private static uint m2 = 3745979853, z2 = 500031303, m3 = 0623716905, z3 = 707339565; private static void init_cmr16(uint seed) { uint[] m = { 3745979853, 4055716687, 3693386591, 3542220329, 1775851103, 1866916287, 4188393139, 4141129223, 1173908643, 3198474053, 11119693, 1282266473, 4076777453, 3908725387, 3293562383, 2492630213, 1818407027, 608828557, 872259061, 2075607481, 1573125557, 2615661665, 1402711077, 3212405133, 680154359, 2023590663, 3458456891, 4184846215, 2408125305, 2558924297, 3008413683, 466035855, 1647905439, 2930730743, 733571709, 3997625831, 1919196763, 3392242035, 100431167, 579587817, 3074845609, 1931914705, 3131462569, 623716905}; seed %= 22; m2 = m[seed]; m3 = m[43 - seed]; z2 = z3 = 1; cmr16(); } private static uint cmr16() // 2^63.99 < period < 2^64 { z2 *= m2; z3 *= m3; z2 = z2 << 16 | z2 >> 16; z3 = z3 << 16 | z3 >> 16; return z2 ^ z3; } }
No comments:
Post a Comment