FLUKA: Question
Dear FLUKA users,
I have got a problem to run an input file.
I don't know why it doesn't run.
I send you my input file and the mgdraw.f file I use for this
calculation
Is there anyone who could tell me what's wrong ?
Thanks for answering
It's quite urgent
Yann FOUCHER
Paul Scherrer Institut
TITLE
Maelle's experiment
*2345678901234567890123456789012345678901234567890123456789012345678901234567890
DEFAULTS CALORIME
* |_________|_________|_________|_________|_________|_________|
BEAM -0.0627 0.0 0.0 +2.0 0.0 -1. NEUTRON
* |_________|_________|_________|_________|_________|_________|
BEAMPOS 0.0 0.0 -1.5e-2
GEOBEGIN COMBINAT
TEST 1
* Bodies
* |_________|_________|_________|_________|_________|_________|
SPH 1 +0.0 +0.0 +0.0 +30.0
* |_________|_________|_________|_________|_________|_________|
SPH 2 +0.0 +0.0 +0.0 +20.0
RPP 3 -3.0 +3.0 -3.0 +3.0 -1.5e-2 +1.5e-2
END
* Regions
* |_________|_________|_________|_________|_________|_________|
*_AAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIII
001 +1 -2
002 +2 -3
003 +3
*_AAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIIIAAIIIII
END
* |_________|_________|_________|_________|_________|_________|
GEOEND
* |_________|_________|_________|_________|_________|_________|
*
MATERIAL 82. 207.2 11.35 26.0 LEAD
* |_________|_________|_________|_________|_________|_________|
ASSIGNMAT 1.0 1.0
ASSIGNMAT 2.0 2.0
ASSIGNMAT 26.0 3.0
* |_________|_________|_________|_________|_________|_________|
LOW-NEUT 72. 22.0 0.0196 1.0 1.0
* |_________|_________|_________|_________|_________|_________|
LOW-MAT 26. 82. -2. 293.
* |_________|_________|_________|_________|_________|_________|
EVENTYPE 10. EVAP
EXTRAWEI 1. 1.0
SCORE 208. 201.
* |_________|_________|_________|_________|_________|_________|
DUMPTHEM 100. 30. 6.0 1.0 dump
RANDOMIZ 1.0
START 5e3 9999999. 80. 1.0 0.0 9999999.
STOP
*$ CREATE MGDRAW.FOR
*COPY MGDRAW
* *
*=== mgdraw ===========================================================*
* *
SUBROUTINE MGDRAW ( ICODE, MREG )
INCLUDE '(DBLPRC)'
INCLUDE '(DIMPAR)'
INCLUDE '(IOUNIT)'
*
*----------------------------------------------------------------------*
* *
* MaGnetic field trajectory DRAWing: actually this entry manages *
* all trajectory dumping for *
* drawing *
* *
* Created by Alfredo Ferrari *
* INFN - Milan *
* last change 19-oct-98 by Alfredo Ferrari *
* INFN - Milan *
* *
*----------------------------------------------------------------------*
*
INCLUDE '(CASLIM)'
INCLUDE '(COMPUT)'
INCLUDE '(EPISOR)'
INCLUDE '(FINUC)'
*----------------------------------------------------------------------*
* *
* Include file: Finuc (new version of old Finuc of FLUKA86) *
* *
* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! *
* !!!! S E E A L S O I N C L U D E F I L E !!!! *
* !!!! F I N U C 2 !!!! *
* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! *
* *
* Created on 20 january 1996 by Alfredo Ferrari & Paola Sala *
* Infn - Milan *
* *
* Last change on 26-jul-97 by Alfredo Ferrari *
* *
* included in the following subroutines or functions: not updated *
* *
* description of the common block(s) and variable(s) *
* *
* /finuc/ is the storage for secondaries created in event *
* np = number of secondaries *
* kpart (ip) = type of the secondary ip *
* cxr (ip) = direction cosine of the secondary ip *
* with respect to x-axis *
* cyr (ip) = direction cosine of the secondary ip *
* with respect to y-axis *
* czr (ip) = direction cosine of the secondary ip *
* with respect to z-axis *
* cxrpol (ip) = direction cosine of the secondary ip polarization *
* with respect to x-axis *
* cyrpol (ip) = direction cosine of the secondary ip polarization *
* with respect to y-axis *
* czrpol (ip) = direction cosine of the secondary ip polarization *
* with respect to z-axis *
* tki (ip) = kinetic energy of secondary ip *
* plr (ip) = momentum of the secondary ip *
* wei (ip) = weight of the secondary ip *
* agesec (ip) = "age" of the secondary ip with respect to the *
* interaction time *
* tv = excitation energy *
* tvcms = actual excitation energy of the residual nucleus *
* tvrecl = recoil kinetic energy of the residual nucleus *
* tvheav = recoil kinetic energies of heavy (2-H, 3-H, 3-He, *
* 4-He) fragments after evaporation *
* tvbind = approximate energy wasted in nuclear binding *
* effects (not yet operational) *
* *
*----------------------------------------------------------------------*
INCLUDE '(MGDDCM)'
INCLUDE '(PAPROP)'
INCLUDE '(STACK)'
INCLUDE '(STARS)'
INCLUDE '(TRACKR)'
*
INCLUDE '(RESNUC)'
*---------------------------------------------------------------------*
* *
* Include file Resnuc *
* *
* Created on 20 april 1990 by Alfredo Ferrari *
* INFN Milan *
* *
* Last change on 09-nov-99 by Alfredo Ferrari, INFN - Milan *
* *
* *
* Description of variables (incomplete): *
* *
* Icres = residual nucleus atomic number *
* Ibres = residual nucleus mass number *
* Istres = residual nucleus stable level index *
* Ismres = residual nucleus isomeric state index *
* Ihyres = residual nucleus hyperon number *
* Amnres = residual nucleus nuclear mass *
* Ammres = residual nucleus atomic mass *
* Eres = residual nucleus total energy *
* Ekres = residual nucleus kinetic energy *
* Px,y,zres = residual nucleus momentum components *
* Ptres2 = residual nucleus squared momentum *
* Angres = residual nucleus angular momentum (GeV/c fm) *
* Anx,y,zres = residual nucleus angular momentum components *
* Khyres(jp) = id of the jp_th hyperon in the residual nucleus *
* Bhyres(jp) = (nuclear) binding energy of the jp_th hyperon *
* in the residual nucleus *
* *
* *
* *
*---------------------------------------------------------------------*
INCLUDE '(FHEAVY)'
*----------------------------------------------------------------------*
* *
* include file: fheavy *
* *
* Created on 5 april 1990 by Alfredo Ferrari, INFN Milan *
* *
* Last change on 26-jul-97 by Alfredo Ferrari, INFN Milan *
* *
* included in the following subroutines or functions: not updated *
* *
* description of the common block(s) and variable(s) *
* *
* /fheavy/ is the storage for heavy secondaries created in the *
* nuclear evaporation *
* npheav = number of secondaries *
* kheavy(ip) = type of the secondary ip *
* ( 3 = deuteron, 4 = 3-H, 5 = 3-He, 6 = 4-He, *
* 7-12 = "Heavy" fragment specified by Ibheav and *
* Icheav ) *
* cxheav(ip) = direction cosine of the secondary ip *
* with respect to x-axis *
* cyheav(ip) = direction cosine of the secondary ip *
* with respect to y-axis *
* czheav(ip) = direction cosine of the secondary ip *
* with respect to z-axis *
* tkheav(ip) = kinetic energy of secondary ip *
* pheavy(ip) = momentum of the secondary ip *
* wheavy(ip) = weight of the secondary ip *
* agheav(ip) = "age" of the secondary ip with respect to the *
* interaction time *
* amheav(kp) = atomic masses of the twelve types of evaporated *
* or fragmented or fissioned particles *
* amnhea(kp) = nuclear masses of the twelve types of evaporated *
* or fragmented or fissioned particles *
* bhheav(jp,kp) = (nuclear) binding energy of the jp_th hyperon of *
* the kp-type heavy particle *
* anheav(kp) = name of the kp-type heavy particle *
* icheav(kp) = charge of the kp-type heavy particle *
* ibheav(kp) = mass number of the kp-type heavy particle *
* imheav(kp) = isomeric state of the kp-type heavy particle *
* ihheav(kp) = number of hyperons of the kp-type heavy particle *
* khheav(jp,kp) = id of the jp_th hyperon of the kp-type heavy *
* particle *
* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! *
* !!! There is now the possibility to produce up to 6 "heavy" !!!! *
* !!! fragments besides the residual nucleus recorded in !!!! *
* !!! Resnuc: they are identified by indeces 7-12, of course !!!! *
* !!! the corresponding physical properties (Z,A,m..) must be !!!! *
* !!! updated every time they are produced !!!! *
* !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! *
*----------------------------------------------------------------------*
*
CHARACTER*20 FILNAM
LOGICAL LFCOPE,LFIRST
SAVE LFCOPE,LFIRST
DATA LFCOPE,LFIRST / .FALSE.,.TRUE. /
*
*----------------------------------------------------------------------*
* *
* Icode = 1: call from Kaskad *
* Icode = 2: call from Emfsco *
* Icode = 3: call from Kasneu *
* Icode = 4: call from Kashea *
* Icode = 5: call from Kasoph *
* *
*----------------------------------------------------------------------*
* *
IF ( .NOT. LFCOPE ) THEN
LFCOPE = .TRUE.
IF ( KOMPUT .EQ. 2 ) THEN
FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
ELSE
FILNAM = CFDRAW
END IF
OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
& 'UNFORMATTED' )
END IF
C WRITE (IODRAW) NTRACK, MTRACK, JTRACK, SNGL (ETRACK),
C & SNGL (WTRACK)
C WRITE (IODRAW) ( SNGL (XTRACK (I)), SNGL (YTRACK (I)),
C & SNGL (ZTRACK (I)), I = 0, NTRACK ),
C & ( SNGL (DTRACK (I)), I = 1,MTRACK ),
C & SNGL (CTRACK)
RETURN
*
*----------------------------------------------------------------------*
* *
* ENergy deposition DRAWing: *
* *
* Icode = 1x: call from Kaskad *
* 10: elastic interaction recoil *
* 11: inelastic interaction recoil *
* 12: stopping particle *
* 13: pseudo-neutron deposition *
* 14: escape *
* Icode = 2x: call from Emfsco *
* 20: local energy deposition (i.e. photoelectric) *
* 21: below threshold, iarg=1 *
* 22: below threshold, iarg=2 *
* 23: escape *
* Icode = 3x: call from Kasneu *
* 30: target recoil *
* 31: below threshold *
* 32: escape *
* Icode = 4x: call from Kashea *
* 40: escape *
* Icode = 5x: call from Kasoph *
* 50: optical photon absorption *
* 51: escape *
* *
*----------------------------------------------------------------------*
* *
ENTRY ENDRAW ( ICODE, MREG, RULL, XSCO, YSCO, ZSCO )
IF ( .NOT. LFCOPE ) THEN
LFCOPE = .TRUE.
IF ( KOMPUT .EQ. 2 ) THEN
FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
ELSE
FILNAM = CFDRAW
END IF
OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
& 'UNFORMATTED' )
END IF
C WRITE (IODRAW) 0, ICODE, JTRACK, SNGL (ETRACK), SNGL (WTRACK)
C WRITE (IODRAW) SNGL (XSCO), SNGL (YSCO), SNGL (ZSCO), SNGL (RULL)
RETURN
*
*======================================================================*
* *
* SOurce particle DRAWing: *
* *
*======================================================================*
*
ENTRY SODRAW
IF ( .NOT. LFCOPE ) THEN
LFCOPE = .TRUE.
IF ( KOMPUT .EQ. 2 ) THEN
FILNAM = '/'//CFDRAW(1:8)//' DUMP A'
ELSE
FILNAM = CFDRAW
END IF
OPEN ( UNIT = IODRAW, FILE = FILNAM, STATUS = 'NEW', FORM =
& 'UNFORMATTED' )
END IF
C WRITE (IODRAW) -NCASE, LSTACK, LSTMAX, SNGL (TKESUM),
C & SNGL (WEIPRI)
C WRITE (IODRAW) ( ILO (I), SNGL (TKE(I)+AM(ILO(I))), SNGL (WT(I)),
C & SNGL (XA(I)), SNGL (YA(I)), SNGL (ZA(I)),
C & SNGL (TX(I)), SNGL (TY(I)), SNGL (TZ(I)),
C & I = 1, LSTACK )
RETURN
*
*======================================================================*
* *
* USer dependent DRAWing: *
* *
* Icode = 10x: call from Kaskad *
* 100: elastic interaction secondaries *
* 101: inelastic interaction secondaries *
* 102: particle decay secondaries *
* 103: delta ray generation secondaries *
* 104: pair production secondaries *
* 105: bremsstrahlung secondaries *
* Icode = 20x: call from Emfsco *
* 208: bremsstrahlung secondaries *
* 210: Moller secondaries *
* 212: Bhabha secondaries *
* 214: in-flight annihilation secondaries *
* 215: annihilation at rest secondaries *
* 217: pair production secondaries *
* 219: Compton scattering secondaries *
* 221: photoelectric secondaries *
* 225: Rayleigh scattering secondaries *
* Icode = 30x: call from Kasneu *
* 300: interaction secondaries *
* Icode = 40x: call from Kashea *
* 400: delta ray generation secondaries *
* For all interactions secondaries are put on FINUC common (kp=1,np) *
* but for KASHEA delta ray generation where only the secondary elec- *
* tron is present and stacked on STACK common for kp=lstack *
* *
*======================================================================*
*
ENTRY USDRAW ( ICODE, MREG, XSCO, YSCO, ZSCO )
*
IF (LFIRST) THEN
OPEN(50,FILE='recoil.dbg',STATUS='UNKNOWN')
OPEN(51,FILE='recoil.out',STATUS='UNKNOWN')
OPEN(52,FILE='neutron.out',STATUS='UNKNOWN')
OPEN(53,FILE='heavy.out',STATUS='UNKNOWN')
OPEN(54,FILE='outpart',STATUS='UNKNOWN')
LFIRST = .FALSE.
ENDIF
*
* zero the momentum-balance accumulators
PXBAL = ZERZER
PYBAL = ZERZER
PZBAL = ZERZER
*
* select elastic and inelastic interactions
IF (ICODE.EQ.101) THEN
c WRITE(50,1000) INT(WEIPRI),ICODE,XSCO,YSCO,ZSCO
c 1000 FORMAT('B',I10,2X,I5,3F8.4)
* outgoing particles
DO 10 IP=1,NP
PXBAL = PXBAL-CXR(IP)*PLR(IP)
PYBAL = PYBAL-CYR(IP)*PLR(IP)
PZBAL = PZBAL-CZR(IP)*PLR(IP)
IF (IP.EQ.1.and.KPART(IP).eq.1) THEN
WRITE(54,1003)XSCO,YSCO,ZSCO,INT(WEIPRI),KPART(IP)
& ,WEI(IP),TKI(IP),CXR(IP),CYR(IP),CZR(IP)
1003 FORMAT('O',3F8.4,I10,I4,2E12.4,2X,3F8.4)
ENDIF
IF (IP.NE.1.and.KPART(IP).eq.1) THEN
WRITE(54,1003)XSCO,YSCO,ZSCO,INT(WEIPRI),KPART(IP)
& ,WEI(IP),TKI(IP),CXR(IP),CYR(IP),CZR(IP)
ENDIF
10 CONTINUE
* outgoing heavy particles
c DO 11 IP=1,NPHEAV
c PXBAL = PXBAL-CXHEAV(IP)*PHEAVY(IP)
c PYBAL = PYBAL-CYHEAV(IP)*PHEAVY(IP)
c PZBAL = PZBAL-CZHEAV(IP)*PHEAVY(IP)
c IF (IP.EQ.1) THEN
c WRITE(50,1005) KHEAVY(IP),WHEAVY(IP),TKHEAV(IP),
c & CXHEAV(IP)*PHEAVY(IP),CYHEAV(IP)*PHEAVY(IP),
c & CZHEAV(IP)*PHEAVY(IP)
c 1005 FORMAT(1X,'Outgoing heavies ',I4,2E12.4,2X,3F8.4)
c ELSE
c WRITE(50,1004) KHEAVY(IP),WHEAVY(IP),TKHEAV(IP),
c & CXHEAV(IP)*PHEAVY(IP),CYHEAV(IP)*PHEAVY(IP),
c & CZHEAV(IP)*PHEAVY(IP)
c ENDIF
* dump heavies
c IF (ICODE.LE.101) THEN
c WRITE(53,1012)
c & INT(WEIPRI),XSCO,YSCO,ZSCO,WHEAVY(IP),TKHEAV(IP),
c & CXHEAV(IP)*PHEAVY(IP),CYHEAV(IP)*PHEAVY(IP),
c & CZHEAV(IP)*PHEAVY(IP),ICHEAV(KHEAVY(IP)),
c & IBHEAV(KHEAVY(IP))
c 1012 FORMAT(I6,3F9.4,1P,5E12.4,2I4)
c ENDIF
c 11 CONTINUE
*
ENDIF
*
RETURN
*=== End of subrutine Mgdraw ==========================================*
END
Partial thread listing: