mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 05:35:49 -04:00
Make Nuclide%reactions allocatable by using associate constructs
This commit is contained in:
Paul Romano
parent
9b87c1b77d
commit
a3199df9d9
6 changed files with 324 additions and 340 deletions
341
src/ace.F90
341
src/ace.F90
|
|
@ -715,7 +715,6 @@ contains
|
|||
integer :: IE ! reaction's starting index on energy grid
|
||||
integer :: NE ! number of energies
|
||||
integer :: NR ! number of interpolation regions
|
||||
type(Reaction), pointer :: rxn
|
||||
type(ListInt) :: MTs
|
||||
|
||||
LMT = JXS(3)
|
||||
|
|
@ -733,14 +732,15 @@ contains
|
|||
! Store elastic scattering cross-section on reaction one -- note that the
|
||||
! sigma array is not allocated or stored for elastic scattering since it is
|
||||
! already stored in nuc % elastic
|
||||
rxn => nuc % reactions(1)
|
||||
rxn % MT = 2
|
||||
rxn % Q_value = ZERO
|
||||
rxn % multiplicity = 1
|
||||
rxn % threshold = 1
|
||||
rxn % scatter_in_cm = .true.
|
||||
rxn % has_angle_dist = .false.
|
||||
rxn % has_energy_dist = .false.
|
||||
associate (rxn => nuc % reactions(1))
|
||||
rxn % MT = 2
|
||||
rxn % Q_value = ZERO
|
||||
rxn % multiplicity = 1
|
||||
rxn % threshold = 1
|
||||
rxn % scatter_in_cm = .true.
|
||||
rxn % has_angle_dist = .false.
|
||||
rxn % has_energy_dist = .false.
|
||||
end associate
|
||||
|
||||
! Add contribution of elastic scattering to total cross section
|
||||
nuc % total = nuc % total + nuc % elastic
|
||||
|
|
@ -753,64 +753,64 @@ contains
|
|||
i_fission = 0
|
||||
|
||||
do i = 1, NMT
|
||||
rxn => nuc % reactions(i+1)
|
||||
associate (rxn => nuc % reactions(i+1))
|
||||
! set defaults
|
||||
rxn % has_angle_dist = .false.
|
||||
rxn % has_energy_dist = .false.
|
||||
|
||||
! set defaults
|
||||
rxn % has_angle_dist = .false.
|
||||
rxn % has_energy_dist = .false.
|
||||
! read MT number, Q-value, and neutrons produced
|
||||
rxn % MT = int(XSS(LMT + i - 1))
|
||||
rxn % Q_value = XSS(JXS4 + i - 1)
|
||||
rxn % multiplicity = abs(nint(XSS(JXS5 + i - 1)))
|
||||
rxn % scatter_in_cm = (nint(XSS(JXS5 + i - 1)) < 0)
|
||||
|
||||
! read MT number, Q-value, and neutrons produced
|
||||
rxn % MT = int(XSS(LMT + i - 1))
|
||||
rxn % Q_value = XSS(JXS4 + i - 1)
|
||||
rxn % multiplicity = abs(nint(XSS(JXS5 + i - 1)))
|
||||
rxn % scatter_in_cm = (nint(XSS(JXS5 + i - 1)) < 0)
|
||||
! Read energy-dependent multiplicities
|
||||
if (rxn % multiplicity > 100) then
|
||||
! Set flag and allocate space for Tab1 to store yield
|
||||
rxn % multiplicity_with_E = .true.
|
||||
allocate(rxn % multiplicity_E)
|
||||
|
||||
! Read energy-dependent multiplicities
|
||||
if (rxn % multiplicity > 100) then
|
||||
! Set flag and allocate space for Tab1 to store yield
|
||||
rxn % multiplicity_with_E = .true.
|
||||
allocate(rxn % multiplicity_E)
|
||||
XSS_index = JXS(11) + rxn % multiplicity - 101
|
||||
NR = nint(XSS(XSS_index))
|
||||
rxn % multiplicity_E % n_regions = NR
|
||||
|
||||
XSS_index = JXS(11) + rxn % multiplicity - 101
|
||||
NR = nint(XSS(XSS_index))
|
||||
rxn % multiplicity_E % n_regions = NR
|
||||
! allocate space for ENDF interpolation parameters
|
||||
if (NR > 0) then
|
||||
allocate(rxn % multiplicity_E % nbt(NR))
|
||||
allocate(rxn % multiplicity_E % int(NR))
|
||||
end if
|
||||
|
||||
! allocate space for ENDF interpolation parameters
|
||||
if (NR > 0) then
|
||||
allocate(rxn % multiplicity_E % nbt(NR))
|
||||
allocate(rxn % multiplicity_E % int(NR))
|
||||
! read ENDF interpolation parameters
|
||||
XSS_index = XSS_index + 1
|
||||
if (NR > 0) then
|
||||
rxn % multiplicity_E % nbt = get_int(NR)
|
||||
rxn % multiplicity_E % int = get_int(NR)
|
||||
end if
|
||||
|
||||
! allocate space for yield data
|
||||
XSS_index = XSS_index + 2*NR
|
||||
NE = nint(XSS(XSS_index))
|
||||
rxn % multiplicity_E % n_pairs = NE
|
||||
allocate(rxn % multiplicity_E % x(NE))
|
||||
allocate(rxn % multiplicity_E % y(NE))
|
||||
|
||||
! read yield data
|
||||
XSS_index = XSS_index + 1
|
||||
rxn % multiplicity_E % x = get_real(NE)
|
||||
rxn % multiplicity_E % y = get_real(NE)
|
||||
end if
|
||||
|
||||
! read ENDF interpolation parameters
|
||||
XSS_index = XSS_index + 1
|
||||
if (NR > 0) then
|
||||
rxn % multiplicity_E % nbt = get_int(NR)
|
||||
rxn % multiplicity_E % int = get_int(NR)
|
||||
end if
|
||||
! read starting energy index
|
||||
LOCA = int(XSS(LXS + i - 1))
|
||||
IE = int(XSS(JXS7 + LOCA - 1))
|
||||
rxn % threshold = IE
|
||||
|
||||
! allocate space for yield data
|
||||
XSS_index = XSS_index + 2*NR
|
||||
NE = nint(XSS(XSS_index))
|
||||
rxn % multiplicity_E % n_pairs = NE
|
||||
allocate(rxn % multiplicity_E % x(NE))
|
||||
allocate(rxn % multiplicity_E % y(NE))
|
||||
|
||||
! read yield data
|
||||
XSS_index = XSS_index + 1
|
||||
rxn % multiplicity_E % x = get_real(NE)
|
||||
rxn % multiplicity_E % y = get_real(NE)
|
||||
end if
|
||||
|
||||
! read starting energy index
|
||||
LOCA = int(XSS(LXS + i - 1))
|
||||
IE = int(XSS(JXS7 + LOCA - 1))
|
||||
rxn % threshold = IE
|
||||
|
||||
! read number of energies cross section values
|
||||
NE = int(XSS(JXS7 + LOCA))
|
||||
allocate(rxn % sigma(NE))
|
||||
XSS_index = JXS7 + LOCA + 1
|
||||
rxn % sigma = get_real(NE)
|
||||
! read number of energies cross section values
|
||||
NE = int(XSS(JXS7 + LOCA))
|
||||
allocate(rxn % sigma(NE))
|
||||
XSS_index = JXS7 + LOCA + 1
|
||||
rxn % sigma = get_real(NE)
|
||||
end associate
|
||||
end do
|
||||
|
||||
! Create set of MT values
|
||||
|
|
@ -821,56 +821,57 @@ contains
|
|||
|
||||
! Create total, absorption, and fission cross sections
|
||||
do i = 2, size(nuc % reactions)
|
||||
rxn => nuc % reactions(i)
|
||||
IE = rxn % threshold
|
||||
NE = size(rxn % sigma)
|
||||
associate (rxn => nuc % reactions(i))
|
||||
IE = rxn % threshold
|
||||
NE = size(rxn % sigma)
|
||||
|
||||
! Skip total inelastic level scattering, gas production cross sections
|
||||
! (MT=200+), etc.
|
||||
if (rxn % MT == N_LEVEL) cycle
|
||||
if (rxn % MT > N_5N2P .and. rxn % MT < N_P0) cycle
|
||||
! Skip total inelastic level scattering, gas production cross sections
|
||||
! (MT=200+), etc.
|
||||
if (rxn % MT == N_LEVEL) cycle
|
||||
if (rxn % MT > N_5N2P .and. rxn % MT < N_P0) cycle
|
||||
|
||||
! Skip level cross sections if total is available
|
||||
if (rxn % MT >= N_P0 .and. rxn % MT <= N_PC .and. MTs % contains(N_P)) cycle
|
||||
if (rxn % MT >= N_D0 .and. rxn % MT <= N_DC .and. MTs % contains(N_D)) cycle
|
||||
if (rxn % MT >= N_T0 .and. rxn % MT <= N_TC .and. MTs % contains(N_T)) cycle
|
||||
if (rxn % MT >= N_3HE0 .and. rxn % MT <= N_3HEC .and. MTs % contains(N_3HE)) cycle
|
||||
if (rxn % MT >= N_A0 .and. rxn % MT <= N_AC .and. MTs % contains(N_A)) cycle
|
||||
if (rxn % MT >= N_2N0 .and. rxn % MT <= N_2NC .and. MTs % contains(N_2N)) cycle
|
||||
! Skip level cross sections if total is available
|
||||
if (rxn % MT >= N_P0 .and. rxn % MT <= N_PC .and. MTs % contains(N_P)) cycle
|
||||
if (rxn % MT >= N_D0 .and. rxn % MT <= N_DC .and. MTs % contains(N_D)) cycle
|
||||
if (rxn % MT >= N_T0 .and. rxn % MT <= N_TC .and. MTs % contains(N_T)) cycle
|
||||
if (rxn % MT >= N_3HE0 .and. rxn % MT <= N_3HEC .and. MTs % contains(N_3HE)) cycle
|
||||
if (rxn % MT >= N_A0 .and. rxn % MT <= N_AC .and. MTs % contains(N_A)) cycle
|
||||
if (rxn % MT >= N_2N0 .and. rxn % MT <= N_2NC .and. MTs % contains(N_2N)) cycle
|
||||
|
||||
! Add contribution to total cross section
|
||||
nuc % total(IE:IE+NE-1) = nuc % total(IE:IE+NE-1) + rxn % sigma
|
||||
! Add contribution to total cross section
|
||||
nuc % total(IE:IE+NE-1) = nuc % total(IE:IE+NE-1) + rxn % sigma
|
||||
|
||||
! Add contribution to absorption cross section
|
||||
if (is_disappearance(rxn % MT)) then
|
||||
nuc % absorption(IE:IE+NE-1) = nuc % absorption(IE:IE+NE-1) + rxn % sigma
|
||||
end if
|
||||
! Add contribution to absorption cross section
|
||||
if (is_disappearance(rxn % MT)) then
|
||||
nuc % absorption(IE:IE+NE-1) = nuc % absorption(IE:IE+NE-1) + rxn % sigma
|
||||
end if
|
||||
|
||||
! Information about fission reactions
|
||||
if (rxn % MT == N_FISSION) then
|
||||
allocate(nuc % index_fission(1))
|
||||
elseif (rxn % MT == N_F) then
|
||||
allocate(nuc % index_fission(PARTIAL_FISSION_MAX))
|
||||
nuc % has_partial_fission = .true.
|
||||
end if
|
||||
! Information about fission reactions
|
||||
if (rxn % MT == N_FISSION) then
|
||||
allocate(nuc % index_fission(1))
|
||||
elseif (rxn % MT == N_F) then
|
||||
allocate(nuc % index_fission(PARTIAL_FISSION_MAX))
|
||||
nuc % has_partial_fission = .true.
|
||||
end if
|
||||
|
||||
! Add contribution to fission cross section
|
||||
if (is_fission(rxn % MT)) then
|
||||
nuc % fissionable = .true.
|
||||
nuc % fission(IE:IE+NE-1) = nuc % fission(IE:IE+NE-1) + rxn % sigma
|
||||
! Add contribution to fission cross section
|
||||
if (is_fission(rxn % MT)) then
|
||||
nuc % fissionable = .true.
|
||||
nuc % fission(IE:IE+NE-1) = nuc % fission(IE:IE+NE-1) + rxn % sigma
|
||||
|
||||
! Also need to add fission cross sections to absorption
|
||||
nuc % absorption(IE:IE+NE-1) = nuc % absorption(IE:IE+NE-1) + rxn % sigma
|
||||
! Also need to add fission cross sections to absorption
|
||||
nuc % absorption(IE:IE+NE-1) = nuc % absorption(IE:IE+NE-1) + rxn % sigma
|
||||
|
||||
! If total fission reaction is present, there's no need to store the
|
||||
! reaction cross-section since it was copied to nuc % fission
|
||||
if (rxn % MT == N_FISSION) deallocate(rxn % sigma)
|
||||
! If total fission reaction is present, there's no need to store the
|
||||
! reaction cross-section since it was copied to nuc % fission
|
||||
if (rxn % MT == N_FISSION) deallocate(rxn % sigma)
|
||||
|
||||
! Keep track of this reaction for easy searching later
|
||||
i_fission = i_fission + 1
|
||||
nuc % index_fission(i_fission) = i
|
||||
nuc % n_fission = nuc % n_fission + 1
|
||||
end if
|
||||
! Keep track of this reaction for easy searching later
|
||||
i_fission = i_fission + 1
|
||||
nuc % index_fission(i_fission) = i
|
||||
nuc % n_fission = nuc % n_fission + 1
|
||||
end if
|
||||
end associate
|
||||
end do
|
||||
|
||||
! Clear MTs set
|
||||
|
|
@ -895,7 +896,6 @@ contains
|
|||
integer :: i ! index in reactions array
|
||||
integer :: j ! index over incoming energies
|
||||
integer :: length ! length of data array to allocate
|
||||
type(Reaction), pointer :: rxn
|
||||
|
||||
JXS8 = JXS(8)
|
||||
JXS9 = JXS(9)
|
||||
|
|
@ -903,71 +903,72 @@ contains
|
|||
! loop over all reactions with secondary neutrons -- NXS(5) does not include
|
||||
! elastic scattering
|
||||
do i = 1, NXS(5) + 1
|
||||
rxn => nuc%reactions(i)
|
||||
associate (rxn => nuc%reactions(i))
|
||||
|
||||
! find location of angular distribution
|
||||
LOCB = int(XSS(JXS8 + i - 1))
|
||||
if (LOCB == -1) then
|
||||
! Angular distribution data are specified through LAWi = 44 in the DLW
|
||||
! block
|
||||
cycle
|
||||
elseif (LOCB == 0) then
|
||||
! No angular distribution data are given for this reaction, isotropic
|
||||
! scattering is asssumed (in CM if TY < 0 and in LAB if TY > 0)
|
||||
cycle
|
||||
end if
|
||||
rxn % has_angle_dist = .true.
|
||||
|
||||
! allocate space for incoming energies and locations
|
||||
NE = int(XSS(JXS9 + LOCB - 1))
|
||||
rxn % adist % n_energy = NE
|
||||
allocate(rxn % adist % energy(NE))
|
||||
allocate(rxn % adist % type(NE))
|
||||
allocate(rxn % adist % location(NE))
|
||||
|
||||
! read incoming energy grid and location of nucs
|
||||
XSS_index = JXS9 + LOCB
|
||||
rxn % adist % energy = get_real(NE)
|
||||
rxn % adist % location = get_int(NE)
|
||||
|
||||
! determine dize of data block
|
||||
length = 0
|
||||
do j = 1, NE
|
||||
LC = rxn % adist % location(j)
|
||||
if (LC == 0) then
|
||||
! isotropic
|
||||
rxn % adist % type(j) = ANGLE_ISOTROPIC
|
||||
elseif (LC > 0) then
|
||||
! 32 equiprobable bins
|
||||
rxn % adist % type(j) = ANGLE_32_EQUI
|
||||
length = length + 33
|
||||
elseif (LC < 0) then
|
||||
! tabular distribution
|
||||
rxn % adist % type(j) = ANGLE_TABULAR
|
||||
NP = int(XSS(JXS9 + abs(LC)))
|
||||
length = length + 2 + 3*NP
|
||||
! find location of angular distribution
|
||||
LOCB = int(XSS(JXS8 + i - 1))
|
||||
if (LOCB == -1) then
|
||||
! Angular distribution data are specified through LAWi = 44 in the DLW
|
||||
! block
|
||||
cycle
|
||||
elseif (LOCB == 0) then
|
||||
! No angular distribution data are given for this reaction, isotropic
|
||||
! scattering is asssumed (in CM if TY < 0 and in LAB if TY > 0)
|
||||
cycle
|
||||
end if
|
||||
end do
|
||||
rxn % has_angle_dist = .true.
|
||||
|
||||
! allocate angular distribution data and read
|
||||
allocate(rxn % adist % data(length))
|
||||
! allocate space for incoming energies and locations
|
||||
NE = int(XSS(JXS9 + LOCB - 1))
|
||||
rxn % adist % n_energy = NE
|
||||
allocate(rxn % adist % energy(NE))
|
||||
allocate(rxn % adist % type(NE))
|
||||
allocate(rxn % adist % location(NE))
|
||||
|
||||
! read angular distribution -- currently this does not actually parse the
|
||||
! angular distribution tables for each incoming energy, that must be done
|
||||
! on-the-fly
|
||||
XSS_index = JXS9 + LOCB + 2 * NE
|
||||
rxn % adist % data = get_real(length)
|
||||
! read incoming energy grid and location of nucs
|
||||
XSS_index = JXS9 + LOCB
|
||||
rxn % adist % energy = get_real(NE)
|
||||
rxn % adist % location = get_int(NE)
|
||||
|
||||
! change location pointers since they are currently relative to JXS(9)
|
||||
LC = LOCB + 2 * NE + 1
|
||||
do j = 1, NE
|
||||
! For consistency, leave location as 0 if type is isotropic.
|
||||
! This is not necessary for current correctness, but can avoid
|
||||
! future issues
|
||||
if (rxn % adist % location(j) /= 0) then
|
||||
rxn % adist % location(j) = abs(rxn % adist % location(j)) - LC
|
||||
end if
|
||||
end do
|
||||
! determine dize of data block
|
||||
length = 0
|
||||
do j = 1, NE
|
||||
LC = rxn % adist % location(j)
|
||||
if (LC == 0) then
|
||||
! isotropic
|
||||
rxn % adist % type(j) = ANGLE_ISOTROPIC
|
||||
elseif (LC > 0) then
|
||||
! 32 equiprobable bins
|
||||
rxn % adist % type(j) = ANGLE_32_EQUI
|
||||
length = length + 33
|
||||
elseif (LC < 0) then
|
||||
! tabular distribution
|
||||
rxn % adist % type(j) = ANGLE_TABULAR
|
||||
NP = int(XSS(JXS9 + abs(LC)))
|
||||
length = length + 2 + 3*NP
|
||||
end if
|
||||
end do
|
||||
|
||||
! allocate angular distribution data and read
|
||||
allocate(rxn % adist % data(length))
|
||||
|
||||
! read angular distribution -- currently this does not actually parse the
|
||||
! angular distribution tables for each incoming energy, that must be done
|
||||
! on-the-fly
|
||||
XSS_index = JXS9 + LOCB + 2 * NE
|
||||
rxn % adist % data = get_real(length)
|
||||
|
||||
! change location pointers since they are currently relative to JXS(9)
|
||||
LC = LOCB + 2 * NE + 1
|
||||
do j = 1, NE
|
||||
! For consistency, leave location as 0 if type is isotropic.
|
||||
! This is not necessary for current correctness, but can avoid
|
||||
! future issues
|
||||
if (rxn % adist % location(j) /= 0) then
|
||||
rxn % adist % location(j) = abs(rxn % adist % location(j)) - LC
|
||||
end if
|
||||
end do
|
||||
end associate
|
||||
end do
|
||||
|
||||
end subroutine read_angular_dist
|
||||
|
|
@ -983,23 +984,23 @@ contains
|
|||
integer :: LED ! location of energy distribution locators
|
||||
integer :: LOCC ! location of energy distributions for given MT
|
||||
integer :: i ! loop index
|
||||
type(Reaction), pointer :: rxn
|
||||
|
||||
LED = JXS(10)
|
||||
|
||||
! Loop over all reactions
|
||||
do i = 1, NXS(5)
|
||||
rxn => nuc % reactions(i+1) ! skip over elastic scattering
|
||||
rxn % has_energy_dist = .true.
|
||||
associate (rxn => nuc % reactions(i+1)) ! skip over elastic scattering
|
||||
rxn % has_energy_dist = .true.
|
||||
|
||||
! find location of energy distribution data
|
||||
LOCC = int(XSS(LED + i - 1))
|
||||
! find location of energy distribution data
|
||||
LOCC = int(XSS(LED + i - 1))
|
||||
|
||||
! allocate energy distribution
|
||||
allocate(rxn % edist)
|
||||
! allocate energy distribution
|
||||
allocate(rxn % edist)
|
||||
|
||||
! read data for energy distribution
|
||||
call get_energy_dist(rxn % edist, LOCC)
|
||||
! read data for energy distribution
|
||||
call get_energy_dist(rxn % edist, LOCC)
|
||||
end associate
|
||||
end do
|
||||
|
||||
end subroutine read_energy_dist
|
||||
|
|
|
|||
|
|
@ -146,7 +146,7 @@ module ace_header
|
|||
|
||||
! Reactions
|
||||
integer :: n_reaction ! # of reactions
|
||||
type(Reaction), pointer :: reactions(:) => null()
|
||||
type(Reaction), allocatable :: reactions(:)
|
||||
type(DictIntInt) :: reaction_index ! map MT values to index in reactions
|
||||
! array; used at tally-time
|
||||
|
||||
|
|
@ -341,11 +341,10 @@ module ace_header
|
|||
end if
|
||||
|
||||
if (associated(this % urr_data)) then
|
||||
call this % urr_data % clear()
|
||||
deallocate(this % urr_data)
|
||||
end if
|
||||
|
||||
if (associated(this % reactions)) then
|
||||
if (allocated(this % reactions)) then
|
||||
do i = 1, size(this % reactions)
|
||||
call this % reactions(i) % clear()
|
||||
end do
|
||||
|
|
|
|||
|
|
@ -379,7 +379,6 @@ contains
|
|||
logical :: same_nuc ! do we know the xs for this nuclide at this energy?
|
||||
type(UrrData), pointer :: urr
|
||||
type(Nuclide), pointer :: nuc
|
||||
type(Reaction), pointer :: rxn
|
||||
|
||||
micro_xs(i_nuclide) % use_ptable = .true.
|
||||
|
||||
|
|
@ -475,18 +474,17 @@ contains
|
|||
! Determine treatment of inelastic scattering
|
||||
inelastic = ZERO
|
||||
if (urr % inelastic_flag > 0) then
|
||||
! Get pointer to inelastic scattering reaction
|
||||
rxn => nuc % reactions(nuc % urr_inelastic)
|
||||
|
||||
! Get index on energy grid and interpolation factor
|
||||
i_energy = micro_xs(i_nuclide) % index_grid
|
||||
f = micro_xs(i_nuclide) % interp_factor
|
||||
|
||||
! Determine inelastic scattering cross section
|
||||
if (i_energy >= rxn % threshold) then
|
||||
inelastic = (ONE - f) * rxn % sigma(i_energy - rxn%threshold + 1) + &
|
||||
f * rxn % sigma(i_energy - rxn%threshold + 2)
|
||||
end if
|
||||
associate (rxn => nuc % reactions(nuc % urr_inelastic))
|
||||
if (i_energy >= rxn % threshold) then
|
||||
inelastic = (ONE - f) * rxn % sigma(i_energy - rxn%threshold + 1) + &
|
||||
f * rxn % sigma(i_energy - rxn%threshold + 2)
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
! Multiply by smooth cross-section if needed
|
||||
|
|
|
|||
|
|
@ -330,7 +330,6 @@ contains
|
|||
integer :: size_energy ! memory used for a energy distributions (bytes)
|
||||
integer :: size_urr ! memory used for probability tables (bytes)
|
||||
character(11) :: law ! secondary energy distribution law
|
||||
type(Reaction), pointer :: rxn
|
||||
type(UrrData), pointer :: urr
|
||||
|
||||
! set default unit for writing information
|
||||
|
|
@ -359,32 +358,32 @@ contains
|
|||
! Information on each reaction
|
||||
write(unit_,*) ' Reaction Q-value COM Law IE size(angle) size(energy)'
|
||||
do i = 1, nuc % n_reaction
|
||||
rxn => nuc % reactions(i)
|
||||
associate (rxn => nuc % reactions(i))
|
||||
! Determine size of angle distribution
|
||||
if (rxn % has_angle_dist) then
|
||||
size_angle = rxn % adist % n_energy * 16 + size(rxn % adist % data) * 8
|
||||
else
|
||||
size_angle = 0
|
||||
end if
|
||||
|
||||
! Determine size of angle distribution
|
||||
if (rxn % has_angle_dist) then
|
||||
size_angle = rxn % adist % n_energy * 16 + size(rxn % adist % data) * 8
|
||||
else
|
||||
size_angle = 0
|
||||
end if
|
||||
! Determine size of energy distribution and law
|
||||
if (rxn % has_energy_dist) then
|
||||
size_energy = size(rxn % edist % data) * 8
|
||||
law = to_str(rxn % edist % law)
|
||||
else
|
||||
size_energy = 0
|
||||
law = 'None'
|
||||
end if
|
||||
|
||||
! Determine size of energy distribution and law
|
||||
if (rxn % has_energy_dist) then
|
||||
size_energy = size(rxn % edist % data) * 8
|
||||
law = to_str(rxn % edist % law)
|
||||
else
|
||||
size_energy = 0
|
||||
law = 'None'
|
||||
end if
|
||||
write(unit_,'(3X,A11,1X,F8.3,3X,L1,3X,A4,1X,I6,1X,I11,1X,I11)') &
|
||||
reaction_name(rxn % MT), rxn % Q_value, rxn % scatter_in_cm, &
|
||||
law(1:4), rxn % threshold, size_angle, size_energy
|
||||
|
||||
write(unit_,'(3X,A11,1X,F8.3,3X,L1,3X,A4,1X,I6,1X,I11,1X,I11)') &
|
||||
reaction_name(rxn % MT), rxn % Q_value, rxn % scatter_in_cm, &
|
||||
law(1:4), rxn % threshold, size_angle, size_energy
|
||||
|
||||
! Accumulate data size
|
||||
size_xs = size_xs + (nuc % n_grid - rxn%threshold + 1) * 8
|
||||
size_angle_total = size_angle_total + size_angle
|
||||
size_energy_total = size_energy_total + size_energy
|
||||
! Accumulate data size
|
||||
size_xs = size_xs + (nuc % n_grid - rxn%threshold + 1) * 8
|
||||
size_angle_total = size_angle_total + size_angle
|
||||
size_energy_total = size_energy_total + size_energy
|
||||
end associate
|
||||
end do
|
||||
|
||||
! Add memory required for summary reactions (total, absorption, fission,
|
||||
|
|
|
|||
|
|
@ -185,7 +185,6 @@ contains
|
|||
!===============================================================================
|
||||
|
||||
subroutine sample_fission(i_nuclide, i_reaction)
|
||||
|
||||
integer, intent(in) :: i_nuclide ! index in nuclides array
|
||||
integer, intent(out) :: i_reaction ! index in nuc % reactions array
|
||||
|
||||
|
|
@ -195,7 +194,6 @@ contains
|
|||
real(8) :: prob
|
||||
real(8) :: cutoff
|
||||
type(Nuclide), pointer :: nuc
|
||||
type(Reaction), pointer :: rxn
|
||||
|
||||
! Get pointer to nuclide
|
||||
nuc => nuclides(i_nuclide)
|
||||
|
|
@ -220,14 +218,15 @@ contains
|
|||
|
||||
FISSION_REACTION_LOOP: do i = 1, nuc % n_fission
|
||||
i_reaction = nuc % index_fission(i)
|
||||
rxn => nuc % reactions(i_reaction)
|
||||
|
||||
! if energy is below threshold for this reaction, skip it
|
||||
if (i_grid < rxn % threshold) cycle
|
||||
associate (rxn => nuc % reactions(i_reaction))
|
||||
! if energy is below threshold for this reaction, skip it
|
||||
if (i_grid < rxn % threshold) cycle
|
||||
|
||||
! add to cumulative probability
|
||||
prob = prob + ((ONE - f)*rxn%sigma(i_grid - rxn%threshold + 1) &
|
||||
+ f*(rxn%sigma(i_grid - rxn%threshold + 2)))
|
||||
! add to cumulative probability
|
||||
prob = prob + ((ONE - f)*rxn%sigma(i_grid - rxn%threshold + 1) &
|
||||
+ f*(rxn%sigma(i_grid - rxn%threshold + 2)))
|
||||
end associate
|
||||
|
||||
! Create fission bank sites if fission occurs
|
||||
if (prob > cutoff) exit FISSION_REACTION_LOOP
|
||||
|
|
@ -312,11 +311,10 @@ contains
|
|||
real(8) :: f
|
||||
real(8) :: prob
|
||||
real(8) :: cutoff
|
||||
type(Nuclide), pointer :: nuc
|
||||
type(Reaction), pointer :: rxn
|
||||
real(8) :: uvw_new(3) ! outgoing uvw for iso-in-lab scattering
|
||||
real(8) :: uvw_old(3) ! incoming uvw for iso-in-lab scattering
|
||||
real(8) :: phi ! azimuthal angle for iso-in-lab scattering
|
||||
type(Nuclide), pointer :: nuc
|
||||
|
||||
! copy incoming direction
|
||||
uvw_old(:) = p % coord(1) % uvw
|
||||
|
|
@ -343,11 +341,8 @@ contains
|
|||
p % E, p % coord(1) % uvw, p % mu)
|
||||
|
||||
else
|
||||
! get pointer to elastic scattering reaction
|
||||
rxn => nuc % reactions(1)
|
||||
|
||||
! Perform collision physics for elastic scattering
|
||||
call elastic_scatter(i_nuclide, rxn, &
|
||||
call elastic_scatter(i_nuclide, nuc % reactions(1), &
|
||||
p % E, p % coord(1) % uvw, p % mu, p % wgt)
|
||||
end if
|
||||
|
||||
|
|
@ -370,28 +365,28 @@ contains
|
|||
&// trim(nuc % name))
|
||||
end if
|
||||
|
||||
rxn => nuc % reactions(i)
|
||||
associate (rxn => nuc % reactions(i))
|
||||
! Skip fission reactions
|
||||
if (rxn % MT == N_FISSION .or. rxn % MT == N_F .or. rxn % MT == N_NF &
|
||||
.or. rxn % MT == N_2NF .or. rxn % MT == N_3NF) cycle
|
||||
|
||||
! Skip fission reactions
|
||||
if (rxn % MT == N_FISSION .or. rxn % MT == N_F .or. rxn % MT == N_NF &
|
||||
.or. rxn % MT == N_2NF .or. rxn % MT == N_3NF) cycle
|
||||
! some materials have gas production cross sections with MT > 200 that
|
||||
! are duplicates. Also MT=4 is total level inelastic scattering which
|
||||
! should be skipped
|
||||
if (rxn % MT >= 200 .or. rxn % MT == N_LEVEL) cycle
|
||||
|
||||
! some materials have gas production cross sections with MT > 200 that
|
||||
! are duplicates. Also MT=4 is total level inelastic scattering which
|
||||
! should be skipped
|
||||
if (rxn % MT >= 200 .or. rxn % MT == N_LEVEL) cycle
|
||||
! if energy is below threshold for this reaction, skip it
|
||||
if (i_grid < rxn % threshold) cycle
|
||||
|
||||
! if energy is below threshold for this reaction, skip it
|
||||
if (i_grid < rxn % threshold) cycle
|
||||
|
||||
! add to cumulative probability
|
||||
prob = prob + ((ONE - f)*rxn%sigma(i_grid - rxn%threshold + 1) &
|
||||
+ f*(rxn%sigma(i_grid - rxn%threshold + 2)))
|
||||
! add to cumulative probability
|
||||
prob = prob + ((ONE - f)*rxn%sigma(i_grid - rxn%threshold + 1) &
|
||||
+ f*(rxn%sigma(i_grid - rxn%threshold + 2)))
|
||||
end associate
|
||||
end do
|
||||
|
||||
! Perform collision physics for inelastic scattering
|
||||
call inelastic_scatter(nuc, rxn, p)
|
||||
p % event_MT = rxn % MT
|
||||
call inelastic_scatter(nuc, nuc%reactions(i), p)
|
||||
p % event_MT = nuc%reactions(i)%MT
|
||||
|
||||
end if
|
||||
|
||||
|
|
@ -1090,11 +1085,9 @@ contains
|
|||
real(8) :: weight ! weight adjustment for ufs method
|
||||
logical :: in_mesh ! source site in ufs mesh?
|
||||
type(Nuclide), pointer :: nuc
|
||||
type(Reaction), pointer :: rxn
|
||||
|
||||
! Get pointers
|
||||
nuc => nuclides(i_nuclide)
|
||||
rxn => nuc % reactions(i_reaction)
|
||||
|
||||
! TODO: Heat generation from fission
|
||||
|
||||
|
|
@ -1165,7 +1158,8 @@ contains
|
|||
|
||||
! Sample secondary energy distribution for fission reaction and set energy
|
||||
! in fission bank
|
||||
fission_bank(i) % E = sample_fission_energy(nuc, rxn, p)
|
||||
fission_bank(i) % E = sample_fission_energy(nuc, nuc%reactions(&
|
||||
i_reaction), p)
|
||||
|
||||
! Set the delayed group of the neutron
|
||||
fission_bank(i) % delayed_group = p % delayed_group
|
||||
|
|
|
|||
197
src/tally.F90
197
src/tally.F90
|
|
@ -66,9 +66,6 @@ contains
|
|||
real(8) :: macro_scatt ! material macro scatt xs
|
||||
real(8) :: uvw(3) ! particle direction
|
||||
real(8) :: E ! particle energy
|
||||
type(Material), pointer :: mat
|
||||
type(Reaction), pointer :: rxn
|
||||
type(Nuclide), pointer :: nuc
|
||||
|
||||
i = 0
|
||||
SCORE_LOOP: do q = 1, t % n_user_score_bins
|
||||
|
|
@ -220,24 +217,20 @@ contains
|
|||
! of one.
|
||||
score = p % last_wgt
|
||||
else
|
||||
do m = 1, nuclides(p % event_nuclide) % n_reaction
|
||||
! Check if this is the desired MT
|
||||
if (p % event_MT == nuclides(p % event_nuclide) % reactions(m) % MT) then
|
||||
! Found the reaction, set our pointer and move on with life
|
||||
rxn => nuclides(p % event_nuclide) % reactions(m)
|
||||
exit
|
||||
end if
|
||||
end do
|
||||
m = nuclides(p%event_nuclide)%reaction_index% &
|
||||
get_key(p % event_MT)
|
||||
|
||||
! Get multiplicity and apply to score
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
associate (rxn => nuclides(p%event_nuclide)%reactions(m))
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
|
||||
|
|
@ -257,24 +250,20 @@ contains
|
|||
! of one.
|
||||
score = p % last_wgt
|
||||
else
|
||||
do m = 1, nuclides(p % event_nuclide) % n_reaction
|
||||
! Check if this is the desired MT
|
||||
if (p % event_MT == nuclides(p % event_nuclide) % reactions(m) % MT) then
|
||||
! Found the reaction, set our pointer and move on with life
|
||||
rxn => nuclides(p % event_nuclide) % reactions(m)
|
||||
exit
|
||||
end if
|
||||
end do
|
||||
m = nuclides(p%event_nuclide)%reaction_index% &
|
||||
get_key(p % event_MT)
|
||||
|
||||
! Get multiplicity and apply to score
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
associate (rxn => nuclides(p%event_nuclide)%reactions(m))
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
|
||||
|
|
@ -294,24 +283,20 @@ contains
|
|||
! of one.
|
||||
score = p % last_wgt
|
||||
else
|
||||
do m = 1, nuclides(p % event_nuclide) % n_reaction
|
||||
! Check if this is the desired MT
|
||||
if (p % event_MT == nuclides(p % event_nuclide) % reactions(m) % MT) then
|
||||
! Found the reaction, set our pointer and move on with life
|
||||
rxn => nuclides(p % event_nuclide) % reactions(m)
|
||||
exit
|
||||
end if
|
||||
end do
|
||||
m = nuclides(p%event_nuclide)%reaction_index% &
|
||||
get_key(p % event_MT)
|
||||
|
||||
! Get multiplicity and apply to score
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
associate (rxn => nuclides(p%event_nuclide)%reactions(m))
|
||||
if (rxn % multiplicity_with_E) then
|
||||
! Then the multiplicity was already incorporated in to p % wgt
|
||||
! per the scattering routine,
|
||||
score = p % wgt
|
||||
else
|
||||
! Grab the multiplicity from the rxn
|
||||
score = p % last_wgt * rxn % multiplicity
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
|
||||
|
|
@ -466,9 +451,6 @@ contains
|
|||
! delayed-nu-fission
|
||||
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
|
||||
|
||||
! Get the event nuclide
|
||||
nuc => nuclides(p % event_nuclide)
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
|
||||
|
|
@ -480,11 +462,11 @@ contains
|
|||
d = t % filters(dg_filter) % int_bins(d_bin)
|
||||
|
||||
! Compute the yield for this delayed group
|
||||
yield = yield_delayed(nuc, E, d)
|
||||
yield = yield_delayed(nuclides(p % event_nuclide), E, d)
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = p % absorb_wgt * yield * micro_xs(p % event_nuclide) &
|
||||
% fission * nu_delayed(nuc, E) / &
|
||||
% fission * nu_delayed(nuclides(p % event_nuclide), E) / &
|
||||
micro_xs(p % event_nuclide) % absorption
|
||||
call score_fission_delayed_dg(t, d_bin, score, score_index)
|
||||
end do
|
||||
|
|
@ -494,7 +476,7 @@ contains
|
|||
! by multiplying the absorbed weight by the fraction of the
|
||||
! delayed-nu-fission xs to the absorption xs
|
||||
score = p % absorb_wgt * micro_xs(p % event_nuclide) &
|
||||
% fission * nu_delayed(nuc, E) / &
|
||||
% fission * nu_delayed(nuclides(p % event_nuclide), E) / &
|
||||
micro_xs(p % event_nuclide) % absorption
|
||||
end if
|
||||
end if
|
||||
|
|
@ -535,9 +517,6 @@ contains
|
|||
! Check if tally is on a single nuclide
|
||||
if (i_nuclide > 0) then
|
||||
|
||||
! Get the nuclide of interest
|
||||
nuc => nuclides(i_nuclide)
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
|
||||
|
|
@ -548,11 +527,19 @@ contains
|
|||
d = t % filters(dg_filter) % int_bins(d_bin)
|
||||
|
||||
! Compute the yield for this delayed group
|
||||
<<<<<<< HEAD
|
||||
yield = yield_delayed(nuc, E, d)
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = micro_xs(i_nuclide) % fission * yield &
|
||||
* nu_delayed(nuc, E) * atom_density * flux
|
||||
=======
|
||||
yield = yield_delayed(nuclides(i_nuclide), p % E, d)
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = micro_xs(i_nuclide) % fission * yield &
|
||||
* nu_delayed(nuclides(i_nuclide), p % E) * atom_density * flux
|
||||
>>>>>>> Make Nuclide%reactions allocatable by using associate constructs
|
||||
call score_fission_delayed_dg(t, d_bin, score, score_index)
|
||||
end do
|
||||
cycle SCORE_LOOP
|
||||
|
|
@ -560,27 +547,29 @@ contains
|
|||
|
||||
! If the delayed group filter is not present, compute the score
|
||||
! by multiplying the delayed-nu-fission macro xs by the flux
|
||||
<<<<<<< HEAD
|
||||
score = micro_xs(i_nuclide) % fission * nu_delayed(nuc, E)&
|
||||
* atom_density * flux
|
||||
=======
|
||||
score = micro_xs(i_nuclide) % fission * &
|
||||
nu_delayed(nuclides(i_nuclide), p % E) * atom_density * flux
|
||||
>>>>>>> Make Nuclide%reactions allocatable by using associate constructs
|
||||
end if
|
||||
|
||||
! Tally is on total nuclides
|
||||
else
|
||||
|
||||
! Get pointer to current material
|
||||
mat => materials(p % material)
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
|
||||
! Loop over all nuclides in the current material
|
||||
do l = 1, mat % n_nuclides
|
||||
do l = 1, materials(p % material) % n_nuclides
|
||||
|
||||
! Get atom density
|
||||
atom_density_ = mat % atom_density(l)
|
||||
atom_density_ = materials(p % material) % atom_density(l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = mat % nuclide(l)
|
||||
i_nuc = materials(p % material) % nuclide(l)
|
||||
|
||||
! Loop over all delayed group bins and tally to them individually
|
||||
do d_bin = 1, t % filters(dg_filter) % n_bins
|
||||
|
|
@ -588,15 +577,20 @@ contains
|
|||
! Get the delayed group for this bin
|
||||
d = t % filters(dg_filter) % int_bins(d_bin)
|
||||
|
||||
! Get the current nuclide
|
||||
nuc => nuclides(i_nuc)
|
||||
|
||||
! Get the yield for the desired nuclide and delayed group
|
||||
<<<<<<< HEAD
|
||||
yield = yield_delayed(nuc, E, d)
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = micro_xs(i_nuc) % fission * yield &
|
||||
* nu_delayed(nuc, E) * atom_density_ * flux
|
||||
=======
|
||||
yield = yield_delayed(nuclides(i_nuc), p % E, d)
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = micro_xs(i_nuc) % fission * yield &
|
||||
* nu_delayed(nuclides(i_nuc), p % E) * atom_density_ * flux
|
||||
>>>>>>> Make Nuclide%reactions allocatable by using associate constructs
|
||||
call score_fission_delayed_dg(t, d_bin, score, score_index)
|
||||
end do
|
||||
end do
|
||||
|
|
@ -606,13 +600,13 @@ contains
|
|||
score = ZERO
|
||||
|
||||
! Loop over all nuclides in the current material
|
||||
do l = 1, mat % n_nuclides
|
||||
do l = 1, materials(p % material) % n_nuclides
|
||||
|
||||
! Get atom density
|
||||
atom_density_ = mat % atom_density(l)
|
||||
atom_density_ = materials(p % material) % atom_density(l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = mat % nuclide(l)
|
||||
i_nuc = materials(p % material) % nuclide(l)
|
||||
|
||||
! Accumulate the contribution from each nuclide
|
||||
score = score + micro_xs(i_nuc) % fission &
|
||||
|
|
@ -693,42 +687,41 @@ contains
|
|||
if (i_nuclide > 0) then
|
||||
if (nuclides(i_nuclide)%reaction_index%has_key(score_bin)) then
|
||||
m = nuclides(i_nuclide)%reaction_index%get_key(score_bin)
|
||||
rxn => nuclides(i_nuclide) % reactions(m)
|
||||
|
||||
! Retrieve index on nuclide energy grid and interpolation
|
||||
! factor
|
||||
i_energy = micro_xs(i_nuclide) % index_grid
|
||||
f = micro_xs(i_nuclide) % interp_factor
|
||||
if (i_energy >= rxn % threshold) then
|
||||
score = ((ONE - f) * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 1) + f * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 2)) * atom_density * flux
|
||||
end if
|
||||
end if
|
||||
|
||||
else
|
||||
! Get pointer to current material
|
||||
mat => materials(p % material)
|
||||
do l = 1, mat % n_nuclides
|
||||
! Get atom density
|
||||
atom_density_ = mat % atom_density(l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = mat % nuclide(l)
|
||||
|
||||
if (nuclides(i_nuc)%reaction_index%has_key(score_bin)) then
|
||||
m = nuclides(i_nuc)%reaction_index%get_key(score_bin)
|
||||
rxn => nuclides(i_nuc) % reactions(m)
|
||||
associate (rxn => nuclides(i_nuclide) % reactions(m))
|
||||
|
||||
! Retrieve index on nuclide energy grid and interpolation
|
||||
! factor
|
||||
i_energy = micro_xs(i_nuc) % index_grid
|
||||
f = micro_xs(i_nuc) % interp_factor
|
||||
i_energy = micro_xs(i_nuclide) % index_grid
|
||||
f = micro_xs(i_nuclide) % interp_factor
|
||||
if (i_energy >= rxn % threshold) then
|
||||
score = score + ((ONE - f) * rxn % sigma(i_energy - &
|
||||
score = ((ONE - f) * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 1) + f * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 2)) * atom_density_ * flux
|
||||
rxn%threshold + 2)) * atom_density * flux
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
else
|
||||
do l = 1, materials(p % material) % n_nuclides
|
||||
! Get atom density
|
||||
atom_density_ = materials(p % material) % atom_density(l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = materials(p % material) % nuclide(l)
|
||||
|
||||
if (nuclides(i_nuc)%reaction_index%has_key(score_bin)) then
|
||||
m = nuclides(i_nuc)%reaction_index%get_key(score_bin)
|
||||
associate (rxn => nuclides(i_nuc) % reactions(m))
|
||||
! Retrieve index on nuclide energy grid and interpolation
|
||||
! factor
|
||||
i_energy = micro_xs(i_nuc) % index_grid
|
||||
f = micro_xs(i_nuc) % interp_factor
|
||||
if (i_energy >= rxn % threshold) then
|
||||
score = score + ((ONE - f) * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 1) + f * rxn % sigma(i_energy - &
|
||||
rxn%threshold + 2)) * atom_density_ * flux
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
end do
|
||||
end if
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue