1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
|
INTEGER FUNCTION ignlgi()
C**********************************************************************
C
C INTEGER FUNCTION IGNLGI()
C GeNerate LarGe Integer
C
C Returns a random integer following a uniform distribution over
C (1, 2147483562) using the current generator.
C
C This is a transcription from Pascal to Fortran of routine
C Random from the paper
C
C L'Ecuyer, P. and Cote, S. "Implementing a Random Number Package
C with Splitting Facilities." ACM Transactions on Mathematical
C Software, 17:98-111 (1991)
C
C**********************************************************************
C .. Parameters ..
INTEGER numg
PARAMETER (numg=32)
C ..
C .. Scalars in Common ..
INTEGER a1,a1vw,a1w,a2,a2vw,a2w,m1,m2
C ..
C .. Arrays in Common ..
INTEGER cg1(numg),cg2(numg),ig1(numg),ig2(numg),lg1(numg),
+ lg2(numg)
LOGICAL qanti(numg)
C ..
C .. Local Scalars ..
INTEGER curntg,k,s1,s2,z
LOGICAL qqssd
C ..
C .. External Functions ..
LOGICAL qrgnin
EXTERNAL qrgnin
C ..
C .. External Subroutines ..
EXTERNAL getcgn,inrgcm,rgnqsd,setall
C ..
C .. Common blocks ..
COMMON /globe/m1,m2,a1,a2,a1w,a2w,a1vw,a2vw,ig1,ig2,lg1,lg2,cg1,
+ cg2,qanti
C ..
C .. Save statement ..
SAVE /globe/
C ..
C .. Executable Statements ..
C
C IF THE RANDOM NUMBER PACKAGE HAS NOT BEEN INITIALIZED YET, DO SO.
C IT CAN BE INITIALIZED IN ONE OF TWO WAYS : 1) THE FIRST CALL TO
C THIS ROUTINE 2) A CALL TO SETALL.
C
IF (.NOT. (qrgnin())) CALL inrgcm()
CALL rgnqsd(qqssd)
IF (.NOT. (qqssd)) CALL setall(1234567890,123456789)
C
C Get Current Generator
C
CALL getcgn(curntg)
s1 = cg1(curntg)
s2 = cg2(curntg)
k = s1/53668
s1 = a1* (s1-k*53668) - k*12211
IF (s1.LT.0) s1 = s1 + m1
k = s2/52774
s2 = a2* (s2-k*52774) - k*3791
IF (s2.LT.0) s2 = s2 + m2
cg1(curntg) = s1
cg2(curntg) = s2
z = s1 - s2
IF (z.LT.1) z = z + m1 - 1
IF (qanti(curntg)) z = m1 - z
ignlgi = z
RETURN
END
|
