mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-28 22:26:08 -04:00
Finish moving score_general_ce to C++
This commit is contained in:
Sterling Harper
parent
355ce75f0a
commit
1c1fa6d753
3 changed files with 462 additions and 704 deletions
|
|
@ -274,7 +274,7 @@ void Nuclide::create_derived()
|
|||
xs_.emplace_back(shape, 0.0);
|
||||
}
|
||||
|
||||
reaction_index_.fill(-1);
|
||||
reaction_index_.fill(C_NONE);
|
||||
for (int i = 0; i < reactions_.size(); ++i) {
|
||||
const auto& rx {reactions_[i]};
|
||||
|
||||
|
|
|
|||
|
|
@ -47,16 +47,16 @@ module tally
|
|||
end interface
|
||||
|
||||
interface
|
||||
subroutine score_general_ce_c(p, i_tally, start_index, filter_index, &
|
||||
subroutine score_general_ce(p, i_tally, start_index, filter_index, &
|
||||
i_nuclide, atom_density, flux) bind(C)
|
||||
import Particle, C_INT, C_DOUBLE
|
||||
type(Particle) :: p
|
||||
integer(C_INT), value :: i_tally
|
||||
integer(C_INT), value :: start_index
|
||||
integer(C_INT), value :: filter_index
|
||||
integer(C_INT), value :: i_nuclide
|
||||
real(C_DOUBLE), value :: atom_density
|
||||
real(C_DOUBLE), value :: flux
|
||||
type(Particle), intent(in) :: p
|
||||
integer(C_INT), intent(in), value :: i_tally
|
||||
integer(C_INT), intent(in), value :: start_index
|
||||
integer(C_INT), intent(in), value :: filter_index
|
||||
integer(C_INT), intent(in), value :: i_nuclide
|
||||
real(C_DOUBLE), intent(in), value :: atom_density
|
||||
real(C_DOUBLE), intent(in), value :: flux
|
||||
end subroutine
|
||||
|
||||
subroutine score_fission_delayed_dg(i_tally, d_bin, score, score_index) bind(C)
|
||||
|
|
@ -135,699 +135,6 @@ contains
|
|||
! analog tallies.
|
||||
!===============================================================================
|
||||
|
||||
subroutine score_general_ce(p, i_tally, start_index, filter_index, i_nuclide, &
|
||||
atom_density, flux) bind(C)
|
||||
type(Particle), intent(in) :: p
|
||||
integer(C_INT), intent(in), value :: i_tally
|
||||
integer(C_INT), intent(in), value :: start_index
|
||||
integer(C_INT), intent(in), value :: i_nuclide
|
||||
integer(C_INT), intent(in), value :: filter_index ! for % results
|
||||
real(C_DOUBLE), intent(in), value :: flux ! flux estimate
|
||||
real(C_DOUBLE), intent(in), value :: atom_density ! atom/b-cm
|
||||
|
||||
integer :: i ! loop index for scoring bins
|
||||
integer :: l ! loop index for nuclides in material
|
||||
integer :: m ! loop index for reactions
|
||||
integer :: i_temp ! temperature index
|
||||
integer :: i_nuc ! index in nuclides array (from material)
|
||||
integer :: i_energy ! index in nuclide energy grid
|
||||
integer :: score_bin ! scoring bin, e.g. SCORE_FLUX
|
||||
integer :: score_index ! scoring bin index
|
||||
integer :: d ! delayed neutron index
|
||||
integer :: g ! delayed neutron index
|
||||
integer :: k ! loop index for bank sites
|
||||
integer :: d_bin ! delayed group bin index
|
||||
integer :: dg_filter ! index of delayed group filter
|
||||
integer :: threshold ! threshold energy index
|
||||
real(8) :: yield ! delayed neutron yield
|
||||
real(8) :: atom_density_ ! atom/b-cm
|
||||
real(8) :: f ! interpolation factor
|
||||
real(8) :: score ! analog tally score
|
||||
real(8) :: E ! particle energy
|
||||
real(8) :: xs ! cross section
|
||||
|
||||
call score_general_ce_c(p, i_tally, start_index, filter_index, &
|
||||
i_nuclide, atom_density, flux)
|
||||
|
||||
associate (t => tallies(i_tally) % obj)
|
||||
|
||||
! Pre-collision energy of particle
|
||||
E = p % last_E
|
||||
|
||||
SCORE_LOOP: do i = 1, t % n_score_bins()
|
||||
|
||||
! determine what type of score bin
|
||||
score_bin = t % score_bins(i)
|
||||
|
||||
! determine scoring bin index
|
||||
score_index = start_index + i
|
||||
|
||||
!#########################################################################
|
||||
! Determine appropirate scoring value.
|
||||
|
||||
select case(score_bin)
|
||||
|
||||
|
||||
case (SCORE_FLUX)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_TOTAL)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_INVERSE_VELOCITY)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_SCATTER)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_NU_SCATTER)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_ABSORPTION)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_FISSION)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_NU_FISSION)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
|
||||
case (SCORE_PROMPT_NU_FISSION)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
case (SCORE_DELAYED_NU_FISSION)
|
||||
cycle SCORE_LOOP
|
||||
|
||||
case (SCORE_DECAY_RATE)
|
||||
|
||||
if (material_xs % absorption == ZERO) cycle SCORE_LOOP
|
||||
|
||||
! Set the delayedgroup filter index
|
||||
dg_filter = t % delayedgroup_filter()
|
||||
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
if (survival_biasing) then
|
||||
! No fission events occur if survival biasing is on -- need to
|
||||
! calculate fraction of absorptions that would have resulted in
|
||||
! delayed-nu-fission
|
||||
if (micro_xs(p % event_nuclide) % absorption > ZERO .and. &
|
||||
nuclides(p % event_nuclide) % fissionable) then
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
select type(filt => filters(t % filter(dg_filter) + 1) % obj)
|
||||
type is (DelayedGroupFilter)
|
||||
|
||||
! Loop over all delayed group bins and tally to them
|
||||
! individually
|
||||
do d_bin = 1, filt % n_bins
|
||||
|
||||
! Get the delayed group for this bin
|
||||
d = filt % groups(d_bin)
|
||||
|
||||
! Compute the yield for this delayed group
|
||||
yield = nuclides(p % event_nuclide) &
|
||||
% nu(E, EMISSION_DELAYED, d)
|
||||
|
||||
associate (rxn => nuclides(p % event_nuclide) % &
|
||||
reactions(nuclides(p % event_nuclide) % &
|
||||
index_fission(1)))
|
||||
|
||||
! Compute the score
|
||||
score = p % absorb_wgt * yield * &
|
||||
micro_xs(p % event_nuclide) % fission &
|
||||
/ micro_xs(p % event_nuclide) % absorption &
|
||||
* rxn % product_decay_rate(1 + d) * flux
|
||||
end associate
|
||||
|
||||
! Tally to bin
|
||||
call score_fission_delayed_dg(i_tally, d_bin, score, score_index)
|
||||
end do
|
||||
cycle SCORE_LOOP
|
||||
end select
|
||||
else
|
||||
|
||||
! If the delayed group filter is not present, compute the score
|
||||
! by accumulating the absorbed weight times the decay rate times
|
||||
! the fraction of the delayed-nu-fission xs to the absorption xs
|
||||
! for all delayed groups.
|
||||
score = ZERO
|
||||
|
||||
associate (rxn => nuclides(p % event_nuclide) % &
|
||||
reactions(nuclides(p % event_nuclide) % index_fission(1)))
|
||||
|
||||
! We need to be careful not to overshoot the number of delayed
|
||||
! groups since this could cause the range of the
|
||||
! rxn % products array to be exceeded. Hence, we use the size
|
||||
! of this array and not the MAX_DELAYED_GROUPS constant for
|
||||
! this loop.
|
||||
do d = 1, rxn % products_size() - 2
|
||||
score = score + rxn % product_decay_rate(1 + d) * &
|
||||
p % absorb_wgt &
|
||||
* micro_xs(p % event_nuclide) % fission &
|
||||
* nuclides(p % event_nuclide) % &
|
||||
nu(E, EMISSION_DELAYED, d) &
|
||||
/ micro_xs(p % event_nuclide) % absorption * flux
|
||||
end do
|
||||
end associate
|
||||
end if
|
||||
end if
|
||||
else
|
||||
|
||||
! Skip any non-fission events
|
||||
if (.not. p % fission) cycle SCORE_LOOP
|
||||
! If there is no outgoing energy filter, than we only need to
|
||||
! score to one bin. For the score to be 'analog', we need to
|
||||
! score the number of particles that were banked in the fission
|
||||
! bank. Since this was weighted by 1/keff, we multiply by keff
|
||||
! to get the proper score. Loop over the neutrons produced from
|
||||
! fission and check which ones are delayed. If a delayed neutron is
|
||||
! encountered, add its contribution to the fission bank to the
|
||||
! score.
|
||||
|
||||
score = ZERO
|
||||
|
||||
! loop over number of particles banked
|
||||
do k = 1, p % n_bank
|
||||
|
||||
! get the delayed group
|
||||
g = fission_bank_delayed_group(n_bank - p % n_bank + k)
|
||||
|
||||
! Case for tallying delayed emissions
|
||||
if (g /= 0) then
|
||||
|
||||
! Accumulate the decay rate times delayed nu fission score
|
||||
associate (rxn => nuclides(p % event_nuclide) % &
|
||||
reactions(nuclides(p % event_nuclide) % index_fission(1)))
|
||||
|
||||
! determine score based on bank site weight and keff.
|
||||
score = score + keff * fission_bank_wgt(n_bank - p % n_bank + k) &
|
||||
* rxn % product_decay_rate(1 + g) * flux
|
||||
end associate
|
||||
|
||||
! if the delayed group filter is present, tally to corresponding
|
||||
! delayed group bin if it exists
|
||||
if (dg_filter > 0) then
|
||||
|
||||
! declare the delayed group filter type
|
||||
select type(filt => filters(t % filter(dg_filter) + 1) % obj)
|
||||
type is (DelayedGroupFilter)
|
||||
|
||||
! loop over delayed group bins until the corresponding bin
|
||||
! is found
|
||||
do d_bin = 1, filt % n_bins
|
||||
d = filt % groups(d_bin)
|
||||
|
||||
! check whether the delayed group of the particle is equal
|
||||
! to the delayed group of this bin
|
||||
if (d == g) then
|
||||
call score_fission_delayed_dg(i_tally, d_bin, score, &
|
||||
score_index)
|
||||
end if
|
||||
end do
|
||||
end select
|
||||
|
||||
! Reset the score to zero
|
||||
score = ZERO
|
||||
end if
|
||||
end if
|
||||
end do
|
||||
end if
|
||||
else
|
||||
|
||||
! Check if tally is on a single nuclide
|
||||
if (i_nuclide >= 0) then
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
select type(filt => filters(t % filter(dg_filter) + 1) % obj)
|
||||
type is (DelayedGroupFilter)
|
||||
|
||||
! Loop over all delayed group bins and tally to them
|
||||
! individually
|
||||
do d_bin = 1, filt % n_bins
|
||||
|
||||
! Get the delayed group for this bin
|
||||
d = filt % groups(d_bin)
|
||||
|
||||
! Compute the yield for this delayed group
|
||||
yield = nuclides(i_nuclide+1) % nu(E, EMISSION_DELAYED, d)
|
||||
|
||||
associate (rxn => nuclides(i_nuclide+1) % &
|
||||
reactions(nuclides(i_nuclide+1) % index_fission(1)))
|
||||
|
||||
! Compute the score and tally to bin
|
||||
score = micro_xs(i_nuclide+1) % fission * yield * flux * &
|
||||
atom_density * rxn % product_decay_rate(1 + d)
|
||||
end associate
|
||||
|
||||
! Tally to bin
|
||||
call score_fission_delayed_dg(i_tally, d_bin, score, score_index)
|
||||
end do
|
||||
cycle SCORE_LOOP
|
||||
end select
|
||||
else
|
||||
|
||||
! If the delayed group filter is not present, compute the score
|
||||
! by accumulating the absorbed weight times the decay rate times
|
||||
! the fraction of the delayed-nu-fission xs to the absorption xs
|
||||
! for all delayed groups.
|
||||
score = ZERO
|
||||
|
||||
associate (rxn => nuclides(p % event_nuclide) % &
|
||||
reactions(nuclides(p % event_nuclide) % index_fission(1)))
|
||||
|
||||
! We need to be careful not to overshoot the number of delayed
|
||||
! groups since this could cause the range of the rxn % products
|
||||
! array to be exceeded. Hence, we use the size of this array
|
||||
! and not the MAX_DELAYED_GROUPS constant for this loop.
|
||||
do d = 1, rxn % products_size() - 2
|
||||
score = score + micro_xs(i_nuclide+1) % fission * flux * &
|
||||
nuclides(i_nuclide+1) % nu(E, EMISSION_DELAYED) * &
|
||||
atom_density * rxn % product_decay_rate(1 + d)
|
||||
end do
|
||||
end associate
|
||||
end if
|
||||
|
||||
! Tally is on total nuclides
|
||||
else
|
||||
|
||||
! Check if the delayed group filter is present
|
||||
if (dg_filter > 0) then
|
||||
select type(filt => filters(t % filter(dg_filter) + 1) % obj)
|
||||
type is (DelayedGroupFilter)
|
||||
|
||||
! Loop over all nuclides in the current material
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
|
||||
! Get atom density
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
|
||||
if (nuclides(i_nuc) % fissionable) then
|
||||
|
||||
! Loop over all delayed group bins and tally to them
|
||||
! individually
|
||||
do d_bin = 1, filt % n_bins
|
||||
|
||||
! Get the delayed group for this bin
|
||||
d = filt % groups(d_bin)
|
||||
|
||||
! Get the yield for the desired nuclide and delayed
|
||||
! group
|
||||
yield = nuclides(i_nuc) % nu(E, EMISSION_DELAYED, d)
|
||||
|
||||
associate (rxn => nuclides(i_nuc) % &
|
||||
reactions(nuclides(i_nuc) % index_fission(1)))
|
||||
|
||||
! Compute the score
|
||||
score = micro_xs(i_nuc) % fission * yield * flux * &
|
||||
atom_density_ &
|
||||
* rxn % product_decay_rate(1 + d)
|
||||
end associate
|
||||
|
||||
! Tally to bin
|
||||
call score_fission_delayed_dg(i_tally, d_bin, score, &
|
||||
score_index)
|
||||
end do
|
||||
end if
|
||||
end do
|
||||
end if
|
||||
cycle SCORE_LOOP
|
||||
end select
|
||||
else
|
||||
|
||||
score = ZERO
|
||||
|
||||
! Loop over all nuclides in the current material
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
|
||||
! Get atom density
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
|
||||
if (nuclides(i_nuc) % fissionable) then
|
||||
|
||||
associate (rxn => nuclides(i_nuc) % &
|
||||
reactions(nuclides(i_nuc) % index_fission(1)))
|
||||
|
||||
! We need to be careful not to overshoot the number of
|
||||
! delayed groups since this could cause the range of the
|
||||
! rxn % products array to be exceeded. Hence, we use the
|
||||
! size of this array and not the MAX_DELAYED_GROUPS
|
||||
! constant for this loop.
|
||||
do d = 1, rxn % products_size() - 2
|
||||
|
||||
! Accumulate the contribution from each nuclide
|
||||
score = score + micro_xs(i_nuc) % fission &
|
||||
* nuclides(i_nuc) % nu(E, EMISSION_DELAYED) &
|
||||
* atom_density_ * flux &
|
||||
* rxn % product_decay_rate(1 + d)
|
||||
end do
|
||||
end associate
|
||||
end if
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
|
||||
case (SCORE_KAPPA_FISSION)
|
||||
|
||||
if (material_xs % absorption == ZERO) cycle SCORE_LOOP
|
||||
|
||||
! Determine kappa-fission cross section on the fly. The ENDF standard
|
||||
! (ENDF-102) states that MT 18 stores the fission energy as the Q_value
|
||||
! (fission(1))
|
||||
|
||||
score = ZERO
|
||||
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
if (survival_biasing) then
|
||||
! No fission events occur if survival biasing is on -- need to
|
||||
! calculate fraction of absorptions that would have resulted in
|
||||
! fission scaled by kappa-fission
|
||||
associate (nuc => nuclides(p % event_nuclide))
|
||||
if (micro_xs(p % event_nuclide) % absorption > ZERO .and. &
|
||||
nuc % fissionable) then
|
||||
score = p % absorb_wgt * &
|
||||
nuc % reactions(nuc % index_fission(1)) % Q_value * &
|
||||
micro_xs(p % event_nuclide) % fission / &
|
||||
micro_xs(p % event_nuclide) % absorption * flux
|
||||
end if
|
||||
end associate
|
||||
else
|
||||
! Skip any non-absorption events
|
||||
if (p % event == EVENT_SCATTER) cycle SCORE_LOOP
|
||||
! All fission events will contribute, so again we can use
|
||||
! particle's weight entering the collision as the estimate for
|
||||
! the fission energy production rate
|
||||
associate (nuc => nuclides(p % event_nuclide))
|
||||
if (nuc % fissionable) then
|
||||
score = p % last_wgt * &
|
||||
nuc % reactions(nuc % index_fission(1)) % Q_value * &
|
||||
micro_xs(p % event_nuclide) % fission / &
|
||||
micro_xs(p % event_nuclide) % absorption * flux
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
else
|
||||
if (i_nuclide >= 0) then
|
||||
associate (nuc => nuclides(i_nuclide+1))
|
||||
if (nuc % fissionable) then
|
||||
score = nuc % reactions(nuc % index_fission(1)) % Q_value * &
|
||||
micro_xs(i_nuclide+1) % fission * atom_density * flux
|
||||
end if
|
||||
end associate
|
||||
else
|
||||
if (p % material == MATERIAL_VOID) then
|
||||
score = ZERO
|
||||
else
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
! Determine atom density and index of nuclide
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
|
||||
! If nuclide is fissionable, accumulate kappa fission
|
||||
associate(nuc => nuclides(i_nuc))
|
||||
if (nuc % fissionable) then
|
||||
score = score + &
|
||||
nuc % reactions(nuc % index_fission(1)) % Q_value * &
|
||||
micro_xs(i_nuc) % fission * atom_density_ * flux
|
||||
end if
|
||||
end associate
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
|
||||
case (SCORE_EVENTS)
|
||||
! Simply count number of scoring events
|
||||
score = ONE
|
||||
|
||||
case (ELASTIC)
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
! Check if event MT matches
|
||||
if (p % event_MT /= ELASTIC) cycle SCORE_LOOP
|
||||
score = p % last_wgt * flux
|
||||
|
||||
else
|
||||
if (i_nuclide >= 0) then
|
||||
if (micro_xs(i_nuclide+1) % elastic == CACHE_INVALID) then
|
||||
call nuclides(i_nuclide+1) % calculate_elastic_xs()
|
||||
end if
|
||||
score = micro_xs(i_nuclide+1) % elastic * atom_density * flux
|
||||
else
|
||||
score = ZERO
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
! Get atom density
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
if (micro_xs(i_nuc) % elastic == CACHE_INVALID) then
|
||||
call nuclides(i_nuc) % calculate_elastic_xs()
|
||||
end if
|
||||
|
||||
score = score + micro_xs(i_nuc) % elastic * atom_density_ * flux
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
|
||||
case (SCORE_FISS_Q_PROMPT, SCORE_FISS_Q_RECOV)
|
||||
|
||||
if (material_xs % absorption == ZERO) cycle SCORE_LOOP
|
||||
|
||||
score = ZERO
|
||||
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
if (survival_biasing) then
|
||||
! No fission events occur if survival biasing is on -- need to
|
||||
! calculate fraction of absorptions that would have resulted in
|
||||
! fission scaled by Q-value
|
||||
associate (nuc => nuclides(p % event_nuclide))
|
||||
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) then
|
||||
xs = nuclide_fission_q_prompt(nuc % ptr, p % last_E)
|
||||
else if (score_bin == SCORE_FISS_Q_RECOV) then
|
||||
xs = nuclide_fission_q_recov(nuc % ptr, p % last_E)
|
||||
end if
|
||||
|
||||
score = p % absorb_wgt * xs * flux &
|
||||
* micro_xs(p % event_nuclide) % fission &
|
||||
/ micro_xs(p % event_nuclide) % absorption
|
||||
end if
|
||||
end associate
|
||||
else
|
||||
! Skip any non-absorption events
|
||||
if (p % event == EVENT_SCATTER) cycle SCORE_LOOP
|
||||
! All fission events will contribute, so again we can use
|
||||
! particle's weight entering the collision as the estimate for
|
||||
! the fission energy production rate
|
||||
associate (nuc => nuclides(p % event_nuclide))
|
||||
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) then
|
||||
xs = nuclide_fission_q_prompt(nuc % ptr, p % last_E)
|
||||
else if (score_bin == SCORE_FISS_Q_RECOV) then
|
||||
xs = nuclide_fission_q_recov(nuc % ptr, p % last_E)
|
||||
end if
|
||||
|
||||
score = p % last_wgt * xs * flux &
|
||||
* micro_xs(p % event_nuclide) % fission &
|
||||
/ micro_xs(p % event_nuclide) % absorption
|
||||
end if
|
||||
end associate
|
||||
end if
|
||||
|
||||
else
|
||||
if (i_nuclide >= 0) then
|
||||
associate (nuc => nuclides(i_nuclide+1))
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) then
|
||||
xs = nuclide_fission_q_prompt(nuc % ptr, E)
|
||||
else if (score_bin == SCORE_FISS_Q_RECOV) then
|
||||
xs = nuclide_fission_q_recov(nuc % ptr, E)
|
||||
end if
|
||||
|
||||
score = micro_xs(i_nuclide+1) % fission * atom_density * flux * xs
|
||||
end associate
|
||||
else
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
|
||||
associate (nuc => nuclides(i_nuc))
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) then
|
||||
xs = nuclide_fission_q_prompt(nuc % ptr, E)
|
||||
else if (score_bin == SCORE_FISS_Q_RECOV) then
|
||||
xs = nuclide_fission_q_recov(nuc % ptr, E)
|
||||
end if
|
||||
|
||||
score = score + micro_xs(i_nuc) % fission * atom_density_ &
|
||||
* flux * xs
|
||||
end associate
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
|
||||
case (N_2N, N_3N, N_4N, N_GAMMA, N_P, N_A)
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
! Check if event MT matches
|
||||
if (p % event_MT /= score_bin) cycle SCORE_LOOP
|
||||
score = p % last_wgt * flux
|
||||
|
||||
else
|
||||
! Determine index in NuclideMicroXS % reaction array
|
||||
select case (score_bin)
|
||||
case (N_GAMMA)
|
||||
m = 1
|
||||
case (N_P)
|
||||
m = 2
|
||||
case (N_A)
|
||||
m = 3
|
||||
case (N_2N)
|
||||
m = 4
|
||||
case (N_3N)
|
||||
m = 5
|
||||
case (N_4N)
|
||||
m = 6
|
||||
end select
|
||||
|
||||
if (i_nuclide >= 0) then
|
||||
score = micro_xs(i_nuclide+1) % reaction(m) * atom_density * flux
|
||||
else
|
||||
score = ZERO
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
score = score + micro_xs(i_nuc) % reaction(m) * atom_density_ * flux
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
end if
|
||||
|
||||
case default
|
||||
if (t % estimator() == ESTIMATOR_ANALOG) then
|
||||
! Any other score is assumed to be a MT number. Thus, we just need
|
||||
! to check if it matches the MT number of the event
|
||||
if (p % event_MT /= score_bin) cycle SCORE_LOOP
|
||||
score = p % last_wgt * flux
|
||||
|
||||
else
|
||||
! Any other cross section has to be calculated on-the-fly. For
|
||||
! cross sections that are used often (e.g. n2n, ngamma, etc. for
|
||||
! depletion), it might make sense to optimize this section or
|
||||
! pre-calculate cross sections
|
||||
if (score_bin > 1) then
|
||||
! Set default score
|
||||
score = ZERO
|
||||
|
||||
if (i_nuclide >= 0) then
|
||||
m = nuclides(i_nuclide+1) % reaction_index(score_bin)
|
||||
if (m /= 0) then
|
||||
! Retrieve temperature and energy grid index and interpolation
|
||||
! factor
|
||||
i_temp = micro_xs(i_nuclide+1) % index_temp + 1
|
||||
if (i_temp > 0) then
|
||||
i_energy = micro_xs(i_nuclide+1) % index_grid
|
||||
f = micro_xs(i_nuclide+1) % interp_factor
|
||||
|
||||
associate (rx => nuclides(i_nuclide+1) % reactions(m))
|
||||
threshold = rx % xs_threshold(i_temp)
|
||||
if (i_energy >= threshold) then
|
||||
score = ((ONE - f) * rx % xs(i_temp, i_energy - &
|
||||
threshold + 1) + f * rx % xs(i_temp, i_energy - &
|
||||
threshold + 2)) * atom_density * flux
|
||||
end if
|
||||
end associate
|
||||
else
|
||||
! This block is reached if multipole is turned on and we're in
|
||||
! the resolved range. Assume xs is zero.
|
||||
score = ZERO
|
||||
end if
|
||||
end if
|
||||
|
||||
else
|
||||
if (p % material /= MATERIAL_VOID) then
|
||||
do l = 1, material_nuclide_size(p % material)
|
||||
! Get atom density
|
||||
atom_density_ = material_atom_density(p % material, l)
|
||||
|
||||
! Get index in nuclides array
|
||||
i_nuc = material_nuclide(p % material, l)
|
||||
|
||||
m = nuclides(i_nuc) % reaction_index(score_bin)
|
||||
if (m /= 0) then
|
||||
! Retrieve temperature and energy grid index and
|
||||
! interpolation factor
|
||||
i_temp = micro_xs(i_nuc) % index_temp + 1
|
||||
if (i_temp > 0) then
|
||||
i_energy = micro_xs(i_nuc) % index_grid
|
||||
f = micro_xs(i_nuc) % interp_factor
|
||||
|
||||
associate (rx => nuclides(i_nuc) % reactions(m))
|
||||
threshold = rx % xs_threshold(i_temp)
|
||||
if (i_energy >= threshold) then
|
||||
score = score + ((ONE - f) * rx % xs(i_temp, i_energy - &
|
||||
threshold + 1) + f * rx % xs(i_temp, i_energy - &
|
||||
threshold + 2)) * atom_density_ * flux
|
||||
end if
|
||||
end associate
|
||||
else
|
||||
! This block is reached if multipole is turned on and
|
||||
! we're in the resolved range. Assume xs is zero.
|
||||
score = ZERO
|
||||
end if
|
||||
end if
|
||||
end do
|
||||
end if
|
||||
end if
|
||||
|
||||
else
|
||||
call fatal_error("Invalid score type on tally " &
|
||||
// to_str(t % id()) // ".")
|
||||
end if
|
||||
end if
|
||||
|
||||
end select
|
||||
|
||||
!#########################################################################
|
||||
! Add derivative information on score for differential tallies.
|
||||
|
||||
if (t % deriv() /= C_NONE) then
|
||||
call apply_derivative_to_score(p, i_tally, i_nuclide, atom_density, &
|
||||
score_bin, score)
|
||||
end if
|
||||
|
||||
!#########################################################################
|
||||
! Expand score if necessary and add to tally results.
|
||||
!$omp atomic
|
||||
t % results(RESULT_VALUE, score_index, filter_index) = &
|
||||
t % results(RESULT_VALUE, score_index, filter_index) + score
|
||||
|
||||
end do SCORE_LOOP
|
||||
end associate
|
||||
end subroutine score_general_ce
|
||||
|
||||
subroutine score_general_mg(p, i_tally, start_index, filter_index, i_nuclide, &
|
||||
atom_density, flux) bind(C)
|
||||
type(Particle), intent(in) :: p
|
||||
|
|
|
|||
|
|
@ -848,7 +848,7 @@ score_fission_eout(Particle* p, int i_tally, int i_score, int score_bin)
|
|||
}
|
||||
|
||||
extern "C" void
|
||||
score_general_ce_c(Particle* p, int i_tally, int start_index, int filter_index,
|
||||
score_general_ce(Particle* p, int i_tally, int start_index, int filter_index,
|
||||
int i_nuclide, double atom_density, double flux)
|
||||
{
|
||||
//TODO: off-by-one
|
||||
|
|
@ -1287,8 +1287,459 @@ score_general_ce_c(Particle* p, int i_tally, int start_index, int filter_index,
|
|||
break;
|
||||
|
||||
|
||||
case SCORE_DECAY_RATE:
|
||||
if (simulation::material_xs.absorption == 0) continue;
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
if (settings::survival_biasing) {
|
||||
// No fission events occur if survival biasing is on -- need to
|
||||
// calculate fraction of absorptions that would have resulted in
|
||||
// delayed-nu-fission
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0
|
||||
&& nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
if (tally.delayedgroup_filter_ > 0) {
|
||||
auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_-1];
|
||||
const DelayedGroupFilter& filt
|
||||
{*dynamic_cast<DelayedGroupFilter*>(
|
||||
model::tally_filters[i_dg_filt].get())};
|
||||
// Tally each delayed group bin individually
|
||||
for (auto d_bin = 0; d_bin < filt.n_bins_; ++d_bin) {
|
||||
auto d = filt.groups_[d_bin];
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
|
||||
auto rate = rxn.products_[d].decay_rate_;
|
||||
score = p->absorb_wgt * yield
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption
|
||||
* rate * flux;
|
||||
score_fission_delayed_dg(i_tally, d_bin+1, score,
|
||||
score_index+1);
|
||||
}
|
||||
continue;
|
||||
} else {
|
||||
// If the delayed group filter is not present, compute the score
|
||||
// by multiplying the absorbed weight by the fraction of the
|
||||
// delayed-nu-fission xs to the absorption xs for all delayed
|
||||
// groups
|
||||
score = 0.;
|
||||
// We need to be careful not to overshoot the number of
|
||||
// delayed groups since this could cause the range of the
|
||||
// rxn.products_ array to be exceeded. Hence, we use the size
|
||||
// of this array and not the MAX_DELAYED_GROUPS constant for
|
||||
// this loop.
|
||||
for (auto d = 0; d < rxn.products_.size() - 2; ++d)
|
||||
{
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
|
||||
auto rate = rxn.products_[d+1].decay_rate_;
|
||||
score += rate * p->absorb_wgt
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission * yield
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
} else {
|
||||
// Skip any non-fission events
|
||||
if (!p->fission) continue;
|
||||
// If there is no outgoing energy filter, than we only need to score
|
||||
// to one bin. For the score to be 'analog', we need to score the
|
||||
// number of particles that were banked in the fission bank. Since
|
||||
// this was weighted by 1/keff, we multiply by keff to get the proper
|
||||
// score. Loop over the neutrons produced from fission and check which
|
||||
// ones are delayed. If a delayed neutron is encountered, add its
|
||||
// contribution to the fission bank to the score.
|
||||
score = 0.;
|
||||
for (auto i = 0; i < p->n_bank; ++i) {
|
||||
auto i_bank = simulation::n_bank - p->n_bank + i;
|
||||
const auto& bank = simulation::fission_bank[i_bank];
|
||||
auto g = bank.delayed_group;
|
||||
if (g != 0) {
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
auto rate = rxn.products_[g].decay_rate_;
|
||||
//score += simulation::keff * bank.wgt * rate * flux;
|
||||
if (tally.delayedgroup_filter_ > 0) {
|
||||
auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_-1];
|
||||
const DelayedGroupFilter& filt
|
||||
{*dynamic_cast<DelayedGroupFilter*>(
|
||||
model::tally_filters[i_dg_filt].get())};
|
||||
// Find the corresponding filter bin and then score
|
||||
for (auto d_bin = 0; d_bin < filt.n_bins_; ++d_bin) {
|
||||
auto d = filt.groups_[d_bin];
|
||||
if (d == g)
|
||||
score_fission_delayed_dg(i_tally, d_bin+1, score,
|
||||
score_index+1);
|
||||
}
|
||||
score = 0.;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
} else {
|
||||
if (i_nuclide >= 0) {
|
||||
const auto& nuc {*data::nuclides[i_nuclide]};
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
if (tally.delayedgroup_filter_ > 0) {
|
||||
auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_-1];
|
||||
const DelayedGroupFilter& filt
|
||||
{*dynamic_cast<DelayedGroupFilter*>(
|
||||
model::tally_filters[i_dg_filt].get())};
|
||||
// Tally each delayed group bin individually
|
||||
for (auto d_bin = 0; d_bin < filt.n_bins_; ++d_bin) {
|
||||
auto d = filt.groups_[d_bin];
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
|
||||
auto rate = rxn.products_[d].decay_rate_;
|
||||
score = simulation::micro_xs[i_nuclide].fission * yield * flux
|
||||
* atom_density * rate;
|
||||
score_fission_delayed_dg(i_tally, d_bin+1, score, score_index+1);
|
||||
}
|
||||
continue;
|
||||
} else {
|
||||
score = 0.;
|
||||
// We need to be careful not to overshoot the number of
|
||||
// delayed groups since this could cause the range of the
|
||||
// rxn.products_ array to be exceeded. Hence, we use the size
|
||||
// of this array and not the MAX_DELAYED_GROUPS constant for
|
||||
// this loop.
|
||||
for (auto d = 0; d < rxn.products_.size() - 2; ++d)
|
||||
{
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
|
||||
auto rate = rxn.products_[d+1].decay_rate_;
|
||||
score += simulation::micro_xs[i_nuclide].fission * flux
|
||||
* yield * atom_density * rate;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
if (tally.delayedgroup_filter_ > 0) {
|
||||
auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_-1];
|
||||
const DelayedGroupFilter& filt
|
||||
{*dynamic_cast<DelayedGroupFilter*>(
|
||||
model::tally_filters[i_dg_filt].get())};
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
const auto& nuc {*data::nuclides[j_nuclide]};
|
||||
if (nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
// Tally each delayed group bin individually
|
||||
for (auto d_bin = 0; d_bin < filt.n_bins_; ++d_bin) {
|
||||
auto d = filt.groups_[d_bin];
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
|
||||
auto rate = rxn.products_[d].decay_rate_;
|
||||
score = simulation::micro_xs[j_nuclide].fission * yield
|
||||
* flux * atom_density * rate;
|
||||
score_fission_delayed_dg(i_tally, d_bin+1, score,
|
||||
score_index+1);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
continue;
|
||||
} else {
|
||||
score = 0.;
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
const auto& nuc {*data::nuclides[j_nuclide]};
|
||||
if (nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
// We need to be careful not to overshoot the number of
|
||||
// delayed groups since this could cause the range of the
|
||||
// rxn.products_ array to be exceeded. Hence, we use the size
|
||||
// of this array and not the MAX_DELAYED_GROUPS constant for
|
||||
// this loop.
|
||||
for (auto d = 0; d < rxn.products_.size() - 2; ++d)
|
||||
{
|
||||
auto yield
|
||||
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
|
||||
auto rate = rxn.products_[d+1].decay_rate_;
|
||||
score += simulation::micro_xs[j_nuclide].fission
|
||||
* yield * atom_density * flux * rate;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
break;
|
||||
|
||||
|
||||
case SCORE_KAPPA_FISSION:
|
||||
if (simulation::material_xs.absorption == 0.) continue;
|
||||
score = 0.;
|
||||
// Kappa-fission values are determined from the Q-value listed for the
|
||||
// fission cross section.
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
if (settings::survival_biasing) {
|
||||
// No fission events occur if survival biasing is on -- need to
|
||||
// calculate fraction of absorptions that would have resulted in
|
||||
// fission scaled by the Q-value
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0
|
||||
&& nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
score = p->absorb_wgt * rxn.q_value_
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
|
||||
}
|
||||
} else {
|
||||
// Skip any non-absorption events
|
||||
if (p->event == EVENT_SCATTER) continue;
|
||||
// All fission events will contribute, so again we can use particle's
|
||||
// weight entering the collision as the estimate for the fission
|
||||
// reaction rate
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0
|
||||
&& nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
score = p->last_wgt * rxn.q_value_
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
if (i_nuclide >= 0) {
|
||||
const auto& nuc {*data::nuclides[i_nuclide]};
|
||||
if (nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
score = rxn.q_value_ * simulation::micro_xs[i_nuclide].fission
|
||||
* atom_density * flux;
|
||||
}
|
||||
} else {
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
const auto& nuc {*data::nuclides[j_nuclide]};
|
||||
if (nuc.fissionable_) {
|
||||
const auto& rxn {*nuc.fission_rx_[0]};
|
||||
score += rxn.q_value_ * simulation::micro_xs[j_nuclide].fission
|
||||
* atom_density * flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
break;
|
||||
|
||||
|
||||
case SCORE_EVENTS:
|
||||
// Simply count the number of scoring events
|
||||
score = 1.;
|
||||
break;
|
||||
|
||||
|
||||
case ELASTIC:
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
// Check if event MT matches
|
||||
if (p->event_MT != ELASTIC) continue;
|
||||
score = p->last_wgt * flux;
|
||||
} else {
|
||||
if (i_nuclide >= 0) {
|
||||
if (simulation::micro_xs[i_nuclide].elastic == CACHE_INVALID)
|
||||
data::nuclides[i_nuclide]->calculate_elastic_xs();
|
||||
score = simulation::micro_xs[i_nuclide].elastic * atom_density * flux;
|
||||
} else {
|
||||
score = 0.;
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
if (simulation::micro_xs[j_nuclide].elastic == CACHE_INVALID)
|
||||
data::nuclides[j_nuclide]->calculate_elastic_xs();
|
||||
score += simulation::micro_xs[j_nuclide].elastic * atom_density
|
||||
* flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
break;
|
||||
|
||||
case SCORE_FISS_Q_PROMPT:
|
||||
case SCORE_FISS_Q_RECOV:
|
||||
//continue;
|
||||
if (simulation::material_xs.absorption == 0.) continue;
|
||||
score = 0.;
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
if (settings::survival_biasing) {
|
||||
// No fission events occur if survival biasing is on -- need to
|
||||
// calculate fraction of absorptions that would have resulted in
|
||||
// fission scaled by the Q-value
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
|
||||
double q_value = 0.;
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) {
|
||||
if (nuc.fission_q_prompt_)
|
||||
q_value = (*nuc.fission_q_prompt_)(p->last_E);
|
||||
} else if (score_bin == SCORE_FISS_Q_RECOV) {
|
||||
if (nuc.fission_q_recov_)
|
||||
q_value = (*nuc.fission_q_recov_)(p->last_E);
|
||||
}
|
||||
score = p->absorb_wgt * q_value
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
|
||||
}
|
||||
} else {
|
||||
// Skip any non-absorption events
|
||||
if (p->event == EVENT_SCATTER) continue;
|
||||
// All fission events will contribute, so again we can use particle's
|
||||
// weight entering the collision as the estimate for the fission
|
||||
// reaction rate
|
||||
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
|
||||
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
|
||||
double q_value = 0.;
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) {
|
||||
if (nuc.fission_q_prompt_)
|
||||
q_value = (*nuc.fission_q_prompt_)(p->last_E);
|
||||
} else if (score_bin == SCORE_FISS_Q_RECOV) {
|
||||
if (nuc.fission_q_recov_)
|
||||
q_value = (*nuc.fission_q_recov_)(p->last_E);
|
||||
}
|
||||
score = p->last_wgt * q_value
|
||||
* simulation::micro_xs[p->event_nuclide-1].fission
|
||||
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
if (i_nuclide >= 0) {
|
||||
const auto& nuc {*data::nuclides[i_nuclide]};
|
||||
double q_value = 0.;
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) {
|
||||
if (nuc.fission_q_prompt_)
|
||||
q_value = (*nuc.fission_q_prompt_)(p->last_E);
|
||||
} else if (score_bin == SCORE_FISS_Q_RECOV) {
|
||||
if (nuc.fission_q_recov_)
|
||||
q_value = (*nuc.fission_q_recov_)(p->last_E);
|
||||
}
|
||||
score = q_value * simulation::micro_xs[i_nuclide].fission
|
||||
* atom_density * flux;
|
||||
} else {
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
const auto& nuc {*data::nuclides[j_nuclide]};
|
||||
double q_value = 0.;
|
||||
if (score_bin == SCORE_FISS_Q_PROMPT) {
|
||||
if (nuc.fission_q_prompt_)
|
||||
q_value = (*nuc.fission_q_prompt_)(p->last_E);
|
||||
} else if (score_bin == SCORE_FISS_Q_RECOV) {
|
||||
if (nuc.fission_q_recov_)
|
||||
q_value = (*nuc.fission_q_recov_)(p->last_E);
|
||||
}
|
||||
score += q_value * simulation::micro_xs[j_nuclide].fission
|
||||
* atom_density * flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
break;
|
||||
|
||||
|
||||
case N_2N:
|
||||
case N_3N:
|
||||
case N_4N:
|
||||
case N_GAMMA:
|
||||
case N_P:
|
||||
case N_A:
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
// Check if the event MT matches
|
||||
if (p->event_MT != score_bin) continue;
|
||||
score = p->last_wgt * flux;
|
||||
} else {
|
||||
int m;
|
||||
switch (score_bin) {
|
||||
case N_GAMMA: m = 0; break;
|
||||
case N_P: m = 1; break;
|
||||
case N_A: m = 2; break;
|
||||
case N_2N: m = 3; break;
|
||||
case N_3N: m = 4; break;
|
||||
case N_4N: m = 5; break;
|
||||
}
|
||||
if (i_nuclide >= 0) {
|
||||
score = simulation::micro_xs[i_nuclide].reaction[m] * atom_density
|
||||
* flux;
|
||||
} else {
|
||||
score = 0.;
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
score += simulation::micro_xs[j_nuclide].reaction[m]
|
||||
* atom_density * flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
break;
|
||||
|
||||
|
||||
default:
|
||||
continue;
|
||||
if (tally.estimator_ == ESTIMATOR_ANALOG) {
|
||||
// Any other score is assumed to be a MT number. Thus, we just need
|
||||
// to check if it matches the MT number of the event
|
||||
if (p->event_MT != score_bin) continue;
|
||||
score = p->last_wgt*flux;
|
||||
} else {
|
||||
// Any other cross section has to be calculated on-the-fly
|
||||
if (score_bin < 2) fatal_error("Invalid score type on tally "
|
||||
+ std::to_string(tally.id_));
|
||||
score = 0.;
|
||||
if (i_nuclide >= 0) {
|
||||
const auto& nuc {*data::nuclides[i_nuclide]};
|
||||
auto m = nuc.reaction_index_[score_bin];
|
||||
if (m == C_NONE) continue;
|
||||
const auto& rxn {*nuc.reactions_[m]};
|
||||
auto i_temp = simulation::micro_xs[i_nuclide].index_temp;
|
||||
if (i_temp >= 0) { // Can be false due to multipole
|
||||
auto i_grid = simulation::micro_xs[i_nuclide].index_grid - 1;
|
||||
auto f = simulation::micro_xs[i_nuclide].interp_factor;
|
||||
const auto& xs {rxn.xs_[i_temp]};
|
||||
auto threshold = xs.threshold - 1;
|
||||
if (i_grid >= xs.threshold) {
|
||||
score = ((1.0 - f) * xs.value[i_grid-threshold]
|
||||
+ f * xs.value[i_grid-threshold+1]) * atom_density * flux;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
if (p->material != MATERIAL_VOID) {
|
||||
const Material& material {*model::materials[p->material-1]};
|
||||
for (auto i = 0; i < material.nuclide_.size(); ++i) {
|
||||
auto j_nuclide = material.nuclide_[i];
|
||||
auto atom_density = material.atom_density_(i);
|
||||
const auto& nuc {*data::nuclides[j_nuclide]};
|
||||
auto m = nuc.reaction_index_[score_bin];
|
||||
if (m == C_NONE) continue;
|
||||
const auto& rxn {*nuc.reactions_[m]};
|
||||
auto i_temp = simulation::micro_xs[j_nuclide].index_temp;
|
||||
if (i_temp >= 0) { // Can be false due to multipole
|
||||
auto i_grid = simulation::micro_xs[j_nuclide].index_grid - 1;
|
||||
auto f = simulation::micro_xs[j_nuclide].interp_factor;
|
||||
const auto& xs {rxn.xs_[i_temp]};
|
||||
auto threshold = xs.threshold - 1;
|
||||
if (i_grid >= threshold) {
|
||||
score += ((1.0 - f) * xs.value[i_grid-threshold]
|
||||
+ f * xs.value[i_grid-threshold+1]) * atom_density
|
||||
* flux;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
// Add derivative information on score for differnetial tallies.
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue