From 524b5344ebe4cf91fea89c08d8a9d5652378b4c9 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Sat, 9 Nov 2019 09:08:56 -0600 Subject: [PATCH 1/9] removed MGXS C interface methods (calculate_xs_c, get_name_c, and get_awr_c) that are no longer necessary now that the code is full C++ --- include/openmc/mgxs_interface.h | 14 ------------ src/mgxs_interface.cpp | 39 +-------------------------------- src/particle.cpp | 14 ++++++------ 3 files changed, 8 insertions(+), 59 deletions(-) diff --git a/include/openmc/mgxs_interface.h b/include/openmc/mgxs_interface.h index afd183e87..2cbe5f21d 100644 --- a/include/openmc/mgxs_interface.h +++ b/include/openmc/mgxs_interface.h @@ -82,10 +82,6 @@ void mark_fissionable_mgxs_materials(); // Mgxs tracking/transport/tallying interface methods //============================================================================== -extern "C" void -calculate_xs_c(int i_mat, int gin, double sqrtkT, Direction u, - double& total_xs, double& abs_xs, double& nu_fiss_xs); - double get_nuclide_xs(int index, int xstype, int gin, const int* gout, const double* mu, const int* dg); @@ -102,15 +98,5 @@ inline double get_macro_xs(int index, int xstype, int gin) {return get_macro_xs(index, xstype, gin, nullptr, nullptr, nullptr);} -//============================================================================== -// General Mgxs methods -//============================================================================== - -extern "C" void -get_name_c(int index, int name_len, char* name); - -extern "C" double -get_awr_c(int index); - } // namespace openmc #endif // OPENMC_MGXS_INTERFACE_H diff --git a/src/mgxs_interface.cpp b/src/mgxs_interface.cpp index 1b0f2675a..0bc61fa00 100644 --- a/src/mgxs_interface.cpp +++ b/src/mgxs_interface.cpp @@ -38,7 +38,7 @@ MgxsInterface::MgxsInterface(const std::string& path_cross_sections, void MgxsInterface::set_nuclides_and_temperatures( std::vector xs_to_read, std::vector> xs_temps) -{ +{ // Check to remove all duplicates xs_to_read_ = xs_to_read; xs_temps_to_read_ = xs_temps; @@ -291,16 +291,6 @@ void mark_fissionable_mgxs_materials() // Mgxs tracking/transport/tallying interface methods //============================================================================== -void -calculate_xs_c(int i_mat, int gin, double sqrtkT, Direction u, - double& total_xs, double& abs_xs, double& nu_fiss_xs) -{ - data::mg.macro_xs_[i_mat].calculate_xs(gin - 1, sqrtkT, u, total_xs, abs_xs, - nu_fiss_xs); -} - -//============================================================================== - double get_nuclide_xs(int index, int xstype, int gin, const int* gout, const double* mu, const int* dg) @@ -333,31 +323,4 @@ get_macro_xs(int index, int xstype, int gin, const int* gout, return data::mg.macro_xs_[index].get_xs(xstype, gin - 1, gout_c_p, mu, dg); } -//============================================================================== -// General Mgxs methods -//============================================================================== - -void -get_name_c(int index, int name_len, char* name) -{ - // First blank out our input string - std::string str(name_len - 1, ' '); - std::strcpy(name, str.c_str()); - - // Now get the data and copy to the C-string - str = data::mg.nuclides_[index - 1].name; - std::strcpy(name, str.c_str()); - - // Finally, remove the null terminator - name[std::strlen(name)] = ' '; -} - -//============================================================================== - -double -get_awr_c(int index) -{ - return data::mg.nuclides_[index - 1].awr; -} - } // namespace openmc diff --git a/src/particle.cpp b/src/particle.cpp index 030c8dc63..3ce3cf84c 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -203,11 +203,11 @@ Particle::transport() } } else { // Get the MG data - calculate_xs_c(material_, g_, sqrtkT_, this->u_local(), - macro_xs_.total, macro_xs_.absorption, macro_xs_.nu_fission); + data::mg.macro_xs_[material_].calculate_xs(g_ - 1, sqrtkT_, + this->u_local(), macro_xs_.total, macro_xs_.absorption, + macro_xs_.nu_fission); - // Finally, update the particle group while we have already checked - // for if multi-group + // Finally, update the particle group since we know we are multi-group g_last_ = g_; } } else { @@ -428,7 +428,7 @@ Particle::cross_surface() } return; - } else if ((surf->bc_ == BC_REFLECT || surf->bc_ == BC_WHITE) + } else if ((surf->bc_ == BC_REFLECT || surf->bc_ == BC_WHITE) && (settings::run_mode != RUN_MODE_PLOTTING)) { // ======================================================================= // PARTICLE REFLECTS FROM SURFACE @@ -458,11 +458,11 @@ Particle::cross_surface() score_surface_tally(this, model::active_meshsurf_tallies); this->r() = r; } - + Direction u = (surf->bc_ == BC_REFLECT) ? surf->reflect(this->r(), this->u()) : surf->diffuse_reflect(this->r(), this->u()); - + // Make sure new particle direction is normalized this->u() = u / u.norm(); From 2903970129d4c643072bd4a3d7d3e6cf2fbda5d3 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Sat, 9 Nov 2019 09:16:37 -0600 Subject: [PATCH 2/9] Added test for making an MGXS lib with nuclides, and then re-running in MG mode. Prior to this we only tested macroscopic MGXS --- .../__init__.py | 0 .../inputs_true.dat | 251 ++++++++++++++++++ .../results_true.dat | 2 + .../mgxs_library_ce_to_mg_nuclides/test.py | 85 ++++++ 4 files changed, 338 insertions(+) create mode 100644 tests/regression_tests/mgxs_library_ce_to_mg_nuclides/__init__.py create mode 100644 tests/regression_tests/mgxs_library_ce_to_mg_nuclides/inputs_true.dat create mode 100644 tests/regression_tests/mgxs_library_ce_to_mg_nuclides/results_true.dat create mode 100644 tests/regression_tests/mgxs_library_ce_to_mg_nuclides/test.py diff --git a/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/__init__.py b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/inputs_true.dat b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/inputs_true.dat new file mode 100644 index 000000000..576f27966 --- /dev/null +++ b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/inputs_true.dat @@ -0,0 +1,251 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 10 + 5 + + + -0.63 -0.63 -1 0.63 0.63 1 + + + + + + + 1 + + + 0.0 0.625 20000000.0 + + + 0.0 0.625 20000000.0 + + + 3 + + + 2 + + + 3 + + + 1 2 + total + flux + tracklength + + + 1 2 + U234 U235 U238 O16 + total + tracklength + + + 1 2 + total + flux + tracklength + + + 1 2 + U234 U235 U238 O16 + absorption + tracklength + + + 1 2 + total + flux + analog + + + 1 2 7 + U234 U235 U238 O16 + nu-fission + analog + + + 1 2 + total + flux + analog + + + 1 2 7 11 + U234 U235 U238 O16 + nu-scatter + analog + + + 1 2 7 + U234 U235 U238 O16 + nu-scatter + analog + + + 1 2 7 + U234 U235 U238 O16 + scatter + analog + + + 15 2 + total + flux + tracklength + + + 15 2 + Zr90 Zr91 Zr92 Zr94 Zr96 + total + tracklength + + + 15 2 + total + flux + tracklength + + + 15 2 + Zr90 Zr91 Zr92 Zr94 Zr96 + absorption + tracklength + + + 15 2 + total + flux + analog + + + 15 2 7 + Zr90 Zr91 Zr92 Zr94 Zr96 + nu-fission + analog + + + 15 2 + total + flux + analog + + + 15 2 7 11 + Zr90 Zr91 Zr92 Zr94 Zr96 + nu-scatter + analog + + + 15 2 7 + Zr90 Zr91 Zr92 Zr94 Zr96 + nu-scatter + analog + + + 15 2 7 + Zr90 Zr91 Zr92 Zr94 Zr96 + scatter + analog + + + 29 2 + total + flux + tracklength + + + 29 2 + H1 O16 B10 B11 + total + tracklength + + + 29 2 + total + flux + tracklength + + + 29 2 + H1 O16 B10 B11 + absorption + tracklength + + + 29 2 + total + flux + analog + + + 29 2 7 + H1 O16 B10 B11 + nu-fission + analog + + + 29 2 + total + flux + analog + + + 29 2 7 11 + H1 O16 B10 B11 + nu-scatter + analog + + + 29 2 7 + H1 O16 B10 B11 + nu-scatter + analog + + + 29 2 7 + H1 O16 B10 B11 + scatter + analog + + diff --git a/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/results_true.dat b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/results_true.dat new file mode 100644 index 000000000..b494a4b4f --- /dev/null +++ b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +5.350442E-01 1.484490E-02 diff --git a/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/test.py b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/test.py new file mode 100644 index 000000000..b624140d4 --- /dev/null +++ b/tests/regression_tests/mgxs_library_ce_to_mg_nuclides/test.py @@ -0,0 +1,85 @@ +import os + +import openmc +import openmc.mgxs +from openmc.examples import pwr_pin_cell + +from tests.testing_harness import PyAPITestHarness +from tests.regression_tests import config + + +class MGXSTestHarness(PyAPITestHarness): + def __init__(self, *args, **kwargs): + # Generate inputs using parent class routine + super().__init__(*args, **kwargs) + + # Initialize a two-group structure + energy_groups = openmc.mgxs.EnergyGroups(group_edges=[0, 0.625, 20.e6]) + + # Initialize MGXS Library for a few cross section types + self.mgxs_lib = openmc.mgxs.Library(self._model.geometry) + self.mgxs_lib.by_nuclide = True + self.mgxs_lib.mgxs_types = ['total', 'absorption', 'nu-fission matrix', + 'nu-scatter matrix', 'multiplicity matrix'] + self.mgxs_lib.energy_groups = energy_groups + self.mgxs_lib.correction = None + self.mgxs_lib.legendre_order = 3 + self.mgxs_lib.domain_type = 'material' + self.mgxs_lib.build_library() + + # Initialize a tallies file + self.mgxs_lib.add_to_tallies_file(self._model.tallies, merge=False) + + def _run_openmc(self): + # Initial run + if config['mpi']: + mpi_args = [config['mpiexec'], '-n', config['mpi_np']] + openmc.run(openmc_exec=config['exe'], mpi_args=mpi_args) + else: + openmc.run(openmc_exec=config['exe']) + + # Build MG Inputs + # Get data needed to execute Library calculations. + sp = openmc.StatePoint(self._sp_name) + self.mgxs_lib.load_from_statepoint(sp) + self._model.mgxs_file, self._model.materials, \ + self._model.geometry = self.mgxs_lib.create_mg_mode() + + # Modify materials and settings so we can run in MG mode + self._model.materials.cross_sections = './mgxs.h5' + self._model.settings.energy_mode = 'multi-group' + + # Write modified input files + self._model.settings.export_to_xml() + self._model.geometry.export_to_xml() + self._model.materials.export_to_xml() + self._model.mgxs_file.export_to_hdf5() + # Dont need tallies.xml, so remove the file + if os.path.exists('tallies.xml'): + os.remove('tallies.xml') + + # Enforce closing statepoint and summary files so HDF5 + # does not throw an error during the next OpenMC execution + sp._f.close() + sp._summary._f.close() + + # Re-run MG mode. + if config['mpi']: + mpi_args = [config['mpiexec'], '-n', config['mpi_np']] + openmc.run(openmc_exec=config['exe'], mpi_args=mpi_args) + else: + openmc.run(openmc_exec=config['exe']) + + def _cleanup(self): + super()._cleanup() + f = 'mgxs.h5' + if os.path.exists(f): + os.remove(f) + + +def test_mgxs_library_ce_to_mg(): + # Set the input set to use the pincell model + model = pwr_pin_cell() + + harness = MGXSTestHarness('statepoint.10.h5', model) + harness.main() From 3ccfa08a115c70833a0ec2a7f4d3995ba0d60372 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Sat, 9 Nov 2019 09:24:30 -0600 Subject: [PATCH 3/9] minor change to add an inline version of Mgxs::get_xs just like @gridley did for mg.get_*_xs --- include/openmc/mgxs.h | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/include/openmc/mgxs.h b/include/openmc/mgxs.h index 313677f6c..71e8c7ec9 100644 --- a/include/openmc/mgxs.h +++ b/include/openmc/mgxs.h @@ -143,6 +143,11 @@ class Mgxs { get_xs(int xstype, int gin, const int* gout, const double* mu, const int* dg); + inline double + get_xs(int xstype, int gin) + {return get_xs(xstype, gin, nullptr, nullptr, nullptr);} + + //! \brief Samples the fission neutron energy and if prompt or delayed. //! //! @param gin Incoming energy group. From 05fac29164411233d423cf3170c1856c3f622646 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Mon, 11 Nov 2019 19:50:52 -0600 Subject: [PATCH 4/9] Revised MG code so that the energy group indices are 0-based instead of 1-based. Next commit will perform some code cleanup enabled by this change. --- src/mgxs_interface.cpp | 8 ++++---- src/particle.cpp | 4 ++-- src/particle_restart.cpp | 2 +- src/physics_mg.cpp | 18 +++++------------- src/source.cpp | 4 ++-- src/tallies/filter_energy.cpp | 6 +++--- src/tallies/tally_scoring.cpp | 4 ++-- 7 files changed, 19 insertions(+), 27 deletions(-) diff --git a/src/mgxs_interface.cpp b/src/mgxs_interface.cpp index 0bc61fa00..3076aedc8 100644 --- a/src/mgxs_interface.cpp +++ b/src/mgxs_interface.cpp @@ -298,12 +298,12 @@ get_nuclide_xs(int index, int xstype, int gin, const int* gout, int gout_c; const int* gout_c_p; if (gout != nullptr) { - gout_c = *gout - 1; + gout_c = *gout; gout_c_p = &gout_c; } else { gout_c_p = gout; } - return data::mg.nuclides_[index].get_xs(xstype, gin - 1, gout_c_p, mu, dg); + return data::mg.nuclides_[index].get_xs(xstype, gin, gout_c_p, mu, dg); } //============================================================================== @@ -315,12 +315,12 @@ get_macro_xs(int index, int xstype, int gin, const int* gout, int gout_c; const int* gout_c_p; if (gout != nullptr) { - gout_c = *gout - 1; + gout_c = *gout; gout_c_p = &gout_c; } else { gout_c_p = gout; } - return data::mg.macro_xs_[index].get_xs(xstype, gin - 1, gout_c_p, mu, dg); + return data::mg.macro_xs_[index].get_xs(xstype, gin, gout_c_p, mu, dg); } } // namespace openmc diff --git a/src/particle.cpp b/src/particle.cpp index 3ce3cf84c..da93c6723 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -123,7 +123,7 @@ Particle::from_source(const Bank* src) } else { g_ = static_cast(src->E); g_last_ = static_cast(src->E); - E_ = data::mg.energy_bin_avg_[g_ - 1]; + E_ = data::mg.energy_bin_avg_[g_]; } E_last_ = E_; } @@ -203,7 +203,7 @@ Particle::transport() } } else { // Get the MG data - data::mg.macro_xs_[material_].calculate_xs(g_ - 1, sqrtkT_, + data::mg.macro_xs_[material_].calculate_xs(g_, sqrtkT_, this->u_local(), macro_xs_.total, macro_xs_.absorption, macro_xs_.nu_fission); diff --git a/src/particle_restart.cpp b/src/particle_restart.cpp index 73e6d5d37..ce6555608 100644 --- a/src/particle_restart.cpp +++ b/src/particle_restart.cpp @@ -52,7 +52,7 @@ void read_particle_restart(Particle& p, int& previous_run_mode) // Set energy group and average energy in multi-group mode if (!settings::run_CE) { p.g_ = p.E_; - p.E_ = data::mg.energy_bin_avg_[p.g_ - 1]; + p.E_ = data::mg.energy_bin_avg_[p.g_]; } // Set particle last attributes diff --git a/src/physics_mg.cpp b/src/physics_mg.cpp index 90e76c67a..924356431 100644 --- a/src/physics_mg.cpp +++ b/src/physics_mg.cpp @@ -77,22 +77,14 @@ sample_reaction(Particle* p) void scatter(Particle* p) { - // Adjust indices for Fortran to C++ indexing - // TODO: Remove when no longer needed - int gin = p->g_last_ - 1; - int gout = p->g_ - 1; - int i_mat = p->material_; - data::mg.macro_xs_[i_mat].sample_scatter(gin, gout, p->mu_, p->wgt_); - - // Adjust return value for fortran indexing - // TODO: Remove when no longer needed - p->g_ = gout + 1; + data::mg.macro_xs_[p->material_].sample_scatter(p->g_last_, p->g_, p->mu_, + p->wgt_); // Rotate the angle p->u() = rotate_angle(p->u(), p->mu_, nullptr); // Update energy value for downstream compatability (in tallying) - p->E_ = data::mg.energy_bin_avg_[gout]; + p->E_ = data::mg.energy_bin_avg_[p->g_]; // Set event component p->event_ = EVENT_SCATTER; @@ -148,8 +140,8 @@ create_fission_sites(Particle* p, std::vector& bank) // the energy in the fission bank int dg; int gout; - data::mg.macro_xs_[p->material_].sample_fission_energy(p->g_ - 1, dg, gout); - site.E = gout + 1; + data::mg.macro_xs_[p->material_].sample_fission_energy(p->g_, dg, gout); + site.E = gout; site.delayed_group = dg + 1; // Set the delayed group on the particle as well diff --git a/src/source.cpp b/src/source.cpp index 6b5cb8aba..1d082240e 100644 --- a/src/source.cpp +++ b/src/source.cpp @@ -309,11 +309,11 @@ Particle::Bank sample_external_source() // Sample source site from i-th source distribution Particle::Bank site {model::external_sources[i].sample()}; - // If running in MG, convert site % E to group + // If running in MG, convert site.E to group if (!settings::run_CE) { site.E = lower_bound_index(data::mg.rev_energy_bins_.begin(), data::mg.rev_energy_bins_.end(), site.E); - site.E = data::mg.num_energy_groups_ - site.E; + site.E = data::mg.num_energy_groups_ - site.E - 1.; } // Set the random number generator back to the tracking stream. diff --git a/src/tallies/filter_energy.cpp b/src/tallies/filter_energy.cpp index 7737e4fac..7c8b5fbee 100644 --- a/src/tallies/filter_energy.cpp +++ b/src/tallies/filter_energy.cpp @@ -61,9 +61,9 @@ const { if (p->g_ != F90_NONE && matches_transport_groups_) { if (estimator == ESTIMATOR_TRACKLENGTH) { - match.bins_.push_back(data::mg.num_energy_groups_ - p->g_); + match.bins_.push_back(data::mg.num_energy_groups_ - p->g_ - 1); } else { - match.bins_.push_back(data::mg.num_energy_groups_ - p->g_last_); + match.bins_.push_back(data::mg.num_energy_groups_ - p->g_last_ - 1); } match.weights_.push_back(1.0); @@ -104,7 +104,7 @@ EnergyoutFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match) const { if (p->g_ != F90_NONE && matches_transport_groups_) { - match.bins_.push_back(data::mg.num_energy_groups_ - p->g_); + match.bins_.push_back(data::mg.num_energy_groups_ - p->g_ - 1); match.weights_.push_back(1.0); } else { diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index a2532c437..95e674735 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -348,9 +348,9 @@ score_fission_eout(const Particle* p, int i_tally, int i_score, int score_bin) // determine outgoing energy group from fission bank auto g_out = static_cast(bank.E); - // modify the value so that g_out = 1 corresponds to the highest energy + // modify the value so that g_out = 0 corresponds to the highest energy // bin - g_out = eo_filt.n_bins() - g_out; + g_out = eo_filt.n_bins() - g_out - 1; // change outgoing energy bin simulation::filter_matches[i_eout_filt].bins_[i_bin] = g_out; From 64dbaff587156803b957debb3dabe9f63992a4bc Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Mon, 11 Nov 2019 20:58:22 -0600 Subject: [PATCH 5/9] Finished cleaning up group indices off-by-ones. Still have to to delayed groups --- include/openmc/mgxs_interface.h | 20 --- src/mgxs_interface.cpp | 36 ---- src/tallies/tally_scoring.cpp | 297 ++++++++++++++------------------ 3 files changed, 131 insertions(+), 222 deletions(-) diff --git a/include/openmc/mgxs_interface.h b/include/openmc/mgxs_interface.h index 2cbe5f21d..6beda5778 100644 --- a/include/openmc/mgxs_interface.h +++ b/include/openmc/mgxs_interface.h @@ -78,25 +78,5 @@ void set_mg_interface_nuclides_and_temps(); // After macro XS have been read, materials can be marked as fissionable void mark_fissionable_mgxs_materials(); -//============================================================================== -// Mgxs tracking/transport/tallying interface methods -//============================================================================== - -double -get_nuclide_xs(int index, int xstype, int gin, const int* gout, - const double* mu, const int* dg); - -inline double -get_nuclide_xs(int index, int xstype, int gin) -{return get_nuclide_xs(index, xstype, gin, nullptr, nullptr, nullptr);} - -double -get_macro_xs(int index, int xstype, int gin, const int* gout, - const double* mu, const int* dg); - -inline double -get_macro_xs(int index, int xstype, int gin) -{return get_macro_xs(index, xstype, gin, nullptr, nullptr, nullptr);} - } // namespace openmc #endif // OPENMC_MGXS_INTERFACE_H diff --git a/src/mgxs_interface.cpp b/src/mgxs_interface.cpp index 3076aedc8..936eeb31f 100644 --- a/src/mgxs_interface.cpp +++ b/src/mgxs_interface.cpp @@ -287,40 +287,4 @@ void mark_fissionable_mgxs_materials() } } -//============================================================================== -// Mgxs tracking/transport/tallying interface methods -//============================================================================== - -double -get_nuclide_xs(int index, int xstype, int gin, const int* gout, - const double* mu, const int* dg) -{ - int gout_c; - const int* gout_c_p; - if (gout != nullptr) { - gout_c = *gout; - gout_c_p = &gout_c; - } else { - gout_c_p = gout; - } - return data::mg.nuclides_[index].get_xs(xstype, gin, gout_c_p, mu, dg); -} - -//============================================================================== - -double -get_macro_xs(int index, int xstype, int gin, const int* gout, - const double* mu, const int* dg) -{ - int gout_c; - const int* gout_c_p; - if (gout != nullptr) { - gout_c = *gout; - gout_c_p = &gout_c; - } else { - gout_c_p = gout; - } - return data::mg.macro_xs_[index].get_xs(xstype, gin, gout_c_p, mu, dg); -} - } // namespace openmc diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index 95e674735..bbbd5411a 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -1374,15 +1374,15 @@ score_general_mg(const Particle* p, int i_tally, int start_index, p_g = p->g_; } - // To significantly reduce de-referencing, point matxs to the macroscopic - // Mgxs for the material of interest - data::mg.macro_xs_[p->material_].set_angle_index(p_u); + // For shorthand, assign pointers to the material and nuclide xs set + Mgxs* nuc_xs = &data::mg.nuclides_[i_nuclide]; + Mgxs* macro_xs = &data::mg.macro_xs_[p->material_]; - // Do same for nucxs, point it to the microscopic nuclide data of interest + // Find the temperature and angle indices of interest + macro_xs->set_angle_index(p_u); if (i_nuclide >= 0) { - // And since we haven't calculated this temperature index yet, do so now - data::mg.nuclides_[i_nuclide].set_temperature_index(p->sqrtkT_); - data::mg.nuclides_[i_nuclide].set_angle_index(p_u); + nuc_xs->set_temperature_index(p->sqrtkT_); + nuc_xs->set_angle_index(p_u); } for (auto i = 0; i < tally.scores_.size(); ++i) { @@ -1426,14 +1426,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } //TODO: should flux be multiplied in above instead of below? if (i_nuclide >= 0) { - score *= flux * atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_TOTAL, p_g) - / get_macro_xs(p->material_, MG_GET_XS_TOTAL, p_g); + score *= flux * atom_density * nuc_xs->get_xs(MG_GET_XS_TOTAL, p_g) + / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); } } else { if (i_nuclide >= 0) { - score = get_nuclide_xs(i_nuclide, MG_GET_XS_TOTAL, p_g) - * atom_density * flux; + score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_TOTAL, p_g); } else { score = p->macro_xs_.total * flux; } @@ -1455,21 +1453,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_; } if (i_nuclide >= 0) { - score *= flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_INVERSE_VELOCITY, p_g) - / get_macro_xs(p->material_, MG_GET_XS_TOTAL, p_g); + score *= flux * nuc_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) + / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); } else { - score *= flux - * get_macro_xs(p->material_, MG_GET_XS_INVERSE_VELOCITY, p_g) - / get_macro_xs(p->material_, MG_GET_XS_TOTAL, p_g); + score *= flux * macro_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) + / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); } } else { if (i_nuclide >= 0) { - score = flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_INVERSE_VELOCITY, p_g); + score = flux * nuc_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); } else { - score = flux - * get_macro_xs(p->material_, MG_GET_XS_INVERSE_VELOCITY, p_g); + score = flux * macro_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); } } break; @@ -1483,19 +1477,18 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // weight entering the collision as the estimator for the reaction rate score = p->wgt_last_ * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_SCATTER_FMU_MULT, - p->g_last_, &p->g_, &p->mu_, nullptr) - / get_macro_xs(p->material_, MG_GET_XS_SCATTER_FMU_MULT, - p->g_last_, &p->g_, &p->mu_, nullptr); + score *= atom_density * nuc_xs->get_xs( + MG_GET_XS_SCATTER_FMU_MULT, p->g_last_, &p->g_, &p->mu_, nullptr) + / macro_xs->get_xs( + MG_GET_XS_SCATTER_FMU_MULT, p->g_last_, &p->g_, &p->mu_, nullptr); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * get_nuclide_xs( - i_nuclide, MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); + score = atom_density * flux * nuc_xs->get_xs( + MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); } else { - score = flux * get_macro_xs( - p->material_, MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); + score = flux * macro_xs->get_xs( + MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); } } break; @@ -1514,17 +1507,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // adjust the score by the actual probability for that nuclide. if (i_nuclide >= 0) { score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_SCATTER_FMU, - p->g_last_, &p->g_, &p->mu_, nullptr) - / get_macro_xs(p->material_, MG_GET_XS_SCATTER_FMU, - p->g_last_, &p->g_, &p->mu_, nullptr); + * nuc_xs->get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, + &p->mu_, nullptr) + / macro_xs->get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, + &p->mu_, nullptr); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * get_nuclide_xs( - i_nuclide, MG_GET_XS_SCATTER, p_g); + score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_SCATTER, p_g); } else { - score = flux * get_macro_xs(p->material_, MG_GET_XS_SCATTER, p_g); + score = flux * macro_xs->get_xs(MG_GET_XS_SCATTER, p_g); } } break; @@ -1544,14 +1536,13 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_ * flux; } if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_ABSORPTION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_ABSORPTION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { score = atom_density * flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_ABSORPTION, p_g); + * nuc_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } else { score = p->macro_xs_.absorption * flux; } @@ -1575,20 +1566,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_ * flux; } if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { - score = get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - * atom_density * flux; + score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g); } else { - score = get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g) * flux; + score = flux * macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } } break; @@ -1610,13 +1598,11 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // nu-fission score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_NU_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_NU_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { // Skip any non-fission events @@ -1628,17 +1614,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // score. score = simulation::keff * p->wgt_bank_ * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } } } else { if (i_nuclide >= 0) { - score = get_nuclide_xs(i_nuclide, MG_GET_XS_NU_FISSION, p_g) - * atom_density * flux; + score = atom_density * flux + * nuc_xs->get_xs(MG_GET_XS_NU_FISSION, p_g); } else { - score = get_macro_xs(p->material_, MG_GET_XS_NU_FISSION, p_g) * flux; + score = flux * macro_xs->get_xs(MG_GET_XS_NU_FISSION, p_g); } } break; @@ -1660,13 +1645,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // prompt-nu-fission score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_PROMPT_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * + nuc_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_PROMPT_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { // Skip any non-fission events @@ -1681,18 +1665,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto prompt_frac = 1. - n_delayed / static_cast(p->n_bank_); score = simulation::keff * p->wgt_bank_ * prompt_frac * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } } } else { if (i_nuclide >= 0) { - score = get_nuclide_xs(i_nuclide, MG_GET_XS_PROMPT_NU_FISSION, p_g) - * atom_density * flux; + score = atom_density * flux * + nuc_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); } else { - score = get_macro_xs(p->material_, MG_GET_XS_PROMPT_NU_FISSION, p_g) - * flux; + score = flux * macro_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); } } break; @@ -1712,7 +1694,8 @@ score_general_mg(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 (get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g) > 0) { + double abs_xs = macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + if (abs_xs > 0.) { if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; const DelayedGroupFilter& filt @@ -1723,18 +1706,15 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= - get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) + / abs_xs; } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) + / abs_xs; } - score_fission_delayed_dg(i_tally, d_bin, score, - score_index); + score_fission_delayed_dg(i_tally, d_bin, score, score_index); } continue; } else { @@ -1743,13 +1723,11 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // delayed-nu-fission xs to the absorption xs score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= - get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) + / abs_xs; } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) + / abs_xs; } } } @@ -1774,9 +1752,8 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = simulation::keff * p->wgt_bank_ / p->n_bank_ * p->n_delayed_bank_[d-1] * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } @@ -1788,9 +1765,8 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = simulation::keff * p->wgt_bank_ / p->n_bank_ * n_delayed * flux; if (i_nuclide >= 0) { - score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } } } @@ -1804,24 +1780,21 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d_bin = 0; d_bin < filt.n_bins(); ++d_bin) { auto d = filt.groups()[d_bin] - 1; if (i_nuclide >= 0) { - score = flux * atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + score = flux * atom_density * nuc_xs->get_xs( + MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } else { - score = flux - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + score = flux * macro_xs->get_xs( + MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } continue; } else { if (i_nuclide >= 0) { - score = flux * atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, p_g); + score = flux * atom_density * + nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); } else { - score = flux - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, p_g); + score = flux * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); } } } @@ -1834,7 +1807,8 @@ score_general_mg(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 (get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g) > 0) { + double abs_xs = macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + if (abs_xs > 0) { if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; const DelayedGroupFilter& filt @@ -1845,22 +1819,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= - get_nuclide_xs(i_nuclide, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, + nullptr, nullptr, &d) + * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } else { - score *= - get_macro_xs(p->material_, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } - score_fission_delayed_dg(i_tally, d_bin, score, - score_index); + score_fission_delayed_dg(i_tally, d_bin, score, score_index); } continue; } else { @@ -1872,18 +1841,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d = 0; d < data::mg.num_delayed_groups_; ++d) { if (i_nuclide >= 0) { score += p->wgt_absorb_ * flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } else { score += p->wgt_absorb_ * flux - * get_macro_xs(p->material_, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } } } @@ -1906,14 +1873,14 @@ score_general_mg(const Particle* p, int i_tally, int start_index, if (g != -1) { if (i_nuclide >= 0) { score += simulation::keff * atom_density * bank.wgt * flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_DECAY_RATE, p_g, - nullptr, nullptr, &g) - * get_nuclide_xs(i_nuclide, MG_GET_XS_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_FISSION, p_g); + * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, + &g) + * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } else { score += simulation::keff * bank.wgt * flux - * get_macro_xs(p->material_, MG_GET_XS_DECAY_RATE, p_g, - nullptr, nullptr, &g); + * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &g); } if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; @@ -1925,7 +1892,7 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin]; if (d == g + 1) score_fission_delayed_dg(i_tally, d_bin, score, - score_index); + score_index); } score = 0.; } @@ -1944,16 +1911,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; if (i_nuclide >= 0) { score += atom_density * flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } else { score += flux - * get_macro_xs(p->material_, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } @@ -1963,16 +1930,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d = 0; d < data::mg.num_delayed_groups_; ++d) { if (i_nuclide >= 0) { score += atom_density * flux - * get_nuclide_xs(i_nuclide, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_nuclide_xs(i_nuclide, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, + &d) + * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } else { score += flux - * get_macro_xs(p->material_, MG_GET_XS_DECAY_RATE, - p_g, nullptr, nullptr, &d) - * get_macro_xs(p->material_, MG_GET_XS_DELAYED_NU_FISSION, - p_g, nullptr, nullptr, &d); + * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } } } @@ -1997,25 +1964,23 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } if (i_nuclide >= 0) { score *= atom_density - * get_nuclide_xs(i_nuclide, MG_GET_XS_KAPPA_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + * nuc_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } else { score *= - get_macro_xs(p->material_, MG_GET_XS_KAPPA_FISSION, p_g) - / get_macro_xs(p->material_, MG_GET_XS_ABSORPTION, p_g); + macro_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g) + / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { - score = get_nuclide_xs(i_nuclide, MG_GET_XS_KAPPA_FISSION, p_g) - * atom_density * flux; + score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_KAPPA_FISSION, + p_g); } else { - score = get_macro_xs(p->material_, MG_GET_XS_KAPPA_FISSION, p_g) - * flux; + score = flux * macro_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g); } } break; - case SCORE_EVENTS: // Simply count the number of scoring events score = 1.; From 12a474d5996d5b5a1873306251a53d5a6e62deae Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Mon, 11 Nov 2019 21:21:01 -0600 Subject: [PATCH 6/9] removed while loops on multigroup outgoing energy sampling, replaced with fors and thus avoding chance for infinite loops --- src/mgxs.cpp | 14 ++++++-------- src/scattdata.cpp | 24 ++++++++++-------------- 2 files changed, 16 insertions(+), 22 deletions(-) diff --git a/src/mgxs.cpp b/src/mgxs.cpp index bc2aa9c0a..53e24c67c 100644 --- a/src/mgxs.cpp +++ b/src/mgxs.cpp @@ -555,11 +555,10 @@ Mgxs::sample_fission_energy(int gin, int& dg, int& gout) dg = -1; // sample the outgoing energy group - gout = 0; - double prob_gout = xs_t->chi_prompt(cache[tid].a, gin, gout); - while (prob_gout < xi_gout) { - gout++; + double prob_gout = 0.; + for (gout = 0; gout < num_groups; ++gout) { prob_gout += xs_t->chi_prompt(cache[tid].a, gin, gout); + if (xi_gout < prob_gout) break; } } else { @@ -575,11 +574,10 @@ Mgxs::sample_fission_energy(int gin, int& dg, int& gout) dg = std::min(dg, num_delayed_groups - 1); // sample the outgoing energy group - gout = 0; - double prob_gout = xs_t->chi_delayed(cache[tid].a, dg, gin, gout); - while (prob_gout < xi_gout) { - gout++; + double prob_gout = 0.; + for (gout = 0; gout < num_groups; ++gout) { prob_gout += xs_t->chi_delayed(cache[tid].a, dg, gin, gout); + if (xi_gout < prob_gout) break; } } } diff --git a/src/scattdata.cpp b/src/scattdata.cpp index a8f63d0f6..5fc19c41e 100644 --- a/src/scattdata.cpp +++ b/src/scattdata.cpp @@ -171,14 +171,12 @@ ScattData::sample_energy(int gin, int& gout, int& i_gout) { // Sample the outgoing group double xi = prn(); - + double prob = 0.; i_gout = 0; - gout = gmin[gin]; - double prob = energy[gin][i_gout]; - while((prob < xi) && (gout < gmax[gin])) { - gout++; - i_gout++; + for (gout = gmin[gin]; gout < gmax[gin]; ++gout) { prob += energy[gin][i_gout]; + if (xi < prob) break; + ++i_gout; } } @@ -358,20 +356,18 @@ ScattDataLegendre::sample(int gin, int& gout, double& mu, double& wgt) // Now we can sample mu using the scattering kernel using rejection // sampling from a rectangular bounding box double M = max_val[gin][i_gout]; - int samples = 0; - - while(true) { + int samples; + for (samples = 0; samples < MAX_SAMPLE; ++samples) { mu = 2. * prn() - 1.; double f = calc_f(gin, gout, mu); if (f > 0.) { double u = prn() * M; if (u <= f) break; } - samples++; - if (samples > MAX_SAMPLE) { - fatal_error("Maximum number of Legendre expansion samples reached"); - } - }; + } + if (samples == MAX_SAMPLE) { + fatal_error("Maximum number of Legendre expansion samples reached!"); + } // Update the weight to reflect neutron multiplicity wgt *= mult[gin][i_gout]; From 7d285b418690e97673033aa9446aa8e975a550aa Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Tue, 12 Nov 2019 04:50:21 -0600 Subject: [PATCH 7/9] Minor comments for delayed group filters --- src/physics_mg.cpp | 6 ++++-- src/tallies/tally_scoring.cpp | 12 ++++++------ 2 files changed, 10 insertions(+), 8 deletions(-) diff --git a/src/physics_mg.cpp b/src/physics_mg.cpp index 924356431..8808adbb2 100644 --- a/src/physics_mg.cpp +++ b/src/physics_mg.cpp @@ -136,12 +136,14 @@ create_fission_sites(Particle* p, std::vector& bank) site.u.y = std::sqrt(1. - mu * mu) * std::cos(phi); site.u.z = std::sqrt(1. - mu * mu) * std::sin(phi); - // Sample secondary energy distribution for the fission reaction and set - // the energy in the fission bank + // Sample secondary energy distribution for the fission reaction int dg; int gout; data::mg.macro_xs_[p->material_].sample_fission_energy(p->g_, dg, gout); + // Store the energy and delayed groups on the fission bank site.E = gout; + // We add 1 to the delayed_group bc in MG, -1 is prompt, but in the rest + // of the code, 0 is prompt. site.delayed_group = dg + 1; // Set the delayed group on the particle as well diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index bbbd5411a..1d23ebb68 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -1869,18 +1869,18 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto i = 0; i < p->n_bank_; ++i) { auto i_bank = simulation::fission_bank.size() - p->n_bank_ + i; const auto& bank = simulation::fission_bank[i_bank]; - auto g = bank.delayed_group - 1; - if (g != -1) { + auto d = bank.delayed_group - 1; + if (d != -1) { if (i_nuclide >= 0) { score += simulation::keff * atom_density * bank.wgt * flux * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, - &g) + &d) * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); } else { score += simulation::keff * bank.wgt * flux * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &g); + nullptr, &d); } if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; @@ -1889,8 +1889,8 @@ score_general_mg(const Particle* p, int i_tally, int start_index, 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 + 1) + auto dg = filt.groups()[d_bin]; + if (dg == d + 1) score_fission_delayed_dg(i_tally, d_bin, score, score_index); } From fe8cbcc487849d4e3cbd6e445cdf534f04faaf29 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Tue, 12 Nov 2019 05:07:43 -0600 Subject: [PATCH 8/9] Updated Mgxs::calculate_xs interface to use the particle, matching interface for continuous energy equivalent --- include/openmc/mgxs.h | 11 +++-------- src/mgxs.cpp | 15 +++++++-------- src/particle.cpp | 10 +++++----- 3 files changed, 15 insertions(+), 21 deletions(-) diff --git a/include/openmc/mgxs.h b/include/openmc/mgxs.h index 71e8c7ec9..3dfd6ba4c 100644 --- a/include/openmc/mgxs.h +++ b/include/openmc/mgxs.h @@ -11,6 +11,7 @@ #include "openmc/constants.h" #include "openmc/hdf5_interface.h" +#include "openmc/particle.h" #include "openmc/xsdata.h" @@ -167,15 +168,9 @@ class Mgxs { //! \brief Calculates cross section quantities needed for tracking. //! - //! @param gin Incoming energy group. - //! @param sqrtkT Temperature of the material. - //! @param u Incoming particle direction. - //! @param total_xs Resultant total cross section. - //! @param abs_xs Resultant absorption cross section. - //! @param nu_fiss_xs Resultant nu-fission cross section. + //! @param p The particle whose attributes set which MGXS to get. void - calculate_xs(int gin, double sqrtkT, Direction u, - double& total_xs, double& abs_xs, double& nu_fiss_xs); + calculate_xs(Particle& p); //! \brief Sets the temperature index in cache given a temperature //! diff --git a/src/mgxs.cpp b/src/mgxs.cpp index 53e24c67c..0dd8c73cc 100644 --- a/src/mgxs.cpp +++ b/src/mgxs.cpp @@ -600,8 +600,7 @@ Mgxs::sample_scatter(int gin, int& gout, double& mu, double& wgt) //============================================================================== void -Mgxs::calculate_xs(int gin, double sqrtkT, Direction u, - double& total_xs, double& abs_xs, double& nu_fiss_xs) +Mgxs::calculate_xs(Particle& p) { // Set our indices #ifdef _OPENMP @@ -609,13 +608,13 @@ Mgxs::calculate_xs(int gin, double sqrtkT, Direction u, #else int tid = 0; #endif - set_temperature_index(sqrtkT); - set_angle_index(u); + set_temperature_index(p.sqrtkT_); + set_angle_index(p.u_local()); XsData* xs_t = &xs[cache[tid].t]; - total_xs = xs_t->total(cache[tid].a, gin); - abs_xs = xs_t->absorption(cache[tid].a, gin); - - nu_fiss_xs = fissionable ? xs_t->nu_fission(cache[tid].a, gin) : 0.; + p.macro_xs_.total = xs_t->total(cache[tid].a, p.g_); + p.macro_xs_.absorption = xs_t->absorption(cache[tid].a, p.g_); + p.macro_xs_.nu_fission = + fissionable ? xs_t->nu_fission(cache[tid].a, p.g_) : 0.; } //============================================================================== diff --git a/src/particle.cpp b/src/particle.cpp index da93c6723..b430d3a06 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -202,12 +202,12 @@ Particle::transport() model::materials[material_]->calculate_xs(*this); } } else { - // Get the MG data - data::mg.macro_xs_[material_].calculate_xs(g_, sqrtkT_, - this->u_local(), macro_xs_.total, macro_xs_.absorption, - macro_xs_.nu_fission); + // Get the MG data; unlike the CE case above, we have to re-calculate + // cross sections for every collision since the cross sections may + // be angle-dependent + data::mg.macro_xs_[material_].calculate_xs(*this); - // Finally, update the particle group since we know we are multi-group + // Update the particle's group while we know we are multi-group g_last_ = g_; } } else { From 4272e9c3f59b950f5b29e8a8ab4f2786c1dbb205 Mon Sep 17 00:00:00 2001 From: Adam G Nelson Date: Fri, 15 Nov 2019 18:50:20 -0600 Subject: [PATCH 9/9] Changed nuc_xs and macro_xs to references vice pointers in score_general_mg --- src/tallies/tally_scoring.cpp | 213 ++++++++++++++++------------------ 1 file changed, 103 insertions(+), 110 deletions(-) diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index 1d23ebb68..7ca279966 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -1375,14 +1375,14 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } // For shorthand, assign pointers to the material and nuclide xs set - Mgxs* nuc_xs = &data::mg.nuclides_[i_nuclide]; - Mgxs* macro_xs = &data::mg.macro_xs_[p->material_]; + auto& nuc_xs = data::mg.nuclides_[i_nuclide]; + auto& macro_xs = data::mg.macro_xs_[p->material_]; // Find the temperature and angle indices of interest - macro_xs->set_angle_index(p_u); + macro_xs.set_angle_index(p_u); if (i_nuclide >= 0) { - nuc_xs->set_temperature_index(p->sqrtkT_); - nuc_xs->set_angle_index(p_u); + nuc_xs.set_temperature_index(p->sqrtkT_); + nuc_xs.set_angle_index(p_u); } for (auto i = 0; i < tally.scores_.size(); ++i) { @@ -1426,12 +1426,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } //TODO: should flux be multiplied in above instead of below? if (i_nuclide >= 0) { - score *= flux * atom_density * nuc_xs->get_xs(MG_GET_XS_TOTAL, p_g) - / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); + score *= flux * atom_density * nuc_xs.get_xs(MG_GET_XS_TOTAL, p_g) + / macro_xs.get_xs(MG_GET_XS_TOTAL, p_g); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_TOTAL, p_g); + score = atom_density * flux * nuc_xs.get_xs(MG_GET_XS_TOTAL, p_g); } else { score = p->macro_xs_.total * flux; } @@ -1453,17 +1453,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_; } if (i_nuclide >= 0) { - score *= flux * nuc_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) - / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); + score *= flux * nuc_xs.get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) + / macro_xs.get_xs(MG_GET_XS_TOTAL, p_g); } else { - score *= flux * macro_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) - / macro_xs->get_xs(MG_GET_XS_TOTAL, p_g); + score *= flux * macro_xs.get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g) + / macro_xs.get_xs(MG_GET_XS_TOTAL, p_g); } } else { if (i_nuclide >= 0) { - score = flux * nuc_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); + score = flux * nuc_xs.get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_INVERSE_VELOCITY, p_g); } } break; @@ -1477,17 +1477,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // weight entering the collision as the estimator for the reaction rate score = p->wgt_last_ * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs( + score *= atom_density * nuc_xs.get_xs( MG_GET_XS_SCATTER_FMU_MULT, p->g_last_, &p->g_, &p->mu_, nullptr) - / macro_xs->get_xs( + / macro_xs.get_xs( MG_GET_XS_SCATTER_FMU_MULT, p->g_last_, &p->g_, &p->mu_, nullptr); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * nuc_xs->get_xs( + score = atom_density * flux * nuc_xs.get_xs( MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); } else { - score = flux * macro_xs->get_xs( + score = flux * macro_xs.get_xs( MG_GET_XS_SCATTER_MULT, p_g, nullptr, &p->mu_, nullptr); } } @@ -1507,16 +1507,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // adjust the score by the actual probability for that nuclide. if (i_nuclide >= 0) { score *= atom_density - * nuc_xs->get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, - &p->mu_, nullptr) - / macro_xs->get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, - &p->mu_, nullptr); + * nuc_xs.get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, + &p->mu_, nullptr) + / macro_xs.get_xs(MG_GET_XS_SCATTER_FMU, p->g_last_, &p->g_, + &p->mu_, nullptr); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_SCATTER, p_g); + score = atom_density * flux * nuc_xs.get_xs(MG_GET_XS_SCATTER, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_SCATTER, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_SCATTER, p_g); } } break; @@ -1536,13 +1536,13 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_ * flux; } if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_ABSORPTION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_ABSORPTION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { score = atom_density * flux - * nuc_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + * nuc_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } else { score = p->macro_xs_.absorption * flux; } @@ -1566,17 +1566,17 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_last_ * flux; } if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= macro_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g); + score = atom_density * flux * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } } break; @@ -1598,11 +1598,11 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // nu-fission score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_NU_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_NU_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= macro_xs->get_xs(MG_GET_XS_NU_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs.get_xs(MG_GET_XS_NU_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { // Skip any non-fission events @@ -1614,16 +1614,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // score. score = simulation::keff * p->wgt_bank_ * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } } } else { if (i_nuclide >= 0) { score = atom_density * flux - * nuc_xs->get_xs(MG_GET_XS_NU_FISSION, p_g); + * nuc_xs.get_xs(MG_GET_XS_NU_FISSION, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_NU_FISSION, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_NU_FISSION, p_g); } } break; @@ -1646,11 +1646,11 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { score *= atom_density * - nuc_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + nuc_xs.get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } else { - score *= macro_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + score *= macro_xs.get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { // Skip any non-fission events @@ -1665,16 +1665,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto prompt_frac = 1. - n_delayed / static_cast(p->n_bank_); score = simulation::keff * p->wgt_bank_ * prompt_frac * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } } } else { if (i_nuclide >= 0) { score = atom_density * flux * - nuc_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); + nuc_xs.get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_PROMPT_NU_FISSION, p_g); } } break; @@ -1694,7 +1694,7 @@ score_general_mg(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 - double abs_xs = macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + double abs_xs = macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); if (abs_xs > 0.) { if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; @@ -1706,12 +1706,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) + score *= nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } else { - score *= macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) + score *= macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } score_fission_delayed_dg(i_tally, d_bin, score, score_index); @@ -1723,10 +1723,10 @@ score_general_mg(const Particle* p, int i_tally, int start_index, // delayed-nu-fission xs to the absorption xs score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) + score *= nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) / abs_xs; } else { - score *= macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) + score *= macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g) / abs_xs; } } @@ -1752,8 +1752,8 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = simulation::keff * p->wgt_bank_ / p->n_bank_ * p->n_delayed_bank_[d-1] * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } @@ -1765,8 +1765,8 @@ score_general_mg(const Particle* p, int i_tally, int start_index, score = simulation::keff * p->wgt_bank_ / p->n_bank_ * n_delayed * flux; if (i_nuclide >= 0) { - score *= atom_density * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + score *= atom_density * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } } } @@ -1780,10 +1780,10 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d_bin = 0; d_bin < filt.n_bins(); ++d_bin) { auto d = filt.groups()[d_bin] - 1; if (i_nuclide >= 0) { - score = flux * atom_density * nuc_xs->get_xs( + score = flux * atom_density * nuc_xs.get_xs( MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } else { - score = flux * macro_xs->get_xs( + score = flux * macro_xs.get_xs( MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); @@ -1792,9 +1792,9 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } else { if (i_nuclide >= 0) { score = flux * atom_density * - nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); + nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g); } } } @@ -1807,7 +1807,7 @@ score_general_mg(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 - double abs_xs = macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + double abs_xs = macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); if (abs_xs > 0) { if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; @@ -1819,15 +1819,15 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; score = p->wgt_absorb_ * flux; if (i_nuclide >= 0) { - score *= nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, - nullptr, nullptr, &d) - * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) / abs_xs; + score *= nuc_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, + nullptr, nullptr, &d) + * nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } else { - score *= macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) / abs_xs; + score *= macro_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } @@ -1841,16 +1841,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d = 0; d < data::mg.num_delayed_groups_; ++d) { if (i_nuclide >= 0) { score += p->wgt_absorb_ * flux - * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) / abs_xs; + * nuc_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } else { score += p->wgt_absorb_ * flux - * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, - nullptr, nullptr, &d) / abs_xs; + * macro_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, + nullptr, nullptr, &d) / abs_xs; } } } @@ -1873,14 +1873,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, if (d != -1) { if (i_nuclide >= 0) { score += simulation::keff * atom_density * bank.wgt * flux - * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, - &d) - * nuc_xs->get_xs(MG_GET_XS_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_FISSION, p_g); + * nuc_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, &d) + * nuc_xs.get_xs(MG_GET_XS_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_FISSION, p_g); } else { score += simulation::keff * bank.wgt * flux - * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d); + * macro_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, &d); } if (tally.delayedgroup_filter_ != C_NONE) { auto i_dg_filt = tally.filters()[tally.delayedgroup_filter_]; @@ -1911,16 +1909,16 @@ score_general_mg(const Particle* p, int i_tally, int start_index, auto d = filt.groups()[d_bin] - 1; if (i_nuclide >= 0) { score += atom_density * flux - * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, - nullptr, &d); + * nuc_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } else { score += flux - * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, - nullptr, &d); + * macro_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, + nullptr, &d) + * macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, + nullptr, &d); } score_fission_delayed_dg(i_tally, d_bin, score, score_index); } @@ -1930,16 +1928,12 @@ score_general_mg(const Particle* p, int i_tally, int start_index, for (auto d = 0; d < data::mg.num_delayed_groups_; ++d) { if (i_nuclide >= 0) { score += atom_density * flux - * nuc_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, - &d) - * nuc_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, - nullptr, &d); + * nuc_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, &d) + * nuc_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } else { score += flux - * macro_xs->get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, - nullptr, &d) - * macro_xs->get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, - nullptr, &d); + * macro_xs.get_xs(MG_GET_XS_DECAY_RATE, p_g, nullptr, nullptr, &d) + * macro_xs.get_xs(MG_GET_XS_DELAYED_NU_FISSION, p_g, nullptr, nullptr, &d); } } } @@ -1964,19 +1958,18 @@ score_general_mg(const Particle* p, int i_tally, int start_index, } if (i_nuclide >= 0) { score *= atom_density - * nuc_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + * nuc_xs.get_xs(MG_GET_XS_KAPPA_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } else { score *= - macro_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g) - / macro_xs->get_xs(MG_GET_XS_ABSORPTION, p_g); + macro_xs.get_xs(MG_GET_XS_KAPPA_FISSION, p_g) + / macro_xs.get_xs(MG_GET_XS_ABSORPTION, p_g); } } else { if (i_nuclide >= 0) { - score = atom_density * flux * nuc_xs->get_xs(MG_GET_XS_KAPPA_FISSION, - p_g); + score = atom_density * flux * nuc_xs.get_xs(MG_GET_XS_KAPPA_FISSION, p_g); } else { - score = flux * macro_xs->get_xs(MG_GET_XS_KAPPA_FISSION, p_g); + score = flux * macro_xs.get_xs(MG_GET_XS_KAPPA_FISSION, p_g); } } break;