1 " ******************************************************
2 " * *
3 " * *
4 " * Copyright c 1972 by Massachusetts Institute of *
5 " * Technology and Honeywell Information Systems, Inc. *
6 " * *
7 " * *
8 " ******************************************************
9
10 """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
11 " This procedure generates pseudo-random numbers using the
12 " Tausworth method. 36 bits are used in the generation.
13 "
14 " There are multiple entry points. For all entry points:
15 " The first argument is a fixed binary input argument,
16 " which is a non-zero integer. This is an optional argument--
17 " if not provided by caller, a value maintained in internal
18 " static is used. This value, from either source,
19 " is the seed for the random number generator. Its value is
20 " modified so that upon return it has the value that should
21 " be used as the seed for the next call.
22 "
23 " There are a set of entry points with two arguments which
24 " are used to generate a single random number. For these:
25 " The second argument is a floating point return argument
26 " that returns the value of the random number generated.
27 "
28 " There are a set of entry points with three arguments which
29 " are used to generate a sequence of random numbers. For these:
30 " The second argument is an array of single precision
31 " floating point numbers. This array returns a sequence of
32 " of random numbers, beginning at the base of the array.
33 " The third argument is a fixed binary17 input
34 " argument the specifies the size of the array.
35 "
36 " Coded 1 January 1970 by Roger R. Schell
37 """"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
38
39 name random_
40
41 " Table of contents
42 entry set_seed
43 entry get_seed
44 entry random_
45 entry uniform
46 entry uniform_ant
47 entry uniform_seq
48 entry uniform_ant_seq
49 entry normal
50 entry normal_ant
51 entry normal_seq
52 entry normal_ant_seq
53 entry exponential
54 entry exponential_seq
55
56
57 equ shift,11
58 equ size,36
59 equ expon,0 exponent to convert integer to floating point
60
61
62 " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " "
63 "
64 " CODING CONVENTIONS
65 "
66 " XR0 used for return address for specific distribution subroutines
67 " XR1 used for return address for generator primitive subroutine
68 " XR2 general purpose register for distribution subroutines
69 " XR3 usedto indicate: 1=> antithetic variable, 0=> usual
70 " XR4 contains the address of the distribution subroutine for this call
71 " XR5 index into return array for the next random number
72 " XR6 count of the number of values to be generated after current one
73 " XR7 general purpose register
74 "
75 " A-reg distribution routine uses to return floating point value
76 " Q-reg always has the seed used by primitive generator
77 "
78 " BP pointer to base of return arguments
79 " AP pointer to the seed
80 "
81 " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " " "
82
83
84
85
86 "
87 " call random_$set_seedseed;
88
89
90 set_seed:
91 ldq ap|2,* qet new seed into Q-reg
92 stq lp|internal_seed save as new value of seed
93 tra return return to caller
94
95
96
97 "
98 " call random_$get_seedseed;
99
100
101 get_seed:
102 ldq lp|internal_seed get current value of seed
103 stq ap|2,* return value to caller
104 tra return return to caller
105
106
107
108 "
109 " call random_$uniformseedrandom_no;
110
111
112 random_:
113 uniform:
114 eax4 uniform_ set XR4 to the address of uniform distribution
115 tra single this entry generates a single random number
116
117
118
119 "
120 " call random_$uniform_antseedrandom_no;
121 " This entry gives negatively correlated value.
122
123
124
125 uniform_ant:
126 eax4 uniform_ant_ set up the proper distribution
127 tra single this entry generates a single random number
128
129
130
131 "
132 " call random_$uniform_seqseedarrayarray_size;
133 " This entry gives an array of return values
134
135
136
137 uniform_seq:
138 eax4 uniform_ we generate sequence from uniform distribution
139 tra sequence this entry gives a sequence of numbers
140
141
142
143 "
144 " call random_$uniform_ant_seqseedarrayarray_size;
145
146
147 uniform_ant_seq:
148 eax4 uniform_ant_ negatively correlated generator
149 tra sequence this entry gives a sequence of numbers
150
151
152
153 "
154 " call random_$normalseedrandom_no;
155
156
157 normal:
158 eax4 normal_ normal distribution
159 eax3 0 not negatively correlated value
160 tra single this entry gives a single number
161
162
163
164 "
165 " call random_$normal_antseedrandom_no;
166
167
168 normal_ant:
169 eax4 normal_ normal distribution
170 eax3 1 negatively correlated
171 tra single this entry gives single number
172
173
174
175 "
176 " call random_$normal_seqseedarrayarray_size;
177
178
179
180 normal_seq:
181 eax4 normal_ normal distribution
182 eax3 0 not negatively correlated value
183 tra sequence this entry gives a sequence of numbers
184
185
186
187 "
188 " call random_$normal_ant_seqseedarrayarray_size;
189
190
191 normal_ant_seq:
192 eax4 normal_ normal distribution
193 eax3 1 negatively correlated
194 tra sequence this entry gives a sequence of numbers
195
196
197
198 "
199 " call random_$exponentialseedrandom_no;
200
201
202 exponential:
203 eax4 exponential_ exponential distribution
204 tra single this entry gives a single value
205
206
207
208 "
209 " call random_$exponential_seqseedarrayarray_size;
210
211
212 exponential_seq:
213 eax4 exponential_ exponential distribution
214 tra sequence this entry gives a sequence of numbers
215
216
217
218 "!!!!!!!!!!--set up the number of values to be generated--!!!!!!!!!!
219
220
221 sequence:
222 ldx7 ap|0 twice number of arguments in XR7
223 lxl6 ap|0,7* length of sequence to XR6
224 eax7 -2,7 subtract two from XR7
225 eax6 -1,6 decrement by one
226 tpl common if positive value, use it
227 tra return if zero or negative, return
228
229 single:
230 ldx7 ap|0 twice number of arguments in XR7
231 eax6 0 use sequence of length one
232
233 common:
234 eppbp ap|0,7* set bp to point to first return value
235 eax5 0 index into array is in XR5
236 eaa -2,7 upper A-reg has offset of seed in arglist
237 ars 19 should be one or zero in A-reg
238 xec set_ap,al set ap to point to the seed
239 szn ap|0 test for a seed of zero
240 tnz loop if non-zero continue
241 eax4 zero_arg if zero, generate zero return values
242
243 loop:
244 tsx0 0,4 go to appropriate generator
245 fst bp|0,5 store value returned by generator
246 eax5 1,5 increment index into array
247 eax6 -1,6 decrement count of remaining
248 tpl loop if not done, loop again
249
250 return:
251 short_return
252
253
254
255 set_ap: "get pointer to seed--from caller or default
256 eppap lp|internal_seed use internal value if not provided in call
257 eppap ap|2,* seed is the first argument if provided
258
259
260
261
262 "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
263 "$ This is the primitive that actually generates the
264 "$ random number in integer form from the seed.
265 "$
266 "$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
267
268
269
270 generate:
271 ldq ap|0 load seed into Q-reg
272 qrl shift shift right the seed
273 ersq ap|0 exclusive or to the seed
274 ldq ap|0 put same value in Q-reg
275 qls size-shift shift the result left
276 ersq ap|0 exclusive or to previous result
277 tra 0,1 return to the caller of primitive
278
279
280
281
282 "!!!!!!!!!!--zero argument generator--!!!!!!!!!!
283
284 zero_arg: "used if input seed is zero
285 fld =0.,du load a floating point zero
286 tra 0,0 return
287
288
289
290
291 "!!!!!!!!!!--uniform generator--!!!!!!!!!!
292
293 uniform_:
294 tsx1 generate generate one random number
295 lda ap|0 load A-reg with integer value
296 arl 1 make it a positive number
297 lde expon,du convert to floating point
298 fad =0.,du normalize
299 tra 0,0 return
300
301
302 uniform_ant_:
303 tsx1 generate generate one random number
304 lda ap|0 load integer value into A-reg
305 arl 1 make it a positive number
306 lde expon,du convert to floating point
307 fneg "take negative value
308 fad =1.,du normalize
309 tra 0,0 return
310
311
312
313
314 "!!!!!!!!!!--exponential generator--!!!!!!!!!!
315
316 exponential_:
317 eax7 -1 count number of 'runs' with XR7
318 outer:
319 eax7 1,7 add one to count of 'runs'
320 eax2 1 use as counter of 'run' length
321 "initialize XR2 with a count of one
322 tsx1 generate go to primitive generator
323 lda ap|0 get seed in A-reg
324 arl 1 make it a positive number
325 lde expon,du convert to floating point
326 fst bp|0,5 store it temporarily in return value
327 inner:
328 lda ap|0 keep value in A-reg for comparison
329 tsx1 generate generate another value
330 eax2 1,2 add one to count of 'run'length
331 cmpa ap|0 compare last number with new one
332 trc inner if still a run down,loop again
333 anx2 =1,du check if 'run' has even length
334 tnz outer if not even, get another run
335 eaa 0,7 no of runs before even length
336 lde =17b25,du convert to floating point
337 fad bp|0,5 add first random number to number of 'runs'
338 tra 0,0 return
339
340
341
342
343 "!!!!!!!!!!--normal distribution generator--!!!!!!!!!!
344
345 normal_:
346 fld =0.,du load a zero
347 eax2 12 use XR2 to count 12 times thru loop
348 n_loop:
349 fst bp|0,5 store the new sum
350 tsx1 generate generate the next random number
351 lda ap|0 load seed into A-reg
352 arl 1 make it a positive number
353 lde expon,du convert to floating point
354 fad bp|0,5 add random number to sum
355 eax2 -1,2 decrement counter by one
356 tnz n_loop accumulate twelve numbers
357 fsb =6.,du give a mean of zero
358 xec n_norm,3 antithetic if appropriate
359 tra 0,0 return
360
361 n_norm:
362 nop "o.k. as is if not antithetic
363 fneg "take negative for antithetic
364
365
366
367
368
369 "
370 " INTERNAL STATIC DATA
371 "
372
373 use .lkstat.
374 join /link/.lkstat.
375
376 internal_seed:
377 dec 4084114320 "initial internal seed for a new process
378
379 use main
380
381
382 end