From: Sebastien WURTH (wurth@ipno.in2p3.fr)
Date: Tue Mar 27 2007 - 11:14:43 CEST
Hello Alfredo, Hello everybody,
I come back two weeks later on that matter.
I'am not able to obtain the same results with the two methods.
Maybe I do something wrong but I cannot find what.
I tried to sample an Am-Be source.
First, I binned the spectrum in the 72 groups, that was a little tricky
but I got excellent results (see attached file).
Then I tried to do it without forcing the energy groups, and my source
routine looks like Sunil's (uploaded earlier today on the fluka discuss)
and my results are not good, the spectrum is not properly rendered (see
file).
I used a USRTRACK detector in a 1 cm radius vaccum sphere centered
around the beam spot (isotropic source point) to score fluence in both
cases presented in the attached file.
The 3 plots in the *.pdf file are :
1) The original spectrum for an Am-Be source
2) The fluence spectrum scored in USRTRACK with the 72 groups and
cumulative spectrum described in source.f
3) The fluence spectrum socred in USRTACK without binning into the 72
groups, I declared 53 energies and 52 groups or intervals in source.f.
So, has anybody some idea of what is happening ? Results are not
identical as Alfredo said, what could lead to such a difference ?
Thank you.
Best Regards.
Sebastien.
Alfredo Ferrari wrote :
> Hi all
>
> you are not forced to bin the spectrum into the 72 groups, however the
> code will do it for you as soon as a neutron is generated by the
> source routine (that is it the corresponding group will be computed
> and used).
> At least you can space some efforts (results will be identical as you
> bin it into the 72 group structure yourself).
>
> Ciao
> Alfredo
>
> On Thu, 15 Mar 2007, Sebastien WURTH wrote:
>
>> Hello Sunil,
>>
>> I did almost the same work with Nisy's basis work on Cf252 (thanks
>> again!).
>> I tried to sample an Am-Be source with respect to the ISO 8529 norm.
>> I respected the 72 neutron groups and to be honnest did some
>> "internal kitchen work" to obtain the spectrum I wished.
>> Indeed the neutron groups intervalls are a bit to large to define
>> properly what really happens between 2 and 11 MeV, but I found my
>> sampling satisfying.
>> I leave your questions to the experts, but I would like to show the
>> differences between the spectrum obtained with the proper definition
>> of Am-Be in terms of ISO 8259 norm and what I managed to sample for
>> FLUKA. See the pdf page attached.
>> In terms of doses obtained, it matched with real measurements with
>> less than 5% error.
>>
>> Regards.
>> Sebastien.
>>
>> Sunil C wrote :
>>
>>> Hi All
>>>
>>> I have the Am-Be neutron spectrum modeled in the source.f attached
>>> here, based on what Nisy had done for Cf-252. I have few questions.
>>> 1. Is it necessary to bin the source spectrum strictly in accordance
>>> to the group structure? What happens if I use the spectrum as given
>>> in IAEA tecdoc 413 (2001), the lethargy converted to fluence and
>>> normalized to one, but energy remaining exactly as in the book?
>>> Such as
>>>
>>> En (MeV) Fluence
>>> 7.94E-02 -
>>> 1.00E-01 -
>>> 1.25E-01 --
>>> 1.58E-01 -
>>> 1.99E-01 -
>>> 2.51E-01 -
>>> 3.16E-01 -
>>> ......
>>> ....
>>> ....
>>> ...
>>>
>>> Should all the 72 groups be listed in the cumulative index? In this
>>> case the numbers are zero below 7e-2 MeV and above 14.9 MeV. The
>>> cumulative index will show lots of zero below 0.07 MeV and lots of
>>> 1.0 above 14.9MeV. In the source attached I have re-binned the
>>> spectrum from the IAEA reference to match the group structure and
>>> have retained only the first groups with 0.0 and 1.0 as cumulative
>>> index.
>>>
>>> 2. Does the neutron group 'NEUGRP' needs to be faithfully passed on
>>> from the source subroutine to the main program? IS the cross section
>>> invoked by this information or is there an internal check?
>>> In the source.f here I am starting from group 3 by setting the
>>> 'NEUGRP' index to 'NEUGRP+2'
>>>
>>> And Thank you very much Nisy.
>>>
>>>
>>> Best Regards
>>>
>>>
>>> Sunil C
>>>
>>>
>>> ------------------------------------------------------------------------
>>>
>>>
>>> *$ CREATE SOURCE.FOR
>>> *COPY SOURCE
>>> *
>>> *=== source
>>> ===========================================================*
>>> *
>>> SUBROUTINE SOURCE ( NOMORE )
>>>
>>> INCLUDE '(DBLPRC)'
>>> INCLUDE '(DIMPAR)'
>>> INCLUDE '(IOUNIT)'
>>> *----------------------------------------------------------------------*
>>>
>>> * Version for isotropic Cf source.
>>> *
>>> c DIMENSION CFCUM(0:72), ENEDGE(73), ANG(0:40), CUMANG(0:40)
>>> DIMENSION CFCUM(0:40), ENEDGE(41), 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 /
>>> c Am-Be neutron spectrum
>>> & 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/ *
>>> * Cumulative Am-Be spectrum in group strcuture
>>> DATA CFCUM / 0.D0,
>>> &1.871941772E-03,4.460515263E-03,
>>> &7.681224425E-03,1.166466537E-02,1.654782161E-02,2.244066112E-02,
>>> &2.947992860E-02,3.777423905E-02,4.725383086E-02,5.774037795E-02,
>>> &6.890897217E-02,8.058700982E-02,9.281815275E-02,1.055467203E-01,
>>> &1.186954005E-01,1.327763145E-01,1.492876327E-01,1.684029556E-01,
>>> &1.901719153E-01,2.149087836E-01,2.455909710E-01,2.844440519E-01,
>>> &3.313109456E-01,3.847880155E-01,4.455136331E-01,5.112235084E-01,
>>> &5.824383965E-01,6.541031298E-01,7.240992794E-01,7.917544886E-01,
>>> &8.561530386E-01,9.152071332E-01,9.564046266E-01,9.777466352E-01,
>>> &9.901464874E-01,9.956112046E-01,9.972629018E-01,9.985685291E-01,
>>> &9.994966257E-01,1.000000000E+00/
>>>
>>> *...+....1....+....2....+....3....+....4....+....5....+....6....+....7..
>>>
>>> * 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, 40
>>> IF(XI .LE. CFCUM(K)) THEN
>>> NEUGRP = K+2
>>> 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
>>>
>>>
>>
>>
>
This archive was generated by hypermail 2.1.6 : Tue Mar 27 2007 - 13:59:21 CEST