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