High+Quality+Random+Number+Generation

//by Richard Russell, March 2010//

Listed below are BBC BASIC implementations of two high performance Pseudo Random Number Generators: the **WELL512A** and **WELL1024A** algorithms developed in 2006 (see this [|Wikipedia article]). They have state sizes of 512 bits and 1024 bits, and therefore sequence lengths of 2^512-1 and 2^1024-1, respectively. They both return 32-bit signed integer values, but these can easily be converted to other formats, for example using the code snippets listed at Notes on the use of RND.

These **Well Equidistributed Long-period Linear** algorithms have long sequences, excellent statistical properties and are computationally simpler than some alternative methods. They are ideally suited to applications such as Monte Carlo simulations.

In the listings below, BBC BASIC's RND function is used to provide the initialisation (seed). This might be adequate for some uses, but for more critical applications a seed with more entropy would be desirable. The 'state' information is contained in the array **State%** and the variable **State_i%** which are passed by reference to the generator function. To 'replay' a sequence it is only necessary to save and restore these values. code format="bb4w" REM BBC BASIC version of the WELL512A Pseudo Random Number Generator REM By Richard Russell, http://www.rtrussell.co.uk/, v1.00, 30-Mar-2010 REM Based on the original C code by François Panneton & Pierre L'Écuyer REM (Université de Montréal) and Makoto Matsumoto (Hiroshima University) REM http://en.wikipedia.org/wiki/Well_Equidistributed_Long-period_Linear DIM State%(15) State_i% = 0 REM Initialisation: FOR j% = 0 TO 15 State%(j%) = RND NEXT j%     FOR j% = 1 TO 32 PRINT RIGHT$("0000000"+STR$~(FN_WELLRNG512a(State_i%, State%)), 8) NEXT j%     END DEF FN_WELLRNG512a(RETURN I%, S%) LOCAL X%, Y%, Z%     X% = S%(I%) EOR S%(I%)<<16 EOR S%(I%+13AND15) EOR S%(I%+13AND15)<<15 Y% = S%(I%+9AND15) EOR S%(I%+9AND15)>>>11 S%(I%) = X% EOR Y%     Z% = S%(I%) EOR S%(I%)<<5 AND &DA442D24 I% = I%+15AND15 S%(I%) = S%(I%) EOR S%(I%)<<2 EOR X% EOR X%<<18 EOR Y%<<28 EOR Z%     = S%(I%) code

code format="bb4w" REM BBC BASIC version of the WELL1024A Pseudo Random Number Generator REM By Richard Russell, http://www.rtrussell.co.uk/, v1.00, 30-Mar-2010 REM Based on the original C code by François Panneton & Pierre L'Écuyer REM (Université de Montréal) and Makoto Matsumoto (Hiroshima University) REM http://en.wikipedia.org/wiki/Well_Equidistributed_Long-period_Linear DIM State%(31) State_i% = 0 REM Initialisation: FOR j% = 0 TO 31 State%(j%) = RND NEXT j%     FOR j% = 1 TO 32 PRINT RIGHT$("0000000"+STR$~(FN_WELLRNG1024a(State_i%, State%)), 8) NEXT j%     END DEF FN_WELLRNG1024a(RETURN I%, S%) LOCAL X%, Y%     X% = S%(I%) EOR S%(I%+3AND31) EOR S%(I%+3AND31)>>>8 Y% = S%(I%+24AND31) EOR S%(I%+24AND31)<<19 EOR S%(I%+10AND31) EOR S%(I%+10AND31)<<14 S%(I%) = X% EOR Y%     I% = I%+31AND31 S%(I%) = S%(I%) EOR S%(I%)<<11 EOR X% EOR X%<<7 EOR Y% EOR Y%<<13 = S%(I%) code