Delayed neutron fraction by isotope with IFP method

docs/source/usersguide/kinetics.rst

@@ -84,6 +84,11 @@ total :math:`\beta_{\text{eff}}` specified by providing a 6-group
     # Add DelayedGroupFilter to enable group-wise tallies
     beta_tally.filters = [openmc.DelayedGroupFilter(list(range(1, 7)))]
 
+A breakdown on isotopes of the :math:`\beta_{\text{eff}}` can be estimated,
+e.g., for :math:`^{U235}U` and :math:`^{U238}U`, by defining::
+
+    beta_tally.nuclides = ['U235', 'U238']
+
 Here is an example showing how to declare the three available IFP scores in a
 single tally::
 
@@ -113,7 +118,7 @@ for ``ifp-denominator``:
     \beta_{\text{eff}} = \frac{S_{\text{ifp-beta-numerator}}}{S_{\text{ifp-denominator}}}
 
 The kinetics parameters can be retrieved directly from a statepoint file using
-the :meth:`openmc.StatePoint.ifp_results` method::
+the :meth:`openmc.StatePoint.get_kinetics_parameters` method::
 
     with openmc.StatePoint(output_path) as sp:
         generation_time, beta_eff = sp.get_kinetics_parameters()

include/openmc/bank.h

@@ -26,10 +26,14 @@ extern SharedArray<SourceSite> fission_bank;
 
 extern vector<vector<int>> ifp_source_delayed_group_bank;
 
+extern vector<vector<int>> ifp_source_ancestor_nuclide_bank;
+
 extern vector<vector<double>> ifp_source_lifetime_bank;
 
 extern vector<vector<int>> ifp_fission_delayed_group_bank;
 
+extern vector<vector<int>> ifp_fission_ancestor_nuclide_bank;
+
 extern vector<vector<double>> ifp_fission_lifetime_bank;
 
 extern vector<int64_t> progeny_per_particle;

include/openmc/ifp.h

@@ -21,13 +21,16 @@ bool is_generation_time_or_both();
 //! Resize IFP vectors
 //!
 //! \param[in,out] delayed_groups List of delayed group numbers
+//! \param[in,out] ancestors List of event nuclides
 //! \param[in,out] lifetimes List of lifetimes
 //! \param[in] n  Dimension to resize vectors
 template<typename T, typename U>
