|
|
1.1 root 1: # Date: 4 Mar 1982 20:54:18-PST
2: # From: cithep!citcsv!kingsley at Berkeley
3: # To: cithep!ucbvax!4bsd-bugs@Berkeley
4: # Subject: random number generation.
5:
6: # I have done a little work with rand supplied with the system and I have
7: # discovered that it is flawed. The manual page claims that it has a
8: # period of 2^32 and returns numbers from 0 to 2^31-1. The code makes it
9: # look like the author thought it was correct, but it is not. Instead of
10: # masking out the most significant (and also most random) bit, you should
11: # do an unsigned shift to throw out the least significant (and least
12: # random) bit. I have also found a multiplier that passes Knuth's
13: # spectral test very well.
14: # I have written a new rand, along with randint(n) which returns 0
15: # to n-1, and flat() which returns 0.0 to <1.0. I did it in assembler
16: # (Mea Maxima Culpa!) to use the extended multiply and some bit fiddling.
17:
18: # Yes, I realize that the bottom bits aren't random. In fact, the bottom
19: # n bits have a period of 2^n. The rng delivered, though, throws out the
20: # most significant bit to produce a 31 bit number, and claims that it has
21: # a period of 2^32.
22:
23: # The actual generator is the routine rand, the global routines just
24: # do range conversion.
25:
26: # Chris Kingsley
27:
28: # Adapted to f77 by David Wasley, U.C.Berkeley
29: # April 1983
30: # @(#)rand_.s 1.1
31:
32: .data
33: .align 2
34: _randx: .long 1 # current value
35: _nitval:.long 1 # seed
36: .text
37: .align 1
38:
39: .globl _isrand_ # set the random seed
40: _isrand_: .word 0 # isrand(seed) int seed;
41: movl _nitval,r0 # return old seed
42: movl *4(ap),_randx
43: movl *4(ap),_nitval
44: ret
45: .align 1
46:
47: .globl _irand_ # give a 31 bit random positive integer
48: _irand_:.word 0 # integer rand(flag) int flag
49: tstl *4(ap) # 0 is normal
50: beql ir1
51: cmpl $1,*4(ap) # if arg is 1, restart
52: bneq ir0
53: movl _nitval,_randx
54: jbr ir1
55: ir0: movl *4(ap),_randx # new seed
56: movl *4(ap),_nitval # new seed
57: ir1: jsb rand
58: bicl2 $1,r0
59: rotl $-1,r0,r0
60: ret
61:
62: .globl _irandn_ # give a random positive integer from 0 to n-1
63: _irandn_:.word 0xc # integer irandn(n) int n;
64: jsb rand
65: emul *4(ap),r0,$0,r2
66: tstl r0
67: jgeq irn1
68: addl3 *4(ap),r3,r0
69: jbr irn2
70: irn1: movl r3,r0
71: irn2: ret
72: .align 1
73:
74: # compute the next 32 bit random number
75: rand: mull3 $505360173,_randx,r0
76: addl2 $907633385,r0
77: movl r0,_randx
78: rsb
79: .align 1
80:
81: .globl _drand_ # give a random double from 0. to <1.
82: _drand_: .word 0xc # double precision drand(flag)
83: dr0: tstl *4(ap) # 0 is normal
84: beql dr2
85: cmpl $1,*4(ap) # if arg is 1, restart
86: bneq dr1
87: movl _nitval,_randx
88: jbr dr2
89: dr1: movl *4(ap),_randx # new seed
90: movl *4(ap),_nitval # new seed
91: dr2: jsb rand
92: movl r0,r2
93: movf $0f1.0,r0
94: extzv $25,$7,r2,r3
95: insv r3,$0,$7,r0
96: extzv $9,$16,r2,r3
97: insv r3,$16,$16,r0
98: extzv $0,$9,r2,r1
99: subd2 $0d1.0,r0
100: ret
101:
102: .globl _rand_ # fake entry for single precision rand
103: _rand_: .word 0xc
104: jbr dr0
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.