MINICERN/packlib/kernlib/kerngen/tcgen/ranmar.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/02/15 17:49:53  mclareni
   6 * Kernlib
   7 *
   8 *
   9 #include "kerngen/pilot.h"
  10       SUBROUTINE RANMAR(RVEC,LENV)
  11 #if defined(CERNLIB_QMCRY)
  12 CDIR$ STACK
  13 #endif
  14 C
  15 C CERN PROGLIB# V113    RANMAR          .VERSION KERNFOR  4.21  890323
  16 C ORIG. 01/03/89 FCA + FJ
  17 C
  18       DIMENSION RVEC(*)
  19 C
  20       COMMON/RANMA1/IJKL,NTOT,NTOT2,I97,J97,C,U(97)
  21       LOGICAL FIRST
  22       PARAMETER (TWOM24=2.**(-24),TWOM48=2.**(-48))
  23       PARAMETER (CD=7654321.*TWOM24,CM=16777213.*TWOM24)
  24       PARAMETER (CINT=362436.*TWOM24,MODCNS=1000000000)
  25       SAVE /RANMA1/, FIRST
  26       DATA FIRST/.TRUE./
  27 C
  28       IF(FIRST) THEN
  29         IJKL = 54217137
  30         NTOT = 0
  31         NTOT2 = 0
  32         GO TO 70
  33       ENDIF
  34 C
  35    80 CONTINUE
  36       DO 100 IVEC= 1, LENV
  37         UNI = U(I97)-U(J97)
  38         IF (UNI .LT. 0.) UNI=UNI+1.
  39         U(I97) = UNI
  40         I97 = I97-1
  41         IF (I97 .EQ. 0)  I97=97
  42         J97 = J97-1
  43         IF (J97 .EQ. 0)  J97=97
  44         C = C - CD
  45         IF (C .LT. 0.)   C=C+CM
  46         UNI = UNI-C
  47         IF (UNI .LT. 0.) UNI=UNI+1.
  48 C
  49 C   Replace exact zeroes by uniform distr. *2**-24
  50 C
  51         IF (UNI .EQ. 0.)  THEN
  52           UNI = TWOM24*U(2)
  53 C
  54 C   An exact zero here is very unlikely, but let's be safe.
  55 C
  56           IF (UNI .EQ. 0.) UNI= TWOM48
  57         ENDIF
  58         RVEC(IVEC) = UNI
  59   100 CONTINUE
  60 C
  61       NTOT = NTOT + LENV
  62       IF (NTOT .GE. MODCNS)  THEN
  63         NTOT2 = NTOT2 + 1
  64         NTOT  = NTOT - MODCNS
  65       ENDIF
  66       RETURN
  67       ENTRY RMARIN(IJKLIN,NTOTIN,NTO2IN)
  68 C
  69       FIRST = .FALSE.
  70       IJKL  = IJKLIN
  71       NTOT  = NTOTIN
  72       NTOT2 = NTO2IN
  73 C
  74    70 CONTINUE
  75       IJ = IJKL/30082
  76       KL = IJKL - 30082*IJ
  77       I = MOD(IJ/177, 177) + 2
  78       J = MOD(IJ, 177)     + 2
  79       K = MOD(KL/169, 178) + 1
  80       L = MOD(KL, 169)
  81       DO 30 II= 1, 97
  82         S = 0.
  83         T = .5
  84         DO 20 JJ= 1, 24
  85           M = MOD(MOD(I*J,179)*K, 179)
  86           I = J
  87           J = K
  88           K = M
  89           L = MOD(53*L+1, 169)
  90           IF (MOD(L*M,64) .GE. 32)  S = S+T
  91           T = 0.5*T
  92   20    CONTINUE
  93         U(II) = S
  94   30  CONTINUE
  95       C   = CINT
  96       I97 = 97
  97       J97 = 33
  98 C       Complete initialization by skipping
  99 C            (NTOT2*MODCNS + NTOT) random numbers
 100       NITER = MODCNS
 101       DO 50 LOOP2= 1, NTOT2+1
 102         IF(LOOP2.GT.NTOT2) NITER=NTOT
 103         DO 40 IDUM = 1, NITER
 104           UNI = U(I97)-U(J97)
 105           IF (UNI .LT. 0.) UNI=UNI+1.
 106           U(I97) = UNI
 107           I97 = I97-1
 108           IF (I97 .EQ. 0)  I97=97
 109           J97 = J97-1
 110           IF (J97 .EQ. 0)  J97=97
 111           C = C - CD
 112           IF (C .LT. 0.)   C=C+CM
 113    40   CONTINUE
 114    50 CONTINUE
 115       NTOT  = 0
 116       NTOT2 = 0
 117       IF(FIRST) THEN
 118         FIRST = .FALSE.
 119         GO TO 80
 120       ENDIF
 121       RETURN
 122       ENTRY RMARUT(IJKLUT,NTOTUT,NTO2UT)
 123 C
 124       NTOTUT = NTOT
 125       NTO2UT = NTOT2
 126       IJKLUT = IJKL
 127 C
 128       END