This code provides support for a 64-bit linear congruential random number generator (in conjunction with the LinearVariate C++ class not shown here). Other than minor formatting for HTML presentation, this code is exactly as found in the original source code.
title Copyright (c) 1993..1996 by Michael Lee Finney. All rights reserved.
subtitle LinearVariate support
include Standard.inc
if nt eq true
public ?advance__13LinearVariateFv
public ?retreat__13LinearVariateFv
else
public advance__13LinearVariateFv
public retreat__13LinearVariateFv
endif
;----------------------------------------------------------------------------*
; *
; This module provides support for the LinearVariate class. *
; *
;----------------------------------------------------------------------------*
;
; This section contains a structure definition for the LinearVariate
; class. The structure elements MUST be defined in the same order and
; be of the same size as those in the C++ definition. These definitions
; are compiler dependent. The IBM VisualAge C++ compilers are assumed.
;
Unsigned_8 struct
;__vptr unsigned4 ?
_lo unsigned4 ?
_hi unsigned4 ?
Unsigned_8 ends
Variate struct
;__vptr unsigned4 ?
_n Unsigned_8 <>
_x extended ?
unsigned2 ?,?,?
Variate ends
LinearVariate struct
__vptr unsigned4 ?
_variate Variate <>
_seed unsigned8 ?
LinearVariate ends
;
; This section contains all data defined in assembly;
;
DATA32 segment PARA USE32 PUBLIC 'DATA'
;
; These constants are the constants required for the linear
; congruential generator defined by:
;
; X = (6364136223846793005 * X + 1442695040888963407) mod 2**64
;
; and its inverse defined by:
;
; X = (13877824140714322085 * X + 11066951453180645397) mod 2**64
;
; Which as measured by Knuth, is a very good generator with a period
; of 2**64.
;
align 4
lineara1 unsigned4 04c957f2dh ; Low 32 bits of LCG multiplier
lineara2 unsigned4 05851f42dh ; Top 32 bits of LCG multiplier
linearc1 unsigned4 0f767814fh ; Low 32 bits of LCG constant
linearc2 unsigned4 014057b7eh ; Top 32 bits of LCG constant
raenila1 unsigned4 0329e28a5h ; Low 32 bits of reverse LCG multiplier
raenila2 unsigned4 0c097ef87h ; Top 32 bits of reverse LCG multiplier
raenilc1 unsigned4 021535015h ; Low 32 bits of reverse LCG constant
raenilc2 unsigned4 09995b5b6h ; Top 32 bits of reverse LCG constant
;
; This value is used to adjust for loading a integral value equal
; or larger than 2**63 into a into a floating point register.
;
two64 real4 18446744073709551616.0
DATA32 ends
;
; This section contains all code defined in assembly;
;
CODE32 segment
;
; void LinearVariate::advance()
;
; The operation implemented by this method is...
;
; _x = 6364136223846793005 * _x + 1442695040888963407;
; _v = _x;
;
; where...
;
; 6364136223846793005 = 0x5851.f42d.4c95.7f2d
; 1442695040888963407 = 0x1405.7b7e.f767.814f
;
align 16
advance__13LinearVariateFv equ $
?advance__13LinearVariateFv equ $
push ebx
mov ebx,eax
mov eax,(LinearVariate ptr [eax])._variate._n._lo
mul lineara2
mov ecx,eax
mov eax,(LinearVariate ptr [ebx])._variate._n._hi
mul lineara1
add ecx,eax
mov eax,(LinearVariate ptr [ebx])._variate._n._lo
mul lineara1
add edx,ecx
add eax,linearc1
adc edx,linearc2
mov (LinearVariate ptr [ebx])._variate._n._lo,eax
mov (LinearVariate ptr [ebx])._variate._n._hi,edx
js normalized1
jz topZero1
mov ecx,3ffeh
normalize1:
dec ecx
add eax,eax
adc edx,edx
jns normalize1
mov unsigned4 ptr (LinearVariate ptr [ebx])._variate._x,eax
mov unsigned4 ptr (LinearVariate ptr [ebx+4])._variate._x,edx
mov unsigned4 ptr (LinearVariate ptr [ebx+8])._variate._x,ecx
pop ebx
ret
normalized1:
mov unsigned4 ptr (LinearVariate ptr [ebx])._variate._x,eax
mov unsigned4 ptr (LinearVariate ptr [ebx+4])._variate._x,edx
mov unsigned4 ptr (LinearVariate ptr [ebx+8])._variate._x,3ffeh
pop ebx
ret
topZero1:
mov ecx,(3ffeh - 32)
or eax,eax
js topNormalized1
jz allZero1
topNormalize1:
dec ecx
add eax,eax
jns topNormalize1
topNormalized1:
mov unsigned4 ptr (LinearVariate ptr [ebx])._variate._x,edx
mov unsigned4 ptr (LinearVariate ptr [ebx+4])._variate._x,eax
mov unsigned4 ptr (LinearVariate ptr [ebx+8])._variate._x,ecx
pop ebx
ret
allZero1:
mov unsigned4 ptr (LinearVariate ptr [ebx])._variate._x,edx
mov unsigned4 ptr (LinearVariate ptr [ebx+4])._variate._x,eax
mov unsigned4 ptr (LinearVariate ptr [ebx+8])._variate._x,edx
pop ebx
ret
;
; void LinearVariate::retreat()
;
; The operation implemented by this method is...
;
; _x = 13877824140714322085 * _x + 11066951453180645397;
; _v = _x;
;
; where...
;
; 13877824140714322085 = 0xc097.ef87.329e.28a5
; 11066951453180645397 = 0x9995.b5b6.2153.5015
;
; Note: This code uses a different algorithm than the advance()
; method for cross checking purposes. The advance() method
; is faster on a 486, but this code may be faster on later
; processors. Since the advance() method is used much more
; frequently than the retreat() method, if the targeted
; processors become Pentium or later then the algorithm
; used here may replace the algorithm used in the advance()
; method after testing the speed differences.
;
align 16
retreat__13LinearVariateFv equ $
?retreat__13LinearVariateFv equ $
push ebx
mov ebx,eax
mov eax,(LinearVariate ptr [ebx])._variate._n._lo
mul raenila2
mov ecx,eax
mov eax,(LinearVariate ptr [ebx])._variate._n._hi
mul raenila1
add ecx,eax
mov eax,(LinearVariate ptr [ebx])._variate._n._lo
mul raenila1
add edx,ecx
add eax,raenilc1
adc edx,raenilc2
mov (LinearVariate ptr [ebx])._variate._n._lo,eax
mov (LinearVariate ptr [ebx])._variate._n._hi,edx
fild unsigned8 ptr (LinearVariate ptr [ebx])._variate._n._lo
jns loadPositive2
fadd two64
loadPositive2:
or edx,eax
mov eax,ebx
pop ebx
fstp (LinearVariate ptr [eax])._variate._x
jz loadZero2
sub unsigned4 ptr (LinearVariate ptr [eax+8])._variate._x,40h
loadZero2:
ret
CODE32 ends
end