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

Get rid of wgt_absorb in ParticleData

Browse source
This commit is contained in:
Paul Romano 2022-01-16 15:13:36 -06:00
parent f741568d06
commit f18caa1712
4 changed files with 37 additions and 33 deletions
Download patch file Download diff file Expand all files Collapse all files

3
include/openmc/particle_data.h
View file

@ -244,7 +244,6 @@ private:
Position r_last_; //!< previous coordinates
Direction u_last_; //!< previous direction coordinates
double wgt_last_ {1.0}; //!< pre-collision particle weight
double wgt_absorb_ {0.0}; //!< weight absorbed for survival biasing
// What event took place
bool fission_ {false}; //!< did particle cause implicit fission
@ -375,8 +374,6 @@ public:
const Position& u_last() const { return u_last_; }
double& wgt_last() { return wgt_last_; }
const double& wgt_last() const { return wgt_last_; }
double& wgt_absorb() { return wgt_absorb_; }
const double& wgt_absorb() const { return wgt_absorb_; }
bool& fission() { return fission_; }
TallyEvent& event() { return event_; }

10
src/physics.cpp
View file

@ -134,8 +134,6 @@ void sample_neutron_reaction(Particle& p)
if (p.neutron_xs(i_nuclide).absorption > 0.0) {
absorption(p, i_nuclide);
} else {
p.wgt_absorb() = 0.0;
}
if (!p.alive())
return;
@ -626,15 +624,15 @@ void absorption(Particle& p, int i_nuclide)
{
if (settings::survival_biasing) {
// Determine weight absorbed in survival biasing
p.wgt_absorb() = p.wgt() * p.neutron_xs(i_nuclide).absorption /
p.neutron_xs(i_nuclide).total;
double 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();
p.wgt() -= wgt_absorb;
// Score implicit absorption estimate of keff
if (settings::run_mode == RunMode::EIGENVALUE) {
p.keff_tally_absorption() += p.wgt_absorb() *
p.keff_tally_absorption() += wgt_absorb *
p.neutron_xs(i_nuclide).nu_fission /
p.neutron_xs(i_nuclide).absorption;
}

8
src/physics_mg.cpp
View file

@ -55,8 +55,6 @@ void sample_reaction(Particle& p)
// weight of the particle. Otherwise, it checks to see if absorption occurs.
if (p.macro_xs().absorption > 0.) {
absorption(p);
} else {
p.wgt_absorb() = 0.;
}
if (!p.alive())
return;
@ -216,14 +214,14 @@ void absorption(Particle& p)
{
if (settings::survival_biasing) {
// Determine weight absorbed in survival biasing
p.wgt_absorb() = p.wgt() * p.macro_xs().absorption / p.macro_xs().total;
double wgt_absorb = p.wgt() * p.macro_xs().absorption / p.macro_xs().total;
// Adjust weight of particle by the probability of absorption
p.wgt() -= p.wgt_absorb();
p.wgt() -= wgt_absorb;
// Score implicit absorpion estimate of keff
p.keff_tally_absorption() +=
p.wgt_absorb() * p.macro_xs().nu_fission / p.macro_xs().absorption;
wgt_absorb * p.macro_xs().nu_fission / p.macro_xs().absorption;
} else {
if (p.macro_xs().absorption > prn(p.current_seed()) * p.macro_xs().total) {
p.keff_tally_absorption() +=

49
src/tallies/tally_scoring.cpp
View file

@ -546,6 +546,14 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
// Get the pre-collision energy of the particle.
auto E = p.E_last();
// Determine how much weight was absorbed due to survival biasing
double wgt_absorb = 0.0;
if (tally.estimator_ == TallyEstimator::ANALOG &&
settings::survival_biasing) {
wgt_absorb = p.wgt_last() * p.neutron_xs(p.event_nuclide()).absorption /
p.neutron_xs(p.event_nuclide()).total;
}
using Type = ParticleType;
for (auto i = 0; i < tally.scores_.size(); ++i) {
@ -614,7 +622,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
continue;
// Since only scattering events make it here, again we can use the
// weight entering the collision as the estimator for the reaction rate
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
} else {
if (i_nuclide >= 0) {
if (p.type() == Type::neutron) {
@ -645,7 +653,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
// For scattering production, we need to use the pre-collision weight
// times the yield as the estimate for the number of neutrons exiting a
// reaction with neutrons in the exit channel
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
// Don't waste time on very common reactions we know have multiplicities
// of one.
@ -667,7 +675,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
if (settings::survival_biasing) {
// No absorption events actually occur if survival biasing is on --
// just use weight absorbed in survival biasing
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
} else {
// Skip any event where the particle wasn't absorbed
if (p.event() == TallyEvent::SCATTER)
@ -1240,7 +1248,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
// Check if event MT matches
if (p.event_mt() != ELASTIC)
continue;
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
} else {
if (i_nuclide >= 0) {
if (p.neutron_xs(i_nuclide).elastic == CACHE_INVALID)
@ -1289,7 +1297,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
// Check if the event MT matches
if (p.event_mt() != score_bin)
continue;
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
} else {
int m;
switch (score_bin) {
@ -1405,7 +1413,7 @@ void score_general_ce(Particle& p, int i_tally, int start_index,
// to check if it matches the MT number of the event
if (p.event_mt() != score_bin)
continue;
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
} else {
// Any other cross section has to be calculated on-the-fly
if (score_bin < 2)
@ -1452,10 +1460,13 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
// Set the direction and group to use with get_xs
Direction p_u;
int p_g;
double wgt_absorb = 0.0;
if (tally.estimator_ == TallyEstimator::ANALOG ||
tally.estimator_ == TallyEstimator::COLLISION) {
if (settings::survival_biasing) {
// Determine weight that was absorbed
wgt_absorb = p.wgt_last() * p.neutron_xs(p.event_nuclide()).absorption /
p.neutron_xs(p.event_nuclide()).total;
// Then we either are alive and had a scatter (and so g changed),
// or are dead and g did not change
@ -1563,7 +1574,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
continue;
// Since only scattering events make it here, again we can use the
// weight entering the collision as the estimator for the reaction rate
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
if (i_nuclide >= 0) {
score *= atom_density *
nuc_xs.get_xs(MgxsType::SCATTER_FMU, p.g_last(), &p.g(),
@ -1591,7 +1602,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
// For scattering production, we need to use the pre-collision weight
// times the multiplicity as the estimate for the number of neutrons
// exiting a reaction with neutrons in the exit channel
score = (p.wgt_last() - p.wgt_absorb()) * flux;
score = (p.wgt_last() - wgt_absorb) * flux;
// Since we transport based on material data, the angle selected
// was not selected from the f(mu) for the nuclide. Therefore
// adjust the score by the actual probability for that nuclide.
@ -1617,7 +1628,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
if (settings::survival_biasing) {
// No absorption events actually occur if survival biasing is on --
// just use weight absorbed in survival biasing
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
} else {
// Skip any event where the particle wasn't absorbed
if (p.event() == TallyEvent::SCATTER)
@ -1646,7 +1657,7 @@ void score_general_mg(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
// fission
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
} else {
// Skip any non-absorption events
if (p.event() == TallyEvent::SCATTER)
@ -1686,7 +1697,7 @@ void score_general_mg(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
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
if (i_nuclide >= 0) {
score *= atom_density * nuc_xs.get_xs(MgxsType::NU_FISSION, p_g) /
macro_xs.get_xs(MgxsType::ABSORPTION, p_g);
@ -1733,7 +1744,7 @@ void score_general_mg(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
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
if (i_nuclide >= 0) {
score *= atom_density *
nuc_xs.get_xs(MgxsType::PROMPT_NU_FISSION, p_g) /
@ -1794,7 +1805,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
// Tally each delayed group bin individually
for (auto d_bin = 0; d_bin < filt.n_bins(); ++d_bin) {
auto d = filt.groups()[d_bin] - 1;
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
if (i_nuclide >= 0) {
score *= nuc_xs.get_xs(MgxsType::DELAYED_NU_FISSION, p_g,
nullptr, nullptr, &d) /
@ -1812,7 +1823,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
// 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
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
if (i_nuclide >= 0) {
score *=
nuc_xs.get_xs(MgxsType::DELAYED_NU_FISSION, p_g) / abs_xs;
@ -1909,7 +1920,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
// Tally each delayed group bin individually
for (auto d_bin = 0; d_bin < filt.n_bins(); ++d_bin) {
auto d = filt.groups()[d_bin] - 1;
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
if (i_nuclide >= 0) {
score *= nuc_xs.get_xs(
MgxsType::DECAY_RATE, p_g, nullptr, nullptr, &d) *
@ -1935,14 +1946,14 @@ void score_general_mg(Particle& p, int i_tally, int start_index,
score = 0.;
for (auto d = 0; d < data::mg.num_delayed_groups_; ++d) {
if (i_nuclide >= 0) {
score += p.wgt_absorb() * flux *
score += wgt_absorb * flux *
nuc_xs.get_xs(
MgxsType::DECAY_RATE, p_g, nullptr, nullptr, &d) *
nuc_xs.get_xs(MgxsType::DELAYED_NU_FISSION, p_g,
nullptr, nullptr, &d) /
abs_xs;
} else {
score += p.wgt_absorb() * flux *
score += wgt_absorb * flux *
macro_xs.get_xs(
MgxsType::DECAY_RATE, p_g, nullptr, nullptr, &d) *
macro_xs.get_xs(MgxsType::DELAYED_NU_FISSION, p_g,
@ -2050,7 +2061,7 @@ void score_general_mg(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
// fission scaled by the Q-value
score = p.wgt_absorb() * flux;
score = wgt_absorb * flux;
} else {
// Skip any non-absorption events
if (p.event() == TallyEvent::SCATTER)