*$ CREATE SOURCE.FOR *COPY SOURCE * *=== source ===========================================================* * SUBROUTINE SOURCE ( NOMORE ) INCLUDE '(DBLPRC)' INCLUDE '(DIMPAR)' INCLUDE '(IOUNIT)' *----------------------------------------------------------------------* * Version for isotropic Cf source. * DIMENSION CFCUM(0:72), ENEDGE(73), ANG(0:40), CUMANG(0:40) *----------------------------------------------------------------------* * INCLUDE '(BEAMCM)' INCLUDE '(CASLIM)' INCLUDE '(FHEAVY)' INCLUDE '(FLKSTK)' INCLUDE '(IOIOCM)' INCLUDE '(LTCLCM)' INCLUDE '(PAPROP)' INCLUDE '(SOURCM)' INCLUDE '(SUMCOU)' * LOGICAL LFIRST SAVE LFIRST DATA LFIRST / .TRUE. / *----------------------------------------------------------------------* * Neutron energy group boundaries DATA ENEDGE / & 1.9600E-02, 1.7500E-02, 1.4918E-02, 1.3499E-02, 1.2214E-02, & 1.1052E-02, 1.0000E-02, 9.0484E-03, 8.1873E-03, 7.4082E-03, & 6.7032E-03, 6.0653E-03, 5.4881E-03, 4.9659E-03, 4.4933E-03, & 4.0657E-03, 3.6788E-03, 3.3287E-03, 3.0119E-03, 2.7253E-03, & 2.4660E-03, 2.2313E-03, 2.0190E-03, 1.8268E-03, 1.6530E-03, & 1.4957E-03, 1.3534E-03, 1.2246E-03, 1.1080E-03, 1.0026E-03, & 9.0718E-04, 8.2085E-04, 7.4274E-04, 6.0810E-04, 4.9787E-04, & 4.0762E-04, 3.3373E-04, 2.7324E-04, 2.2371E-04, 1.8316E-04, & 1.4996E-04, 1.2277E-04, 8.6517E-05, 5.2475E-05, 3.1828E-05, & 2.1852E-05, 1.5034E-05, 1.0332E-05, 7.1018E-06, 4.8809E-06, & 3.3546E-06, 2.3054E-06, 1.5846E-06, 1.0446E-06, 6.8871E-07, & 4.5400E-07, 2.7537E-07, 1.6702E-07, 1.0130E-07, 6.1442E-08, & 3.7267E-08, 2.2603E-08, 1.5535E-08, 1.0677E-08, 7.3375E-09, & 5.0435E-09, 3.4662E-09, 2.3824E-09, 1.6374E-09, 1.1254E-09, & 6.8257E-10, 4.1400E-10, 1.0000E-14/ * Cumulative Cf spectrum DATA CFCUM / 0.D0, *...+....1....+....2....+....3....+....4....+....5....+....6....+....7.. & 1.919922147E-05,1.076421528E-04,2.985396423E-04,6.750280793E-04, & 1.460097913E-03,2.919855883E-03,5.440437967E-03,9.517761940E-03, & 1.572873215E-02,2.550738486E-02,3.904978251E-02,5.635569494E-02, & 7.365912720E-02,9.966894989E-02,1.279101323E-01,1.650685484E-01, & 2.022289923E-01,2.393916485E-01,2.836610364E-01,3.279253392E-01, & 3.733559139E-01,4.193675039E-01,4.654041752E-01,5.064859139E-01, & 5.475773136E-01,5.886592151E-01,6.297539805E-01,6.708701087E-01, & 6.974131472E-01,7.239697428E-01,7.505235450E-01,7.770757115E-01, & 8.301848328E-01,8.661373176E-01,8.849342051E-01,9.037310414E-01, & 9.225258728E-01,9.413224153E-01,9.601180671E-01,9.756795138E-01, & 9.815263394E-01,9.879075224E-01,9.948263962E-01,9.971427050E-01, & 9.985910051E-01,9.990648305E-01,9.995399469E-01,9.996648243E-01, & 9.997897420E-01,9.999146525E-01,9.999395164E-01,9.999643698E-01, & 9.999838018E-01,9.999909490E-01,9.999946297E-01,9.999962702E-01, & 9.999979108E-01,9.999995515E-01,9.999997034E-01,9.999998553E-01, & 9.999999446E-01,9.999999646E-01,9.999999846E-01,9.999999885E-01, & 9.999999924E-01,9.999999963E-01,9.999999981E-01,9.999999988E-01, & 9.999999995E-01,9.999999997E-01,9.999999999E-01,1.000000000E+00/ * Angle cosines DATA ANG/ -1.0D0, -0.95D0, -0.9D0, -0.85D0, -0.8D0, -0.75D0, & -0.7D0, -0.65D0, -0.6D0, -0.55D0, -0.5D0, -0.45D0, & -0.4D0, -0.35D0, -0.3D0, -0.25D0, -0.2D0, -0.15D0, & -0.1D0, -5.D-02, 0.0D0, 5.D-02, 0.1D0, 0.15D0, & 0.2D0, 0.25D0, 0.3D0, 0.35D0, 0.4D0, 0.45D0, & 0.5D0, 0.55D0, 0.6D0, 0.65D0, 0.7D0, 0.75D0, & 0.8D0, 0.85D0, 0.9D0, 0.95D0, 1.0D0 / * Cumulative angular distribution DATA CUMANG / 0.0D0, 4.3787904871D-4, 9.309108266D-4, & 1.5080318037D-3, 2.1835838416D-3, 2.9743548826D-3, & 3.8909576476D-3, 4.8923191479D-3, 5.9862770681D-3, & 7.1813938468D-3, 8.4870236975D-3, 9.9133858209D-3, & 1.1471644393D-2, 1.3097707545D-2, 1.4776914040D-2, & 1.6511000723D-2, 1.8301761200D-2, 2.0151047694D-2, & 2.2060772967D-2, 2.4032912292D-2, 2.6069505497D-2, & 2.9084203728D-2, 3.3546756876D-2, 4.0152519402D-2, & 4.9930800378D-2, 6.4405250301D-2, 8.5831276094D-2, & 0.1175474784D0, 0.1569148680D0, 0.1981130258D0, & 0.2412270912D0, 0.2863461629D0, 0.3335634830D0, & 0.3829766299D0, 0.4346877199D0, 0.4919529577D0, & 0.5559660349D0, 0.6275220820D0, 0.7075099237D0, & 0.8189694928D0, 1.0D0 / *----------------------------------------------------------------------* *...+....1....+....2....+....3....+....4....+....5....+....6....+....7.. *----------------------------------------------------------------------* NOMORE = 0 * +-------------------------------------------------------------------* * | First call initializations: IF ( LFIRST ) THEN * | *** The following 3 cards are mandatory *** TKESUM = ZERZER LFIRST = .FALSE. LUSSRC = .TRUE. * | *** User initialization *** write(lunout,*) write(lunout,'(a,132a)') ("*",i=1,132) write(lunout,*) write(lunout,*) write(lunout,'(a)') "Californium source. Version corrected" write(lunout,*) write(lunout,*) write(lunout,'(a,132a)') ("*",i=1,132) write(lunout,*) END IF * | * +-------------------------------------------------------------------* * Sample the energy group XI = FLRNDM(DUMMY) DO 500 K = 1, 72 IF(XI .LE. CFCUM(K)) THEN NEUGRP = K ENERGY = ENEDGE(K) - & (XI-CFCUM(K-1))*(ENEDGE(K)-ENEDGE(K+1))/(CFCUM(K)-CFCUM(K-1)) GO TO 600 END IF 500 CONTINUE STOP ' Failed to sample the energy group' 600 CONTINUE kount=kount+1 * * Sample the cosine of the polar angle XI = FLRNDM(XI) DO 300 K = 1, 40 IF(XI .LE. CUMANG(K)) THEN COSTHE = ANG(K-1) + & (XI - CUMANG(K-1))*(ANG(K)-ANG(K-1))/(CUMANG(K)-CUMANG(K-1)) GO TO 400 END IF 300 CONTINUE STOP ' Failed to sample the polar angle cosine' 400 CONTINUE SINTHE = SQRT(ONEONE - COSTHE**2) IF(FLRNDM(COSTHE) .LE. HLFHLF) SINTHE = -SINTHE * Samples the azimuthal angle PHI = TWOPIP * FLRNDM(COSTHE) COSPHI = COS(PHI) * Npflka is the stack counter: of course any time source is called it * must be =0 NPFLKA = NPFLKA + 1 * Wt is the weight of the particle WTFLK (NPFLKA) = ONEONE WEIPRI = WEIPRI + WTFLK (NPFLKA) * Particle type (1=proton.....). Ijbeam is the type set by the BEAM * card * +-------------------------------------------------------------------* * | Heavy ion: IF ( IJBEAM .EQ. -2 ) THEN IJHION = IPROZ * 1000 + IPROA IJHION = IJHION * 100 + KXHEAV IONID = IJHION CALL DCDION ( IONID ) CALL SETION ( IONID ) ILOFLK (NPFLKA) = IJHION * | * +-------------------------------------------------------------------* * | Normal hadron: ELSE IONID = IJBEAM ILOFLK (NPFLKA) = IJBEAM END IF * | * +-------------------------------------------------------------------* * From this point ..... * Particle generation (1 for primaries) LOFLK (NPFLKA) = 1 * User dependent flag: LOUSE (NPFLKA) = 0 * User dependent spare variables: DO 100 ISPR = 1, MKBMX1 SPAREK (ISPR,NPFLKA) = ZERZER 100 CONTINUE * User dependent spare flags: DO 200 ISPR = 1, MKBMX2 ISPARK (ISPR,NPFLKA) = 0 200 CONTINUE * Save the track number of the stack particle: ISPARK (MKBMX2,NPFLKA) = NPFLKA NPARMA = NPARMA + 1 NUMPAR (NPFLKA) = NPARMA NEVENT (NPFLKA) = 0 DFNEAR (NPFLKA) = +ZERZER * ... to this point: don't change anything * Particle age (s) AGESTK (NPFLKA) = +ZERZER AKNSHR (NPFLKA) = -TWOTWO * Group number for "low" energy neutrons IGROUP (NPFLKA) = NEUGRP * Kinetic energy of the particle (GeV) TKEFLK (NPFLKA) = ENERGY * Particle momentum PMOFLK (NPFLKA) = SQRT ( TKEFLK (NPFLKA) * ( TKEFLK (NPFLKA) * + TWOTWO * AM (ILOFLK(NPFLKA)) ) ) * Cosines (tx,ty,tz) (make sure they are properly normalized) TZFLK (NPFLKA) = COSTHE TYFLK (NPFLKA) = SINTHE * COSPHI TXFLK (NPFLKA) = SQRT(ONEONE-TZFLK(NPFLKA)**2-TYFLK(NPFLKA)**2) IF(FLRNDM(COSPHI) .LE. HLFHLF) TXFLK (NPFLKA) = -TXFLK (NPFLKA) * Polarization cosines: TXPOL (NPFLKA) = -TWOTWO TYPOL (NPFLKA) = +ZERZER TZPOL (NPFLKA) = +ZERZER * Particle coordinates XFLK (NPFLKA) = XBEAM YFLK (NPFLKA) = YBEAM ZFLK (NPFLKA) = ZBEAM * Calculate the total kinetic energy of the primaries: don't change IF ( ILOFLK (NPFLKA) .EQ. -2 .OR. ILOFLK (NPFLKA) .GT. 100000 ) & THEN TKESUM = TKESUM + TKEFLK (NPFLKA) * WTFLK (NPFLKA) ELSE IF ( ILOFLK (NPFLKA) .NE. 0 ) THEN TKESUM = TKESUM + ( TKEFLK (NPFLKA) + AMDISC (ILOFLK(NPFLKA)) ) & * WTFLK (NPFLKA) ELSE TKESUM = TKESUM + TKEFLK (NPFLKA) * WTFLK (NPFLKA) END IF * Flag this is prompt radiation LRADDC (NPFLKA) = .FALSE. RADDLY (NPFLKA) = ZERZER * Here we ask for the region number of the hitting point. * NREG (NPFLKA) = ... * The following line makes the starting region search much more * robust if particles are starting very close to a boundary: CALL GEOCRS ( TXFLK (NPFLKA), TYFLK (NPFLKA), TZFLK (NPFLKA) ) CALL GEOREG ( XFLK (NPFLKA), YFLK (NPFLKA), ZFLK (NPFLKA), & NRGFLK(NPFLKA), IDISC ) * Do not change these cards: CALL GEOHSM ( NHSPNT (NPFLKA), 1, -11, MLATTC ) NLATTC (NPFLKA) = MLATTC CMPATH (NPFLKA) = ZERZER CALL SOEVSV RETURN *=== End of subroutine Source =========================================* END