10   '
20   '  CRITICAL MASS FOR A FISSION CHAIN REACTION
30   '
40   ' ----------------------- 1/15/88 ---------------------
50   '
60   '    This program computes the survival fraction, F,
     '    for a rectangular slab of fissionable material.  F
     '    is defined as the number of induced fissions per
     '    spontaneous fission.  For values of F greater than
     '    1.0, a chain reaction occurs.  The slab has
     '    "critical mass" if it has F=1.  We assume that two
     '    neutrons are emitted per fission.
120  '
130  '           --> LIST OF PROCEDURES <--
140  '
150        GOSUB 270   'INITIALIZATIONS
160        GOSUB 400   'INPUT
170        GOSUB 680   'LOOP OVER # OF SPONTANEOUS FISSIONS
180        GOSUB 1340  'OUTPUT
190  '
200  '
210        PRINT
220        INPUT "WANT TO RUN ANOTHER (Y/N)? ", AGAIN$
230        IF AGAIN$ = "Y" OR AGAIN$ = "y" THEN 160
240  '
250        END
260  '--------------------------------------------
270  'INITIALIZATIONS
280  '
290  DIM R(9)
300  PI = 3.14159265#
310  DEFDBL I
320  '
330  ' -- IX is set to an odd integer to initialize the
     '    random number generator:
350  '
360  IX = 13: IN = 32767: IC = 899
370  '
380  RETURN
390  '--------------------------------------------
400  'INPUT
410  '
420  CLS
430  PRINT "INPUT VALUES FOR THE FOLLOWING:"
440  PRINT "  THE MASS OF THE RECTANGULAR SLAB (M)"
450  PRINT "  THE RATIO OF LENGTH TO THICKNESS (S)"
460  PRINT "  THE NUMBER OF RANDOMLY GENERATED FISSIONS (N)"
470  PRINT
480  PRINT "SEPARATE INPUT VALUES WITH COMMAS, PRESS RETURN"
490  INPUT M, S, N
500  '
510  ' --  Find the dimensions of the slab from M and S:
520  '
530  A = (M * S) ^ (1 / 3)
540  B = (M / S ^ 2) ^ (1 / 3)
550  L1 = A / 2
560  L2 = A / 2
570  L3 = B / 2
580  '
590  ' --Set NIN, the number of cases in which an emitted
     '   neutron stops within the boundaries of the slab,
     '   thus producing an induced fission, equal to zero
     '   initially.
620  '
630  NIN = 0
640  N = INT(N)
650  '
660  RETURN
670  '---------------------------------------------
680  'LOOP OVER NUMBER OF SPONTANEOUS FISSIONS
690  '
700  PRINT
710  FOR J = 1 TO N
720    '
730    PRINT " > We are at fission #"; J;
740    LOCATE , 1
750    '
760    ' --For each spontaneous fission we need 9 random
       '   numbers
770    '
780    GOSUB 1240 'GET 9 RANDOM NUMBERS R(N) IN (0,1)
790    '
800    ' -- Calculate the X, Y, and Z coordinates of the
       '    nucleus undergoing spontaneous emission.  They
       '    are random numbers inside the boundary of the
       '    block.
830    '    Note: X,Y are in the interval (-A/2,+A/2)
840    '            Z  is in the interval (-B/2,+B/2)
850    '
860    X0 = A * (R(1) - .5)
870    Y0 = A * (R(2) - .5)
880    Z0 = B * (R(3) - .5)
890    '
900    ' -- Begin loop over each of the 2 neutrons; 1 at a
       '    time
910    '
920    FOR K = 1 TO 2
930      '
940      '-- Get two angles that define a random direction
950      '      in space for the neutron:
960      '
970      PHI = 2 * PI * R(2 * K + 2)
980      CSN = 2 * (R(2 * K + 3) - .5)
990      SNE = SQR(1 - CSN ^ 2)
1000      '
1010      ' -- Get the distance traveled by the neutron, 
1020      '    assumed to be a random number in the 
1030      '    interval (0,1):
1040      D = R(K + 7)
1050     '
1060     ' -- Calculate the coordinates of the end point 
1070     '    (interaction point) for the neutron
1080     '
1090     X1 = X0 + D * SNE * COS(PHI)
1100     Y1 = Y0 + D * SNE * SIN(PHI)
1110     Z1 = Z0 + D * CSN
1120     '
1130     ' -- If the end point is inside the block, the
         '    neutron produces an induced fission so add 1
         '    to NIN:
1150     '
1160     IF ABS(X1)<= L1 AND ABS(Y1)<= L2 AND_
             ABS(Z1)<= L3 THEN 
1170     NIN = NIN + 1
1180     NEXT K  'End of 1-neutron loop for this fission.
1190     '
1200   NEXT J  'End of 1-fission loop.
1210   '
1220 RETURN
1230 '--------------------------------------------
1240 'RANDOM NUMBER GENERATOR
1250 '
1260 FOR K = 1 TO 9
1270    IX = ABS(IC * IX)
1280    IX = IX - (IN + 1) * INT(IX / (IN + 1))
1290    R(K) = IX / IN
1300    NEXT K
1310 '
1320 RETURN
1330 '--------------------------------------------
1340 'OUTPUT
1350 '
1360 PRINT
1370 PRINT
1380 INPUT "PRINTER (Y/N)? ", ANS$
1390 IF ANS$ = "y" OR ANS$ = "Y" THEN 
1400 OPEN "LPT1:" FOR OUTPUT AS #1
1410 IF ANS$ <> "y" AND ANS$ <> "Y" THEN 
1415 OPEN "SCRN:" FOR OUTPUT AS #1
1420 PRINT #1,
1430 PRINT #1, "             MASS= "; M
1440 PRINT #1, " LENGTH/THICKNESS= "; S
1450 PRINT #1, "# RANDOM FISSIONS= "; N
1460 '
1470 ' --  The survival fraction is NIN/N:
1480 '
1490 F = NIN / N
1500 '
1510 PRINT #1,
1520 PRINT #1, "THE SURVIVAL FRACTION=   "; F
1530 '
1540 CLOSE #1
1550 '
1560 RETURN
1570 '------------------ END OF PROGRAM LINES ------------------