Linear congruential generator
You are encouraged to solve this task according to the task description, using any language you may know.
The linear congruential generator is a very simple example of a random number generator. All linear congruential generators use this formula:
Where:
- is a seed.
- , , , ..., are the random numbers.
- , , are constants.
If one chooses the values of , and with care, then the generator produces a uniform distribution of integers from to .
LCG numbers have poor quality. and are not independent, as true random numbers would be. Anyone who knows can predict , therefore LCG is not cryptographically secure. The LCG is still good enough for simple tasks like Miller-Rabin primality test, or FreeCell deals. Among the benefits of the LCG, one can easily reproduce a sequence of numbers, from the same . One can also reproduce such sequence with a different programming language, because the formula is so simple.
The task is to replicate two historic random number generators. One is the rand()
function from BSD libc, and the other is the rand()
function from the Microsoft C Runtime (MSCVRT.DLL). Each replica must yield the same sequence of integers as the original generator, when starting from the same seed.
In these formulas, the seed becomes . The random sequence is , and so on.
BSD formula:
- is in range 0 to 2147483647.
Microsoft formula:
- is in range 0 to 32767.
The BSD formula was so awful that FreeBSD switched to a different formula. More info is at Random number generator (included)#C.
Ada
We first specify a generic package LCG:
<lang Ada>generic
type Base_Type is mod <>; Multiplyer, Adder: Base_Type; Output_Divisor: Base_Type := 1;
package LCG is
procedure Initialize(Seed: Base_Type); function Random return Base_Type; -- changes the state and outputs the result
end LCG;</lang>
Then we provide a generic implementation:
<lang Ada>package body LCG is
State: Base_Type := Base_Type'First;
procedure Initialize(Seed: Base_Type) is begin State := Seed; end Initialize;
function Random return Base_Type is begin State := State * Multiplyer + Adder; return State / Output_Divisor; end Random;
end LCG;</lang>
Next, we define the MS- and BSD-instantiations of the generic package:
<lang Ada>with Ada.Text_IO, LCG;
procedure Run_LCGs is
type M31 is mod 2**31;
package BSD_Rand is new LCG(Base_Type => M31, Multiplyer => 1103515245, Adder => 12345);
package MS_Rand is new LCG(Base_Type => M31, Multiplyer => 214013, Adder => 2531011, Output_Divisor => 2**16);
begin
for I in 1 .. 10 loop Ada.Text_IO.Put_Line(M31'Image(BSD_Rand.Random)); end loop; for I in 1 .. 10 loop Ada.Text_IO.Put_Line(M31'Image(MS_Rand.Random)); end loop;
end Run_LCGs;</lang>
Finally, we run the program, which generates the following output (note that the first ten lines are from the BSD generator, the next ten from the MS generator):
12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310 38 7719 21238 2437 8855 11797 8365 32285 10450 30612
bc
As with dc, bc has no bitwise operators. <lang bc>/* BSD rand */
define rand() { randseed = (randseed * 1103515245 + 12345) % 2147483648 return randseed }
randseed = 1 rand(); rand(); rand(); print "\n"
/* Microsoft rand */
define rand() { randseed = (randseed * 214013 + 2531011) % 2147483648 return randseed / 65536 }
randseed = 1 rand(); rand(); rand(); print "\n"</lang>
C
In a pretended lib style, this code produces a rand() function depends on compiler macro: gcc -DMS_RAND
uses MS style, otherwise it's BSD rand by default.
<lang C>#include <stdio.h>
/* always assuming int is at least 32 bits */ int rand(); int rseed = 0;
inline void srand(int x) { rseed = x; }
- ifndef MS_RAND
- define RAND_MAX ((1U << 31) - 1)
inline int rand() { return rseed = (rseed * 1103515245 + 12345) & RAND_MAX; }
- else /* MS rand */
- define RAND_MAX_32 ((1U << 31) - 1)
- define RAND_MAX ((1U << 15) - 1)
inline int rand() { return (rseed = (rseed * 214013 + 2531011) & RAND_MAX_32) >> 16; }
- endif/* MS_RAND */
int main() { int i; printf("rand max is %d\n", RAND_MAX);
for (i = 0; i < 100; i++) printf("%d\n", rand());
return 0; }</lang>
Common Lisp
<lang lisp>(defun make-rng (&key (seed 0) (mode nil))
"returns an RNG according to :seed and :mode keywords default mode: bsd default seed: 0 (should be 1 actually)" (if (eql mode 'ms) #'(lambda ()
(ash (setf seed (mod (+ (* 214013 seed) 2531011) (expt 2 31))) -16))
#'(lambda () (setf seed (mod (+ (* seed 1103515245) 12345) (expt 2 31))))))
(let ((rng (make-rng)))
(dotimes (x 10) (format t "BSD: ~d~%" (funcall rng))))
(let ((rng (make-rng :mode 'ms :seed 1)))
(dotimes (x 10) (format t "MS: ~d~%" (funcall rng))))</lang>
D
<lang d>import std.stdio;
struct LinearCongruentialGenerator {
enum uint RAND_MAX = (1U << 31) - 1; uint seed = 0;
uint randBSD() pure nothrow { seed = (seed * 1_103_515_245 + 12_345) & RAND_MAX; return seed; }
uint randMS() pure nothrow { seed = (seed * 214_013 + 2_531_011) & RAND_MAX; return seed >> 16; }
}
void main() {
LinearCongruentialGenerator rnd;
foreach (i; 0 .. 10) writeln(rnd.randBSD()); writeln();
rnd.seed = 0; foreach (i; 0 .. 10) writeln(rnd.randMS());
}</lang> Output:
12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310 38 7719 21238 2437 8855 11797 8365 32285 10450 30612
dc
dc has no bitwise operations, so this program uses the modulus operator (2147483648 %
) and division (65536 /
). Fortunately, dc numbers cannot overflow to negative, so the modulus calculation involves only non-negative integers.
For BSD rand(): <lang dc>[*
* lrx -- (random number from 0 to 2147483647) * * Returns a number from the BSD rand() sequence. * Seeded by storing a seed in register R. *]sz
[lR 1103515245 * 12345 + 2147483648 % d sR]sr
[* Set seed to 1, then print the first 3 random numbers. *]sz 1 sR lrx psz lrx psz lrx psz</lang>
1103527590 377401575 662824084
For Microsoft rand(): <lang dc>[*
* lrx -- (random number from 0 to 32767) * * Returns a number from the Microsoft rand() sequence. * Seeded by storing a seed in register R. *]sz
[lR 214013 * 2531011 + 2147483648 % d sR 65536 /]sr
[* Set seed to 1, then print the first 3 random numbers. *]sz 1 sR lrx psz lrx psz lrx psz</lang>
41 18467 6334
F#
<lang fsharp>module lcg =
let bsd seed = let state = ref seed (fun (_:unit) -> state := (1103515245 * !state + 12345) % (1<<<31) !state)
let ms seed = let state = ref seed (fun (_:unit) -> state := (214013 * !state + 2531011) % (1<<<31) !state / (1<<<16))
</lang>
let rnd = lcg.bsd 0; [for n in [0 .. 9] -> rnd()];; BSD: [12345; -740551042; -1492899873; -698016724; 229283573; -1038148470; 1051550459; -853684456; -1353011855; 551188310] MS: [38; 7719; -11529; 2437; -23912; 11797; 8365; 32285; -22317; 30612]
Forth
<lang forth>1 31 lshift 1- constant MAX-RAND-BSD 1 15 lshift 1- constant MAX-RAND-MS
variable seed \ seed variable
- (random) seed @ * + dup seed ! ; ( -- n)
- BSDrandom MAX-RAND-BSD 12345 1103515245 (random) and ;
- MSrandom MAX-RAND-MS 2531011 214013 (random) 16 rshift and ;
- test-random
1 seed ! cr ." BSD (seed=1)" cr 5 0 do BSDrandom . cr loop 1 seed ! cr ." MS (seed=1)" cr 5 0 do MSrandom . cr loop
test-random</lang>
Output:
BSD (seed=1) 1103527590 377401575 662824084 1147902781 2035015474 MS (seed=1) 41 18467 6334 26500 19169
Fortran
<lang fortran>module lcgs
implicit none
integer, parameter :: i64 = selected_int_kind(18) integer, parameter :: a1 = 1103515245, a2 = 214013 integer, parameter :: c1 = 12345, c2 = 2531011 integer, parameter :: div = 65536 integer(i64), parameter :: m = 2147483648_i64 ! need to go to 64 bits because ! of the use of signed integers
contains
function bsdrand(seed)
integer :: bsdrand integer, optional, intent(in) :: seed integer(i64) :: x = 0 if(present(seed)) x = seed x = mod(a1 * x + c1, m) bsdrand = x
end function
function msrand(seed)
integer :: msrand integer, optional, intent(in) :: seed integer(i64) :: x = 0 if(present(seed)) x = seed x = mod(a2 * x + c2, m) msrand = x / div
end function end module
program lcgtest
use lcgs implicit none integer :: i write(*, "(a)") " BSD MS" do i = 1, 10 write(*, "(2i12)") bsdrand(), msrand() end do
end program</lang> Output
BSD MS 12345 38 1406932606 7719 654583775 21238 1449466924 2437 229283573 8855 1109335178 11797 1051550459 8365 1293799192 32285 794471793 10450 551188310 30612
Go
<lang go>package main
import "fmt"
// basic linear congruential generator func lcg(a, c, m, seed uint32) func() uint32 {
r := seed return func() uint32 { r = (a*r + c) % m return r }
}
// microsoft generator has extra division step func msg(seed uint32) func() uint32 {
g := lcg(214013, 2531011, 1<<31, seed) return func() uint32 { return g() / (1 << 16) }
}
func example(seed uint32) {
fmt.Printf("\nWith seed = %d\n", seed) bsd := lcg(1103515245, 12345, 1<<31, seed) msf := msg(seed) fmt.Println(" BSD Microsoft") for i := 0; i < 5; i++ { fmt.Printf("%10d %5d\n", bsd(), msf()) }
}
func main() {
example(0) example(1)
}</lang> Output:
With seed = 0 BSD Microsoft 12345 38 1406932606 7719 654583775 21238 1449466924 2437 229283573 8855 With seed = 1 BSD Microsoft 1103527590 41 377401575 18467 662824084 6334 1147902781 26500 2035015474 19169
Icon and Unicon
The following LCRNG's behave in the same way maintaining the state (seed) from round to round. There is an srand procedure for each lcrng that maintains the seed state and allows the user to assign a new state. <lang Icon>link printf
procedure main()
printf(" BSD MS\n") every 1 to 10 do printf("%10s %10s\n",rand_BSD(),rand_MS())
end
procedure srand_BSD(x) #: seed random static seed
return seed := \x | \seed | 0 # parm or seed or zero if none
end
procedure rand_BSD() #: lcrng
return srand_BSD((1103515245 * srand_BSD() + 12345) % 2147483648)
end
procedure srand_MS(x) #: seed random static seed
return seed := \x | \seed | 0 # parm or seed or zero if none
end
procedure rand_MS() #: lcrng
return ishift(srand_MS((214013 * srand_MS() + 2531011) % 2147483648),-16)
end</lang>
J
Solution: <lang j>lcg=: adverb define
0 m lcg y
'a c mod'=. m }. (mod | c + a * ])^:(<y+1) x
)
rand_bsd=: (1103515245 12345 , <.2^31) lcg rand_ms=: (2^16) <.@:%~ (214013 2531011 , <.2^31) lcg</lang> Example Use: <lang j> rand_bsd 10 12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310
654583775 rand_bsd 4
1449466924 229283573 1109335178 1051550459
rand_ms 10
38 7719 21238 2437 8855 11797 8365 32285 10450 30612</lang>
K
<lang K> bsd:{1_ y{((1103515245*x)+12345)!(_2^31)}\x}
ms:{1_(y{_(((214013*x)+2531011)!(_2^31))}\x)%(_2^16)}
bsd[0;10]
12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310
ms[0;10]
38 7719 21238 2437 8855 11797 8365 32285 10450 30612</lang>
Liberty BASIC
<lang lb> 'by default these are 0 global BSDState global MSState
for i = 1 to 10
print randBSD()
next i
for i = 1 to 10
print randMS()
next i
function randBSD()
randBSD = (1103515245 * BSDState + 12345) mod (2 ^ 31) BSDState = randBSD
end function
function randMS()
MSState = (214013 * MSState + 2531011) mod (2 ^ 31) randMS = int(MSState / 2 ^ 16)
end function </lang>
PARI/GP
Note that up to PARI/GP version 2.3.0, random()
used a linear congruential generator.
<lang parigp>BSDseed=Mod(1,1<<31);
MSFTseed=Mod(1,1<<31);
BSD()=BSDseed=1103515245*BSDseed+12345;lift(BSDseed);
MSFT()=MSFTseed=214013*MSFTseed+2531011;lift(MSFTseed)%(1<<31);</lang>
Pascal
<lang pascal>Program LinearCongruentialGenerator(output);
var
x1, x2: int64;
function bsdrand: longint;
const a = 1103515245; c = 12345; m = 2147483648; begin x1 := (a * x1 + c) mod m; bsdrand := x1; end;
function msrand: longint;
const a = 214013; c = 2531011; m = 2147483648; begin x2 := (a * x2 + c) mod m; msrand := x2 div 65536; end;
var
i: longint;
begin
writeln(' BSD MS'); x1 := 0; x2 := 0; for i := 1 to 10 do writeln(bsdrand:12, msrand:12);
end.</lang> Output:
BSD MS 12345 7584 1124652145 3277 1499545833 3067 1558406049 31446 696007321 13069 56579025 17343 1312705865 2510 811881729 5264 1301653753 21298 1318262577 27689
Perl
Creates a magic scalar whose value is next in the LCG sequence when read.<lang perl>use strict; package LCG;
use overload '0+' => \&get;
use integer; sub gen_bsd { (1103515245 * shift() + 12345) % (1 << 31) }
sub gen_ms { my $s = (214013 * shift() + 2531011) % (1 << 31); $s, $s / (1 << 16) }
sub set { $_[0]->{seed} = $_[1] } # srand sub get { my $o = shift; ($o->{seed}, my $r) = $o->{meth}->($o->{seed}); $r //= $o->{seed} }
sub new { my $cls = shift; my %opts = @_; bless { seed => $opts{seed}, meth => $opts{meth} eq 'MS' ? \&gen_ms : \&gen_bsd, }, ref $cls || $cls; }
package main;
my $rand = LCG->new;
print "BSD:\n"; print "$rand\n" for 1 .. 10;
$rand = LCG->new(meth => 'MS');
print "\nMS:\n"; print "$rand\n" for 1 .. 10;</lang>output<lang>BSD: 12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310
MS: 38 7719 21238 2437 8855 11797 8365 32285 10450 30612</lang>
Perl 6
Define subroutines implementing the LCG algorithm for each version then use those to generate lazy infinite lists of values and return the first 10 values from each. <lang perl6>my $mod = 2**31; sub bsd ($seed) { ( 1103515245 * $seed + 12345 ) % $mod }; sub ms ($seed) { ( 214013 * $seed + 2531011 ) % $mod };
say 'BSD LCG first 10 values:'; .say for ( 0.&bsd, -> $seed { $seed.&bsd } ... * )[^10];
say "\nMS LCG first 10 values:"; ($_ +> 16).say for ( 0.&ms, -> $seed { $seed.&ms } ... * )[^10];</lang>
BSD LCG first 10 values: 12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 1293799192 794471793 551188310 MS LCG first 10 values: 38 7719 21238 2437 8855 11797 8365 32285 10450 30612
PicoLisp
<lang PicoLisp>(zero *BsdSeed *MsSeed)
(de bsdRand ()
(setq *BsdSeed (& (+ 12345 (* 1103515245 *BsdSeed)) `(dec (** 2 31))) ) )
(de msRand ()
(>> 16 (setq *MsSeed (& (+ 2531011 (* 214013 *MsSeed)) `(dec (** 2 31))) ) ) )</lang>
Output:
: (do 7 (printsp (bsdRand))) 12345 1406932606 654583775 1449466924 229283573 1109335178 1051550459 -> 1051550459 : (do 12 (printsp (msRand))) 38 7719 21238 2437 8855 11797 8365 32285 10450 30612 5853 28100 -> 28100
PL/I
<lang> (nofixedoverflow, nosize): LCG: procedure options (main);
declare i fixed binary;
put skip list ('BSD', 'MS'); do i = 1 to 20; put skip list (BSD(), MS()); end;
bsd: procedure returns (fixed binary (31));
declare const fixed binary static initial (12345); declare s fixed binary (31) static initial (123456789);
s = s * 1103515245 + const; s = isrl(isll(s,1), 1); return (s);
end bsd; ms: procedure returns (fixed binary (15));
declare const fixed binary (31) static initial (2531011); declare s fixed binary (31) static initial (123456789);
s = s * 214013 + const; s = isrl(isll(s,1), 1); return (isrl(s,16));
end ms;
end LCG; </lang> OUTPUT:
BSD MS 231794730 13259 1126946331 26974 1757975480 13551 850994577 30354 1634557174 18709 707246327 15861 1397699428 16906 1035569613 21981 1904890498 8603 1335160211 12911 1434329552 18110 1273099721 3228 1250890958 27918 1016516591 17989 1097566972 22768 436938117 23599 1175171034 7712 1059748875 15601 308566760 7038 534615297 21512
PureBasic
<lang purebasic>Procedure ms_LCG(seed.q = -1)
Static state.q If seed >= 0 state = seed Else state = (state * 214013 + 2531011) % (1 << 31) ProcedureReturn state >> 16 EndIf
EndProcedure
Procedure.q bsd_LCG(seed.q = -1)
Static state.q If seed >= 0 state = seed Else state = (state * 1103515245 + 12345) % (1 << 31) ProcedureReturn state EndIf
EndProcedure
If OpenConsole()
Define i PrintN("BSD (seed = 1)") bsd_LCG(1) For i = 1 To 5 PrintN(Str(bsd_LCG())) Next PrintN(#CRLF$ + "MS (seed = 1)") ms_LCG(1) For i = 1 To 5 PrintN(Str(ms_LCG())) Next Print(#CRLF$ + #CRLF$ + "Press ENTER to exit"): Input() CloseConsole()
EndIf</lang> Sample output:
BSD (seed = 1) 1103527590 377401575 662824084 1147902781 2035015474 MS (seed = 1) 41 18467 6334 26500 19169
Ruby
You can create multiple instances of LCG::Berkeley or LCG::Microsoft. Each instance privately keeps the original seed in @seed, and the current state in @r. Each class resembles the core Random class, but with fewer features. The .new method takes a seed. The #rand method returns the next random number. The #seed method returns the original seed.
<lang ruby>module LCG
module Common # The original seed of this generator. attr_reader :seed
# Creates a linear congruential generator with the given _seed_. def initialize(seed) @seed = @r = seed end end
# LCG::Berkeley generates 31-bit integers using the same formula # as BSD rand(). class Berkeley include Common def rand @r = (1103515245 * @r + 12345) & 0x7fff_ffff end end
# LCG::Microsoft generates 15-bit integers using the same formula # as rand() from the Microsoft C Runtime. class Microsoft include Common def rand @r = (214013 * @r + 2531011) & 0x7fff_ffff @r >> 16 end end
end</lang>
The next example sets the seed to 1, and prints the first 5 random numbers.
<lang ruby>lcg = LCG::Berkeley.new(1) p (1..5).map {lcg.rand}
- prints [1103527590, 377401575, 662824084, 1147902781, 2035015474]
lcg = LCG::Microsoft.new(1) p (1..5).map {lcg.rand}
- prints [41, 18467, 6334, 26500, 19169]</lang>
Tcl
Using an object-oriented solution, inspired by (but not a translation of) the Ruby solution above. <lang tcl>package require Tcl 8.6
- General form of a linear-congruential RNG
oo::class create LCRNG {
variable seed A B C D constructor {init a b c d} {
if {$init < 1} {set init [clock clicks]} variable seed $init A $a B $b C $c D $d
} method rand {} {
set seed [expr {($A * $seed + $B) % $C}] return [expr {$seed / $D}]
} method srand x {
set seed $x
}
}
- Subclass to introduce constants
oo::class create BSDRNG {
superclass LCRNG constructor Template:InitialSeed -1 {
next $initialSeed 1103515245 12345 [expr {2**31}] 1
}
} oo::class create MSRNG {
superclass LCRNG constructor Template:InitialSeed -1 {
next $initialSeed 214013 2531011 [expr {2**31}] [expr {2**16}]
}
}</lang> Demo code: <lang tcl>proc sample rng {foreach - {1 2 3 4 5} {lappend r [$rng rand]}; join $r ", "} puts BSD:\t\[[sample [BSDRNG new 1]]\] puts MS:\t\[[sample [MSRNG new 1]]\]</lang> Output:
BSD: [1103527590, 377401575, 662824084, 1147902781, 2035015474] MS: [41, 18467, 6334, 26500, 19169]