crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 06:05:58 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Rename cross section caches

Browse source
This commit is contained in:
Paul Romano 2019-03-18 09:18:54 -05:00
parent 8a8b7b2c0a
commit 52fc93e6e0
9 changed files with 239 additions and 238 deletions
Download patch file Download diff file Expand all files Collapse all files

11
include/openmc/particle.h
View file

@ -111,11 +111,11 @@ struct ElementMicroXS {
};
//==============================================================================
// MATERIALMACROXS contains cached macroscopic cross sections for the material a
// MACROXS contains cached macroscopic cross sections for the material a
// particle is traveling through
//==============================================================================
struct MaterialMacroXS {
struct MacroXS {
double total; //!< macroscopic total xs
double absorption; //!< macroscopic absorption xs
double fission; //!< macroscopic fission xs
@ -219,9 +219,10 @@ public:
//==========================================================================
// Data members
std::vector<NuclideMicroXS> micro_xs_;
std::vector<ElementMicroXS> micro_photon_xs_;
MaterialMacroXS material_xs_;
// Cross section caches
std::vector<NuclideMicroXS> neutron_xs_; //!< Microscopic neutron cross sections
std::vector<ElementMicroXS> photon_xs_; //!< Microscopic photon cross sections
MacroXS macro_xs_; //!< Macroscopic cross sections
int64_t id_; //!< Unique ID
Type type_ {Type::neutron}; //!< Particle type (n, p, e, etc.)

38
src/material.cpp
View file

@ -732,10 +732,10 @@ void Material::init_nuclide_index()
void Material::calculate_xs(Particle& p) const
{
// Set all material macroscopic cross sections to zero
p.material_xs_.total = 0.0;
p.material_xs_.absorption = 0.0;
p.material_xs_.fission = 0.0;
p.material_xs_.nu_fission = 0.0;
p.macro_xs_.total = 0.0;
p.macro_xs_.absorption = 0.0;
p.macro_xs_.fission = 0.0;
p.macro_xs_.nu_fission = 0.0;
if (p.type_ == Particle::Type::neutron) {
this->calculate_neutron_xs(p);
@ -792,7 +792,7 @@ void Material::calculate_neutron_xs(Particle& p) const
int i_nuclide = nuclide_[i];
// Calculate microscopic cross section for this nuclide
const auto& micro {p.micro_xs_[i_nuclide]};
const auto& micro {p.neutron_xs_[i_nuclide]};
if (p.E_ != micro.last_E
|| p.sqrtkT_ != micro.last_sqrtkT
|| i_sab != micro.index_sab
@ -807,19 +807,19 @@ void Material::calculate_neutron_xs(Particle& p) const
double atom_density = atom_density_(i);
// Add contributions to cross sections
p.material_xs_.total += atom_density * micro.total;
p.material_xs_.absorption += atom_density * micro.absorption;
p.material_xs_.fission += atom_density * micro.fission;
p.material_xs_.nu_fission += atom_density * micro.nu_fission;
p.macro_xs_.total += atom_density * micro.total;
p.macro_xs_.absorption += atom_density * micro.absorption;
p.macro_xs_.fission += atom_density * micro.fission;
p.macro_xs_.nu_fission += atom_density * micro.nu_fission;
}
}
void Material::calculate_photon_xs(Particle& p) const
{
p.material_xs_.coherent = 0.0;
p.material_xs_.incoherent = 0.0;
p.material_xs_.photoelectric = 0.0;
p.material_xs_.pair_production = 0.0;
p.macro_xs_.coherent = 0.0;
p.macro_xs_.incoherent = 0.0;
p.macro_xs_.photoelectric = 0.0;
p.macro_xs_.pair_production = 0.0;
// Add contribution from each nuclide in material
for (int i = 0; i < nuclide_.size(); ++i) {
@ -830,7 +830,7 @@ void Material::calculate_photon_xs(Particle& p) const
int i_element = element_[i];
// Calculate microscopic cross section for this nuclide
const auto& micro {p.micro_photon_xs_[i_element]};
const auto& micro {p.photon_xs_[i_element]};
if (p.E_ != micro.last_E) {
data::elements[i_element].calculate_xs(p);
}
@ -842,11 +842,11 @@ void Material::calculate_photon_xs(Particle& p) const
double atom_density = atom_density_(i);
// Add contributions to material macroscopic cross sections
p.material_xs_.total += atom_density * micro.total;
p.material_xs_.coherent += atom_density * micro.coherent;
p.material_xs_.incoherent += atom_density * micro.incoherent;
p.material_xs_.photoelectric += atom_density * micro.photoelectric;
p.material_xs_.pair_production += atom_density * micro.pair_production;
p.macro_xs_.total += atom_density * micro.total;
p.macro_xs_.coherent += atom_density * micro.coherent;
p.macro_xs_.incoherent += atom_density * micro.incoherent;
p.macro_xs_.photoelectric += atom_density * micro.photoelectric;
p.macro_xs_.pair_production += atom_density * micro.pair_production;
}
}

8
src/nuclide.cpp
View file

@ -459,7 +459,7 @@ double Nuclide::nu(double E, EmissionMode mode, int group) const
void Nuclide::calculate_elastic_xs(Particle& p) const
{
// Get temperature index, grid index, and interpolation factor
auto& micro {p.micro_xs_[i_nuclide_]};
auto& micro {p.neutron_xs_[i_nuclide_]};
int i_temp = micro.index_temp;
int i_grid = micro.index_grid;
double f = micro.interp_factor;
@ -495,7 +495,7 @@ double Nuclide::elastic_xs_0K(double E) const
void Nuclide::calculate_xs(int i_sab, int i_log_union, double sab_frac, Particle& p)
{
auto& micro {p.micro_xs_[i_nuclide_]};
auto& micro {p.neutron_xs_[i_nuclide_]};
// Initialize cached cross sections to zero
micro.elastic = CACHE_INVALID;
@ -702,7 +702,7 @@ void Nuclide::calculate_xs(int i_sab, int i_log_union, double sab_frac, Particle
void Nuclide::calculate_sab_xs(int i_sab, double sab_frac, Particle& p)
{
auto& micro {p.micro_xs_[i_nuclide_]};
auto& micro {p.neutron_xs_[i_nuclide_]};
// Set flag that S(a,b) treatment should be used for scattering
micro.index_sab = i_sab;
@ -731,7 +731,7 @@ void Nuclide::calculate_sab_xs(int i_sab, double sab_frac, Particle& p)
void Nuclide::calculate_urr_xs(int i_temp, Particle& p) const
{
auto& micro = p.micro_xs_[i_nuclide_];
auto& micro = p.neutron_xs_[i_nuclide_];
micro.use_ptable = true;
// Create a shorthand for the URR data

26
src/particle.cpp
View file

@ -58,8 +58,8 @@ Particle::Particle()
}
// Create microscopic cross section caches
micro_xs_.resize(data::nuclides.size());
micro_photon_xs_.resize(data::elements.size());
neutron_xs_.resize(data::nuclides.size());
photon_xs_.resize(data::elements.size());
}
void
@ -135,7 +135,7 @@ Particle::transport()
// Force calculation of cross-sections by setting last energy to zero
if (settings::run_CE) {
for (auto& micro : micro_xs_) micro.last_E = 0.0;
for (auto& micro : neutron_xs_) micro.last_E = 0.0;
}
// Prepare to write out particle track.
@ -189,17 +189,17 @@ Particle::transport()
} else {
// Get the MG data
calculate_xs_c(material_, g_, sqrtkT_, this->u_local(),
material_xs_.total, material_xs_.absorption, material_xs_.nu_fission);
macro_xs_.total, macro_xs_.absorption, macro_xs_.nu_fission);
// Finally, update the particle group while we have already checked
// for if multi-group
g_last_ = g_;
}
} else {
material_xs_.total = 0.0;
material_xs_.absorption = 0.0;
material_xs_.fission = 0.0;
material_xs_.nu_fission = 0.0;
macro_xs_.total = 0.0;
macro_xs_.absorption = 0.0;
macro_xs_.fission = 0.0;
macro_xs_.nu_fission = 0.0;
}
// Find the distance to the nearest boundary
@ -215,10 +215,10 @@ Particle::transport()
if (type_ == Particle::Type::electron ||
type_ == Particle::Type::positron) {
d_collision = 0.0;
} else if (material_xs_.total == 0.0) {
} else if (macro_xs_.total == 0.0) {
d_collision = INFINITY;
} else {
d_collision = -std::log(prn()) / material_xs_.total;
d_collision = -std::log(prn()) / macro_xs_.total;
}
// Select smaller of the two distances
@ -237,7 +237,7 @@ Particle::transport()
// Score track-length estimate of k-eff
if (settings::run_mode == RUN_MODE_EIGENVALUE &&
type_ == Particle::Type::neutron) {
global_tally_tracklength += wgt_ * distance * material_xs_.nu_fission;
global_tally_tracklength += wgt_ * distance * macro_xs_.nu_fission;
}
// Score flux derivative accumulators for differential tallies.
@ -280,8 +280,8 @@ Particle::transport()
// Score collision estimate of keff
if (settings::run_mode == RUN_MODE_EIGENVALUE &&
type_ == Particle::Type::neutron) {
global_tally_collision += wgt_ * material_xs_.nu_fission
/ material_xs_.total;
global_tally_collision += wgt_ * macro_xs_.nu_fission
/ macro_xs_.total;
}
// Score surface current tallies -- this has to be done before the collision

2
src/photon.cpp
View file

@ -452,7 +452,7 @@ void PhotonInteraction::calculate_xs(Particle& p) const
// calculate interpolation factor
double f = (log_E - energy_(i_grid)) / (energy_(i_grid+1) - energy_(i_grid));
auto& xs {p.micro_photon_xs_[i_element_]};
auto& xs {p.photon_xs_[i_element_]};
xs.index_grid = i_grid;
xs.interp_factor = f;

62
src/physics.cpp
View file

@ -110,7 +110,7 @@ void sample_neutron_reaction(Particle* p)
// If survival biasing is being used, the following subroutine adjusts the
// weight of the particle. Otherwise, it checks to see if absorption occurs
if (p->micro_xs_[i_nuclide].absorption > 0.0) {
if (p->neutron_xs_[i_nuclide].absorption > 0.0) {
absorption(p, i_nuclide);
} else {
p->wgt_absorb_ = 0.0;
@ -146,8 +146,8 @@ create_fission_sites(Particle* p, int i_nuclide, const Reaction* rx,
double weight = settings::ufs_on ? ufs_get_weight(p) : 1.0;
// Determine the expected number of neutrons produced
double nu_t = p->wgt_ / simulation::keff * weight * p->micro_xs_[
i_nuclide].nu_fission / p->micro_xs_[i_nuclide].total;
double nu_t = p->wgt_ / simulation::keff * weight * p->neutron_xs_[
i_nuclide].nu_fission / p->neutron_xs_[i_nuclide].total;
// Sample the number of neutrons produced
int nu = static_cast<int>(nu_t);
@ -227,7 +227,7 @@ void sample_photon_reaction(Particle* p)
// Sample element within material
int i_element = sample_element(p);
p->event_nuclide_ = i_element;
const auto& micro {p->micro_photon_xs_[i_element]};
const auto& micro {p->photon_xs_[i_element]};
const auto& element {data::elements[i_element]};
// Calculate photon energy over electron rest mass equivalent
@ -408,7 +408,7 @@ void sample_positron_reaction(Particle* p)
int sample_nuclide(const Particle* p)
{
// Sample cumulative distribution function
double cutoff = prn() * p->material_xs_.total;
double cutoff = prn() * p->macro_xs_.total;
// Get pointers to nuclide/density arrays
const auto& mat {model::materials[p->material_]};
@ -421,7 +421,7 @@ int sample_nuclide(const Particle* p)
double atom_density = mat->atom_density_[i];
// Increment probability to compare to cutoff
prob += atom_density * p->micro_xs_[i_nuclide].total;
prob += atom_density * p->neutron_xs_[i_nuclide].total;
if (prob >= cutoff) return i_nuclide;
}
@ -433,7 +433,7 @@ int sample_nuclide(const Particle* p)
int sample_element(Particle* p)
{
// Sample cumulative distribution function
double cutoff = prn() * p->material_xs_.total;
double cutoff = prn() * p->macro_xs_.total;
// Get pointers to elements, densities
const auto& mat {model::materials[p->material_]};
@ -454,7 +454,7 @@ int sample_element(Particle* p)
double atom_density = mat->atom_density_[i];
// Determine microscopic cross section
double sigma = atom_density * p->micro_photon_xs_[i_element].total;
double sigma = atom_density * p->photon_xs_[i_element].total;
// Increment probability to compare to cutoff
prob += sigma;
@ -472,7 +472,7 @@ Reaction* sample_fission(int i_nuclide, const Particle* p)
// If we're in the URR, by default use the first fission reaction. We also
// default to the first reaction if we know that there are no partial fission
// reactions
if (p->micro_xs_[i_nuclide].use_ptable || !nuc->has_partial_fission_) {
if (p->neutron_xs_[i_nuclide].use_ptable || !nuc->has_partial_fission_) {
return nuc->fission_rx_[0];
}
@ -485,10 +485,10 @@ Reaction* sample_fission(int i_nuclide, const Particle* p)
}
// Get grid index and interpolatoin factor and sample fission cdf
int i_temp = p->micro_xs_[i_nuclide].index_temp;
int i_grid = p->micro_xs_[i_nuclide].index_grid;
double f = p->micro_xs_[i_nuclide].interp_factor;
double cutoff = prn() * p->micro_xs_[i_nuclide].fission;
int i_temp = p->neutron_xs_[i_nuclide].index_temp;
int i_grid = p->neutron_xs_[i_nuclide].index_grid;
double f = p->neutron_xs_[i_nuclide].interp_factor;
double cutoff = prn() * p->neutron_xs_[i_nuclide].fission;
double prob = 0.0;
// Loop through each partial fission reaction type
@ -512,10 +512,10 @@ Reaction* sample_fission(int i_nuclide, const Particle* p)
void sample_photon_product(int i_nuclide, const Particle* p, int* i_rx, int* i_product)
{
// Get grid index and interpolation factor and sample photon production cdf
int i_temp = p->micro_xs_[i_nuclide].index_temp;
int i_grid = p->micro_xs_[i_nuclide].index_grid;
double f = p->micro_xs_[i_nuclide].interp_factor;
double cutoff = prn() * p->micro_xs_[i_nuclide].photon_prod;
int i_temp = p->neutron_xs_[i_nuclide].index_temp;
int i_grid = p->neutron_xs_[i_nuclide].index_grid;
double f = p->neutron_xs_[i_nuclide].interp_factor;
double cutoff = prn() * p->neutron_xs_[i_nuclide].photon_prod;
double prob = 0.0;
// Loop through each reaction type
@ -548,8 +548,8 @@ void absorption(Particle* p, int i_nuclide)
{
if (settings::survival_biasing) {
// Determine weight absorbed in survival biasing
p->wgt_absorb_ = p->wgt_ * p->micro_xs_[i_nuclide].absorption /
p->micro_xs_[i_nuclide].total;
p->wgt_absorb_ = p->wgt_ * p->neutron_xs_[i_nuclide].absorption /
p->neutron_xs_[i_nuclide].total;
// Adjust weight of particle by probability of absorption
p->wgt_ -= p->wgt_absorb_;
@ -557,17 +557,17 @@ void absorption(Particle* p, int i_nuclide)
// Score implicit absorption estimate of keff
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
global_tally_absorption += p->wgt_absorb_ * p->micro_xs_[
i_nuclide].nu_fission / p->micro_xs_[i_nuclide].absorption;
global_tally_absorption += p->wgt_absorb_ * p->neutron_xs_[
i_nuclide].nu_fission / p->neutron_xs_[i_nuclide].absorption;
}
} else {
// See if disappearance reaction happens
if (p->micro_xs_[i_nuclide].absorption >
prn() * p->micro_xs_[i_nuclide].total) {
if (p->neutron_xs_[i_nuclide].absorption >
prn() * p->neutron_xs_[i_nuclide].total) {
// Score absorption estimate of keff
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
global_tally_absorption += p->wgt_ * p->micro_xs_[
i_nuclide].nu_fission / p->micro_xs_[i_nuclide].absorption;
global_tally_absorption += p->wgt_ * p->neutron_xs_[
i_nuclide].nu_fission / p->neutron_xs_[i_nuclide].absorption;
}
p->alive_ = false;
@ -584,7 +584,7 @@ void scatter(Particle* p, int i_nuclide)
// Get pointer to nuclide and grid index/interpolation factor
const auto& nuc {data::nuclides[i_nuclide]};
const auto& micro {p->micro_xs_[i_nuclide]};
const auto& micro {p->neutron_xs_[i_nuclide]};
int i_temp = micro.index_temp;
int i_grid = micro.index_grid;
double f = micro.interp_factor;
@ -691,9 +691,9 @@ void elastic_scatter(int i_nuclide, const Reaction& rx, double kT,
// Sample velocity of target nucleus
Direction v_t {};
if (!p->micro_xs_[i_nuclide].use_ptable) {
if (!p->neutron_xs_[i_nuclide].use_ptable) {
v_t = sample_target_velocity(nuc.get(), p->E_, p->u(), v_n,
p->micro_xs_[i_nuclide].elastic, kT);
p->neutron_xs_[i_nuclide].elastic, kT);
}
// Velocity of center-of-mass
@ -746,7 +746,7 @@ void elastic_scatter(int i_nuclide, const Reaction& rx, double kT,
void sab_scatter(int i_nuclide, int i_sab, Particle* p)
{
// Determine temperature index
const auto& micro {p->micro_xs_[i_nuclide]};
const auto& micro {p->neutron_xs_[i_nuclide]};
int i_temp = micro.index_temp_sab;
// Sample energy and angle
@ -1104,8 +1104,8 @@ void inelastic_scatter(const Nuclide* nuc, const Reaction* rx, Particle* p)
void sample_secondary_photons(Particle* p, int i_nuclide)
{
// Sample the number of photons produced
double y_t = p->wgt_ * p->micro_xs_[i_nuclide].photon_prod /
p->micro_xs_[i_nuclide].total;
double y_t = p->wgt_ * p->neutron_xs_[i_nuclide].photon_prod /
p->neutron_xs_[i_nuclide].total;
int y = static_cast<int>(y_t);
if (prn() <= y_t - y) ++y;

18
src/physics_mg.cpp
View file

@ -60,7 +60,7 @@ sample_reaction(Particle* p)
// If survival biasing is being used, the following subroutine adjusts the
// weight of the particle. Otherwise, it checks to see if absorption occurs.
if (p->material_xs_.absorption > 0.) {
if (p->macro_xs_.absorption > 0.) {
absorption(p);
} else {
p->wgt_absorb_ = 0.;
@ -113,7 +113,7 @@ create_fission_sites(Particle* p, Particle::Bank* bank_array, int64_t* size_bank
// Determine the expected number of neutrons produced
double nu_t = p->wgt_ / simulation::keff * weight *
p->material_xs_.nu_fission / p->material_xs_.total;
p->macro_xs_.nu_fission / p->macro_xs_.total;
// Sample the number of neutrons produced
int nu = static_cast<int>(nu_t);
@ -204,7 +204,7 @@ absorption(Particle* p)
if (settings::survival_biasing) {
// Determine weight absorbed in survival biasing
p->wgt_absorb_ = p->wgt_ *
p->material_xs_.absorption / p->material_xs_.total;
p->macro_xs_.absorption / p->macro_xs_.total;
// Adjust weight of particle by the probability of absorption
p->wgt_ -= p->wgt_absorb_;
@ -213,14 +213,14 @@ absorption(Particle* p)
// Score implicit absorpion estimate of keff
#pragma omp atomic
global_tally_absorption += p->wgt_absorb_ *
p->material_xs_.nu_fission /
p->material_xs_.absorption;
p->macro_xs_.nu_fission /
p->macro_xs_.absorption;
} else {
if (p->material_xs_.absorption >
prn() * p->material_xs_.total) {
if (p->macro_xs_.absorption >
prn() * p->macro_xs_.total) {
#pragma omp atomic
global_tally_absorption += p->wgt_ * p->material_xs_.nu_fission /
p->material_xs_.absorption;
global_tally_absorption += p->wgt_ * p->macro_xs_.nu_fission /
p->macro_xs_.absorption;
p->alive_ = false;
p->event_ = EVENT_ABSORB;
}

142
src/tallies/derivative.cpp
View file

@ -221,12 +221,12 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
switch (score_bin) {
case SCORE_TOTAL:
if (i_nuclide == -1 && p->material_xs_.total > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.total > 0.0) {
score *= flux_deriv
+ p->micro_xs_[deriv.diff_nuclide].total
/ p->material_xs_.total;
+ p->neutron_xs_[deriv.diff_nuclide].total
/ p->macro_xs_.total;
} else if (i_nuclide == deriv.diff_nuclide
&& p->micro_xs_[i_nuclide].total) {
&& p->neutron_xs_[i_nuclide].total) {
score *= flux_deriv + 1. / atom_density;
} else {
score *= flux_deriv;
@ -234,13 +234,13 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_SCATTER:
if (i_nuclide == -1 && (p->material_xs_.total
- p->material_xs_.absorption) > 0.0) {
if (i_nuclide == -1 && (p->macro_xs_.total
- p->macro_xs_.absorption) > 0.0) {
score *= flux_deriv
+ (p->micro_xs_[deriv.diff_nuclide].total
- p->micro_xs_[deriv.diff_nuclide].absorption)
/ (p->material_xs_.total
- p->material_xs_.absorption);
+ (p->neutron_xs_[deriv.diff_nuclide].total
- p->neutron_xs_[deriv.diff_nuclide].absorption)
/ (p->macro_xs_.total
- p->macro_xs_.absorption);
} else if (i_nuclide == deriv.diff_nuclide) {
score *= flux_deriv + 1. / atom_density;
} else {
@ -249,12 +249,12 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_ABSORPTION:
if (i_nuclide == -1 && p->material_xs_.absorption > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.absorption > 0.0) {
score *= flux_deriv
+ p->micro_xs_[deriv.diff_nuclide].absorption
/ p->material_xs_.absorption;
+ p->neutron_xs_[deriv.diff_nuclide].absorption
/ p->macro_xs_.absorption;
} else if (i_nuclide == deriv.diff_nuclide
&& p->micro_xs_[i_nuclide].absorption) {
&& p->neutron_xs_[i_nuclide].absorption) {
score *= flux_deriv + 1. / atom_density;
} else {
score *= flux_deriv;
@ -262,12 +262,12 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_FISSION:
if (i_nuclide == -1 && p->material_xs_.fission > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.fission > 0.0) {
score *= flux_deriv
+ p->micro_xs_[deriv.diff_nuclide].fission
/ p->material_xs_.fission;
+ p->neutron_xs_[deriv.diff_nuclide].fission
/ p->macro_xs_.fission;
} else if (i_nuclide == deriv.diff_nuclide
&& p->micro_xs_[i_nuclide].fission) {
&& p->neutron_xs_[i_nuclide].fission) {
score *= flux_deriv + 1. / atom_density;
} else {
score *= flux_deriv;
@ -275,12 +275,12 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_NU_FISSION:
if (i_nuclide == -1 && p->material_xs_.nu_fission > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.nu_fission > 0.0) {
score *= flux_deriv
+ p->micro_xs_[deriv.diff_nuclide].nu_fission
/ p->material_xs_.nu_fission;
+ p->neutron_xs_[deriv.diff_nuclide].nu_fission
/ p->macro_xs_.nu_fission;
} else if (i_nuclide == deriv.diff_nuclide
&& p->micro_xs_[i_nuclide].nu_fission) {
&& p->neutron_xs_[i_nuclide].nu_fission) {
score *= flux_deriv + 1. / atom_density;
} else {
score *= flux_deriv;
@ -331,64 +331,64 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
switch (score_bin) {
case SCORE_TOTAL:
if (p->micro_xs_[p->event_nuclide_].total) {
if (p->neutron_xs_[p->event_nuclide_].total) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + (dsig_s + dsig_a) * material.atom_density_(i)
/ p->material_xs_.total;
/ p->macro_xs_.total;
} else {
score *= flux_deriv;
}
break;
case SCORE_SCATTER:
if (p->micro_xs_[p->event_nuclide_].total
- p->micro_xs_[p->event_nuclide_].absorption) {
if (p->neutron_xs_[p->event_nuclide_].total
- p->neutron_xs_[p->event_nuclide_].absorption) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + dsig_s * material.atom_density_(i)
/ (p->material_xs_.total
- p->material_xs_.absorption);
/ (p->macro_xs_.total
- p->macro_xs_.absorption);
} else {
score *= flux_deriv;
}
break;
case SCORE_ABSORPTION:
if (p->micro_xs_[p->event_nuclide_].absorption) {
if (p->neutron_xs_[p->event_nuclide_].absorption) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + dsig_a * material.atom_density_(i)
/ p->material_xs_.absorption;
/ p->macro_xs_.absorption;
} else {
score *= flux_deriv;
}
break;
case SCORE_FISSION:
if (p->micro_xs_[p->event_nuclide_].fission) {
if (p->neutron_xs_[p->event_nuclide_].fission) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + dsig_f * material.atom_density_(i)
/ p->material_xs_.fission;
/ p->macro_xs_.fission;
} else {
score *= flux_deriv;
}
break;
case SCORE_NU_FISSION:
if (p->micro_xs_[p->event_nuclide_].fission) {
double nu = p->micro_xs_[p->event_nuclide_].nu_fission
/ p->micro_xs_[p->event_nuclide_].fission;
if (p->neutron_xs_[p->event_nuclide_].fission) {
double nu = p->neutron_xs_[p->event_nuclide_].nu_fission
/ p->neutron_xs_[p->event_nuclide_].fission;
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + nu * dsig_f * material.atom_density_(i)
/ p->material_xs_.nu_fission;
/ p->macro_xs_.nu_fission;
} else {
score *= flux_deriv;
}
@ -413,141 +413,141 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
switch (score_bin) {
case SCORE_TOTAL:
if (i_nuclide == -1 && p->material_xs_.total > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.total > 0.0) {
double cum_dsig = 0;
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto i_nuc = material.nuclide_[i];
const auto& nuc {*data::nuclides[i_nuc]};
if (multipole_in_range(&nuc, p->E_last_)
&& p->micro_xs_[i_nuc].total) {
&& p->neutron_xs_[i_nuc].total) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
cum_dsig += (dsig_s + dsig_a) * material.atom_density_(i);
}
}
score *= flux_deriv + cum_dsig / p->material_xs_.total;
} else if (p->micro_xs_[i_nuclide].total) {
score *= flux_deriv + cum_dsig / p->macro_xs_.total;
} else if (p->neutron_xs_[i_nuclide].total) {
const auto& nuc {*data::nuclides[i_nuclide]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv
+ (dsig_s + dsig_a) / p->micro_xs_[i_nuclide].total;
+ (dsig_s + dsig_a) / p->neutron_xs_[i_nuclide].total;
} else {
score *= flux_deriv;
}
break;
case SCORE_SCATTER:
if (i_nuclide == -1 && (p->material_xs_.total
- p->material_xs_.absorption)) {
if (i_nuclide == -1 && (p->macro_xs_.total
- p->macro_xs_.absorption)) {
double cum_dsig = 0;
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto i_nuc = material.nuclide_[i];
const auto& nuc {*data::nuclides[i_nuc]};
if (multipole_in_range(&nuc, p->E_last_)
&& (p->micro_xs_[i_nuc].total
- p->micro_xs_[i_nuc].absorption)) {
&& (p->neutron_xs_[i_nuc].total
- p->neutron_xs_[i_nuc].absorption)) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
cum_dsig += dsig_s * material.atom_density_(i);
}
}
score *= flux_deriv + cum_dsig / (p->material_xs_.total
- p->material_xs_.absorption);
} else if (p->micro_xs_[i_nuclide].total
- p->micro_xs_[i_nuclide].absorption) {
score *= flux_deriv + cum_dsig / (p->macro_xs_.total
- p->macro_xs_.absorption);
} else if (p->neutron_xs_[i_nuclide].total
- p->neutron_xs_[i_nuclide].absorption) {
const auto& nuc {*data::nuclides[i_nuclide]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv + dsig_s / (p->micro_xs_[i_nuclide].total
- p->micro_xs_[i_nuclide].absorption);
score *= flux_deriv + dsig_s / (p->neutron_xs_[i_nuclide].total
- p->neutron_xs_[i_nuclide].absorption);
} else {
score *= flux_deriv;
}
break;
case SCORE_ABSORPTION:
if (i_nuclide == -1 && p->material_xs_.absorption > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.absorption > 0.0) {
double cum_dsig = 0;
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto i_nuc = material.nuclide_[i];
const auto& nuc {*data::nuclides[i_nuc]};
if (multipole_in_range(&nuc, p->E_last_)
&& p->micro_xs_[i_nuc].absorption) {
&& p->neutron_xs_[i_nuc].absorption) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
cum_dsig += dsig_a * material.atom_density_(i);
}
}
score *= flux_deriv + cum_dsig / p->material_xs_.absorption;
} else if (p->micro_xs_[i_nuclide].absorption) {
score *= flux_deriv + cum_dsig / p->macro_xs_.absorption;
} else if (p->neutron_xs_[i_nuclide].absorption) {
const auto& nuc {*data::nuclides[i_nuclide]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv
+ dsig_a / p->micro_xs_[i_nuclide].absorption;
+ dsig_a / p->neutron_xs_[i_nuclide].absorption;
} else {
score *= flux_deriv;
}
break;
case SCORE_FISSION:
if (i_nuclide == -1 && p->material_xs_.fission > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.fission > 0.0) {
double cum_dsig = 0;
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto i_nuc = material.nuclide_[i];
const auto& nuc {*data::nuclides[i_nuc]};
if (multipole_in_range(&nuc, p->E_last_)
&& p->micro_xs_[i_nuc].fission) {
&& p->neutron_xs_[i_nuc].fission) {
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
cum_dsig += dsig_f * material.atom_density_(i);
}
}
score *= flux_deriv + cum_dsig / p->material_xs_.fission;
} else if (p->micro_xs_[i_nuclide].fission) {
score *= flux_deriv + cum_dsig / p->macro_xs_.fission;
} else if (p->neutron_xs_[i_nuclide].fission) {
const auto& nuc {*data::nuclides[i_nuclide]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv
+ dsig_f / p->micro_xs_[i_nuclide].fission;
+ dsig_f / p->neutron_xs_[i_nuclide].fission;
} else {
score *= flux_deriv;
}
break;
case SCORE_NU_FISSION:
if (i_nuclide == -1 && p->material_xs_.nu_fission > 0.0) {
if (i_nuclide == -1 && p->macro_xs_.nu_fission > 0.0) {
double cum_dsig = 0;
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto i_nuc = material.nuclide_[i];
const auto& nuc {*data::nuclides[i_nuc]};
if (multipole_in_range(&nuc, p->E_last_)
&& p->micro_xs_[i_nuc].fission) {
double nu = p->micro_xs_[i_nuc].nu_fission
/ p->micro_xs_[i_nuc].fission;
&& p->neutron_xs_[i_nuc].fission) {
double nu = p->neutron_xs_[i_nuc].nu_fission
/ p->neutron_xs_[i_nuc].fission;
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
cum_dsig += nu * dsig_f * material.atom_density_(i);
}
}
score *= flux_deriv + cum_dsig / p->material_xs_.nu_fission;
} else if (p->micro_xs_[i_nuclide].fission) {
score *= flux_deriv + cum_dsig / p->macro_xs_.nu_fission;
} else if (p->neutron_xs_[i_nuclide].fission) {
const auto& nuc {*data::nuclides[i_nuclide]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);
score *= flux_deriv
+ dsig_f / p->micro_xs_[i_nuclide].fission;
+ dsig_f / p->neutron_xs_[i_nuclide].fission;
} else {
score *= flux_deriv;
}
@ -582,7 +582,7 @@ score_track_derivative(const Particle* p, double distance)
// phi is proportional to e^(-Sigma_tot * dist)
// (1 / phi) * (d_phi / d_rho) = - (d_Sigma_tot / d_rho) * dist
// (1 / phi) * (d_phi / d_rho) = - Sigma_tot / rho * dist
deriv.flux_deriv -= distance * p->material_xs_.total
deriv.flux_deriv -= distance * p->macro_xs_.total
/ material.density_gpcc_;
break;
@ -591,7 +591,7 @@ score_track_derivative(const Particle* p, double distance)
// (1 / phi) * (d_phi / d_N) = - (d_Sigma_tot / d_N) * dist
// (1 / phi) * (d_phi / d_N) = - sigma_tot * dist
deriv.flux_deriv -= distance
* p->micro_xs_[deriv.diff_nuclide].total;
* p->neutron_xs_[deriv.diff_nuclide].total;
break;
case DIFF_TEMPERATURE:
@ -660,7 +660,7 @@ void score_collision_derivative(const Particle* p)
// phi is proportional to Sigma_s
// (1 / phi) * (d_phi / d_T) = (d_Sigma_s / d_T) / Sigma_s
// (1 / phi) * (d_phi / d_T) = (d_sigma_s / d_T) / sigma_s
const auto& micro_xs {p->micro_xs_[i_nuc]};
const auto& micro_xs {p->neutron_xs_[i_nuc]};
double dsig_s, dsig_a, dsig_f;
std::tie(dsig_s, dsig_a, dsig_f)
= nuc.multipole_->evaluate_deriv(p->E_last_, p->sqrtkT_);

170
src/tallies/tally_scoring.cpp
View file

@ -349,7 +349,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
if (p->type_ == Particle::Type::neutron ||
p->type_ == Particle::Type::photon) {
score *= flux / p->material_xs_.total;
score *= flux / p->macro_xs_.total;
} else {
score = 0.;
}
@ -373,9 +373,9 @@ score_general_ce(Particle* p, int i_tally, int start_index,
} else {
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].total * atom_density * flux;
score = p->neutron_xs_[i_nuclide].total * atom_density * flux;
} else {
score = p->material_xs_.total * flux;
score = p->macro_xs_.total * flux;
}
}
break;
@ -393,7 +393,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
} else {
score = p->wgt_last_;
}
score *= flux / p->material_xs_.total;
score *= flux / p->macro_xs_.total;
} else {
score = flux;
}
@ -411,11 +411,11 @@ score_general_ce(Particle* p, int i_tally, int start_index,
score = p->wgt_last_ * flux;
} else {
if (i_nuclide >= 0) {
score = (p->micro_xs_[i_nuclide].total
- p->micro_xs_[i_nuclide].absorption) * atom_density * flux;
score = (p->neutron_xs_[i_nuclide].total
- p->neutron_xs_[i_nuclide].absorption) * atom_density * flux;
} else {
score = (p->material_xs_.total
- p->material_xs_.absorption) * flux;
score = (p->macro_xs_.total
- p->macro_xs_.absorption) * flux;
}
}
break;
@ -458,26 +458,26 @@ score_general_ce(Particle* p, int i_tally, int start_index,
}
} else {
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].absorption * atom_density
score = p->neutron_xs_[i_nuclide].absorption * atom_density
* flux;
} else {
score = p->material_xs_.absorption * flux;
score = p->macro_xs_.absorption * flux;
}
}
break;
case SCORE_FISSION:
if (p->material_xs_.absorption == 0) continue;
if (p->macro_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
// fission
if (p->micro_xs_[p->event_nuclide_].absorption > 0) {
if (p->neutron_xs_[p->event_nuclide_].absorption > 0) {
score = p->wgt_absorb_
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
} else {
score = 0.;
}
@ -488,21 +488,21 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// weight entering the collision as the estimate for the fission
// reaction rate
score = p->wgt_last_
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
} else {
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].fission * atom_density * flux;
score = p->neutron_xs_[i_nuclide].fission * atom_density * flux;
} else {
score = p->material_xs_.fission * flux;
score = p->macro_xs_.fission * flux;
}
}
break;
case SCORE_NU_FISSION:
if (p->material_xs_.absorption == 0) continue;
if (p->macro_xs_.absorption == 0) continue;
if (tally.estimator_ == ESTIMATOR_ANALOG) {
if (settings::survival_biasing || p->fission_) {
if (tally.energyout_filter_ != C_NONE) {
@ -516,10 +516,10 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// No fission events occur if survival biasing is on -- need to
// calculate fraction of absorptions that would have resulted in
// nu-fission
if (p->micro_xs_[p->event_nuclide_].absorption > 0) {
if (p->neutron_xs_[p->event_nuclide_].absorption > 0) {
score = p->wgt_absorb_
* p->micro_xs_[p->event_nuclide_].nu_fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].nu_fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
} else {
score = 0.;
}
@ -535,17 +535,17 @@ score_general_ce(Particle* p, int i_tally, int start_index,
}
} else {
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].nu_fission * atom_density
score = p->neutron_xs_[i_nuclide].nu_fission * atom_density
* flux;
} else {
score = p->material_xs_.nu_fission * flux;
score = p->macro_xs_.nu_fission * flux;
}
}
break;
case SCORE_PROMPT_NU_FISSION:
if (p->material_xs_.absorption == 0) continue;
if (p->macro_xs_.absorption == 0) continue;
if (tally.estimator_ == ESTIMATOR_ANALOG) {
if (settings::survival_biasing || p->fission_) {
if (tally.energyout_filter_ != C_NONE) {
@ -559,12 +559,12 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// No fission events occur if survival biasing is on -- need to
// calculate fraction of absorptions that would have resulted in
// prompt-nu-fission
if (p->micro_xs_[p->event_nuclide_].absorption > 0) {
if (p->neutron_xs_[p->event_nuclide_].absorption > 0) {
score = p->wgt_absorb_
* p->micro_xs_[p->event_nuclide_].fission
* p->neutron_xs_[p->event_nuclide_].fission
* data::nuclides[p->event_nuclide_]
->nu(E, ReactionProduct::EmissionMode::prompt)
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
} else {
score = 0.;
}
@ -583,7 +583,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
}
} else {
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].fission
score = p->neutron_xs_[i_nuclide].fission
* data::nuclides[i_nuclide]
->nu(E, ReactionProduct::EmissionMode::prompt)
* atom_density * flux;
@ -595,7 +595,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto j_nuclide = material.nuclide_[i];
auto atom_density = material.atom_density_(i);
score += p->micro_xs_[j_nuclide].fission
score += p->neutron_xs_[j_nuclide].fission
* data::nuclides[j_nuclide]
->nu(E, ReactionProduct::EmissionMode::prompt)
* atom_density * flux;
@ -607,7 +607,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
case SCORE_DELAYED_NU_FISSION:
if (p->material_xs_.absorption == 0) continue;
if (p->macro_xs_.absorption == 0) continue;
if (tally.estimator_ == ESTIMATOR_ANALOG) {
if (settings::survival_biasing || p->fission_) {
if (tally.energyout_filter_ != C_NONE) {
@ -621,7 +621,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// No fission events occur if survival biasing is on -- need to
// calculate fraction of absorptions that would have resulted in
// delayed-nu-fission
if (p->micro_xs_[p->event_nuclide_].absorption > 0
if (p->neutron_xs_[p->event_nuclide_].absorption > 0
&& data::nuclides[p->event_nuclide_]->fissionable_) {
if (tally.delayedgroup_filter_ != C_NONE) {
auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_];
@ -634,8 +634,8 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto yield = data::nuclides[p->event_nuclide_]
->nu(E, ReactionProduct::EmissionMode::delayed, d);
score = p->wgt_absorb_ * yield
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
score_fission_delayed_dg(i_tally, d_bin, score,
score_index);
}
@ -645,10 +645,10 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// by multiplying the absorbed weight by the fraction of the
// delayed-nu-fission xs to the absorption xs
score = p->wgt_absorb_
* p->micro_xs_[p->event_nuclide_].fission
* p->neutron_xs_[p->event_nuclide_].fission
* data::nuclides[p->event_nuclide_]
->nu(E, ReactionProduct::EmissionMode::delayed)
/ p->micro_xs_[p->event_nuclide_].absorption *flux;
/ p->neutron_xs_[p->event_nuclide_].absorption *flux;
}
}
} else {
@ -694,7 +694,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto d = filt.groups_[d_bin];
auto yield = data::nuclides[i_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed, d);
score = p->micro_xs_[i_nuclide].fission * yield
score = p->neutron_xs_[i_nuclide].fission * yield
* atom_density * flux;
score_fission_delayed_dg(i_tally, d_bin, score, score_index);
}
@ -702,7 +702,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
} else {
// If the delayed group filter is not present, compute the score
// by multiplying the delayed-nu-fission macro xs by the flux
score = p->micro_xs_[i_nuclide].fission
score = p->neutron_xs_[i_nuclide].fission
* data::nuclides[i_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed)
* atom_density * flux;
@ -723,7 +723,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto d = filt.groups_[d_bin];
auto yield = data::nuclides[j_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed, d);
score = p->micro_xs_[j_nuclide].fission * yield
score = p->neutron_xs_[j_nuclide].fission * yield
* atom_density * flux;
score_fission_delayed_dg(i_tally, d_bin, score,
score_index);
@ -738,7 +738,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto j_nuclide = material.nuclide_[i];
auto atom_density = material.atom_density_(i);
score += p->micro_xs_[j_nuclide].fission
score += p->neutron_xs_[j_nuclide].fission
* data::nuclides[j_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed)
* atom_density * flux;
@ -751,14 +751,14 @@ score_general_ce(Particle* p, int i_tally, int start_index,
case SCORE_DECAY_RATE:
if (p->material_xs_.absorption == 0) continue;
if (p->macro_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_]};
if (p->micro_xs_[p->event_nuclide_].absorption > 0
if (p->neutron_xs_[p->event_nuclide_].absorption > 0
&& nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
if (tally.delayedgroup_filter_ != C_NONE) {
@ -773,8 +773,8 @@ score_general_ce(Particle* p, int i_tally, int start_index,
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
auto rate = rxn.products_[d].decay_rate_;
score = p->wgt_absorb_ * yield
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption
* rate * flux;
score_fission_delayed_dg(i_tally, d_bin, score,
score_index);
@ -796,8 +796,8 @@ score_general_ce(Particle* p, int i_tally, int start_index,
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
auto rate = rxn.products_[d+1].decay_rate_;
score += rate * p->wgt_absorb_
* p->micro_xs_[p->event_nuclide_].fission * yield
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission * yield
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
}
}
@ -854,7 +854,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto yield
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
auto rate = rxn.products_[d].decay_rate_;
score = p->micro_xs_[i_nuclide].fission * yield * flux
score = p->neutron_xs_[i_nuclide].fission * yield * flux
* atom_density * rate;
score_fission_delayed_dg(i_tally, d_bin, score, score_index);
}
@ -870,7 +870,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto yield
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
auto rate = rxn.products_[d+1].decay_rate_;
score += p->micro_xs_[i_nuclide].fission * flux
score += p->neutron_xs_[i_nuclide].fission * flux
* yield * atom_density * rate;
}
}
@ -894,7 +894,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto yield
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d);
auto rate = rxn.products_[d].decay_rate_;
score = p->micro_xs_[j_nuclide].fission * yield
score = p->neutron_xs_[j_nuclide].fission * yield
* flux * atom_density * rate;
score_fission_delayed_dg(i_tally, d_bin, score,
score_index);
@ -922,7 +922,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto yield
= nuc.nu(E, ReactionProduct::EmissionMode::delayed, d+1);
auto rate = rxn.products_[d+1].decay_rate_;
score += p->micro_xs_[j_nuclide].fission
score += p->neutron_xs_[j_nuclide].fission
* yield * atom_density * flux * rate;
}
}
@ -935,7 +935,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
case SCORE_KAPPA_FISSION:
if (p->material_xs_.absorption == 0.) continue;
if (p->macro_xs_.absorption == 0.) continue;
score = 0.;
// Kappa-fission values are determined from the Q-value listed for the
// fission cross section.
@ -945,12 +945,12 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// calculate fraction of absorptions that would have resulted in
// fission scaled by the Q-value
const auto& nuc {*data::nuclides[p->event_nuclide_]};
if (p->micro_xs_[p->event_nuclide_].absorption > 0
if (p->neutron_xs_[p->event_nuclide_].absorption > 0
&& nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
score = p->wgt_absorb_ * rxn.q_value_
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
} else {
// Skip any non-absorption events
@ -959,12 +959,12 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// weight entering the collision as the estimate for the fission
// reaction rate
const auto& nuc {*data::nuclides[p->event_nuclide_]};
if (p->micro_xs_[p->event_nuclide_].absorption > 0
if (p->neutron_xs_[p->event_nuclide_].absorption > 0
&& nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
score = p->wgt_last_ * rxn.q_value_
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
}
} else {
@ -972,7 +972,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
const auto& nuc {*data::nuclides[i_nuclide]};
if (nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
score = rxn.q_value_ * p->micro_xs_[i_nuclide].fission
score = rxn.q_value_ * p->neutron_xs_[i_nuclide].fission
* atom_density * flux;
}
} else {
@ -984,7 +984,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
const auto& nuc {*data::nuclides[j_nuclide]};
if (nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
score += rxn.q_value_ * p->micro_xs_[j_nuclide].fission
score += rxn.q_value_ * p->neutron_xs_[j_nuclide].fission
* atom_density * flux;
}
}
@ -1007,9 +1007,9 @@ score_general_ce(Particle* p, int i_tally, int start_index,
score = p->wgt_last_ * flux;
} else {
if (i_nuclide >= 0) {
if (p->micro_xs_[i_nuclide].elastic == CACHE_INVALID)
if (p->neutron_xs_[i_nuclide].elastic == CACHE_INVALID)
data::nuclides[i_nuclide]->calculate_elastic_xs(*p);
score = p->micro_xs_[i_nuclide].elastic * atom_density * flux;
score = p->neutron_xs_[i_nuclide].elastic * atom_density * flux;
} else {
score = 0.;
if (p->material_ != MATERIAL_VOID) {
@ -1017,9 +1017,9 @@ score_general_ce(Particle* p, int i_tally, int start_index,
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto j_nuclide = material.nuclide_[i];
auto atom_density = material.atom_density_(i);
if (p->micro_xs_[j_nuclide].elastic == CACHE_INVALID)
if (p->neutron_xs_[j_nuclide].elastic == CACHE_INVALID)
data::nuclides[j_nuclide]->calculate_elastic_xs(*p);
score += p->micro_xs_[j_nuclide].elastic * atom_density
score += p->neutron_xs_[j_nuclide].elastic * atom_density
* flux;
}
}
@ -1030,7 +1030,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
case SCORE_FISS_Q_PROMPT:
case SCORE_FISS_Q_RECOV:
//continue;
if (p->material_xs_.absorption == 0.) continue;
if (p->macro_xs_.absorption == 0.) continue;
score = 0.;
if (tally.estimator_ == ESTIMATOR_ANALOG) {
if (settings::survival_biasing) {
@ -1038,7 +1038,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// calculate fraction of absorptions that would have resulted in
// fission scaled by the Q-value
const auto& nuc {*data::nuclides[p->event_nuclide_]};
if (p->micro_xs_[p->event_nuclide_].absorption > 0) {
if (p->neutron_xs_[p->event_nuclide_].absorption > 0) {
double q_value = 0.;
if (score_bin == SCORE_FISS_Q_PROMPT) {
if (nuc.fission_q_prompt_)
@ -1048,8 +1048,8 @@ score_general_ce(Particle* p, int i_tally, int start_index,
q_value = (*nuc.fission_q_recov_)(p->E_last_);
}
score = p->wgt_absorb_ * q_value
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
} else {
// Skip any non-absorption events
@ -1058,7 +1058,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
// weight entering the collision as the estimate for the fission
// reaction rate
const auto& nuc {*data::nuclides[p->event_nuclide_]};
if (p->micro_xs_[p->event_nuclide_].absorption > 0) {
if (p->neutron_xs_[p->event_nuclide_].absorption > 0) {
double q_value = 0.;
if (score_bin == SCORE_FISS_Q_PROMPT) {
if (nuc.fission_q_prompt_)
@ -1068,8 +1068,8 @@ score_general_ce(Particle* p, int i_tally, int start_index,
q_value = (*nuc.fission_q_recov_)(p->E_last_);
}
score = p->wgt_last_ * q_value
* p->micro_xs_[p->event_nuclide_].fission
/ p->micro_xs_[p->event_nuclide_].absorption * flux;
* p->neutron_xs_[p->event_nuclide_].fission
/ p->neutron_xs_[p->event_nuclide_].absorption * flux;
}
}
} else {
@ -1083,7 +1083,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
if (nuc.fission_q_recov_)
q_value = (*nuc.fission_q_recov_)(p->E_last_);
}
score = q_value * p->micro_xs_[i_nuclide].fission
score = q_value * p->neutron_xs_[i_nuclide].fission
* atom_density * flux;
} else {
if (p->material_ != MATERIAL_VOID) {
@ -1100,7 +1100,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
if (nuc.fission_q_recov_)
q_value = (*nuc.fission_q_recov_)(p->E_last_);
}
score += q_value * p->micro_xs_[j_nuclide].fission
score += q_value * p->neutron_xs_[j_nuclide].fission
* atom_density * flux;
}
}
@ -1130,7 +1130,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
case N_4N: m = 5; break;
}
if (i_nuclide >= 0) {
score = p->micro_xs_[i_nuclide].reaction[m] * atom_density
score = p->neutron_xs_[i_nuclide].reaction[m] * atom_density
* flux;
} else {
score = 0.;
@ -1139,7 +1139,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
for (auto i = 0; i < material.nuclide_.size(); ++i) {
auto j_nuclide = material.nuclide_[i];
auto atom_density = material.atom_density_(i);
score += p->micro_xs_[j_nuclide].reaction[m]
score += p->neutron_xs_[j_nuclide].reaction[m]
* atom_density * flux;
}
}
@ -1164,10 +1164,10 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto m = nuc.reaction_index_[score_bin];
if (m == C_NONE) continue;
const auto& rxn {*nuc.reactions_[m]};
auto i_temp = p->micro_xs_[i_nuclide].index_temp;
auto i_temp = p->neutron_xs_[i_nuclide].index_temp;
if (i_temp >= 0) { // Can be false due to multipole
auto i_grid = p->micro_xs_[i_nuclide].index_grid;
auto f = p->micro_xs_[i_nuclide].interp_factor;
auto i_grid = p->neutron_xs_[i_nuclide].index_grid;
auto f = p->neutron_xs_[i_nuclide].interp_factor;
const auto& xs {rxn.xs_[i_temp]};
if (i_grid >= xs.threshold) {
score = ((1.0 - f) * xs.value[i_grid-xs.threshold]
@ -1184,10 +1184,10 @@ score_general_ce(Particle* p, int i_tally, int start_index,
auto m = nuc.reaction_index_[score_bin];
if (m == C_NONE) continue;
const auto& rxn {*nuc.reactions_[m]};
auto i_temp = p->micro_xs_[j_nuclide].index_temp;
auto i_temp = p->neutron_xs_[j_nuclide].index_temp;
if (i_temp >= 0) { // Can be false due to multipole
auto i_grid = p->micro_xs_[j_nuclide].index_grid;
auto f = p->micro_xs_[j_nuclide].interp_factor;
auto i_grid = p->neutron_xs_[j_nuclide].index_grid;
auto f = p->neutron_xs_[j_nuclide].interp_factor;
const auto& xs {rxn.xs_[i_temp]};
if (i_grid >= xs.threshold) {
score += ((1.0 - f) * xs.value[i_grid-xs.threshold]
@ -1291,7 +1291,7 @@ score_general_mg(const Particle* p, int i_tally, int start_index,
} else {
score = p->wgt_last_;
}
score *= flux / p->material_xs_.total;
score *= flux / p->macro_xs_.total;
} else {
score = flux;
}
@ -1320,7 +1320,7 @@ score_general_mg(const Particle* p, int i_tally, int start_index,
score = get_nuclide_xs(i_nuclide, MG_GET_XS_TOTAL, p_g)
* atom_density * flux;
} else {
score = p->material_xs_.total * flux;
score = p->macro_xs_.total * flux;
}
}
break;
@ -1438,7 +1438,7 @@ score_general_mg(const Particle* p, int i_tally, int start_index,
score = atom_density * flux
* get_nuclide_xs(i_nuclide, MG_GET_XS_ABSORPTION, p_g);
} else {
score = p->material_xs_.absorption * flux;
score = p->macro_xs_.absorption * flux;
}
}
break;
@ -2127,9 +2127,9 @@ void score_collision_tally(Particle* p)
// Determine the collision estimate of the flux
double flux;
if (!settings::survival_biasing) {
flux = p->wgt_last_ / p->material_xs_.total;
flux = p->wgt_last_ / p->macro_xs_.total;
} else {
flux = (p->wgt_last_ + p->wgt_absorb_) / p->material_xs_.total;
flux = (p->wgt_last_ + p->wgt_absorb_) / p->macro_xs_.total;
}
for (auto i_tally : model::active_collision_tallies) {