Temperature feedback support in the random ray solver.

docs/source/usersguide/geometry.rst

@@ -248,6 +248,28 @@ The classes :class:`Halfspace`, :class:`Intersection`, :class:`Union`, and
 :class:`Complement` and all instances of :class:`openmc.Region` and can be
 assigned to the :attr:`Cell.region` attribute.
 
+Cells also contain :attr:`Cell.temperature` and :attr:`Cell.density`
+attributes which override the temperature and density of the fill. These can
+be quite useful when temperatures and densities are spatially varying, as the
+alternative would be to add a unique :class:`Material` for each permutation of
+temperature, density, and composition. You can set the temperature (K) and
+density (g/cc) of a cell like so::
+
+  fuel.temperature = 800.0
+  fuel.density = 10.0
+
+The real utility of cell temperatures and densities occurs when a cell is
+replicated across the geometry, such as when a cell is the root geometric element
+in a replicated :ref:`universe<usersguide_universes>` or :ref:`lattice
+<usersguide_lattices>`. In those cases, you can provide a list of temperatures
+and densities to apply a temperature/density field to all of the distributed cells::
+
+  fuel.temperature = [800.0, 900.0, 800.0, 900.0]
+  fuel.density = [10.0, 9.0, 10.0, 9.0]
+
+In this example, the fuel cell is distributed four times in the geometry. Each
+distributed instance then receives its own temperature and density.
+
 .. _usersguide_universes:
 
 ---------

docs/source/usersguide/random_ray.rst

@@ -645,7 +645,9 @@ model to use these multigroup cross sections. An example is given below::
       overwrite_mgxs_library=False,
       mgxs_path="mgxs.h5",
       correction=None,
-      source_energy=None
+      source_energy=None,
+      temperatures=None,
+      temperature_settings=None
   )
 
 The most important parameter to set is the ``method`` parameter, which can be
@@ -732,6 +734,18 @@ distribution for MGXS generation as::
 
   source_energy = openmc.stats.delta_function(2.45e6)
 
+The ``temperatures`` parameter can be provided if temperature-dependent
+multi-group cross sections are desired for multi-physics simulations. An
+individual cross section generation calculation is run for each temperature
+provided, where the materials in the model are set to the temperature. The
+temperature settings used during cross section generation can be specified with the
+``temperature_settings`` parameter. If no ``temperature_settings`` are provided,
+the settings contained in the model will be used. This approach yields isothermal
+cross section interpolation tables, which can be inaccurate for systems with
+large differences between temperatures in each material (often the case in
+fission reactors). If a more sophisticated temperature-dependence is required,
+we recommend generating cross sections manually.
+
 Ultimately, the methods described above are all just approximations.
 Approximations in the generated MGXS data will fundamentally limit the potential
 accuracy of the random ray solver. However, the methods described above are all

include/openmc/mgxs.h

@@ -113,6 +113,11 @@ class Mgxs {
     const vector<Mgxs*>& micros, const vector<double>& atom_densities,
     int num_group, int num_delay);
 
+  //! \brief Get the number of temperature data points.
+  //!
+  //! @return The number of temperature data points for this MGXS
+  inline int n_temperature_points() { return kTs.size(); }
+
   //! \brief Provides a cross section value given certain parameters
   //!
   //! @param xstype Type of cross section requested, according to the

include/openmc/random_ray/flat_source_domain.h

