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