      PROGRAM SHIELD
************************************************************************
* This program uses the random number generator  RAND  to simulate     *
* particles entering a shield and to determine what percentage reaches *
* the outside.  The particles are assumed to move forward, backward,   *
* left, and right with equal likelihood and to die within the shield   *
* if a certain number of collisions (changes of direction) have        *
* occurred. Variables used are:                                        *
*     THICK  : thickness of shield                                     *
*     LIMCOL : limit on # of collisions (before energy dissipated)     *
*     DIR    : a random integer 1, 2, 3, or 4 representing direction   *
*     OLDDIR : previous direction of particle                          *
*     COLLIS : # of collisions (changes of direction)                  *
*     FORW   : net units forward traveled                              *
*     N      : # of particles simulated                                *
*     COUNT  : # of particles reaching outside of shield               *
*     SEED   : a seed for RAND                                         *
*     I      : index variable                                          *
*                                                                      *
* Input:   THICK, LIMCOL, N, and SEED                                  *
* Output:  Percentage of particles that reach the outside              *
************************************************************************

      INTEGER THICK, LIMCOL, DIR, OLDDIR, COLLIS, FORW, N, COUNT,
     +        SEED, I
      REAL RAND

      PRINT *, 'ENTER THICKNESS OF SHIELD, LIMIT ON # OF COLLISIONS,'
      PRINT *, 'AND THE NUMBER OF PARTICLES TO SIMULATE'
      READ *, THICK, LIMCOL, N
      PRINT *, 'SEED FOR RANDOM NUMBER GENERATOR'
      READ *, SEED
      COUNT = 0

* Begin the simulation

      DO 20 I = 1, N
          FORW = 1
          OLDDIR = 0
          COLLIS = 0

*         Repeat the following until particle reaches outside of
*         shield, returns inside reactor, or dies within shield

10        CONTINUE
              DIR = 1 + 4 * RAND(SEED)
              IF (DIR .NE. OLDDIR) THEN
                  COLLIS = COLLIS + 1
                  OLDDIR = DIR
              END IF
              IF (DIR .EQ. 1) THEN
                  FORW = FORW + 1
              ELSE IF (DIR .EQ. 2) THEN
                  FORW = FORW - 1
              END IF
          IF ((FORW .LT. THICK) .AND. (FORW .GT. 0) .AND.
     +        (COLLIS .LT. LIMCOL)) GO TO 10

          IF (FORW .EQ. THICK) THEN
              COUNT = COUNT + 1
          END IF
20    CONTINUE

      PRINT 30, 100 * COUNT / REAL(N)
30    FORMAT(1X, F6.2, '% OF THE PARTICLES ESCAPED')
      END


**RAND******************************************************************
* This function generates a random real number in the interval from    *
* 0 to 1.  The variable M is initially the seed supplied by the        *
* user;  thereafter, it is the random integer generated on the pre-    *
* ceding call to the function and saved using the SAVE statement.      *
* NOTE:  The constants 2147483647 and .4656613E-9 used in this sub-    *
* program are appropriate when it is executed on a machine having      *
* 32-bit memory words.  For a machine having M-bit words, these two    *
* constants should be replaced by the values of  2**M - 1  and         *
* 1/(2**M - 1) , respectively.                                         *
************************************************************************

      FUNCTION RAND(K)

      INTEGER K, M, CONST1
      REAL RAND, CONST2
      PARAMETER (CONST1 = 2147483647, CONST2 = .4656613E-9)
      SAVE
      DATA M /0/

      IF (M .EQ. 0) M = K
      M = M * 65539
      IF (M .LT. 0) M = (M + 1) + CONST1
      RAND = M * CONST2
      END

