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

Fix off-by-one on event_nuclide and diff_nuclide

Browse source
This commit is contained in:
Paul Romano 2019-02-25 10:07:04 -06:00
parent 6ce85a9f0c
commit a2e058c8cf
5 changed files with 80 additions and 92 deletions
Download patch file Download diff file Expand all files Collapse all files

3
src/output.cpp
View file

@ -657,8 +657,7 @@ write_tallies()
case DIFF_NUCLIDE_DENSITY:
tallies_out << " Nuclide density derivative Material "
<< std::to_string(deriv.diff_material) << " Nuclide "
// TODO: off-by-one
<< data::nuclides[deriv.diff_nuclide-1]->name_ << "\n";
<< data::nuclides[deriv.diff_nuclide]->name_ << "\n";
break;
case DIFF_TEMPERATURE:
tallies_out << " Temperature derivative Material "

10
src/physics.cpp
View file

@ -63,10 +63,10 @@ void collision(Particle* p)
std::stringstream msg;
if (static_cast<ParticleType>(p->type) == ParticleType::neutron) {
msg << " " << reaction_name(p->event_MT) << " with " <<
data::nuclides[p->event_nuclide-1]->name_ << ". Energy = " << p->E << " eV.";
data::nuclides[p->event_nuclide]->name_ << ". Energy = " << p->E << " eV.";
} else {
msg << " " << reaction_name(p->event_MT) << " with " <<
data::elements[p->event_nuclide-1].name_ << ". Energy = " << p->E << " eV.";
data::elements[p->event_nuclide].name_ << ". Energy = " << p->E << " eV.";
}
write_message(msg, 1);
}
@ -78,8 +78,7 @@ void sample_neutron_reaction(Particle* p)
int i_nuclide = sample_nuclide(p);
// Save which nuclide particle had collision with
// TODO: off-by-one
p->event_nuclide = i_nuclide + 1;
p->event_nuclide = i_nuclide;
// Create fission bank sites. Note that while a fission reaction is sampled,
// it never actually "happens", i.e. the weight of the particle does not
@ -228,8 +227,7 @@ void sample_photon_reaction(Particle* p)
// Sample element within material
int i_element = sample_element(p);
// TODO: off-by-one
p->event_nuclide = i_element + 1;
p->event_nuclide = i_element;
const auto& micro {simulation::micro_photon_xs[i_element]};
const auto& element {data::elements[i_element]};

2
src/state_point.cpp
View file