@@ -100,16 +100,17 @@ class FlatSourceDomain {
   // in model::cells
   vector<int64_t> source_region_offsets_;
 
-  // 2D arrays stored in 1D representing values for all materials x energy
-  // groups
+  // 3D arrays stored in 1D representing values for all materials x temperature
+  // points x energy groups
   int n_materials_;
+  int ntemperature_;
   vector<double> sigma_t_;
   vector<double> nu_sigma_f_;
   vector<double> sigma_f_;
   vector<double> chi_;
 
-  // 3D arrays stored in 1D representing values for all materials x energy
-  // groups x energy groups
+  // 4D arrays stored in 1D representing values for all materials x temperature
+  // points x energy groups x energy groups
   vector<double> sigma_s_;
 
   // The abstract container holding all source region-specific data

include/openmc/random_ray/random_ray.h

@@ -65,6 +65,7 @@ class RandomRay : public Particle {
   vector<double> mesh_fractional_lengths_;
 
   int negroups_;
+  int ntemperature_;
   FlatSourceDomain* domain_ {nullptr}; // pointer to domain that has flat source
                                        // data needed for ray transport
   double distance_travelled_ {0};

include/openmc/random_ray/source_region.h

@@ -146,6 +146,7 @@ class SourceRegionHandle {
 
   // Scalar fields
   int* material_;
+  int* temperature_idx_;
   double* density_mult_;
   int* is_small_;
   int* n_hits_;
@@ -199,6 +200,9 @@ class SourceRegionHandle {
   double& density_mult() { return *density_mult_; }
   const double density_mult() const { return *density_mult_; }
 
+  int& temperature_idx() { return *temperature_idx_; }
+  const int temperature_idx() const { return *temperature_idx_; }
+
   int& is_small() { return *is_small_; }
   const int is_small() const { return *is_small_; }
 
@@ -319,8 +323,9 @@ class SourceRegion {
 
   //---------------------------------------
   // Scalar fields
-
-  int material_ {0};          //!< Index in openmc::model::materials array
+  int material_ {0}; //!< Index in openmc::model::materials array
+  int temperature_idx_ {
+    0}; //!< Index into the MGXS array representing temperature
   double density_mult_ {1.0}; //!< A density multiplier queried from the cell
                               //!< corresponding to the source region.
   OpenMPMutex lock_;
@@ -400,6 +405,9 @@ class SourceRegionContainer {
   int& material(int64_t sr) { return material_[sr]; }
   const int material(int64_t sr) const { return material_[sr]; }
 
+  int& temperature_idx(int64_t sr) { return temperature_idx_[sr]; }
+  const int temperature_idx(int64_t sr) const { return temperature_idx_[sr]; }
+
   double& density_mult(int64_t sr) { return density_mult_[sr]; }
   const double density_mult(int64_t sr) const { return density_mult_[sr]; }
 
@@ -634,6 +642,7 @@ class SourceRegionContainer {
 
   // SoA storage for scalar fields (one item per source region)
   vector<int> material_;
+  vector<int> temperature_idx_;
   vector<double> density_mult_;
   vector<int> is_small_;
   vector<int> n_hits_;

openmc/examples.py

@@ -656,13 +656,22 @@ def slab_mg(num_regions=1, mat_names=None, mgxslib_name='2g.h5') -> openmc.Model
     return model
 
 
-def random_ray_lattice() -> openmc.Model:
-    """Create a 2x2 PWR pincell asymmetrical lattic eexample.
+def random_ray_lattice(second_temp: bool = False) -> openmc.Model:
+    """Create a 2x2 PWR pincell asymmetrical lattice example.
 
     This model is a 2x2 reflective lattice of fuel pins with one of the lattice
     locations having just moderator instead of a fuel pin. It uses 7 group
     cross section data.
 
+    Parameters
+    ----------
+    second_temp : bool, optional
+        Whether or not the cross sections should contain two temperature datapoints.
+        The first data point is the C5G7 cross sections, which corresponds to a temperature
+        of 294 K. The second data point is the C5G7 cross sections multiplied by 1/2,
+        which corresponds to a temperature of 394 K. This temperature dependence is
+        ficticious; it is used for testing temperature feedback in the random ray solver.
+
     Returns
     -------
     model : openmc.Model
@@ -675,63 +684,75 @@ def random_ray_lattice() -> openmc.Model:
     # Create MGXS data for the problem
 
     # Instantiate the energy group data
+    # MGXS for the UO2 pins.
     group_edges = [1e-5, 0.0635, 10.0, 1.0e2, 1.0e3, 0.5e6, 1.0e6, 20.0e6]
     groups = openmc.mgxs.EnergyGroups(group_edges)
 
-    # Instantiate the 7-group (C5G7) cross section data
-    uo2_xsdata = openmc.XSdata('UO2', groups)
-    uo2_xsdata.order = 0
-    uo2_xsdata.set_total(
-        [0.1779492, 0.3298048, 0.4803882, 0.5543674, 0.3118013, 0.3951678,
-         0.5644058])
-    uo2_xsdata.set_absorption([8.0248e-03, 3.7174e-03, 2.6769e-02, 9.6236e-02,
-                               3.0020e-02, 1.1126e-01, 2.8278e-01])
-    scatter_matrix = np.array(
+    uo2_total = np.array([0.1779492, 0.3298048, 0.4803882, 0.5543674, 0.3118013, 0.3951678,
+                          0.5644058])
+    uo2_abs = np.array([8.0248e-03, 3.7174e-03, 2.6769e-02, 9.6236e-02, 3.0020e-02,
+                        1.1126e-01, 2.8278e-01])
+    uo2_scatter_matrix = np.array(
         [[[0.1275370, 0.0423780, 0.0000094, 0.0000000, 0.0000000, 0.0000000, 0.0000000],
-          [0.0000000, 0.3244560, 0.0016314, 0.0000000,
-              0.0000000, 0.0000000, 0.0000000],
-          [0.0000000, 0.0000000, 0.4509400, 0.0026792,
-              0.0000000, 0.0000000, 0.0000000],
-          [0.0000000, 0.0000000, 0.0000000, 0.4525650,
-              0.0055664, 0.0000000, 0.0000000],
-          [0.0000000, 0.0000000, 0.0000000, 0.0001253,
-              0.2714010, 0.0102550, 0.0000000],
-          [0.0000000, 0.0000000, 0.0000000, 0.0000000,
-              0.0012968, 0.2658020, 0.0168090],
+          [0.0000000, 0.3244560, 0.0016314, 0.0000000, 0.0000000, 0.0000000, 0.0000000],
+          [0.0000000, 0.0000000, 0.4509400, 0.0026792, 0.0000000, 0.0000000, 0.0000000],
+          [0.0000000, 0.0000000, 0.0000000, 0.4525650, 0.0055664, 0.0000000, 0.0000000],
+          [0.0000000, 0.0000000, 0.0000000, 0.0001253, 0.2714010, 0.0102550, 0.0000000],
+          [0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0012968, 0.2658020, 0.0168090],
           [0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0085458, 0.2730800]]])
-    scatter_matrix = np.rollaxis(scatter_matrix, 0, 3)
-    uo2_xsdata.set_scatter_matrix(scatter_matrix)
-    uo2_xsdata.set_fission([7.21206e-03, 8.19301e-04, 6.45320e-03,
-                            1.85648e-02, 1.78084e-02, 8.30348e-02,
-                            2.16004e-01])
-    uo2_xsdata.set_nu_fission([2.005998e-02, 2.027303e-03, 1.570599e-02,
-                               4.518301e-02, 4.334208e-02, 2.020901e-01,
-                               5.257105e-01])
-    uo2_xsdata.set_chi([5.8791e-01, 4.1176e-01, 3.3906e-04, 1.1761e-07, 0.0000e+00,
+    uo2_scatter_matrix = np.rollaxis(uo2_scatter_matrix, 0, 3)
+    uo2_fission = np.array([7.21206e-03, 8.19301e-04, 6.45320e-03, 1.85648e-02, 1.78084e-02,
+                            8.30348e-02, 2.16004e-01])
+    uo2_nu_fission = np.array([2.005998e-02, 2.027303e-03, 1.570599e-02, 4.518301e-02,
+                               4.334208e-02, 2.020901e-01, 5.257105e-01])
+    uo2_chi = np.array([5.8791e-01, 4.1176e-01, 3.3906e-04, 1.1761e-07, 0.0000e+00,
                         0.0000e+00, 0.0000e+00])
 
-    h2o_xsdata = openmc.XSdata('LWTR', groups)
-    h2o_xsdata.order = 0
-    h2o_xsdata.set_total([0.15920605, 0.412969593, 0.59030986, 0.58435,
-                          0.718, 1.2544497, 2.650379])
-    h2o_xsdata.set_absorption([6.0105e-04, 1.5793e-05, 3.3716e-04,
-                               1.9406e-03, 5.7416e-03, 1.5001e-02,
-                               3.7239e-02])
-    scatter_matrix = np.array(
+    # MGXS for the H2O moderator.
+    h2o_total = np.array([0.15920605, 0.412969593, 0.59030986, 0.58435, 0.718, 1.2544497,
+                          2.650379])
+    h2o_abs = np.array([6.0105e-04, 1.5793e-05, 3.3716e-04, 1.9406e-03, 5.7416e-03,
+                        1.5001e-02, 3.7239e-02])
+    h2o_scatter_matrix = np.array(
         [[[0.0444777, 0.1134000, 0.0007235, 0.0000037, 0.0000001, 0.0000000, 0.0000000],
-          [0.0000000, 0.2823340, 0.1299400, 0.0006234,
-              0.0000480, 0.0000074, 0.0000010],
-          [0.0000000, 0.0000000, 0.3452560, 0.2245700,
-              0.0169990, 0.0026443, 0.0005034],
-          [0.0000000, 0.0000000, 0.0000000, 0.0910284,
-              0.4155100, 0.0637320, 0.0121390],
-          [0.0000000, 0.0000000, 0.0000000, 0.0000714,
-              0.1391380, 0.5118200, 0.0612290],
-          [0.0000000, 0.0000000, 0.0000000, 0.0000000,
-              0.0022157, 0.6999130, 0.5373200],
+          [0.0000000, 0.2823340, 0.1299400, 0.0006234, 0.0000480, 0.0000074, 0.0000010],
+          [0.0000000, 0.0000000, 0.3452560, 0.2245700, 0.0169990, 0.0026443, 0.0005034],
+          [0.0000000, 0.0000000, 0.0000000, 0.0910284, 0.4155100, 0.0637320, 0.0121390],
+          [0.0000000, 0.0000000, 0.0000000, 0.0000714, 0.1391380, 0.5118200, 0.0612290],
+          [0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0022157, 0.6999130, 0.5373200],
           [0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.1324400, 2.4807000]]])
-    scatter_matrix = np.rollaxis(scatter_matrix, 0, 3)
-    h2o_xsdata.set_scatter_matrix(scatter_matrix)
+    h2o_scatter_matrix = np.rollaxis(h2o_scatter_matrix, 0, 3)
+
+    # Instantiate the 7-group (C5G7) cross section data
+    uo2_xsdata = openmc.XSdata('UO2', groups)
+    uo2_xsdata.order = 0
+    uo2_xsdata.set_total(uo2_total, temperature=294.0)
+    uo2_xsdata.set_absorption(uo2_abs, temperature=294.0)
+    uo2_xsdata.set_scatter_matrix(uo2_scatter_matrix, temperature=294.0)
+    uo2_xsdata.set_fission(uo2_fission, temperature=294.0)
+    uo2_xsdata.set_nu_fission(uo2_nu_fission, temperature=294.0)
+    uo2_xsdata.set_chi(uo2_chi, temperature=294.0)
+
+    h2o_xsdata = openmc.XSdata('LWTR', groups)
+    h2o_xsdata.order = 0
+    h2o_xsdata.set_total(h2o_total, temperature=294.0)
+    h2o_xsdata.set_absorption(h2o_abs, temperature=294.0)
+    h2o_xsdata.set_scatter_matrix(h2o_scatter_matrix, temperature=294.0)
+
+    # Add the second temperature data point if requested.
+    if second_temp:
+        uo2_xsdata.add_temperature(394.0)
+        uo2_xsdata.set_total(0.5 * uo2_total, temperature=394.0)
+        uo2_xsdata.set_absorption(0.5 * uo2_abs, temperature=394.0)
+        uo2_xsdata.set_scatter_matrix(0.5 * uo2_scatter_matrix, temperature=394.0)
+        uo2_xsdata.set_fission(0.5 * uo2_fission, temperature=394.0)
+        uo2_xsdata.set_nu_fission(0.5 * uo2_nu_fission, temperature=394.0)
+        uo2_xsdata.set_chi(uo2_chi, temperature=394.0)
+
+        h2o_xsdata.add_temperature(394.0)
+        h2o_xsdata.set_total(0.5 * h2o_total, temperature=394.0)
+        h2o_xsdata.set_absorption(0.5 * h2o_abs, temperature=394.0)
+        h2o_xsdata.set_scatter_matrix(0.5 * h2o_scatter_matrix, temperature=394.0)
 
     mg_cross_sections = openmc.MGXSLibrary(groups)
     mg_cross_sections.add_xsdatas([uo2_xsdata, h2o_xsdata])