-void resize_ifp_data(vector<T>& delayed_groups, vector<U>& lifetimes, int64_t n)
+void resize_ifp_data(vector<T>& delayed_groups, vector<T>& ancestors,
+  vector<U>& lifetimes, int64_t n)
 {
   if (is_beta_effective_or_both()) {
     delayed_groups.resize(n);
+    ancestors.resize(n);
   }
   if (is_generation_time_or_both()) {
     lifetimes.resize(n);
@@ -84,9 +87,10 @@ void resize_simulation_ifp_banks();
 //!
 //! \param[in] i_bank Index in the fission banks
 //! \param[in,out] delayed_groups Delayed group numbers
+//! \param[in,out] ancestors List of Ancestor nuclides
 //! \param[in,out] lifetimes Lifetimes lists
-void copy_ifp_data_from_fission_banks(
-  int i_bank, vector<int>& delayed_groups, vector<double>& lifetimes);
+void copy_ifp_data_from_fission_banks(int i_bank, vector<int>& delayed_groups,
+  vector<int>& ancestors, vector<double>& lifetimes);
 
 #ifdef OPENMC_MPI
 
@@ -101,9 +105,11 @@ struct DeserializationInfo {
 //!
 //! \param[in] n_generation Number of generations
 //! \param[in] delayed_groups List of delayed group numbers lists
+//! \param[in] ancestors List of Ancestor nuclide
 //! \param[in] lifetimes List of lifetimes lists
 void broadcast_ifp_n_generation(int& n_generation,
   const vector<vector<int>>& delayed_groups,
+  const vector<vector<int>>& ancestors,
   const vector<vector<double>>& lifetimes);
 
 //! Send IFP data using MPI.
@@ -115,11 +121,14 @@ void broadcast_ifp_n_generation(int& n_generation,
 //! \param[in] requests MPI requests
 //! \param[in] delayed_groups List of delayed group numbers lists
 //! \param[out] send_delayed_groups Delayed group numbers buffer
+//! \param[in] ancestors List of event nuclide lists
+//! \param[out] send_ancestors Event nuclide buffer
 //! \param[in] lifetimes List of lifetimes lists
 //! \param[out] send_lifetimes Lifetimes buffer
 void send_ifp_info(int64_t idx, int64_t n, int n_generation, int neighbor,
   vector<MPI_Request>& requests, const vector<vector<int>>& delayed_groups,
-  vector<int>& send_delayed_groups, const vector<vector<double>>& lifetimes,
+  vector<int>& send_delayed_groups, const vector<vector<int>>& ancestors,
+  vector<int>& send_ancestors, const vector<vector<double>>& lifetimes,
   vector<double>& send_lifetimes);
 
 //! Receive IFP data using MPI.
@@ -130,21 +139,25 @@ void send_ifp_info(int64_t idx, int64_t n, int n_generation, int neighbor,
 //! \param[in] neighbor Index of the neighboring processor
 //! \param[in] requests MPI requests
 //! \param[in] delayed_groups List of delayed group numbers
+//! \param[in] ancestors List of event nuclide
 //! \param[in] lifetimes List of lifetimes
 //! \param[out] deserialization Information to deserialize the received data
 void receive_ifp_data(int64_t idx, int64_t n, int n_generation, int neighbor,
   vector<MPI_Request>& requests, vector<int>& delayed_groups,
-  vector<double>& lifetimes, vector<DeserializationInfo>& deserialization);
+  vector<int>& ancestors, vector<double>& lifetimes,
+  vector<DeserializationInfo>& deserialization);
 
 //! Copy partial IFP data from local lists to source banks.
 //!
 //! \param[in] idx Index of the first site
 //! \param[in] n Number of sites to copy
 //! \param[in] i_bank Index in the IFP source banks
 //! \param[in] delayed_groups List of delayed group numbers lists
+//! \param[in] ancestors List of event nuclide
 //! \param[in] lifetimes List of lifetimes lists
 void copy_partial_ifp_data_to_source_banks(int64_t idx, int n, int64_t i_bank,
   const vector<vector<int>>& delayed_groups,
+  const vector<vector<int>>& ancestors,
   const vector<vector<double>>& lifetimes);
 
 //! Deserialize IFP information received using MPI and store it in
@@ -153,34 +166,40 @@ void copy_partial_ifp_data_to_source_banks(int64_t idx, int n, int64_t i_bank,
 //! \param[in] n_generation Number of generations
 //! \param[out] deserialization Information to deserialize the received data
 //! \param[in] delayed_groups List of delayed group numbers
+//! \param[in] ancestors List of event nuclide
 //! \param[in] lifetimes List of lifetimes
 void deserialize_ifp_info(int n_generation,
   const vector<DeserializationInfo>& deserialization,
-  const vector<int>& delayed_groups, const vector<double>& lifetimes);
+  const vector<int>& delayed_groups, const vector<int>& ancestors,
+  const vector<double>& lifetimes);
 
 #endif
 
 //! Copy IFP temporary vectors to source banks.
 //!
 //! \param[in] delayed_groups List of delayed group numbers lists
+//! \param[in] ancestors List of event nuclide lists
 //! \param[in] lifetimes List of lifetimes lists
 void copy_complete_ifp_data_to_source_banks(
   const vector<vector<int>>& delayed_groups,
+  const vector<vector<int>>& ancestors,
   const vector<vector<double>>& lifetimes);
 
 //! Allocate temporary vectors for IFP data.
 //!
 //! \param[in,out] delayed_groups List of delayed group numbers lists
+//! \param[in,out] ancestors List of event nuclide lists
 //! \param[in,out] lifetimes List of delayed group numbers lists
-void allocate_temporary_vector_ifp(
-  vector<vector<int>>& delayed_groups, vector<vector<double>>& lifetimes);
+void allocate_temporary_vector_ifp(vector<vector<int>>& delayed_groups,
+  vector<vector<int>>& ancestors, vector<vector<double>>& lifetimes);
 
 //! Copy local IFP data to IFP fission banks.
 //!
 //! \param[in] delayed_groups_ptr Pointer to delayed group numbers
+//! \param[in] ancestors_ptr Pointer to delayed group numbers
 //! \param[in] lifetimes_ptr Pointer to lifetimes
-void copy_ifp_data_to_fission_banks(
-  const vector<int>* delayed_groups_ptr, const vector<double>* lifetimes_ptr);
+void copy_ifp_data_to_fission_banks(const vector<int>* delayed_groups_ptr,
+  const vector<int>* ancestors_ptr, const vector<double>* lifetimes_ptr);
 
 } // namespace openmc
 

include/openmc/particle_data.h

@@ -52,6 +52,7 @@ struct SourceSite {
 
   // Extra attributes that don't show up in source written to file
   int parent_nuclide {-1};
+  int ancestor_nuclide {-1};
   int64_t parent_id;
   int64_t progeny_id;
 };
@@ -65,6 +66,7 @@ struct CollisionTrackSite {
   double wgt {1.0};
   int event_mt {0};
   int delayed_group {0};
+  int ancestor_nuclide {0};
   int cell_id {0};
   int nuclide_id;
   int material_id {0};
@@ -517,6 +519,7 @@ class ParticleData : public GeometryState {
   int event_mt_;
   int delayed_group_ {0};
   int parent_nuclide_ {-1};
+  int ancestor_nuclide_ {-1};
 
   int n_bank_ {0};
   double bank_second_E_ {0.0};
@@ -652,8 +655,10 @@ class ParticleData : public GeometryState {
   const int& event_mt() const { return event_mt_; }
   int& delayed_group() { return delayed_group_; } // delayed group
   const int& delayed_group() const { return delayed_group_; }
-  const int& parent_nuclide() const { return parent_nuclide_; }
   int& parent_nuclide() { return parent_nuclide_; } // Parent nuclide
+  const int& parent_nuclide() const { return parent_nuclide_; }
+  int& ancestor_nuclide() { return ancestor_nuclide_; } // Ancestor nuclide
+  const int& ancestor_nuclide() const { return ancestor_nuclide_; }
 
   // Post-collision data
   double& bank_second_E()

include/openmc/settings.h

@@ -71,6 +71,7 @@ extern "C" bool entropy_on; //!< calculate Shannon entropy?
 extern "C" bool
   event_based;      //!< use event-based mode (instead of history-based)
 extern bool ifp_on; //!< Use IFP for kinetics parameters?
+extern bool ifp_beta_nuclide;    //!< Use IFP for beta calculation by nuclide?
 extern bool legendre_to_tabular; //!< convert Legendre distributions to tabular?
 extern bool material_cell_offsets;   //!< create material cells offsets?
 extern "C" bool output_summary;      //!< write summary.h5?

openmc/lib/core.py

@@ -30,6 +30,7 @@ class _SourceSite(Structure):
                 ('surf_id', c_int),
                 ('particle', c_int),
                 ('parent_nuclide', c_int),
+                ('ancestor_nuclide', c_int),
                 ('parent_id', c_int64),
                 ('progeny_id', c_int64)]
 

openmc/model/model.py

@@ -235,20 +235,23 @@ def _get_all_materials(self) -> dict[int, openmc.Material]:
 
         return materials
 
-    def add_kinetics_parameters_tallies(self, num_groups: int | None = None):
+    def add_kinetics_parameters_tallies(self, num_groups: int | None = None, nuclides: Sequence[str] | None = None):
         """Add tallies for calculating kinetics parameters using the IFP method.
 
         This method adds tallies to the model for calculating two kinetics
         parameters, the generation time and the effective delayed neutron
         fraction (beta effective). After a model is run, these parameters can be
-        determined through the :meth:`openmc.StatePoint.ifp_results` method.
+        determined through the :meth:`openmc.StatePoint.get_kinetics_parameters` method.
 
         Parameters
         ----------
         num_groups : int, optional
             Number of precursor groups to filter the delayed neutron fraction.
             If None, only the total effective delayed neutron fraction is
             tallied.
+        nuclides : int, optional 
+            Nuclides to calculate separate kinetic parameters for.
+            If None, do not separate kinetic parameters per nuclide.
 
         """
         if not any('ifp-time-numerator' in t.scores for t in self.tallies):
@@ -260,6 +263,8 @@ def add_kinetics_parameters_tallies(self, num_groups: int | None = None):
             beta_tally.scores = ['ifp-beta-numerator']
             if num_groups is not None:
                 beta_tally.filters = [openmc.DelayedGroupFilter(list(range(1, num_groups + 1)))]
+            if nuclides:
+                beta_tally.nuclides = list(nuclides)
             self.tallies.append(beta_tally)
         if not any('ifp-denominator' in t.scores for t in self.tallies):
             denom_tally = openmc.Tally(name='IFP denominator')

openmc/statepoint.py

@@ -761,21 +761,17 @@ def get_ufloat(tally, score):
             return uarray(tally.get_values(scores=[score]),
                           tally.get_values(scores=[score], value='std_dev'))
 
-        denom_values = get_ufloat(denom_tally, 'ifp-denominator')
+        denom_value = get_ufloat(denom_tally, 'ifp-denominator').squeeze()
         if gen_time_tally is None:
             generation_time = None
         else:
-            gen_time_values = get_ufloat(gen_time_tally, 'ifp-time-numerator')
-            gen_time_values /= denom_values*self.keff
-            generation_time = gen_time_values.flatten()[0]
+            generation_time = get_ufloat(gen_time_tally, 'ifp-time-numerator').squeeze()
+            generation_time /= denom_value*self.keff
 
         if beta_tally is None:
             beta_effective = None
         else:
-            beta_values = get_ufloat(beta_tally, 'ifp-beta-numerator')
-            beta_values /= denom_values
-            beta_effective = beta_values.flatten()
-            if beta_effective.size == 1:
-                beta_effective = beta_effective[0]
+            beta_effective = get_ufloat(beta_tally, 'ifp-beta-numerator').squeeze()
+            beta_effective /= denom_value
 
         return KineticsParameters(generation_time, beta_effective)

src/bank.cpp

@@ -31,10 +31,14 @@ SharedArray<SourceSite> fission_bank;
 
 vector<vector<int>> ifp_source_delayed_group_bank;
 
+vector<vector<int>> ifp_source_ancestor_nuclide_bank;
+
 vector<vector<double>> ifp_source_lifetime_bank;
 
 vector<vector<int>> ifp_fission_delayed_group_bank;
 
+vector<vector<int>> ifp_fission_ancestor_nuclide_bank;
+
 vector<vector<double>> ifp_fission_lifetime_bank;
 
 // Each entry in this vector corresponds to the number of progeny produced
@@ -56,8 +60,10 @@ void free_memory_bank()
   simulation::fission_bank.clear();
   simulation::progeny_per_particle.clear();
   simulation::ifp_source_delayed_group_bank.clear();
+  simulation::ifp_source_ancestor_nuclide_bank.clear();
   simulation::ifp_source_lifetime_bank.clear();
   simulation::ifp_fission_delayed_group_bank.clear();
+  simulation::ifp_fission_ancestor_nuclide_bank.clear();
   simulation::ifp_fission_lifetime_bank.clear();
 }
 
@@ -92,6 +98,7 @@ void sort_fission_bank()
   SourceSite* sorted_bank;
   vector<SourceSite> sorted_bank_holder;
   vector<vector<int>> sorted_ifp_delayed_group_bank;
+  vector<vector<int>> sorted_ifp_ancestor_nuclide_bank;
   vector<vector<double>> sorted_ifp_lifetime_bank;
 
   // If there is not enough space, allocate a temporary vector and point to it
@@ -104,8 +111,8 @@ void sort_fission_bank()
   }
 
   if (settings::ifp_on) {
-    allocate_temporary_vector_ifp(
-      sorted_ifp_delayed_group_bank, sorted_ifp_lifetime_bank);
+    allocate_temporary_vector_ifp(sorted_ifp_delayed_group_bank,
+      sorted_ifp_ancestor_nuclide_bank, sorted_ifp_lifetime_bank);
   }
 
   // Use parent and progeny indices to sort fission bank
@@ -119,17 +126,17 @@ void sort_fission_bank()
     }
     sorted_bank[idx] = site;
     if (settings::ifp_on) {
-      copy_ifp_data_from_fission_banks(
-        i, sorted_ifp_delayed_group_bank[idx], sorted_ifp_lifetime_bank[idx]);
+      copy_ifp_data_from_fission_banks(i, sorted_ifp_delayed_group_bank[idx],
+        sorted_ifp_ancestor_nuclide_bank[idx], sorted_ifp_lifetime_bank[idx]);
     }
   }
 
   // Copy sorted bank into the fission bank
   std::copy(sorted_bank, sorted_bank + simulation::fission_bank.size(),
     simulation::fission_bank.data());
   if (settings::ifp_on) {
-    copy_ifp_data_to_fission_banks(
-      sorted_ifp_delayed_group_bank.data(), sorted_ifp_lifetime_bank.data());
+    copy_ifp_data_to_fission_banks(sorted_ifp_delayed_group_bank.data(),
+      sorted_ifp_ancestor_nuclide_bank.data(), sorted_ifp_lifetime_bank.data());
   }
 }