@ -125,7 +125,7 @@ openmc_statepoint_write(const char* filename, bool* write_source)
write_dataset(deriv_group, "independent variable", "nuclide_density");
//TODO: off-by-one
write_dataset(deriv_group, "nuclide",
data::nuclides[deriv.diff_nuclide-1]->name_);
data::nuclides[deriv.diff_nuclide]->name_);
} else if (deriv.variable == DIFF_TEMPERATURE) {
write_dataset(deriv_group, "independent variable", "temperature");
} else {

61
src/tallies/derivative.cpp
View file

@ -51,8 +51,7 @@ TallyDerivative::TallyDerivative(pugi::xml_node node)
for (auto i = 0; i < data::nuclides.size(); ++i) {
if (data::nuclides[i]->name_ == nuclide_name) {
found = true;
//TODO: off-by-one
diff_nuclide = i + 1;
diff_nuclide = i;
}
}
if (!found) {
@ -188,7 +187,6 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
// where i is the perturbed nuclide.
case DIFF_NUCLIDE_DENSITY:
//TODO: off-by-one throughout on diff_nuclide
switch (tally.estimator_) {
case ESTIMATOR_ANALOG:
@ -208,7 +206,7 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
// Find the index of the perturbed nuclide.
int i;
for (i = 0; i < material.nuclide_.size(); ++i)
if (material.nuclide_[i] == deriv.diff_nuclide - 1) break;
if (material.nuclide_[i] == deriv.diff_nuclide) break;
score *= flux_deriv + 1. / material.atom_density_(i);
}
break;
@ -225,9 +223,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
case SCORE_TOTAL:
if (i_nuclide == -1 && simulation::material_xs.total > 0.0) {
score *= flux_deriv
+ simulation::micro_xs[deriv.diff_nuclide-1].total
+ simulation::micro_xs[deriv.diff_nuclide].total
/ simulation::material_xs.total;
} else if (i_nuclide == deriv.diff_nuclide-1
} else if (i_nuclide == deriv.diff_nuclide
&& simulation::micro_xs[i_nuclide].total) {
score *= flux_deriv + 1. / atom_density;
} else {
@ -239,11 +237,11 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
if (i_nuclide == -1 && (simulation::material_xs.total
- simulation::material_xs.absorption) > 0.0) {
score *= flux_deriv
+ (simulation::micro_xs[deriv.diff_nuclide-1].total
- simulation::micro_xs[deriv.diff_nuclide-1].absorption)
+ (simulation::micro_xs[deriv.diff_nuclide].total
- simulation::micro_xs[deriv.diff_nuclide].absorption)
/ (simulation::material_xs.total
- simulation::material_xs.absorption);
} else if (i_nuclide == deriv.diff_nuclide-1) {
} else if (i_nuclide == deriv.diff_nuclide) {
score *= flux_deriv + 1. / atom_density;
} else {
score *= flux_deriv;
@ -253,9 +251,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
case SCORE_ABSORPTION:
if (i_nuclide == -1 && simulation::material_xs.absorption > 0.0) {
score *= flux_deriv
+ simulation::micro_xs[deriv.diff_nuclide-1].absorption
+ simulation::micro_xs[deriv.diff_nuclide].absorption
/ simulation::material_xs.absorption;
} else if (i_nuclide == deriv.diff_nuclide-1
} else if (i_nuclide == deriv.diff_nuclide
&& simulation::micro_xs[i_nuclide].absorption) {
score *= flux_deriv + 1. / atom_density;
} else {
@ -266,9 +264,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
case SCORE_FISSION:
if (i_nuclide == -1 && simulation::material_xs.fission > 0.0) {
score *= flux_deriv
+ simulation::micro_xs[deriv.diff_nuclide-1].fission
+ simulation::micro_xs[deriv.diff_nuclide].fission
/ simulation::material_xs.fission;
} else if (i_nuclide == deriv.diff_nuclide-1
} else if (i_nuclide == deriv.diff_nuclide
&& simulation::micro_xs[i_nuclide].fission) {
score *= flux_deriv + 1. / atom_density;
} else {
@ -279,9 +277,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
case SCORE_NU_FISSION:
if (i_nuclide == -1 && simulation::material_xs.nu_fission > 0.0) {
score *= flux_deriv
+ simulation::micro_xs[deriv.diff_nuclide-1].nu_fission
+ simulation::micro_xs[deriv.diff_nuclide].nu_fission
/ simulation::material_xs.nu_fission;
} else if (i_nuclide == deriv.diff_nuclide-1
} else if (i_nuclide == deriv.diff_nuclide
&& simulation::micro_xs[i_nuclide].nu_fission) {
score *= flux_deriv + 1. / atom_density;
} else {
@ -322,9 +320,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
// Find the index of the event nuclide.
int i;
for (i = 0; i < material.nuclide_.size(); ++i)
if (material.nuclide_[i] == p->event_nuclide-1) break;
if (material.nuclide_[i] == p->event_nuclide) break;
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (!multipole_in_range(&nuc, p->last_E)) {
score *= flux_deriv;
break;
@ -333,7 +331,7 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
switch (score_bin) {
case SCORE_TOTAL:
if (simulation::micro_xs[p->event_nuclide-1].total) {
if (simulation::micro_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->last_E, p->sqrtkT);
@ -345,8 +343,8 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_SCATTER:
if (simulation::micro_xs[p->event_nuclide-1].total
- simulation::micro_xs[p->event_nuclide-1].absorption) {
if (simulation::micro_xs[p->event_nuclide].total
- simulation::micro_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->last_E, p->sqrtkT);
@ -359,7 +357,7 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_ABSORPTION:
if (simulation::micro_xs[p->event_nuclide-1].absorption) {
if (simulation::micro_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->last_E, p->sqrtkT);
@ -371,7 +369,7 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_FISSION:
if (simulation::micro_xs[p->event_nuclide-1].fission) {
if (simulation::micro_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->last_E, p->sqrtkT);
@ -383,9 +381,9 @@ apply_derivative_to_score(const Particle* p, int i_tally, int i_nuclide,
break;
case SCORE_NU_FISSION:
if (simulation::micro_xs[p->event_nuclide-1].fission) {
double nu = simulation::micro_xs[p->event_nuclide-1].nu_fission
/ simulation::micro_xs[p->event_nuclide-1].fission;
if (simulation::micro_xs[p->event_nuclide].fission) {
double nu = simulation::micro_xs[p->event_nuclide].nu_fission
/ simulation::micro_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->last_E, p->sqrtkT);
@ -592,9 +590,8 @@ score_track_derivative(const Particle* p, double distance)
// phi is proportional to e^(-Sigma_tot * dist)
// (1 / phi) * (d_phi / d_N) = - (d_Sigma_tot / d_N) * dist
// (1 / phi) * (d_phi / d_N) = - sigma_tot * dist
//TODO: off-by-one
deriv.flux_deriv -= distance
* simulation::micro_xs[deriv.diff_nuclide-1].total;
* simulation::micro_xs[deriv.diff_nuclide].total;
break;
case DIFF_TEMPERATURE:
@ -635,17 +632,16 @@ void score_collision_derivative(const Particle* p)
break;
case DIFF_NUCLIDE_DENSITY:
//TODO: off-by-one throughout on diff_nuclide
if (p->event_nuclide != deriv.diff_nuclide) continue;
// Find the index in this material for the diff_nuclide.
int i;
for (i = 0; i < material.nuclide_.size(); ++i)
if (material.nuclide_[i] == deriv.diff_nuclide - 1) break;
if (material.nuclide_[i] == deriv.diff_nuclide) break;
// Make sure we found the nuclide.
if (material.nuclide_[i] != deriv.diff_nuclide - 1) {
if (material.nuclide_[i] != deriv.diff_nuclide) {
std::stringstream err_msg;
err_msg << "Could not find nuclide "
<< data::nuclides[deriv.diff_nuclide-1]->name_ << " in material "
<< data::nuclides[deriv.diff_nuclide]->name_ << " in material "
<< material.id_ << " for tally derivative " << deriv.id;
fatal_error(err_msg);
}
@ -661,8 +657,7 @@ void score_collision_derivative(const Particle* p)
for (auto i_nuc : material.nuclide_) {
const auto& nuc {*data::nuclides[i_nuc]};
//TODO: off-by-one
if (i_nuc == p->event_nuclide - 1
&& multipole_in_range(&nuc, p->last_E)) {
if (i_nuc == p->event_nuclide && multipole_in_range(&nuc, p->last_E)) {
// 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

96
src/tallies/tally_scoring.cpp
View file

@ -337,9 +337,6 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
double score;
//TODO: off-by-one throughout on p->event_nuclide
//TODO: off-by-one throughout on p->material
switch (score_bin) {
@ -445,8 +442,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
} else {
// Get yield and apply to score
auto m =
data::nuclides[p->event_nuclide-1]->reaction_index_[p->event_MT];
const auto& rxn {*data::nuclides[p->event_nuclide-1]->reactions_[m]};
data::nuclides[p->event_nuclide]->reaction_index_[p->event_MT];
const auto& rxn {*data::nuclides[p->event_nuclide]->reactions_[m]};
score = p->last_wgt * flux * (*rxn.products_[0].yield_)(E);
}
break;
@ -483,10 +480,10 @@ score_general_ce(const 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
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
if (simulation::micro_xs[p->event_nuclide].absorption > 0) {
score = p->absorb_wgt
* simulation::micro_xs[p->event_nuclide-1].fission
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
} else {
score = 0.;
}
@ -497,8 +494,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
// weight entering the collision as the estimate for the fission
// reaction rate
score = p->last_wgt
* simulation::micro_xs[p->event_nuclide-1].fission
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
} else {
if (i_nuclide >= 0) {
@ -525,10 +522,10 @@ score_general_ce(const 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 (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
if (simulation::micro_xs[p->event_nuclide].absorption > 0) {
score = p->absorb_wgt
* simulation::micro_xs[p->event_nuclide-1].nu_fission
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
* simulation::micro_xs[p->event_nuclide].nu_fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
} else {
score = 0.;
}
@ -568,12 +565,12 @@ score_general_ce(const 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 (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
if (simulation::micro_xs[p->event_nuclide].absorption > 0) {
score = p->absorb_wgt
* simulation::micro_xs[p->event_nuclide-1].fission
* data::nuclides[p->event_nuclide-1]
* simulation::micro_xs[p->event_nuclide].fission
* data::nuclides[p->event_nuclide]
->nu(E, ReactionProduct::EmissionMode::prompt)
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
} else {
score = 0.;
}
@ -630,8 +627,8 @@ score_general_ce(const 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 (simulation::micro_xs[p->event_nuclide-1].absorption > 0
&& data::nuclides[p->event_nuclide-1]->fissionable_) {
if (simulation::micro_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_];
const DelayedGroupFilter& filt
@ -640,11 +637,11 @@ score_general_ce(const 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];
auto yield = data::nuclides[p->event_nuclide-1]
auto yield = data::nuclides[p->event_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed, d);
score = p->absorb_wgt * yield
* simulation::micro_xs[p->event_nuclide-1].fission
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
score_fission_delayed_dg(i_tally, d_bin+1, score,
score_index);
}
@ -654,10 +651,10 @@ score_general_ce(const 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->absorb_wgt
* simulation::micro_xs[p->event_nuclide-1].fission
* data::nuclides[p->event_nuclide-1]
* simulation::micro_xs[p->event_nuclide].fission
* data::nuclides[p->event_nuclide]
->nu(E, ReactionProduct::EmissionMode::delayed)
/ simulation::micro_xs[p->event_nuclide-1].absorption *flux;
/ simulation::micro_xs[p->event_nuclide].absorption *flux;
}
}
} else {
@ -766,8 +763,8 @@ score_general_ce(const 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
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (simulation::micro_xs[p->event_nuclide].absorption > 0
&& nuc.fissionable_) {
const auto& rxn {*nuc.fission_rx_[0]};
if (tally.delayedgroup_filter_ != C_NONE) {
@ -782,8 +779,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
= 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
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption
* rate * flux;
score_fission_delayed_dg(i_tally, d_bin+1, score,
score_index);
@ -805,8 +802,8 @@ score_general_ce(const 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->absorb_wgt
* simulation::micro_xs[p->event_nuclide-1].fission * yield
/ simulation::micro_xs[p->event_nuclide-1].absorption * flux;
* simulation::micro_xs[p->event_nuclide].fission * yield
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
}
}
@ -826,7 +823,7 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
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& nuc {*data::nuclides[p->event_nuclide]};
const auto& rxn {*nuc.fission_rx_[0]};
auto rate = rxn.products_[g].decay_rate_;
score += simulation::keff * bank.wgt * rate * flux;
@ -953,13 +950,13 @@ score_general_ce(const 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
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (simulation::micro_xs[p->event_nuclide].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;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
} else {
// Skip any non-absorption events
@ -967,13 +964,13 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
// 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
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (simulation::micro_xs[p->event_nuclide].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;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
}
} else {
@ -1046,8 +1043,8 @@ score_general_ce(const 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
const auto& nuc {*data::nuclides[p->event_nuclide-1]};
if (simulation::micro_xs[p->event_nuclide-1].absorption > 0) {
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (simulation::micro_xs[p->event_nuclide].absorption > 0) {
double q_value = 0.;
if (score_bin == SCORE_FISS_Q_PROMPT) {
if (nuc.fission_q_prompt_)
@ -1057,8 +1054,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
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;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
} else {
// Skip any non-absorption events
@ -1066,8 +1063,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
// 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) {
const auto& nuc {*data::nuclides[p->event_nuclide]};
if (simulation::micro_xs[p->event_nuclide].absorption > 0) {
double q_value = 0.;
if (score_bin == SCORE_FISS_Q_PROMPT) {
if (nuc.fission_q_prompt_)
@ -1077,8 +1074,8 @@ score_general_ce(const Particle* p, int i_tally, int start_index,
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;
* simulation::micro_xs[p->event_nuclide].fission
/ simulation::micro_xs[p->event_nuclide].absorption * flux;
}
}
} else {
@ -1992,8 +1989,7 @@ void score_analog_tally_ce(const Particle* p)
// the event nuclide or the total material. Note that the i_nuclide
// and flux arguments for score_general are not used for analog
// tallies.
//TODO: off-by-one
if (i_nuclide == p->event_nuclide-1 || i_nuclide == -1)
if (i_nuclide == p->event_nuclide || i_nuclide == -1)
score_general_ce(p, i_tally, i*tally.scores_.size(), filter_index,
-1, -1., filter_weight);
}