
      DOUBLE PRECISION FUNCTION NORMAL(DUMMY)
C     ---------------------------------------

      DOUBLE PRECISION     DUMMY

C     ----------------------------------------------------------------------
C     |  Returns a random REAL value, normally distributed with standard   |
C     |   deviation = SQRT (2).                                            |
C     |  The algorithm actually generates 2 independent normals at once -  |
C     |   so the second is saved and returned at the next call.            |
C     ----------------------------------------------------------------------

      INTEGER  SEED
      DOUBLE PRECISION     N2,  RAND
      LOGICAL  STORED
      COMMON / RANDOM / N2, SEED, STORED

      DOUBLE PRECISION        TWOPI
      PARAMETER  (TWOPI = 6.28318531D0)
      DOUBLE PRECISION     LLOG, THETA, N1, X, Y

      IF  (STORED)  THEN
          NORMAL = N2
          STORED = .FALSE.

C------... else generate another 2 normals, 'n1' and 'n2'  (and return 'n1')

        ELSE
  10      CONTINUE
            X = RAND (SEED)
            Y = RAND (SEED)
          IF  (X  .EQ.  0.0D0) GO TO 10
          LLOG    = SQRT (- LOG (X))
          THETA  = TWOPI * Y
          N1     = LLOG * SIN (THETA)
          N2     = LLOG * COS (THETA)
          NORMAL = N1
          STORED = .TRUE.
      ENDIF
      END
