/* $NetBSD: mertwist.c,v 1.7 2008/02/17 22:49:11 matt Exp $ */ /*- * Copyright (c) 2008 The NetBSD Foundation, Inc. * All rights reserved. * * This code is derived from software contributed to The NetBSD Foundation * by Matt Thomas . * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions * are met: * 1. Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * 2. Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the distribution. * 3. All advertising materials mentioning features or use of this software * must display the following acknowledgement: * This product includes software developed by the NetBSD * Foundation, Inc. and its contributors. * 4. Neither the name of The NetBSD Foundation nor the names of its * contributors may be used to endorse or promote products derived * from this software without specific prior written permission. * * THIS SOFTWARE IS PROVIDED BY THE NETBSD FOUNDATION, INC. AND CONTRIBUTORS * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED * TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR * PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR CONTRIBUTORS * BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE * POSSIBILITY OF SUCH DAMAGE. */ #if defined(_KERNEL) || defined(_STANDALONE) #include #include #include #include #else #include #include #include #include #define KASSERT(x) assert(x) #endif /* * Mersenne Twister. Produces identical output compared to mt19937ar.c * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html */ #define MATRIX_A(a) (((a) & 1) ? 0x9908b0df : 0) #define TEMPERING_MASK_B 0x9d2c5680 #define TEMPERING_MASK_C 0xefc60000 #define UPPER_MASK 0x80000000 #define LOWER_MASK 0x7fffffff #define MIX(u,l) (((u) & UPPER_MASK) | ((l) & LOWER_MASK)) #define KNUTH_MULTIPLIER 0x6c078965 #ifndef MTPRNG_RLEN #define MTPRNG_RLEN 624 #endif #define MTPRNG_POS1 397 static void mtprng_refresh(struct mtprng_state *); /* * Initialize the generator from a seed */ void mtprng_init32(struct mtprng_state *mt, uint32_t seed) { size_t i; /* * Use Knuth's algorithm for expanding this seed over its * portion of the key space. */ mt->mt_elem[0] = seed; for (i = 1; i < MTPRNG_RLEN; i++) { mt->mt_elem[i] = KNUTH_MULTIPLIER * (mt->mt_elem[i-1] ^ (mt->mt_elem[i-1] >> 30)) + i; } mtprng_refresh(mt); } void mtprng_initarray(struct mtprng_state *mt, const uint32_t *key, size_t keylen) { uint32_t *mp; size_t i, j, k; /* * Use Knuth's algorithm for expanding this seed over its * portion of the key space. */ mt->mt_elem[0] = 19650218UL; for (i = 1; i < MTPRNG_RLEN; i++) { mt->mt_elem[i] = KNUTH_MULTIPLIER * (mt->mt_elem[i-1] ^ (mt->mt_elem[i-1] >> 30)) + i; } KASSERT(keylen > 0); i = 1; j = 0; k = (keylen < MTPRNG_RLEN ? MTPRNG_RLEN : keylen); mp = &mt->mt_elem[1]; for (; k-- > 0; mp++) { mp[0] ^= (mp[-1] ^ (mp[-1] >> 30)) * 1664525UL; mp[0] += key[j] + j; if (++i == MTPRNG_RLEN) { KASSERT(mp == mt->mt_elem + MTPRNG_RLEN - 1); mt->mt_elem[0] = mp[0]; i = 1; mp = mt->mt_elem; } if (++j == keylen) j = 0; } for (j = MTPRNG_RLEN; --j > 0; mp++) { mp[0] ^= (mp[-1] ^ (mp[-1] >> 30)) * 1566083941UL; mp[0] -= i; if (++i == MTPRNG_RLEN) { KASSERT(mp == mt->mt_elem + MTPRNG_RLEN - 1); mt->mt_elem[0] = mp[0]; i = 1; mp = mt->mt_elem; } } mt->mt_elem[0] = 0x80000000; mtprng_refresh(mt); } /* * Generate an array of 624 untempered numbers */ void mtprng_refresh(struct mtprng_state *mt) { uint32_t y; size_t i, j; /* * The following has been refactored to avoid the need for 'mod 624' */ for (i = 0, j = MTPRNG_POS1; j < MTPRNG_RLEN; i++, j++) { y = MIX(mt->mt_elem[i], mt->mt_elem[i+1]); mt->mt_elem[i] = mt->mt_elem[j] ^ (y >> 1) ^ MATRIX_A(y); } for (j = 0; i < MTPRNG_RLEN - 1; i++, j++) { y = MIX(mt->mt_elem[i], mt->mt_elem[i+1]); mt->mt_elem[i] = mt->mt_elem[j] ^ (y >> 1) ^ MATRIX_A(y); } y = MIX(mt->mt_elem[MTPRNG_RLEN - 1], mt->mt_elem[0]); mt->mt_elem[MTPRNG_RLEN - 1] = mt->mt_elem[MTPRNG_POS1 - 1] ^ (y >> 1) ^ MATRIX_A(y); } /* * Extract a tempered PRN based on the current index. Then recompute a * new value for the index. This avoids having to regenerate the array * every 624 iterations. We can do this since recomputing only the next * element and the [(i + 397) % 624] one. */ uint32_t mtprng_rawrandom(struct mtprng_state *mt) { uint32_t x, y; const size_t i = mt->mt_idx; size_t j; /* * First generate the random value for the current position. */ x = mt->mt_elem[i]; x ^= x >> 11; x ^= (x << 7) & TEMPERING_MASK_B; x ^= (x << 15) & TEMPERING_MASK_C; x ^= x >> 18; /* * Next recalculate the next sequence for the current position. */ y = mt->mt_elem[i]; if (__predict_true(i < MTPRNG_RLEN - 1)) { /* * Avoid doing % since it can be expensive. * j = (i + MTPRNG_POS1) % MTPRNG_RLEN; */ j = i + MTPRNG_POS1; if (j >= MTPRNG_RLEN) j -= MTPRNG_RLEN; mt->mt_idx++; } else { j = MTPRNG_POS1 - 1; mt->mt_idx = 0; } y = MIX(y, mt->mt_elem[mt->mt_idx]); mt->mt_elem[i] = mt->mt_elem[j] ^ (y >> 1) ^ MATRIX_A(y); /* * Return the value calculated in the first step. */ return x; } /* * This is a non-standard routine which attempts to return a cryptographically * strong random number by collapsing 2 32bit values outputed by the twister * into one 32bit value. */ uint32_t mtprng_random(struct mtprng_state *mt) { uint32_t a; mt->mt_count = (mt->mt_count + 13) & 31; a = mtprng_rawrandom(mt); a = (a << mt->mt_count) | (a >> (32 - mt->mt_count)); return a + mtprng_rawrandom(mt); }