Re: FLUKA: Production of Hydrogen
Hi Laurent,
please find attached:
- the include file for the TRACKR COMMON, with comments
- Chap. 22 of the manual, with instructions on how to write a
user routine.
- the include file for the FINUC COMMON, with comments
- a description of option DUMPTHEM, from the manual, describing how
to activate calls to various entries of the MGDRAW routine.
Look at the instructions for routine MGDRAW (and its entries ENDRAW,
USDRAW)
Stopping hadrons and muons are accessible trough calls to ENTRY ENDRAW
with ICODE=12. Protons are identified by JTRACK=1.
Secondary particles from a low-energy neutron interactions are
accessible in COMMON FINUC after a call to ENTRY USDRAW with ICODE=300.
Best regards,
Alberto
On Wed, 12 Sep 2001, Laurent Aphecetche wrote:
> Alberto Fasso' wrote:
> >
> > ...
> >
> > Anyway, one way to attack your problem could be to score
> > separately "stopping protons" (protons reaching energy cutoff). This is
> > possible using the MGDRAW user routine. Unfortunately it will not be a
> > complete solution, because you will still be missing protons produced
> > in low energy neutron reactions, which FLUKA does not generate
> > explicitely. (But it will catch the hyrogen recoils, although I am not
> > sure that it is a desirable feature for you)
> >
>
> Hi Alberto,
>
> I'm also interested of trapping those "stopping protons". Unfortunately
> I do have a FLUKA version where there's no comment in the include files,
> thus it's difficult to have an idea about which variable is relevant to
> this problem. Could you provide an example on how to store the particles
> reaching energy cut-off ?
>
> Thanks,
>
>
--
Alberto Fassò
CERN-EP/AIP, CH-1211 Geneve 23 (Switzerland)
Phone: (41 22) 767 2398 Fax: (41 22) 767 9480 Alberto.Fasso@cern.ch
*$ CREATE TRACKR.ADD
*COPY TRACKR
* *
*=== trackr ===========================================================*
* *
*----------------------------------------------------------------------*
* *
* TRACKs Recording by Alfredo Ferrari, INFN - Milan *
* *
* last change 31 january 2001 by Alfredo Ferrari *
* *
* included in : *
* electr *
* emfsco *
* kaskad (new version) *
* kashea *
* kasneu *
* geoden (new version) *
* mageas *
* magmov *
* magnew *
* move *
* photon *
* usrsco *
* *
* Ntrack = number of track segments *
* Mtrack = number of energy deposition events along the track *
* 0 < i < Ntrack *
* Xtrack = end x-point of the ith track segment *
* Ytrack = end y-point of the ith track segment *
* Ztrack = end z-point of the ith track segment *
* 1 < i < Ntrack *
* Ttrack = length of the ith track segment *
* 1 < j < Mtrack *
* Dtrack = energy deposition of the jth deposition event *
* *
* Jtrack = identity number of the particle *
* Etrack = total energy of the particle *
* Ptrack = momentum of the particle (not always defined, if *
* < 0 must be obtained from Etrack) *
* Cx,y,ztrck = direction cosines of the current particle *
* Cx,y,ztrpl = polarization cosines of the current particle *
* Wtrack = weight of the particle *
* Wscrng = scoring weight: it can differ from Wtrack if some *
* biasing techniques are used (for example inelastic *
* interaction length biasing) *
* Ctrack = total curved path *
* Zfftrk = <Z_eff> of the particle *
* Zfrttk = actual Z_eff of the particle *
* Atrack = age of the particle *
* Akshrt = Kshrt amplitude for K0/K0bar *
* Aklong = Klong amplitude for K0/K0bar *
* Wninou = neutron algebraic balance of interactions (both *
* for "high" energy particles and "low" energy *
* neutrons) *
* Spausr = user defined spare variables for the current *
* particle *
* Sttrck = macroscopic total cross section for low energy *
* neutron collisions *
* Satrck = macroscopic absorption cross section for low energy*
* neutron collisions (it can be negative for Pnab>1) *
* Ktrack = if > 0 neutron group of the particle (neutron) *
* *
* Ntrack > 0, Mtrack > 0 : energy loss distributed along the *
* track *
* Ntrack > 0, Mtrack = 0 : no energy loss along the track *
* Ntrack = 0, Mtrack = 0 : local energy deposition (the *
* value and the point are not re- *
* corded in Trackr) *
* Mmtrck = flag recording the material index for low energy *
* neutron collisions *
* Lt1trk = initial lattice cell of the current track *
* (or lattice cell for a point energy deposition) *
* Lt2trk = final lattice cell of the current track *
* Ihspnt = current geometry history pointer (not set if -1) *
* Ltrack = flag recording the generation number *
* Llouse = user defined flag for the current particle *
* Ispusr = user defined spare flags for the current particle *
* Lfsssc = logical flag for inelastic interactions ending with*
* fission (used also for low energy neutrons) *
* *
*----------------------------------------------------------------------*
* *
PARAMETER ( MXTRCK = 2500 )
LOGICAL LFSSSC
COMMON / TRACKR / XTRACK ( 0:MXTRCK ), YTRACK ( 0:MXTRCK ),
& ZTRACK ( 0:MXTRCK ), TTRACK ( MXTRCK ),
& DTRACK ( MXTRCK ), ETRACK, PTRACK, CXTRCK,
& CYTRCK, CZTRCK, WTRACK, CXTRPL, CYTRPL, CZTRPL,
& ZFFTRK, ZFRTTK, ATRACK, CTRACK, AKSHRT, AKLONG,
& WSCRNG, WNINOU, SPAUSR(MKBMX1), STTRCK, SATRCK,
& NTRACK, MTRACK, JTRACK, KTRACK, MMTRCK, LT1TRK,
& LT2TRK, IHSPNT, LTRACK, LLOUSE, ISPUSR(MKBMX2),
& LFSSSC
EQUIVALENCE ( SPAUSE, SPAUSR (1) )
EQUIVALENCE ( ISPUSE, ISPUSR (1) )
1***********************************************************************
22} How to write, compile and link a user routine
***********************************************************************
Unlike some other Monte Carlo particle transport codes, FLUKA gets its input
mainly from a simple file. It offers a rich choice of options for scoring
most
quantities of possible interest and for applying different variance
reduction
techniques, without requiring the user to write a single line of code.
However, although normally there is no need for any "user code", there are
special cases where this is unavoidable, either because of the complexity of
the problem, or because the desired information is too unusual or too
problem-specific to be offered as a standard option.
And on the other hand, even when this is not strictly necessary, experienced
programmers may like to create customized input/output interfaces.
A number of user routines (available on LINUX and UNIX platforms in
directory
usermvax) allow to define non-standard input and output, and in some cases
even to modify to a limited extent the normal particle transport. Most of
them
are already present in the FLUKA library as dummy or template routines, and
require a special command in the standard input file to be activated. The
user
can modify any one of these routines, and even insert into them further
calls
to his own private ones, or to external packages (at his own risk!). This
increased flexibility must be balanced against the advantage of using as far
as
possible the FLUKA standard facilities, which are known to be reliable and
well
tested.
To implement his own code, the user must perform the following steps:
1) make a modified copy of one or more of these routines. It is recommended
that each modified routine should always print an informative message
when
called for the first time, to confirm that it has been successfully
activated, for future documentation, and to avoid misinterpretations of
the
standard output. It is important to remember that when calling modified
user routines, the units, titles etc. reported in the normal FLUKA output
become often meaningless.
A typical way to do this is:
..............
LOGICAL LFIRST
SAVE LFIRST
DATA LFIRST /.TRUE./
* return message from first call
IF (LFIRST) THEN
WRITE(LUNOUT,*) 'Version xxx of Routine yyy called'
LFIRST = .FALSE.
ENDIF
..............
IMPORTANT: The user should not modify the value of any argument in a
routine
calling list, except when marked as "returned" in the description of the
routine here below.
Similarly, no variable contained in COMMON blocks should be overwritten
unless explicitely indicated.
2) compile the modified routines (with the fff script on LINUX/UNIX):
fff yyy.f (produces a new file yyy.o)
3) link them (with the lfluka script on LINUX/UNIX) to the FLUKA
library and any additional library of interest (for instance
CERNLIB):
lfluka -o myfluka -m fluka yyy.o
This will produce a new executable (indicated here as myfluka).
To run the new executable, launch the usual rfluka script with the option
-e myfluka.
INCLUDE files
-------------
It is recommended that at least the following lines be present at the
beginning
of each routine:
INCLUDE '(DBLPRC)'
INCLUDE '(DIMPAR)'
INCLUDE '(IOUNIT)'
Each INCLUDE contains a COMMON block, plus related constants.
Additional INCLUDEs may be useful, in particular BEAM, CASLIM, EMFSTK,
EPISOR,
EVTFLG, FHEAVY, FINUC, LTCLCM, MAPA, RESNUCLEi, SCOHLP, SOUEVT, STACK,
STARS,
TRACKR, USRBIN, USRBDX, USRTRC, USRYLD.
Files flukaadd,add and emfadd.add contain a full documentation about the
meaning of the variables of these INCLUDE files.
Here is a summary of their content:
DBLPRC: included in ALL routines of FLUKA, contains the declaration
IMPLICIT DOUBLE PRECISION (A-H,O-Z)
and sets many mathematical and physical constants. Users are
strongly
encouraged to adhere to "FLUKA style" by using systematically double
precision (except for very good reasons such as calling external
single
precision scoring packages), and to use constants defined in this
file
for maximum accuracy.
DIMPAR: dimensions of the most important arrays
IOUNIT: logical input and output unit numbers
BEAM: properties of primary particles as defined by options BEAM
and BEAMPOS
CASLIM: number of primary particles followed (needed for normalization)
EMFSTK: particle stack for electrons and photons
EPISOR: user variables and information for a user-written source
EVTFLG: event flags (still undocumented)
FHEAVY: stack of heavy secondaries created in nuclear evaporation
FINUC: properties of each secondary created in a hadronic event
LTCLCM: LaTtice CeLl CoMmon (needed when writing symmetry transformations
for
Lattice Geometry)
MAPA: material properties
RESNUC: properties of the current residual nucleus
SCOHLP: SCOring HeLP (information on current estimator type). It contains a
flag ISCRNG, indicating the quantity being scored by the current
binning or estimator, and one JSCRNG corresponding to the
binning/estimator number. Binnings and estimators are sequentially
numbered according to their order of appearance in standard input.
Note that several estimators can have the same JSCRNG number (for
instance Binning N. 3 and Tracklength estimator N. 3). They can be
distinguished based on the value of ISCRNG. However, note that the
same value of ISCRNG may have different meanings in functions FLUSCW
and
COMSCW.
SOUEVT: SOUrce EVenT (useful when source particles are obtained from an
external event generator)
STACK: main FLUKA particle stack
STARS: total numbers and total weights relative to many physical and Monte
Carlo events (needed for normalization, energy balance etc.)
TRACKR: TRACKs Recording (properties of the currently transported particle
and
its path)
USRBIN, USRBDX, USRSNC, USRTRC, USRYLD: parameters of the requested
estimators
User routines
-------------
===============================================================================
--------
comscw.f: weighting deposited energy or stars
--------
Argument list (all variables are input only):
IJ : particle type (1 = proton, 8 = neutron, etc.: see code in
4})
XA,YA,ZA : current particle position
MREG : current geometry region
RULL : amount to be deposited (unweighted)
LLO : particle generation
ICALL : internal code calling flag (not for general use)
This routine is activated by option EXTRAWEIgh, with WHAT(6) > 0.0.
Energy and star densities obtained via SCORE and USRBIN, and energy
and stars obtained via EVENTBIN are multiplied by the value
returned
by this function. The user can implement any desired logic to
differentiate the returned value according to any information
contained in the argument list (particle type, position, region,
amount deposited, particle generation), or information available in
the SCOHLP COMMON block (binning number, type of scored quantity).
The scored quantity is given by the flag ISCRNG (in SCOHLP):
ISCRNG = 1 --> Energy density binning
ISCRNG = 2 --> Star density binning
ISCRNG = 3 --> Residual nuclei scoring
The binning/estimator number is given by JSCRNG (in SCOHLP) and is
printed in output:
Res. nuclei n. 3 "any-name" , "high" energy products, region n.
4
R-Phi-Z binning n. 5 "other-name" , generalized particle n.
1
. Note that an estimator of residual nuclei can have the same JSCRNG
number as a binning (use the value of ISCRNG to discriminate).
Further information can be obtained including COMMON TRACKR (for
instance particle's total energy, direction cosines, age). TRACKR
contains also special user variables (both integer and in double
precision) which can be used to save information about particles
which
have undergone some particular event. If data concerning the
current
material are needed, it can be accessed as MEDIUM(MREG) if file
MAPA
is included. Indeed, a common simple application of COMSCW is to
score
dose according to the local density (especially useful to get the
correct average dose in bins straddling a boundary between two
different media):
..................
INCLUDE '(MAPA)'
INCLUDE '(SCOHLP)'
..................
* ========== In order to compute doses ========= *
* Medium(n) is the material number of region n
* Rho(m) is the density of material m (in g/cm3)
* Iscrng = 1 means we are depositing energy (not stars)
IF ( ISCRNG .EQ. 1 ) THEN
* to get dose in Gy (elcmks is the electron charge in C)
COMSCW = ELCMKS * 1.D12 / RHO (MEDIUM(MREG))
ELSE
* oneone is defined as 1.D0 in include DBLPRC
COMSCW = ONEONE
ENDIF
..................
Note: setting the variable LSCZER = .TRUE. before RETURN (LSCZER is
in
COMMON SCOHLP), will cause zero scoring whatever the value returned
by COMSCW. This is more efficient than returning a zero value.
===============================================================================
--------
fluscw.f: weighting fluence, current and yield
--------
Argument list (all variables are input only):
IJ : particle type
PLA: particle momentum (if > 0), or -PLA = kinetic energy
(if <0)
TXX,TYY,TZZ: particle current direction cosines
WEE: particle weight
XX,YY,ZZ: particle position
NREG: current region (after boundary crossing)
IOLREG: previous region (before boundary crossing). Useful only
whith boundary crossing estimators; for other estimators it
has no meaming.
LLO: particle generation
ICALL: internal code calling flag (not for general use)
Function FLUSCW is activated by option EXTRAWEIgh, with WHAT(3) >
0.0.
Yields obtained via USRYIELD, fluences calculated with USRBDX,
USRTRACK, USRCOLL, USRBIN, and currents calculated with USRBDX are
multiplied by the value returned by this function. The user can
implement any desired logic to differentiate the returned value
according to any information contained in the argument list
(particle
type, energy, direction, weight, position, region, boundary,
particle
generation), or information available in the SCOHLP COMMON block
(binning or estimator number, estimator type).
The scored quantity is given by the flag ISCRNG (in SCOHLP):
ISCRNG = 1 --> Boundary crossing estimator
ISCRNG = 2 --> Track length binning
ISCRNG = 3 --> Track length estimator
ISCRNG = 4 --> Collision density estimator
ISCRNG = 5 --> Yield estimator
The binning/estimator number is given by JSCRNG (in SCOHLP) and is
printed in output:
Bdrx n. 2 "bdxname" , generalized particle n. 8,
from region n. 22 to region n.
78
Track n. 6 "trkname" , generalized particle n. 14, region n. 9
. Note that a track-length estimator can have the same JSCRNG number
as
a boundary crossing one or a binning etc. (use the value of ISCRNG
to
discriminate).
Further information can be obtained including COMMON TRACKR (for
instance particle age). TRACKR contains also special user variables
(both integer and in double precision) which can be used to save
information about particles which have undergone some particular
event.
Function FLUSCW has many applications. A common one in shielding
calculations is to multiply selected scored fluences by
particle/energy-dependent fluence-to-dose equivalent conversion
coefficients, or by some instrument response, radiation damage
curve,
etc. Another application is conditional scoring (score only if
within a certain distance from a point, etc.): for instance it is
possible to implement a sort of 2-dimensional fluence binning on
a plane boundary.
Other interesting applications are based on the fact that
FLUSCW is called at every boundary crossing, provided that at least
one USRBDX estimator has been requested. Although the function has
been designed mainly to weight scored quantities, it can be
"cheated"
to do all sorts of side things, even not directly connected with
scoring.
Note: setting the variable LSCZER (in SCOHLP) = .TRUE. before
RETURN, will cause zero scoring whatever the value returned
by FLUSCW. This is more efficient than returning a zero value.
===============================================================================
--------
frghns.f material roughness (for optical photons)
--------
Argument list (all variables are input only):
TXX, TYY, TZZ: particle direction cosines
UXSRFC, UYSRFC, UZSRFC: direction of the normal to the surface
MREG: region the particle is coming from
NEWREG: region the particle is going to
MMAT: material the particle is coming from
MMATNW: material the particle is going to
Function FRGHNS can be used to return a non-zero value for the
roughness of a boundary between two materials, relevant for
optical photon transport (default roughness is zero for all
boundaries). Meaningful only if options OPT-PROP or OPT-PROD have
been requested.
===============================================================================
--------
musrbr.f 1st user defined quantity for special binning
lusrbl.f 2nd user defined quantity for special binning
fusrbv.f 3rd user defined quantity for special binning
--------
These three functions are used to define 3-dimensional
distributions
to be calculated by special user-defined binnings (see USRBIN with
WHAT(1) = 8.0 in the first card).
MUSRBR defines a discrete (integer) variable (by default: region
number).
Argument list (all variables are input only):
IJ : particle type
PCONTR: particle momentum
XA,YA,ZA: particle position
MREG: current region
LATCLL: current lattice cell
ICALL: internal code calling flag (not for general use)
LUSRBL defines another discrete (integer) variable (by default:
lattice number)
Argument list: same as for MUSRBR above
FUSRBV defines a continuous (double precision) variable (by
default:
pseudorapidity with respect to the Z axis)
Argument list: same as for MUSRBR above (except LATCLL)
The 3 functions are called at tracklength events. What is scored is
the particle track-length multiplied by the particle's weight,
possibly modified by a user-written FLUSCW (see above).
===============================================================================
--------
geoden.f FLUKA general scoring routine
--------
Subroutine GEODEN is the highest level scoring routine in FLUKA,
which is called at all events where something happens which could
contribute to an estimator response (energy deposition, end of a
step, crossing of a region boundary, elastic or inelastic
interaction).
It is not really a user routine, but it has been included in the
list
of user routines because it can be modified to replace standard
FLUKA
scoring with alternative histogramming packages.
Each type of scorable event corresponds to a different ENTRY with
its
own argument list, which acts as a "switch" to call one or more
corresponding specific lower-level subroutines. These calls can be
replaced by calls to other packages desired by the user. Please
remember that the scored amount (variable RULL), as all other FLUKA
variables, is in double precision and may need to be converted to
SNGL(RULL) if transmitted to single precision histogramming
packages.
Entries:
GEODEN (energy deposition events).
Estimators: SCORE, USRBIN, EVENTBIN, DETECT
Argument list (all variables are input only):
IJ : particle or generalized particle type
XA,YA,ZA: particle position
MREG: current region
RULL: amount of energy to be deposited (multiplied by
particle's
weight)
ICALL: internal code calling flag (not for general use)
LLO: particle generation
ECONTR: kinetic energy of contributing particle
ICNTR: flag = 1 for particle escaping from system, discarded, or
in general depositing energy in the blackhole; otherwise
= 0
GEODEN handles directly the energy deposition scoring requested
by
the SCORE command. Energy deposition binning (requested by
the USRBIN or EVENTBIN command) and event by event energy
deposition scoring (requested by the DETECT command) are
handled indirectly via calls to USRSCO and to DETECT,
respectively.
Note that energy weighting, via user routine COMSCW or
Birks
quenching, is done only inside USRSCO: if the latter is
bypassed by the user, no weighting or quenching can be
applied, even if EXTRAWEIgh or TCQUENCH options have been
requested.
GEOSTR (spallations and fission events).
Estimators: SCORE, USRBIN, EVENTDAT
Argument list (all variables are input only):
IJ, XA, YA, ZA, MREG: see GEODEN above
RULL: amount to be scored (weight of the particle producing a
star
or a fission; or neutron balance multiplied by the
particle's
weight)
ICALL, LLO, ECONTR: see GEODEN above
GEOSTR handles directly star/fission/neutron balance scoring
requested by the SCORE command.
Star/fission/neutron balance binning (requested by the
USRBIN or EVENTBIN command) is handled indirectly via a
call
to USRSCS.
Note that star weighting via user routine COMSCW is done
only inside USRSCS: if the latter is bypassed by the user,
no weighting can be applied, even if EXTRAWEIgh has been
requested.
GEOTRK (particle step event)
Estimators: USRTRACK, USRBIN
Argument list (all variables are input only):
IJ, XA, YA, ZA, MREG: see GEODEN above
RULL: amount to be scored (particle track length multiplied by
particle's weight)
ICALL, LLO, ECONTR: see GEODEN above
GEOSTR handles indirectly track length fluence (requested by the
USRBIN/EVENTBIN or USRTRACK commands) via calls to USRSTB
and USRSCT, respectively.
Note that fluence weighting via user routine FLUSCW is
done
only inside USRSTB and USRSCT: if the latter are bypassed
by
the user, no weighting can be applied, even if EXTRAWEIgh
has been requested.
GEOBDX (boundary crossing event)
Estimator: USRBDX
Argument list (all variables are input only):
IJ, XA, YA, ZA: see GEODEN above
TXX, TYY, TZZ: direction cosines of the particle
MREG: region before boundary crossing
NWREG: region after boundary crossing
RULL: amount to be scored (weight of the particle if current is
being scored; weight of the particle divided by the cosine
between its trajectory and the normal to the surface if
fluence is being scored))
ICALL, LLO, ECONTR: see GEODEN above
GEOBDX handles indirectly boundary crossing fluence (requested by
the USRBDX commands) via a call to USRSCX.
Note that fluence weighting via user routine FLUSCW is
done
only inside USRSCX: if the latter are bypassed by the
user,
no weighting can be applied, even if EXTRAWEIgh has been
requested.
GEOCLL (collision events)
Estimator: USRCOLL
Argument list (all variables are input only):
IJ, XA, YA, ZA, MREG: see GEODEN above
RULL: amount to be scored (weight of the colliding particle
divided
by its total macroscopic cross section)
ICALL, LLO, ECONTR: see GEODEN above
GEOCLL handles indirectly collision fluence (requested by the
USRCOLL commands) via a call to USRSCC.
Note that fluence weighting via user routine FLUSCW is
done
only inside USRSCC: if the latter are bypassed by the
user,
no weighting can be applied, even if EXTRAWEIgh has been
requested.
GEORSN (production of residual nuclei)
Estimator: RESNUCLEi
Argument list (all variables are input only):
IJ, XA, YA, ZA, MREG: see GEODEN above
RULL: amount to be scored (weight of the particle producing the
residual nucleus)
ICALL, LLO, ECONTR: see GEODEN above
GEORSN handles indirectly residual nuclei production (requested
by
the RESNUCLEi commands) via a call to USRSRN.
The properties (atomic number, mass number, energy etc.)
of
the main residual nucleus produced in a collision are
available in COMMON RESNUC. However, treatment of
isomers and identification of further possible nuclei
created in the same collision is done only inside USRSRN.
Therefore, overriding of the standard FLUKA treatment by
the
user is not recommended for this particular estimator.
===============================================================================
--------
lattic.f (symmetry transformation for lattice geometry)
--------
Subroutine LATTIC is activated by one or more LATTICE cards in the
geometry input (see 10}). It is expected to transform coordinates
and
direction cosines from any lattice cell (defined by card LATTICE)
to
the reference system in which the basic structure has been
described.
The user is expected to provide a transformation of coordinates and
vector direction cosines from each lattice cell to the
corresponding
basic structure (in entry LATTIC) and of direction cosines from the
basic structure to each corresponding lattice cell (in entry
LATNOR).
Entries:
LATTIC (position and direction symmetry transformation from
lattice
cell to prototype structure)
Argument list:
XB(1), XB(2), XB(3) : actual physical position coordinates in
IRLTGG
lattice cell
WB(1), WB(2), WB(3) : actual physical direction cosines in IRLTGG
lattice cell
DIST : reserved variable???
SB(1), SB(2), SB(3) : transformed coordinates in prototype cell
UB(1), UB(2), UB(3) : transformed cosines in prototype cell
IR : region number corresponding to cell IRLTGG
IRLTGG : lattice cell number
IFLAG : reserved variable
LATTIC returns the tracking point coordinates (SB) and direction
cosines(UB) in the reference prototype geometrical structure,
corresponding to real position/direction XB, WB in the actual cell
IRLTGG (defined as input region IR by a LATTICE card).
When the lattice option is activated, the tracking proceeds in two
different systems: the "real" one, and that of the basic symmetry
unit.
Particle positions and directions are swapped from their real
values
to their symmetric ones in the basic cell, to perform the physical
transport in the regions and materials that form the prototype
geometrical structure and back again to the real world. The
correspondence between "real" and "basic" position/direction
depends
on the symmetry transformation and on the lattice cell number.
LATNOR (LATtice cell NORmal transformation from prototype structure
to
lattice cell)
Argument list:
UN(1), UN(2), UN(3) : direction cosines of the vector normal to
the
surface, in the prototype cell (entry
values)
and in the lattice cell (returned values)
IRLTNO : present lattice cell number
Entry LATNOR transforms the direction cosines stored in the vector
UN(3) from the system of the basic prototype unit to that of the
real
world in lattice cell number IRLTNO. Therefore, this cosine
transformation must be the inverse of that performed on the cosines
by
the LATTIC entry: but while LATTIC maps vector UB to a different
vector WB, LATNOR maps the UN vector to itself.
Note that if the transformation implies a rotation, it is necessary
to
save first the incoming UN cosines to local variables, to avoid
overwriting the vector before all transformation statements are
executed.
NOTES:
Different symmetry transformations can of course be implemented
in the same LATTIC routine (each being activated by a different
cell
number or range of cell numbers).
The advantage of the lattice geometry is to avoid describing in
detail
the geometry of repetitive multi-modular structures. It must be
realized, however, that a penalty is generally paid in computer
efficiency.
Also, a region contained in the prototype cell and all those
"mirrored"
to it inside lattice cells are treated by the program as if they
were connected with "non-overlapping ORs" into a single region.
Therefore, any region-based scoring (options SCORE, USRTRACK, etc.)
can only provide quantities averaged over the whole structure. More
detailed information must be obtained by region-independent options
such as USRBIN or by user-written routines (MGDRAW).
===============================================================================
--------
magfld.f definition of a magnetic field
--------
Argument list:
X,Y,Z: current position (input only)
BTX,BTY,BTZ: direction cosines of the magnetic field vector
(returned).
B: magnetic field intensity in Tesla (returned)
NREG: current region (input only)
IDISC: if returned = 1, the particle will be discarded
MAGFLD is activated by option MGNFIELD with WHAT(4-6) = 0.0 and
is
used to return intensity and direction of a magnetic
field
based on the current position and region. It is called
only
if the current region has been flagged as having a
non-zero
magnetic field by option ASSIGNMAt, with WHAT(5) = 1.0.
The magnetic field spatial distribution is often read and
interpolated from an external field map. Note that in any
case the direction cosines MUST be properly normalized in
double precision (e.g. BTX = SQRT(ONEONE - BTY**2 -
BTZ**2)),
even if B = 0.0.
Please read carefully the notes on option MGNFIELD.
===============================================================================
--------
mdstck.f management of the stack of secondaries
--------
Argument list:
IFLAG: type of nuclear interaction which has produced
secondaries:
1: inelastic
2: elastic
3: low-energy neutron
NUMSEC: number of secondary particles produced in the interaction
MDSTCK is called after a nuclear interaction in which at least
one
secondary particle has been produced, before any biasing
is applied to decide which secondary will be loaded in
the main stack for further transport. The properties of
the secondaries are stored in the secondary stack (COMMON
FINUC). With MDSTCK, the user can analyze those
secondaries,
write them to a file, or even modify the content of FINUC
(for instance applying his own biasing). In the latter
case,
however, it is his responsibility to make sure that energy
is
conserved, the various physical quantities are still
consistent, etc.
===============================================================================
--------
mgdraw.f general event interface
--------
Subroutine MGDRAW, activated by option DUMPTHEM with WHAT(1) >=
100.,
usually writes a "collision tape", i.e. a file where all or
selected
transport events are recorded. The default version (unmodified by
the
user) offers several possibilities, selected by WHAT(3) in
DUMPTHEM.
Details are given in 17}.
Additional flexibility is offered by a user entry USDRAW,
interfaced
with the most important physical events happening during particle
transport.
The user can modify of course also any other entry of this
subroutine (MGDRAW for trajectory drawing, ENDRAW for recording of
energy deposition events, and SODRAW for recording of source
events):
for instance the format of the output file can be changed, and
different combinations of events can be written to file.
But the most interesting aspect of the routine is that the four
entries (all of which, if desired, can be activated at the same
time
by setting DUMPTHEM with WHAT(3) = 0.0 and WHAT(4) >= 1.0)
constitute
a complete interface to the whole FLUKA transport. Therefore,
MGDRAW
can be used not only to write a collision tape, but to do any kind
of complex analysis (for instance studying correlations between
events).
Entries:
MGDRAW (trajectory dumping for drawing)
Argument list (all variables are input only):
ICODE: FLUKA physical compartment originating the call
= 1: call from Kaskad (hadrons and muons)
= 2: call from Emfsco (electrons, positrons and photons)
= 3: call from Kasneu (low-energy neutrons)
= 4: call from Kashea (heavy ions)
= 5: call from Kasoph (optical photons)
MREG: current region
MGDRAW writes by default, for each trajectory, the following
variables (contained in COMMON TRACKR):
NTRACK: number of track segments
MTRACK: number of energy deposition events along the track
JTRACK: type of particle
ETRACK: total energy of the particle
WTRACK: weight of the particle
Ntrack values of XTRACK, YTRACK, ZTRACK: end of each track
segment
Mtrack values of DTRACK: energy deposited at each
deposition
event
CTRACK: total length of the curved path
Other variables are available in TRACKR (but not written
by
MGDRAW unless the latter is modified by the user: particle
momentum, direction cosines, cosines of the polarization
vector, age, generation, etc. (see a full list in the
comment
in the INCLUDE file).
ENDRAW (energy deposition dumping)
Argument list (all variables are input only):
ICODE: type of event originating energy deposition
1x: call from Kaskad (hadrons and muons)
10: elastic interaction recoil
11: inelastic interaction recoil
12: stopping particle
13: pseudo-neutron deposition
14: particle escaping (energy deposited in blackhole)
2x: call from Emfsco (electrons, positrons and photons)
20: local energy deposition (i.e. photoelectric)
21 or 22: particle below threshold
23: particle escaping (energy deposited in blackhole)
3x: call from Kasneu (low-energy neutrons)
30: target recoil
31: neutron below threshold
32: neutron escaping (energy deposited in blackhole)
4x: call from Kashea (heavy ions)
40: ion escaping (energy deposited in blackhole)
5x: call from Kasoph (optical photons)
50: optical photon absorption
51: optical photon escaping (energy deposited in
blackhole)
MREG: current region
RULL: energy amount deposited
XSCO, YSCO, ZSCO: point where energy is deposited
ENDRAW writes by default, for each energy deposition point:
0: flag identifying Endraw output from that of other
entries
ICODE: see argument list
JTRACK, ETRACK, WTRACK: see MGDRAW above
XSCO, YSCO, ZSCO, RULL: see argument list
SODRAW (source particle dumping)
No arguments
SODRAW writes by default, for each source or beam particle:
-NCASE (in COMMON CASLIM, with a minus sign to identify
Sodraw output): number of primaries followed so far
LSTACK (in COMMON STACK): stack pointer
LSTMAX (in COMMON STACK): highest value of the stack
pointer
encountered so far
TKESUM (in COMMON EPISOR): total kinetic energy of the
primaries of a user written source (see source.f
here below), if applicable. Otherwise = 0.0
WEIPRI (in COMMON STARS): total weight of the primaries
handled so far
Lstack times: ILO: type of source particle
TKE+AM: total particle energy (kinetic+mass)
WT: source particle weight
XA, YA, ZA: source particle position
TX, TY, TZ: source particle direction
cosines
USDRAW
Argument list (all variables are input only):
ICODE: 10x: call from Kaskad (hadron and muon interactions)
100: elastic interaction secondaries
101: inelastic interaction secondaries
102: particle decay secondaries
103: delta ray generation secondaries
104: pair production secondaries
105: bremsstrahlung secondaries
20x: call from Emfsco (electron, positron and photon
interactions)
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
30x: call from Kasneu (low-energy neutron interactions)
300: interaction secondaries
40x: call from Kashea (heavy ion interactions)
400: delta ray generation secondaries
MREG: current region
XSCO, YSCO, ZSCO: interaction point
USDRAW is called after each particle interaction (if requested by
the user with DUMPTHEM, WHAT(4) >= 1.0). There is no default
output:
any output must be supplied by the user.
Information about the secondary particles produced is available
in COMMON FINUC, except that concerning delta rays produced by
heavy
ions (in which case the properties of the single electron
produced
are available in COMMON EMFSTK, with index NP).
Information about the interacting particle and its trajectory can
be
found in COMMON TRACKR (see description under the MGDRAW entry
above). In TRACKR there are also some spare variables at the
user's
disposal: LLOUSE (integer), ISPUSR (integer array) and SPAUSR
(double precision array). Like many other TRACKR variables, each
of
them has a correspondent in the particle stacks, i.e. the COMMONs
from which the particles are unloaded at the beginning of their
transport: STACK, EMFSTK and OPPHST (respectively, the stack of
hadrons/muons, electrons/photons, and optical photons). The
correspondence with TRACKR is shown below under STUPRF/STUPRE.
When a particle is generated, its properties (weight, momentum,
energy, coordinates etc., as well as the values of the user
flags)
are loaded into one of the stacks. The user can write a STUPRF or
STUPRE subroutine (see description below) to change anyone of
such
flags just before it is saved in stack.
When a particle starts to be transported, its stack variables are
copied to the corresponding TRACKR ones. Unlike the other TRACKR
variables, which in general become modified during transport due
to
energy loss, scattering etc., the user flags keep their original
value copied from stack until they are changed by the user
himself
(generally in USDRAW).
One common application is the following: after an interaction
which
has produced sencondaries, make USDRAW copy some properties of
the
interacting particle into the TRACKR user variables. When STUPRF
is
called next to load the secondaries into stack, by default copies
the TRACKR user variables to the stack ones. In this way,
information about the parent can be still carried by its
daughters
(and possibly by further descendants). This technique is
sometimes
referred to as "latching".
===============================================================================
--------
soevsv.f SOurce EVent SaVing
--------
No arguments
Subroutine SOEVSV is always called after a beam particle is loaded
into STACK, but a call to SOEVS can be inserted by the user
anywhere
in a user routine.
SOEVSV copies the whole STACK to another COMMON, SOUEVT, which can
be included in other user routines. In other words, this routine is
used to "take a snapshot" of the particle bank at a particular time
for further use (interfacing to independent generators, etc.)
===============================================================================
--------
source.f user-written source
--------
Argument list:
NOMORE: if set = 1, no more calls will occur (the run will be
terminated after exhausting the primary particles loaded
into STACK in the present call). The history number limit
set with option START will be overridden.
Subroutine SOURCE is probably the most frequently used user
routine.
It is activated by option SOURCE and is used to sample primary
particle properties from distributions (in space, energy, time,
direction or mixture of particles) too complicated to be described
with the BEAM and BEAMPOS cards alone. For each phase-space
variable,
a value must be loaded into the STACK COMMON (particle bank) before
returning control. These values can be read from a file, generated
by
some sampling algorithm, or just assigned.
Reading from a file:
--------------------
Reading from a file is needed, for instance, when the particle data
are taken from a collision file, writtem by FLUKA or by another
program. The user must open the file with a unit number > 20 (unit
numbers lower than 20 are reserved), in one of the following ways:
1) Using option OPEN, with SDUM = OLD
2) In a user subroutine USRINI (see below), with a Fortran OPEN
statement. Option USRICALL is needed to activate the call to
the routine.
3) With an OPEN statement in the initialization part of
subroutine
SOURCE itself.
Then, a READ statement in SOURCE can be used to get the data to
load
in stack, for instance:
READ(21,*) IPART, X, Y, Z, COSX, COSY, COSZ, ENERGY, WEIGHT
ILO (LSTACK) = IPART
XA (LSTACK) = X
YA (LSTACK) = Y
ZA (LSTACK) = Z
TX (LSTACK) = COSX
...etc...
(LSTACK is the current stack index).
Direct assignment:
------------------
Direct assignment can be done explicitly, for instance:
PMOM (LSTACK) = 305.2D0
or implicitely, leaving unmodified values input with BEAM or
BEAMPOS:
PMOM (LSTACK) = PBEAM
(PBEAM is the momentum value input as WHAT(1) in option BEAM).
A set of direct assignments, one for each of several different
stack
entries, can be useful, for example, to define a series of RAYs
through
the geometry (see 23}):
DO 10 I = 1, 20
LSTACK = LSTACK + 1
ILO (LSTACK) = 0 (0 is the RAY particle id
number)
XA (LSTACK) = 500.D0 + DBLE(I) * 40.D0
YA (LSTACK) = 200.D0
...etc...
10 CONTINUE
Sampling from a uniform distribution:
-------------------------------------
To sample from a uniform distribution, the user must use the
function FLRNDM(DUMMY), which returns a double precision
pseudo-random number uniformly distributed between 0 (included)
and 1 (not included). Actually, DUMMY can be any variable name.
A simple example of sampling from a uniform distribution is that of
a linear source along the Z axis, between Z = 10 and Z = 80:
Z1 = 10.D0
Z2 = 80.D0
ZA (LSTACK) = 10.D0 + (Z2 - Z1) * FLRNDM(XXX)
Sampling from a generic distribution:
-------------------------------------
One way to sample a value XX from a generic distribution f(x) is
the following.
First integrate the distribution function, analytically or
numerically, and normalize to 1 the obtained cumulative
distribution:
/x /xmax
F(x) = | f(x)dx / | f(x)dx
/xmin /xmin
Then, sample a uniform pseudo-random number t using FLRNDM and
-1
get the desired result by finding the inverse value XX = F (t)
(analytically or most often by interpolation).
A FLUKA subroutine is available to sample directly from a Gaussian
distribution:
CALL FLNRRN (RGAUSS)
or, if two independent Gaussian distributed numbers are needed:
CALL FLNRR2 (RGAUS1, RGAUS2)
(faster than calling FLNRRN twice).
Sampling from a biased distribution:
------------------------------------
The sampling technique 1) described above can be extended to modify
the probability of sampling in different parts of the interval
(importance sampling). We replace f(x) by a weighted function
g(x) = f(x) * h(x), where h(x) is any appropriate function of x we
like to choose. We normalize g(x) in the same way as f(x) before:
/x /xmax /x
G(x) = | g(x)dx / | g(x)dx = | f(x)*h(x)dx / B
/xmin /xmin /xmin
and we need also the integral of f(x) over the whole interval:
/xmax
A = | f(x)dx
/xmin
All the sampling is done using the biased cumulative normalized
function G instead of the original unbiased F: we sample a uniform
pseudo-random number t as in 1), and we get the sampled value XX by
inverting G(x):
-1
XX = G (t)
The particle is assigned a weight B/(A * h(XX))
A special case of importance sampling is when the biasing function
chosen is the inverse of the unbiased distribution function:
h(x) = 1/f(x)
g(x) = f(x) * h(x) = 1
G(x) = (x - xmin) / (xmax - xmin)
In this case we sample a uniform pseudo-random number t using
FLRNDM as shown above. The sampled value XX is simply given by:
XX = xmin + (xmax - xmin)*t
and the particle is assigned a weight = f(XX)*(xmax - xmin).
Because XX is sampled with the same probability over all possible
values of x, independently of the value f(XX) of the function,
this technique is used to ensure that sampling is done uniformly
over the whole interval, even though f(x) might have very small
values somewhere. For instance it may be important to avoid
undersampling in the high-energy tail of a spectrum, steeply
falling
with energy but more penetrating, such as that of cosmic rays or
synchrotron radiation.
Option SOURCE allows the user to input up to 12 numerical values
(WHASOU(1),(2)...(12)) and one 8-character string (SDUSOU) which
can
be accessed by the subroutine by including the following lines:
INCLUDE '(EPISOR)'
INCLUDE '(CHEPSR)'
These values can be used as parameters or switches for a
multi-source
routine capable to handle several cases, or to identify an external
file to be read, etc., without having to compile and link again the
routine.
In the SOURCE routine there are a number of mandatory statements,
(clearly marked as such in accompanying comments) which must not be
removed or modified. The following IF block intializes the total
kinetic energy of the primary particles and sets two flags: the
first to skip the IF block in all next calls, and the second to
remind the program, when writing the final output, that a user
source
has been used:
* +-------------------------------------------------------------------*
* | First call initializations:
IF ( LFIRST ) THEN
* | *** The following 3 cards are mandatory ***
TKESUM = ZERZER
LFIRST = .FALSE.
LUSSRC = .TRUE.
* | *** User initialization ***
END IF
* |
* +-------------------------------------------------------------------*
The user can insert into the above IF block any other
initialization
needed, for instance the preparation of a cumulative spectrum
array from which to sample the energy of the source particles.
Note that user initialization can take place also in routine USRINI
(activated at input time by input option USRICALL) and USREIN
(called
before unloading from stack the first source particle of an event,
i.e., just after the call to SOURCE)
At the time SOURCE is called, the particle bank STACK is always
empty
and the stack pointer LSTACK has value 0.
The user can load into STACK one or more source particles at each
call: for each particle loaded the pointer must be increased by 1.
The template version of SOURCE loads only one particle: if several
are loaded the following sequence, until the statement CALL SOEVSV
not
included, must be repeated once for each particle, possibly inside
a
DO loop:
LSTACK = LSTACK + 1 increases the pointer
The following statements assign a value to each of the STACK
variables
concerning the particle being loaded.
WT (LSTACK) = ONEONE sets the weight of the particle =
1.
This must be changed if the sampling of one or more of the
particle properties are biased. In that case, generally the
weight
must be set after the sampling, and its value depends on the
sampling outcome.
WEIPRI = WEIPRI + WT (LSTACK) updates the total weight of the
primaries (don't change)
ILO (LSTACK) = IJBEAM by default sets the type of
particle
equal to the one defined by the BEAM card (if no BEAM card is
given in input, IJBEAM is = 1 (proton)).
The above statement is followed by several others that must not be
changed or removed. In the template routine, they are encompassed
by
the comment lines:
From this point .... / ... to this point: don't change anything
These statements are:
* From this point ....
LO (LSTACK) = 1 Generation is 1 for source
particles
LOUSE (LSTACK) = 0 User variables: the user can set
DO 100 ISPR = 1, MKBMX1 different values in the STUPRF or
SPAREK (1,LSTACK) = ZERZER STUPRE routine, but it is better
100 CONTINUE not to do it here
DO 200 ISPR = 1, MKBMX2
ISPARK (ISPR,LSTACK) = 0 More user variables (integer)
200 CONTINUE
NPARMA = NPARMA + 1 Updating the maximum particle
number
NUMPAR (LSTACK) = NPARMA Setting the current particle
number
NEVENT (LSTACK) = 0 Resetting the current event
number
DFNEAR (LSTACK) = +ZERZER Resetting the distance to the
nearest boundary
* ... to this point: don't change anything
The following statements can be overridden or rewritten by the
user,
assigning new values or sampling them from problem-dependent
distributions.
First three statements which are rarely modified:
AGESTK (LSTACK) = +ZERZER Particle age is zero by default
AKNSHR (LSTACK) = -TWOTWO Resets the Kshort component of
K0/K0bar. Usually not to be
changed.
IGROUP (LSTACK) = 0 Group number for low energy
neutrons: if set to 0, the
program
derives it from the kinetic
energy
Then the most frequently changed lines: both energy and momentum of
the particle must be loaded onto STACK, but the two cannot be
defined
independently. Appropriate kinematical (relativistic) relations
must
be applied to derive one from the other.
In the template routine, the momentum is assumed to be assigned by
BEAM option (its value, PBEAM, is taken from COMMON BEAM, which
contains all values defined by options BEAM and BEAMPOS).
PMOM (LSTACK) = PBEAM
Therefore, the kinetic energy (in GeV) must be derived:
TKE (LSTACK) = SQRT ( PBEAM**2 + AM (IJBEAM)**2 ) - AM (IJBEAM)
(where AM is the rest mass, in COMMON PAPROP, and IJBEAM is the
particle type, in COMMON BEAM)
If instead the energy had been sampled first from some spectrum,
and
ENSAMP would be the sampled value, the two statements above would
become:
TKE (LSTACK) = ENSAMP
PMOM (LSTACK) = SQRT(ENSAMP * (ENSAMP + TWOTWO * AM(IJBEAM)))
The direction cosines are loaded next:
TX (LSTACK) = TINX Assumed here to be the same as
TY (LSTACK) = TINY defined by option BEAMPOS. TINX,
TZ (LSTACK) = TINZ TINY, TINZ are some among the
beam
properties in COMMON BEAM.
(If BEAMPOS is not given, by default TINX = TINY = 0.0, TINZ = 1.0)
Remember to make sure that the cosines are normalized! One could
replace the last statement by:
TZ (LSTACK) = SQRT ( ONEONE - TX(LSTACK)**2 - TY(LSTACK)**2 )
The polarization cosines are not set by default:
TXPOL (LSTACK) = -TWOTWO -2 is a flag for "no
polarization"
TYPOL (LSTACK) = +ZERZER
TZPOL (LSTACK) = +ZERZER
but appropriate values need to be given in some cases, for instance
in
synchrotron radiation shielding problems.
Finally the particle coordinates, set again by default equal to
those
input with BEAMPOS:
XA (LSTACK) = XINA Assumed here to be the same as
YA (LSTACK) = YINA defined by option BEAMPOS. XINA,
ZA (LSTACK) = ZINA YINA, ZINA are also in COMMON
BEAM.
(If BEAMPOS is not given, by default XINA = YINA = ZINA = 0.0)
If for example our problem required instead a linear source
uniformly
distributed along Z between Z1 and Z2, we could replace the last
statement by:
ZA (LSTACK) = Z1 + FLRNDM(UGH) * (Z2 - Z1)
The following lines in the template SOURCE routine should never be
changed. They calculate the total energy of the primary particles,
define the remaining properties of the particles (starting region
and
lattice cell) and do some geometry initialization.
The last line calls the SOEVSV user routine (see description above)
to save the stack for possible further use.
===============================================================================
--------
stupre.f SeT User PRoperties for EMF particles
stuprf.f SeT User PRoperties for FLUKA particles
--------
These two functions are used to assign a value to one or more stack
user variables when the corresponding particle is loaded into one
of
the stacks (STACK for hadrons/muons, EMFSTK for electrons/photons,
OPPHST for optical photon).
In each of these stacks the user has access to one integer
variable,
one integer array and one double precision array. Each of them is
copied to a correspondent variable or array in COMMON TRACKR at the
beginning of transport:
Correspondence: STACK EMFSTK OPPHST TRACKR
-------------- ----- ------ ------ ------
integer variable: LOUSE LOUEMF LOUOPP --> LLOUSE
integer array: ISPARK IESPAK ISPORK --> ISPUSR
double precision array: SPAREK ESPARK SPAROK --> SPAUSR
The user can access and modify the TRACKR variables via subroutine
MGDRAW and its entries ENDRAW, SODRAW and especially USDRAW (see
description above). STUPRF and STUPRE can be used to do the
reverse,
namely to copy TRACKR user variables to those of the relevant stack
(see USDRAW above).
STUPRE is called before loading into stack electrons, positrons and
photons.
No arguments
The default version does nothing (the user variables of the
parent
particle are already set equal to the original projectile by the
various electromagnetic interaction routines. Also the
region/position etc. are already set inside the stack arrays.
STUPRF is called before loading into stack hadrons, muons,
neutrinos
low-energy neutrons and optical photons
Argument list:
IJ : type of the parent particle
MREG : current region
XX, YY, ZZ: particle position
NUMSEC : index in the FINUC COMMON of the secondary being loaded
onto stack
NPPRMR : if > 0, the secondary being loaded is actually still the
interacting particle (it can happen in some biasing
situations)
The default version copies to stack the user flags of the parent.
===============================================================================
--------
ubsset.f User BiaSing SETting
--------
Argument list:
IR : region number
RRHADR : multiplicity biasing factor to be applied to the
secondaries from hadronic interactions in region IR
(WHAT(2) of card BIASING)
HMPHAD : Importance of region IR for hadrons and muons (WHAT(3)
of
card BIASING, with WHAT(1) = 0 or 1).
Actually the routine argument is an integer, IMPHAD,
equal
to importance multiplied by 10000, but the user should
consider only the double precision version HMPHAD
(a conversion from and to the integer version is
provided
at the beginning and at the end of the routine, and
should
not be changed).
HMPLOW : Importance of region IR for low energy neutrons (WHAT(3)
of card BIASING, with WHAT(1) = 0 or 3).
Actually the routine argument is an integer, IMPLOW,
equal
to importance multiplied by 10000, but the user should
consider only the double precision version HMPLOW
(a conversion from and to the integer version is
provided
at the beginning and at the end of the routine, and
should
not be changed).
HMPEMF : Importance of region IR for electrons and photons
(WHAT(3)
of card BIASING, with WHAT(1) = 0 or 2).
Actually the routine argument is an integer, IMPEMF,
equal
to importance multiplied by 10000, but the user should
consider only the double precision version HMPEMF
(a conversion from and to the integer version is
provided
at the beginning and at the end of the routine, and
should
not be changed).
IGCUTO : Cut-off group index for low energy neutrons in region IR
(WHAT(1) in card LOW-BIAS)
IGNONA : Non-analog absorption group limit for low energy
neutrons
in region IR (WHAT(2) in card LOW-BIAS)
PNONAN : Non-analog survival probability for low energy neutrons
in
region IR (WHAT(3) in card LOW-BIAS)
IGDWSC : Group limit for biased downscattering for low energy
neutrons in region IR (WHAT(1) in card LOW-DOWN)
FDOWSC : Biased downscattering factor for low energy neutrons in
region IR (WHAT(2) in card LOW-DOWN)
JWSHPP : Weight-Window/importance profile index for low energy
neutrons in region IR (SDUM in WW-FACTO)
WWLOW : Weight-Window lower level in region IR (WHAT(1) in card
WW-FACTO, possibly modified by WHAT(4) in WW-THRES or
WHAT(2) in WW-PROFI)
WWHIG : Weight-Window upper level in region IR (WHAT(2) in card
WW-FACTO, possibly modified by WHAT(4) in WW-THRES or
WHAT(2) in WW-PROFI)
WWMUL : Weight-Window multiplicative factor applied to the two
energy thresholds defined with WW-THRES, for region IR
(WHAT(3) in card WW-FACTO)
EXPTR : Exponential transform parameter for region IR (WHAT(2)
in
card EXPTRANS) (not implemented yet!!!!!!!!)
ELECUT : e+,e- cut-off in region IR (WHAT(1) in card EMFCUT)
GAMCUT : Photon cut-off in region IR (WHAT(2) in card EMFCUT)
LPEMF : Leading Particle Biasing flag in region IR (SDUM =
LPBEMF
in card EMF-BIAS, or WHAT(3) in card EMFCUT)
ELPEMF : Maximum e+/e- energy for applying Leading Particle
Biasing
(WHAT(2) in card EMF-BIAS with SDUM = LPBEMF)
PLPEMF : Maximum photon energy for applying leading particle
biasing (WHAT(3) in card EMF-BIAS with SDUM = LPBEMF)
Subroutine UBSSET does not require a special command to be
activated:
is always called several times for each region: (once for every
biasing option or suboption) after the end of input reading and
before starting the calculations. The default version is a dummy
and
does nothing. The user can replace it to override any biasing
parameters specified in input.
The UBSSET subroutine is used especially in cases with a large
number
of regions, because it allows to derive the biasing parameters from
simple algorithms instead of entering each input value by hand.
Choosing an appropriate numbering scheme for the geometry regions
can
often facilitate the task.
For instance, assuming a simple slab geometry with an expected
exponential hadron attenuation from region 3 to region 20, each
region being one half-value-layer thick, one could write the
following in order to set importances that would keep the hadron
number about constant in all regions:
IF(IR .GE. 3 .AND. IR .LE. 20) HMPHAD = ONEONE * TWOTWO**(IR-3)
It is important, however, not to do recursive assignments of the
type:
GAMCUT(5) = GAMCUT(5) * HLFHLF (HLFHLF is a FLUKA constant for
0.5D0)
because that would halve the value of the photon cut-off for
Region 5 AT EVERY CALL, and the number of calls is not known to the
user.
===============================================================================
--------
udcdrl.f User defined DeCay DiRection biasing and Lambda (for neutrinos
only)
--------
Argument list:
IJ : type of decaying particle
KPB : outgoing neutrino (???)
NDCY : (???)
UDCDRB, VDCDRB, WDCDRB : cosines of the preferential outgoing
direction for the neutrino (???)
Function UDCDRL is used to bias the direction of a neutrino emitted
by a decaying particle of type IJ. The value returned by UDCDRL is
the Lambda for direction biasing (1-cos(theta)) (???)
===============================================================================
--------
usrein.f USeR Event INitialization
--------
No arguments
Subroutine USREIN is called just before the first source particle
of
an event is unloaded from stack and begins to be transported.
An event is the full history of a group of related particles and
their
descendants. If primaries are loaded into stack by the input option
BEAM, there is only one source particle per event; but there can be
more if the user routine SOURCE is used to load particles into
stack.
USREIN does not need any special command to be activated, but the
default version of USREIN does nothing: the user can write here any
kind of initialization.
===============================================================================
--------
usreou.f USeR Event OUtput (called at the end of each event)
--------
No arguments
Subroutine USREOU is called at the end of each event, namely after
all
event primary particles and their descendants have been
transported.
(See USREIN above for a definition of an event).
USREOU does not need any special command to be activated, but the
default version of USREOU does nothing: the user can write here
any
kind of event analysis, output, etc.
===============================================================================
--------
usrini.f USeR INItialisation
--------
Argument list:
WHAT(1), (2), (3), (4), (5), (6) : user-given numerical
parameters
SDUM : user-given character string (8 characters)
Subroutine USRINI is called every time a USRICALL card is read in
input. It can be used to do any kind of initialization: reading
and manipulating data from one or more files, calling other private
routines, etc.
The calling parameters can carry any kind of useful information or
can
be used as flags to choose between different possible actions to be
performed before any particle transport takes place.
===============================================================================
--------
usrmed.f USeR MEDium dependent directives
--------
Argument list:
IJ : particle type
EKSCO : particle kinetic energy (GeV)
PLA : particle momentum (GeV/c)
WEE : particle weight
MREG : previous region number
NEWREG : current region number
XX, YY, ZZ : particle position
TXX, TYY, TZZ : particle direction
Subroutine USRMED is activated by option MAT-PROP with SDUM =
USERDIRE,
for one or more materials indicated by the user.
It is called every time a particle is going to be transported in
one of
the user-flagged materials.
Two cases are possible:
1) MREG = NEWREG: the particle is going to move from a point inside
the medium. The user is allowed to change only the particle
weight.
Typical application: simulating attenuation of optical photons
in an
absorbing medium by reducing the photon weight.
2) MREG not = NEWREG: the particle is going to move from a point on
a
boundary between two different regions. The user may change any
of
the following: particle weight, current region number, direction
cosines.
Typical applications:
- simulating refraction, by changing the direction cosines so
that
the particle is still inside the same region. To do this, one
generally needs the direction cosines of the normal to the
surface: TXNOR(LSTACK), TYNOR(LSTACK), TZNOR(LSTACK) (COMMON
STACK must be included).
- simulating reflection (albedo) at a boundary. The direction
cosines must be modified according to some reflection law or
albedo angular distribution, and NEWREG must be set = MREG.
In both cases the weight can also be reduced to account for
surface
reflectivity or similar (if the particle is an optical photon,
the
FRGHNS user function can be called to establish a surface
roughness)
===============================================================================
--------
usrout.f USeR OUTput
--------
Argument list:
WHAT(1), (2), (3), (4), (5), (6) : user-given numerical
parameters
SDUM : user-given character string (8 characters)
Subroutine USROUT is called every time a USROCALL card is read in
input. It is used to print special user-written output in addition
to the standard one provided by default.
The calling parameters can carry any kind of useful information or
can
be used as flags to choose between different possible actions to be
performed after all particle transport has taken place.
*$ CREATE FINUC.ADD
*COPY FINUC
*
*=== 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) *
* *
*----------------------------------------------------------------------*
*
PARAMETER (MXP=999)
*
COMMON / FINUC / CXR (MXP), CYR (MXP), CZR (MXP),
& CXRPOL (MXP), CYRPOL (MXP), CZRPOL (MXP),
& TKI (MXP), PLR (MXP), WEI (MXP),
& AGESEC (MXP), TV, TVCMS, TVRECL, TVHEAV, TVBIND,
& NP0, NP, KPART (MXP)
DUMPTHEM}
defines a collision "tape" (or better a collision file) to be
written.
WHAT(1) >= 100 : A complete dump is written in unformatted
(binary) form of one or more of the following:
- all source particles
- all particle trajectories
- all local (pointlike) energy deposition events
- all continuous energy deposition events
- user-defined events
= 0.0 : ignored
=< 0.0 : the default is reset, i.e. no dump is written
> 0.0 and < 100: not allowed! Originally used to request
another form of collision tape. Presently only
the "new" form of collision tape is possible
(the old one being incompatible with the present version
of FLUKA).
WHAT(2) = number of the unformatted output unit.
Values of WHAT(2) < 21 must be avoided because of possible
conflicts with FLUKA pre-defined units.
Default: 49
WHAT(3) <= 0: source particles, trajectories, continuous and local
energy losses are all dumped
1: only source particles are dumped
2: only trajectories and continuous energy losses are
dumped
3: only local energy losses are dumped (e.g., heavy recoil
kerma, cut-off energy). Proton recoils are not included
(since recoil protons are transported by FLUKA)
4: source particles, trajectories and continuous
energy losses are dumped
5: source particles and local energy losses are
dumped
6: trajectories and all energy losses (continuous and
local) are dumped
>= 7: source particles, trajectories, continuous and local
energy losses are not dumped (but user-defined dumps
required by WHAT(4) are unaffected)
Default: source particles, trajectories, continuous and
local energy losses are all dumped, provided WHAT(1) >
1.0
WHAT(4)>= 1: user dependent dumps after collisions are activated
= 0: ignored
< 0: resets to default (user dependent dumps after
collisions are de-activated)
WHAT(5...6) = not used
SDUM = filename of the output file (max. 10 characters)
Default (option DUMPTHEM not given): no dump will be written.
Note: The format of the collision tape, the code number of the
events and the variables written for the different type of
events, are described in 17}.
Be careful, the amount of output can be enormous.
Example:
*...+....1....+....2....+....3....+....4....+....5....+....6....+....7...+...8
DUMPTHEM 200. 37.0 2.0
TRAKFILE
* It is requested to write a file TRAKFILE, containing all trajectories and
* continuous energy losses, and pre-connected to the logical output unit 37.
Partial thread listing:
- Re: FLUKA: Production of Hydrogen, (continued)
foucher
FLUKA: Help me...plz...,
=?ks_c_5601-1987?B?wK+828Dn?=
FLUKA: Help me...
=?ks_c_5601-1987?B?wK+828Dn?=