diff --git a/docs/source/usersguide/geometry.rst b/docs/source/usersguide/geometry.rst index 3224d5fff51..8f83b9b084e 100644 --- a/docs/source/usersguide/geometry.rst +++ b/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` or :ref:`lattice +`. 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: --------- diff --git a/docs/source/usersguide/random_ray.rst b/docs/source/usersguide/random_ray.rst index 64ef7de6d54..0c9a0402814 100644 --- a/docs/source/usersguide/random_ray.rst +++ b/docs/source/usersguide/random_ray.rst @@ -646,7 +646,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 @@ -733,6 +735,20 @@ 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. The valid keys and values in the +``temperature_settings`` dictionary are identical to +:attr:`openmc.Settings.temperature_settings`; more information can be found in +:class:`openmc.Settings` . 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 diff --git a/include/openmc/mgxs.h b/include/openmc/mgxs.h index 9b1602f299a..43505f432d6 100644 --- a/include/openmc/mgxs.h +++ b/include/openmc/mgxs.h @@ -113,6 +113,11 @@ class Mgxs { const vector& micros, const vector& 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 diff --git a/include/openmc/random_ray/flat_source_domain.h b/include/openmc/random_ray/flat_source_domain.h index 0d4086966dd..693093683a9 100644 --- a/include/openmc/random_ray/flat_source_domain.h +++ b/include/openmc/random_ray/flat_source_domain.h @@ -100,17 +100,18 @@ class FlatSourceDomain { // in model::cells vector 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 sigma_t_; vector nu_sigma_f_; vector sigma_f_; vector chi_; vector kappa_fission_; - // 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 sigma_s_; // The abstract container holding all source region-specific data diff --git a/include/openmc/random_ray/random_ray.h b/include/openmc/random_ray/random_ray.h index 40c67ef9549..7fd68b83e05 100644 --- a/include/openmc/random_ray/random_ray.h +++ b/include/openmc/random_ray/random_ray.h @@ -65,6 +65,7 @@ class RandomRay : public Particle { vector 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}; diff --git a/include/openmc/random_ray/source_region.h b/include/openmc/random_ray/source_region.h index c20d46abe9c..65f2e6bc41c 100644 --- a/include/openmc/random_ray/source_region.h +++ b/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 material_; + vector temperature_idx_; vector density_mult_; vector is_small_; vector n_hits_; diff --git a/openmc/examples.py b/openmc/examples.py index 350a4d24d5e..8eaeb29fe24 100644 --- a/openmc/examples.py +++ b/openmc/examples.py @@ -656,10 +656,19 @@ def slab_mg(num_regions=1, mat_names=None, mgxslib_name='2g.h5') -> openmc.Model return model -def _generate_c5g7_materials() -> openmc.Materials: +def _generate_c5g7_materials(second_temp = False) -> openmc.Materials: """Generate materials utilizing multi-group cross sections based on the the C5G7 Benchmark. + 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 ------- materials : openmc.Materials @@ -672,9 +681,45 @@ def _generate_c5g7_materials() -> openmc.Materials: assembly transport calculations without spatial homogenization" """ # 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) + 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.0000000, 0.0000000, 0.0000000, 0.0000000, 0.0085458, 0.2730800]]]) + 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]) + + # 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.0000000, 0.0000000, 0.0000000, 0.0000000, 0.1324400, 2.4807000]]]) + 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 @@ -707,29 +752,33 @@ def _generate_c5g7_materials() -> openmc.Materials: uo2_xsdata.set_nu_fission(nu_fission) uo2_xsdata.set_chi([5.8791e-01, 4.1176e-01, 3.3906e-04, 1.1761e-07, 0.0000e+00, 0.0000e+00, 0.0000e+00]) + 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([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( - [[[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.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_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]) @@ -823,10 +872,19 @@ def _generate_subdivided_pin_cell(uo2, water) -> openmc.Universe: return pincell -def random_ray_pin_cell() -> openmc.Model: +def random_ray_pin_cell(second_temp = False) -> openmc.Model: """Create a PWR pin cell example using C5G7 cross section data. 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 3934 K. This temperature dependence is + ficticious; it is used for testing temperature feedback in the random ray solver. + Returns ------- model : openmc.Model @@ -837,7 +895,7 @@ def random_ray_pin_cell() -> openmc.Model: ########################################################################### # Create Materials for the problem - materials = _generate_c5g7_materials() + materials = _generate_c5g7_materials(second_temp) uo2 = materials[0] water = materials[1] @@ -897,13 +955,22 @@ def random_ray_pin_cell() -> openmc.Model: return model -def random_ray_lattice() -> openmc.Model: +def random_ray_lattice(second_temp = False) -> openmc.Model: """Create a 2x2 PWR pin cell 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 C5G7 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 3934 K. This temperature dependence is + ficticious; it is used for testing temperature feedback in the random ray solver. + Returns ------- model : openmc.Model @@ -914,7 +981,7 @@ def random_ray_lattice() -> openmc.Model: ########################################################################### # Create Materials for the problem - materials = _generate_c5g7_materials() + materials = _generate_c5g7_materials(second_temp) uo2 = materials[0] water = materials[1] diff --git a/openmc/mgxs/library.py b/openmc/mgxs/library.py index b476de90202..78ab3ca1989 100644 --- a/openmc/mgxs/library.py +++ b/openmc/mgxs/library.py @@ -13,6 +13,7 @@ from openmc.checkvalue import PathLike from ..tallies import ESTIMATOR_TYPES +ROOM_TEMPERATURE_KELVIN = 294.0 class Library: """A multi-energy-group and multi-delayed-group cross section library for @@ -954,7 +955,7 @@ def load_from_file(filename='mgxs', directory='mgxs'): return pickle.load(f) def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', - subdomain=None, apply_domain_chi=False): + subdomain=None, apply_domain_chi=False, temperature=ROOM_TEMPERATURE_KELVIN): """Generates an openmc.XSdata object describing a multi-group cross section dataset for writing to an openmc.MGXSLibrary object. @@ -990,6 +991,9 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', downstream multigroup solvers that precompute a material-specific chi before the transport solve provides group-wise fluxes. Defaults to False. + temperature : float, optional + The temperature to set in the XSdata object. Defaults to 294 K + (room temperature). Returns ------- @@ -1036,6 +1040,7 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', else: representation = 'isotropic' xsdata = openmc.XSdata(name, self.energy_groups, + temperatures=[temperature], representation=representation) xsdata.num_delayed_groups = self.num_delayed_groups if self.num_polar > 1 or self.num_azimuthal > 1: @@ -1053,45 +1058,61 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', # Now get xs data itself if 'nu-transport' in self.mgxs_types and self.correction == 'P0': mymgxs = self.get_mgxs(domain, 'nu-transport') - xsdata.set_total_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide], + xsdata.set_total_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuclide], subdomain=subdomain) elif 'transport' in self.mgxs_types and self.correction == 'P0': mymgxs = self.get_mgxs(domain, 'transport') - xsdata.set_total_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide], + xsdata.set_total_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuclide], subdomain=subdomain) elif 'total' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'total') - xsdata.set_total_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide], + xsdata.set_total_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuclide], subdomain=subdomain) if 'absorption' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'absorption') - xsdata.set_absorption_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_absorption_mgxs(mymgxs, + temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'fission' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'fission') - xsdata.set_fission_mgxs(mymgxs, xs_type=xs_type, - nuclide=[nuclide], subdomain=subdomain) + xsdata.set_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuclide], + subdomain=subdomain) if 'kappa-fission' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'kappa-fission') - xsdata.set_kappa_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_kappa_fission_mgxs(mymgxs, + temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'inverse-velocity' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'inverse-velocity') - xsdata.set_inverse_velocity_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_inverse_velocity_mgxs(mymgxs, + temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'nu-fission matrix' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'nu-fission matrix') - xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_nu_fission_mgxs(mymgxs, + temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) @@ -1101,7 +1122,9 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuc = "sum" else: nuc = nuclide - xsdata.set_chi_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuc], + xsdata.set_chi_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuc], subdomain=subdomain) if 'chi-prompt' in self.mgxs_types: @@ -1110,8 +1133,10 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuc = "sum" else: nuc = nuclide - xsdata.set_chi_prompt_mgxs(mymgxs, xs_type=xs_type, - nuclide=[nuc], subdomain=subdomain) + xsdata.set_chi_prompt_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuc], + subdomain=subdomain) if 'chi-delayed' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'chi-delayed') @@ -1119,53 +1144,61 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuc = "sum" else: nuc = nuclide - xsdata.set_chi_delayed_mgxs(mymgxs, xs_type=xs_type, - nuclide=[nuc], subdomain=subdomain) + xsdata.set_chi_delayed_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, + nuclide=[nuc], + subdomain=subdomain) if 'nu-fission' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'nu-fission') - xsdata.set_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_nu_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'prompt-nu-fission' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'prompt-nu-fission') - xsdata.set_prompt_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_prompt_nu_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'prompt-nu-fission matrix' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'prompt-nu-fission matrix') - xsdata.set_prompt_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_prompt_nu_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'delayed-nu-fission' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'delayed-nu-fission') - xsdata.set_delayed_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_delayed_nu_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'delayed-nu-fission matrix' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'delayed-nu-fission matrix') - xsdata.set_delayed_nu_fission_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_delayed_nu_fission_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) if 'beta' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'beta') - xsdata.set_beta_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide], - subdomain=subdomain) + xsdata.set_beta_mgxs(mymgxs, temperature=temperature, xs_type=xs_type, + nuclide=[nuclide], subdomain=subdomain) if 'decay-rate' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'decay-rate') - xsdata.set_decay_rate_mgxs(mymgxs, xs_type=xs_type, nuclide=[nuclide], - subdomain=subdomain) + xsdata.set_decay_rate_mgxs(mymgxs, temperature=temperature, xs_type=xs_type, + nuclide=[nuclide], subdomain=subdomain) # If multiplicity matrix is available, prefer that if 'multiplicity matrix' in self.mgxs_types: mymgxs = self.get_mgxs(domain, 'multiplicity matrix') - xsdata.set_multiplicity_matrix_mgxs(mymgxs, xs_type=xs_type, + xsdata.set_multiplicity_matrix_mgxs(mymgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) using_multiplicity = True @@ -1176,6 +1209,7 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', scatt_mgxs = self.get_mgxs(domain, 'scatter matrix') nuscatt_mgxs = self.get_mgxs(domain, 'nu-scatter matrix') xsdata.set_multiplicity_matrix_mgxs(nuscatt_mgxs, scatt_mgxs, + temperature=temperature, xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) @@ -1188,6 +1222,7 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', nuscatt_mgxs = \ self.get_mgxs(domain, 'consistent nu-scatter matrix') xsdata.set_multiplicity_matrix_mgxs(nuscatt_mgxs, scatt_mgxs, + temperature=temperature, xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) @@ -1202,7 +1237,8 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', else: nuscatt_mgxs = \ self.get_mgxs(domain, 'consistent nu-scatter matrix') - xsdata.set_scatter_matrix_mgxs(nuscatt_mgxs, xs_type=xs_type, + xsdata.set_scatter_matrix_mgxs(nuscatt_mgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) else: @@ -1213,7 +1249,8 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', else: nuscatt_mgxs = \ self.get_mgxs(domain, 'consistent nu-scatter matrix') - xsdata.set_scatter_matrix_mgxs(nuscatt_mgxs, xs_type=xs_type, + xsdata.set_scatter_matrix_mgxs(nuscatt_mgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) @@ -1253,14 +1290,15 @@ def get_xsdata(self, domain, xsdata_name, nuclide='total', xs_type='macro', 'are ignored since multiplicity or nu-scatter matrices '\ 'were not tallied for ' + xsdata_name warn(msg, RuntimeWarning) - xsdata.set_scatter_matrix_mgxs(scatt_mgxs, xs_type=xs_type, + xsdata.set_scatter_matrix_mgxs(scatt_mgxs, temperature=temperature, + xs_type=xs_type, nuclide=[nuclide], subdomain=subdomain) return xsdata def create_mg_library(self, xs_type='macro', xsdata_names=None, - apply_domain_chi=False): + apply_domain_chi=False, temperature=ROOM_TEMPERATURE_KELVIN): """Creates an openmc.MGXSLibrary object to contain the MGXS data for the Multi-Group mode of OpenMC. @@ -1286,6 +1324,9 @@ def create_mg_library(self, xs_type='macro', xsdata_names=None, downstream multigroup solvers that precompute a material-specific chi before the transport solve provides group-wise fluxes. Defaults to False. + temperature : float, optional + The temperature to set in the MGXSLibrary object. Defaults to 294 K + (room temperature). Returns ------- diff --git a/openmc/mgxs_library.py b/openmc/mgxs_library.py index c1a15998e19..ab9b58b7a3a 100644 --- a/openmc/mgxs_library.py +++ b/openmc/mgxs_library.py @@ -515,6 +515,75 @@ def _set_fissionable(self, array): if np.sum(array) > 0: self._fissionable = True + def add_temperature_data(self, other): + """This method adds temperature-dependent cross section + values from another XSdata object to this XSdata object. + Note: if a temperature datapoint from 'other' already exists in this + object, it will be overridden. + + Parameters + ---------- + other: openmc.XSdata + The other XSdata object to fetch data from + """ + + # Sanity check to make sure they have the same name, energy group structure, + # and delayed group structure + check_value('name', other.name, self.name) + check_value('energy_groups', other.energy_groups, [self.energy_groups]) + check_value('delayed_groups', other.num_delayed_groups, [self.num_delayed_groups]) + + # Add the temperature data. + for temp in other.temperatures: + if temp not in self.temperatures: + self.add_temperature(temp) + + if np.all(other.absorption[other._temperature_index(temp)] != None): + self.set_absorption(other.absorption[other._temperature_index(temp)], temp) + + if np.all(other.beta[other._temperature_index(temp)] != None): + self.set_beta(other.beta[other._temperature_index(temp)], temp) + + if np.all(other.chi[other._temperature_index(temp)] != None): + self.set_chi(other.chi[other._temperature_index(temp)], temp) + + if np.all(other.chi_delayed[other._temperature_index(temp)] != None): + self.set_chi_delayed(other.chi_delayed[other._temperature_index(temp)], temp) + + if np.all(other.chi_prompt[other._temperature_index(temp)] != None): + self.set_chi_prompt(other.chi_prompt[other._temperature_index(temp)], temp) + + if np.all(other.decay_rate[other._temperature_index(temp)] != None): + self.set_decay_rate(other.decay_rate[other._temperature_index(temp)], temp) + + if np.all(other.delayed_nu_fission[other._temperature_index(temp)] != None): + self.set_delayed_nu_fission(other.delayed_nu_fission[other._temperature_index(temp)], temp) + + if np.all(other.fission[other._temperature_index(temp)] != None): + self.set_fission(other.fission[other._temperature_index(temp)], temp) + + if np.all(other.inverse_velocity[other._temperature_index(temp)] != None): + self.set_inverse_velocity(other.inverse_velocity[other._temperature_index(temp)], temp) + + if np.all(other.kappa_fission[other._temperature_index(temp)] != None): + self.set_kappa_fission(other.kappa_fission[other._temperature_index(temp)], temp) + + if np.all(other.multiplicity_matrix[other._temperature_index(temp)] != None): + self.set_multiplicity_matrix(other.multiplicity_matrix[other._temperature_index(temp)], temp) + + if np.all(other.nu_fission[other._temperature_index(temp)] != None): + self.set_nu_fission(other.nu_fission[other._temperature_index(temp)], temp) + + if np.all(other.prompt_nu_fission[other._temperature_index(temp)] != None): + self.set_prompt_nu_fission(other.prompt_nu_fission[other._temperature_index(temp)], temp) + + if np.all(other.scatter_matrix[other._temperature_index(temp)] != None): + self.set_scatter_matrix(other.scatter_matrix[other._temperature_index(temp)], temp) + + if np.all(other.fission[other._temperature_index(temp)] != None): + self.set_total(other.total[other._temperature_index(temp)], temp) + + def set_total(self, total, temperature=ROOM_TEMPERATURE_KELVIN): """This method sets the cross section for this XSdata object at the provided temperature. diff --git a/openmc/model/model.py b/openmc/model/model.py index 90b43f9f74f..ecdddb2c9c6 100644 --- a/openmc/model/model.py +++ b/openmc/model/model.py @@ -1685,8 +1685,8 @@ def differentiate_mats(self, diff_volume_method: str = None, depletable_only: bo self.geometry.get_all_materials().values() ) + @staticmethod def _auto_generate_mgxs_lib( - self, model: openmc.model.model, groups: openmc.mgxs.EnergyGroups, correction: str | none, @@ -1851,6 +1851,97 @@ def _create_mgxs_sources( return sources + @staticmethod + def _isothermal_infinite_media_mgxs( + material: openmc.Material, + groups: openmc.mgxs.EnergyGroups, + nparticles: int, + correction: str | None, + directory: PathLike, + source: openmc.IndependentSource, + temperature_settings: dict, + temperature: float | None = None, + ) -> openmc.XSdata: + """Generate a single MGXS set for one material, where the geometry is an + infinite medium composed of that material at an isothermal temperature value. + + Parameters + ---------- + material : openmc.Material + The material to generate MGXS for + groups : openmc.mgxs.EnergyGroups + Energy group structure for the MGXS. + nparticles : int + Number of particles to simulate per batch when generating MGXS. + correction : str + Transport correction to apply to the MGXS. Options are None and + "P0". + directory : str + Directory to run the simulation in, so as to contain XML files. + source : openmc.IndependentSource + Source to use when generating MGXS. + temperature_settings : dict + A dictionary of temperature settings to use when generating MGXS. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. + temperature : float, optional + The isothermal temperature value to apply to the material. If not specified, + defaults to the temperature in the material. + + Returns + ------- + data : openmc.XSdata + The material MGXS for the given temperature isotherm. + """ + model = openmc.Model() + + # Set materials on the model + model.materials = [material] + if temperature != None: + model.materials[-1].temperature = temperature + + # Settings + model.settings.batches = 100 + model.settings.particles = nparticles + + model.settings.source = source + + model.settings.run_mode = 'fixed source' + model.settings.create_fission_neutrons = False + + model.settings.output = {'summary': True, 'tallies': False} + model.settings.temperature = temperature_settings + + # Geometry + box = openmc.model.RectangularPrism( + 100000.0, 100000.0, boundary_type='reflective') + name = material.name + infinite_cell = openmc.Cell(name=name, fill=model.materials[-1], region=-box) + infinite_universe = openmc.Universe(name=name, cells=[infinite_cell]) + model.geometry.root_universe = infinite_universe + + # Generate MGXS + mgxs_lib = Model._auto_generate_mgxs_lib( + model, groups, correction, directory) + + if temperature != None: + return mgxs_lib.get_xsdata(domain=material, xsdata_name=name, + temperature=temperature) + else: + return mgxs_lib.get_xsdata(domain=material, xsdata_name=name) + def _generate_infinite_medium_mgxs( self, groups: openmc.mgxs.EnergyGroups, @@ -1859,13 +1950,17 @@ def _generate_infinite_medium_mgxs( correction: str | None, directory: PathLike, source_energy: openmc.stats.Univariate | None = None, - ): + temperatures: Sequence[float] | None = None, + temperature_settings: dict | None = None, + ) -> None: """Generate a MGXS library by running multiple OpenMC simulations, each representing an infinite medium simulation of a single isolated material. A discrete source is used to sample particles, with an equal strength spread across each of the energy groups. This is a highly naive method that ignores all spatial self shielding effects and all resonance - shielding effects between materials. + shielding effects between materials. If temperature data points are provided, + isothermal cross sections are generated at each temperature point for + each material to build a temperature interpolation table. Note that in all cases, a discrete source that is uniform over all energy groups is created (strength = 0.01) to ensure that total cross @@ -1897,50 +1992,91 @@ def _generate_infinite_medium_mgxs( source_energy : openmc.stats.Univariate, optional Energy distribution to use when generating MGXS data, replacing any existing sources in the model. + temperatures : Sequence[float], optional + A list of temperatures to generate MGXS at. Each infinite material region + is isothermal at a given temperature data point for cross + section generation. + temperature_settings : dict, optional + A dictionary of temperature settings to use when generating MGXS. If not + provided, the settings stored in the model will be used. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. """ - mgxs_sets = [] - for material in self.materials: - model = openmc.Model() - # Set materials on the model - model.materials = [material] - - # Settings - model.settings.batches = 100 - model.settings.particles = nparticles - - model.settings.source = self._create_mgxs_sources( - groups, - spatial_dist=openmc.stats.Point(), - source_energy=source_energy - ) - - model.settings.run_mode = 'fixed source' - model.settings.create_fission_neutrons = False - - model.settings.output = {'summary': True, 'tallies': False} - - # Geometry - box = openmc.model.RectangularPrism( - 100000.0, 100000.0, boundary_type='reflective') - name = material.name - infinite_cell = openmc.Cell(name=name, fill=material, region=-box) - infinite_universe = openmc.Universe(name=name, cells=[infinite_cell]) - model.geometry.root_universe = infinite_universe + src = self._create_mgxs_sources( + groups, + spatial_dist=openmc.stats.Point(), + source_energy=source_energy + ) - # Add MGXS Tallies - mgxs_lib = self._auto_generate_mgxs_lib( - model, groups, correction, directory) + temp_settings = {} + if temperature_settings == None: + temp_settings = self.settings.temperature + else: + temp_settings = temperature_settings - # Create a MGXS File which can then be written to disk - mgxs_set = mgxs_lib.get_xsdata(domain=material, xsdata_name=name) - mgxs_sets.append(mgxs_set) + if temperatures == None: + mgxs_sets = [] + for material in self.materials: + xs_data = Model._isothermal_infinite_media_mgxs( + material, + groups, + nparticles, + correction, + directory, + src, + temp_settings + ) + mgxs_sets.append(xs_data) - # Write the file to disk - mgxs_file = openmc.MGXSLibrary(energy_groups=groups) - for mgxs_set in mgxs_sets: - mgxs_file.add_xsdata(mgxs_set) - mgxs_file.export_to_hdf5(mgxs_path) + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) + else: + # Build a series of XSData objects, one for each isothermal temperature value. + raw_mgxs_sets = {} + for temperature in temperatures: + raw_mgxs_sets[temperature] = [] + for material in self.materials: + xs_data = Model._isothermal_infinite_media_mgxs( + material, + groups, + nparticles, + correction, + directory, + src, + temp_settings, + temperature + ) + raw_mgxs_sets[temperature].append(xs_data) + + # Unpack the isothermal XSData objects and build a single XSData object per material. + mgxs_sets = [] + for m in range(len(self.materials)): + mgxs_sets.append(openmc.XSdata(self.materials[m].name, groups)) + mgxs_sets[-1].order = 0 + for temperature in temperatures: + mgxs_sets[-1].add_temperature_data(raw_mgxs_sets[temperature][m]) + + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) @staticmethod def _create_stochastic_slab_geometry( @@ -2016,6 +2152,101 @@ def _create_stochastic_slab_geometry( return geometry, box + @staticmethod + def _isothermal_stochastic_slab_mgxs( + stoch_geom: openmc.Geometry, + groups: openmc.mgxs.EnergyGroups, + nparticles: int, + correction: str | None, + directory: PathLike, + source: openmc.IndependentSource, + temperature_settings: dict, + temperature: float | None = None, + ) -> dict[str, openmc.XSdata]: + """Generate MGXS assuming a stochastic "sandwich" of materials in a layered + slab geometry. If a temperature is specified, all materials in the slab have + their temperatures set to be isothermal at this temperature. + + Parameters + ---------- + stoch_geom : openmc.Geometry + The stochastic slab geometry. + groups : openmc.mgxs.EnergyGroups + Energy group structure for the MGXS. + nparticles : int + Number of particles to simulate per batch when generating MGXS. + correction : str + Transport correction to apply to the MGXS. Options are None and + "P0". + directory : str + Directory to run the simulation in, so as to contain XML files. + source : openmc.IndependentSource + Source to use when generating MGXS. + temperature_settings : dict + A dictionary of temperature settings to use when generating MGXS. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. + temperature : float, optional + The isothermal temperature value to apply to the materials in the + slab. If not specified, defaults to the temperature in the materials. + + Returns + ------- + data : dict[str, openmc.XSdata] + A dictionary where the key is the name of the material and the value is the isothermal MGXS. + """ + + model = openmc.Model() + model.geometry = stoch_geom + + if temperature != None: + for material in model.geometry.get_all_materials().values(): + material.temperature = temperature + + # Settings + model.settings.batches = 200 + model.settings.inactive = 100 + model.settings.particles = nparticles + model.settings.output = {'summary': True, 'tallies': False} + model.settings.temperature = temperature_settings + + # Define the sources + model.settings.source = source + + model.settings.run_mode = 'fixed source' + model.settings.create_fission_neutrons = False + + model.settings.output = {'summary': True, 'tallies': False} + + # Generate MGXS + mgxs_lib = Model._auto_generate_mgxs_lib( + model, groups, correction, directory) + + # Fetch all of the isothermal results. + if temperature != None: + return { + mat.name : mgxs_lib.get_xsdata(domain=mat, xsdata_name=mat.name, + temperature=temperature) + for mat in mgxs_lib.domains + } + else: + return { + mat.name : mgxs_lib.get_xsdata(domain=mat, xsdata_name=mat.name) + for mat in mgxs_lib.domains + } + def _generate_stochastic_slab_mgxs( self, groups: openmc.mgxs.EnergyGroups, @@ -2024,6 +2255,8 @@ def _generate_stochastic_slab_mgxs( correction: str | None, directory: PathLike, source_energy: openmc.stats.Univariate | None = None, + temperatures: Sequence[float] | None = None, + temperature_settings: dict | None = None, ) -> None: """Generate MGXS assuming a stochastic "sandwich" of materials in a layered slab geometry. While geometry-specific spatial shielding effects are not @@ -2033,7 +2266,9 @@ def _generate_stochastic_slab_mgxs( will generate cross sections for all materials in the problem regardless of type. If this is a fixed source problem, a discrete source is used to sample particles, with an equal strength spread across each of the - energy groups. + energy groups. If temperature data points are provided, + isothermal cross sections are generated at each temperature point for + the stochastic slab to build a temperature interpolation table. Parameters ---------- @@ -2065,41 +2300,177 @@ def _generate_stochastic_slab_mgxs( no sources are defined on the model and the run mode is 'eigenvalue', then a default Watt spectrum source (strength = 0.99) is added. + temperatures : Sequence[float], optional + A list of temperatures to generate MGXS at. Each infinite material region + is isothermal at a given temperature data point for cross + section generation. + temperature_settings : dict, optional + A dictionary of temperature settings to use when generating MGXS. If not + provided, the settings stored in the model will be used. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. """ - model = openmc.Model() - model.materials = self.materials - - # Settings - model.settings.batches = 200 - model.settings.inactive = 100 - model.settings.particles = nparticles - model.settings.output = {'summary': True, 'tallies': False} # Stochastic slab geometry - model.geometry, spatial_distribution = Model._create_stochastic_slab_geometry( - model.materials) + geo, spatial_distribution = Model._create_stochastic_slab_geometry( + self.materials) - # Define the sources - model.settings.source = self._create_mgxs_sources( + src = self._create_mgxs_sources( groups, spatial_dist=spatial_distribution, source_energy=source_energy ) - model.settings.run_mode = 'fixed source' - model.settings.create_fission_neutrons = False + temp_settings = {} + if temperature_settings == None: + temp_settings = self.settings.temperature + else: + temp_settings = temperature_settings + if temperatures == None: + mgxs_sets = Model._isothermal_stochastic_slab_mgxs( + geo, + groups, + nparticles, + correction, + directory, + src, + temp_settings + ).values() + + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) + else: + # Build a series of XSData objects, one for each isothermal temperature value. + raw_mgxs_sets = {} + for temperature in temperatures: + raw_mgxs_sets[temperature] = Model._isothermal_stochastic_slab_mgxs( + geo, + groups, + nparticles, + correction, + directory, + src, + temp_settings, + temperature + ) + + # Unpack the isothermal XSData objects and build a single XSData object per material. + mgxs_sets = [] + for mat in self.materials: + mgxs_sets.append(openmc.XSdata(mat.name, groups)) + mgxs_sets[-1].order = 0 + for temperature in temperatures: + mgxs_sets[-1].add_temperature_data(raw_mgxs_sets[temperature][mat.name]) + + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) + + @staticmethod + def _isothermal_materialwise_mgxs( + input_model: openmc.Model, + groups: openmc.mgxs.EnergyGroups, + nparticles: int, + correction: str | None, + directory: PathLike, + temperature_settings: dict, + temperature: float | None = None, + ) -> dict[str, openmc.XSdata]: + """Generate a material-wise MGXS library for the model by running the + original continuous energy OpenMC simulation. If a temperature is + specified, each material in the input model is set to that temperature. + Otherwise, the original material temperatures are used. If temperature + data points are provided, isothermal cross sections are generated at + each temperature point for the whole model to build a temperature + interpolation table. + + Parameters + ---------- + input_model : openmc.Model + The model to use when computing material-wise MGXS. + groups : openmc.mgxs.EnergyGroups + Energy group structure for the MGXS. + nparticles : int + Number of particles to simulate per batch when generating MGXS. + correction : str + Transport correction to apply to the MGXS. Options are None and + "P0". + directory : str + Directory to run the simulation in, so as to contain XML files. + temperature_settings : dict + A dictionary of temperature settings to use when generating MGXS. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. + temperature : float, optional + The isothermal temperature value to apply to the materials in the + input model. If not specified, defaults to the temperatures in the + materials. + + Returns + ------- + data : dict[str, openmc.XSdata] + A dictionary where the key is the name of the material and the value is the isothermal MGXS. + """ + model = copy.deepcopy(input_model) + model.tallies = openmc.Tallies() + + if temperature != None: + for material in model.geometry.get_all_materials().values(): + material.temperature = temperature + + # Settings + model.settings.batches = 200 + model.settings.inactive = 100 + model.settings.particles = nparticles model.settings.output = {'summary': True, 'tallies': False} + model.settings.temperature = temperature_settings - # Add MGXS Tallies - mgxs_lib = self._auto_generate_mgxs_lib( + # Generate MGXS + mgxs_lib = Model._auto_generate_mgxs_lib( model, groups, correction, directory) - names = [mat.name for mat in mgxs_lib.domains] - - # Create a MGXS File which can then be written to disk - mgxs_file = mgxs_lib.create_mg_library(xs_type='macro', xsdata_names=names) - mgxs_file.export_to_hdf5(mgxs_path) + # Fetch all of the isothermal results. + if temperature != None: + return { + mat.name : mgxs_lib.get_xsdata(domain=mat, xsdata_name=mat.name, + temperature=temperature) + for mat in mgxs_lib.domains + } + else: + return { + mat.name : mgxs_lib.get_xsdata(domain=mat, xsdata_name=mat.name) + for mat in mgxs_lib.domains + } def _generate_material_wise_mgxs( self, @@ -2108,6 +2479,8 @@ def _generate_material_wise_mgxs( mgxs_path: PathLike, correction: str | None, directory: PathLike, + temperatures: Sequence[float] | None = None, + temperature_settings: dict | None = None, ) -> None: """Generate a material-wise MGXS library for the model by running the original continuous energy OpenMC simulation of the full material @@ -2132,26 +2505,76 @@ def _generate_material_wise_mgxs( "P0". directory : PathLike Directory to run the simulation in, so as to contain XML files. + temperatures : Sequence[float], optional + A list of temperatures to generate MGXS at. Each infinite material region + is isothermal at a given temperature data point for cross + section generation. + temperature_settings : dict, optional + A dictionary of temperature settings to use when generating MGXS. If not + provided, the settings stored in the model will be used. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. """ - model = copy.deepcopy(self) - model.tallies = openmc.Tallies() - - # Settings - model.settings.batches = 200 - model.settings.inactive = 100 - model.settings.particles = nparticles - model.settings.output = {'summary': True, 'tallies': False} + temp_settings = {} + if temperature_settings == None: + temp_settings = self.settings.temperature + else: + temp_settings = temperature_settings - # Add MGXS Tallies - mgxs_lib = self._auto_generate_mgxs_lib( - model, groups, correction, directory) + if temperatures == None: + mgxs_sets = Model._isothermal_materialwise_mgxs( + self, + groups, + nparticles, + correction, + directory, + temp_settings + ).values() + + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) + else: + # Build a series of XSData objects, one for each isothermal temperature value. + raw_mgxs_sets = {} + for temperature in temperatures: + raw_mgxs_sets[temperature] = Model._isothermal_materialwise_mgxs( + self, + groups, + nparticles, + correction, + directory, + temp_settings, + temperature + ) - names = [mat.name for mat in mgxs_lib.domains] + # Unpack the isothermal XSData objects and build a single XSData object per material. + mgxs_sets = [] + for mat in self.materials: + mgxs_sets.append(openmc.XSdata(mat.name, groups)) + mgxs_sets[-1].order = 0 + for temperature in temperatures: + mgxs_sets[-1].add_temperature_data(raw_mgxs_sets[temperature][mat.name]) - # Create a MGXS File which can then be written to disk - mgxs_file = mgxs_lib.create_mg_library( - xs_type='macro', xsdata_names=names) - mgxs_file.export_to_hdf5(mgxs_path) + # Write the file to disk. + mgxs_file = openmc.MGXSLibrary(energy_groups=groups) + for mgxs_set in mgxs_sets: + mgxs_file.add_xsdata(mgxs_set) + mgxs_file.export_to_hdf5(mgxs_path) def convert_to_multigroup( self, @@ -2162,6 +2585,8 @@ def convert_to_multigroup( mgxs_path: PathLike = "mgxs.h5", correction: str | None = None, source_energy: openmc.stats.Univariate | None = None, + temperatures: Sequence[float] | None = None, + temperature_settings: dict | None = None, ): """Convert all materials from continuous energy to multigroup. @@ -2202,6 +2627,27 @@ def convert_to_multigroup( 'eigenvalue', then a default Watt spectrum source (strength = 0.99) is added. Note that this argument is only used when using the "stochastic_slab" or "infinite_medium" MGXS generation methods. + temperatures : Sequence[float], optional + A list of temperatures to generate MGXS at. Each infinite material region + is isothermal at a given temperature data point for cross + section generation. + temperature_settings : dict, optional + A dictionary of temperature settings to use when generating MGXS. If not + provided, the settings stored in the model will be used. + Valid entries for temperature_settings are the same as the valid + entries in openmc.Settings.temperature_settings. Accepted keys are + 'default', 'method', 'range', 'tolerance', and 'multipole'. The value + for 'default' should be a float representing the default temperature + in Kelvin. The value for 'method' should be 'nearest' or 'interpolation'. + If the method is 'nearest', 'tolerance' indicates a range of temperature + within which cross sections may be used. If the method is 'interpolation', + 'tolerance' indicates the range of temperatures outside of the available + cross section temperatures where cross sections will evaluate to the nearer + bound. The value for 'range' should be a pair of minimum and maximum + temperatures which are used to indicate that cross sections be loaded + at all temperatures within the range. 'multipole' is a boolean indicating + whether or not the windowed multipole method should be used to evaluate + resolved resonance cross sections. """ if isinstance(groups, str): groups = openmc.mgxs.EnergyGroups(groups) @@ -2231,13 +2677,16 @@ def convert_to_multigroup( if not Path(mgxs_path).is_file() or overwrite_mgxs_library: if method == "infinite_medium": self._generate_infinite_medium_mgxs( - groups, nparticles, mgxs_path, correction, tmpdir, source_energy) + groups, nparticles, mgxs_path, correction, tmpdir, source_energy, + temperatures, temperature_settings) elif method == "material_wise": self._generate_material_wise_mgxs( - groups, nparticles, mgxs_path, correction, tmpdir) + groups, nparticles, mgxs_path, correction, tmpdir, + temperatures, temperature_settings) elif method == "stochastic_slab": self._generate_stochastic_slab_mgxs( - groups, nparticles, mgxs_path, correction, tmpdir, source_energy) + groups, nparticles, mgxs_path, correction, tmpdir, source_energy, + temperatures, temperature_settings) else: raise ValueError( f'MGXS generation method "{method}" not recognized') diff --git a/src/mgxs_interface.cpp b/src/mgxs_interface.cpp index ed734d401ea..d31f5a0b0bb 100644 --- a/src/mgxs_interface.cpp +++ b/src/mgxs_interface.cpp @@ -169,13 +169,13 @@ vector> MgxsInterface::get_mat_kTs() continue; // Get temperature of cell (rounding to nearest integer) - double sqrtkT = - cell->sqrtkT_.size() == 1 ? cell->sqrtkT_[j] : cell->sqrtkT_[0]; - double kT = sqrtkT * sqrtkT; + for (int k = 0; k < cell->sqrtkT_.size(); ++k) { + double kT = cell->sqrtkT_[k] * cell->sqrtkT_[k]; - // Add temperature if it hasn't already been added - if (!contains(kTs[i_material], kT)) { - kTs[i_material].push_back(kT); + // Add temperature if it hasn't already been added + if (!contains(kTs[i_material], kT)) { + kTs[i_material].push_back(kT); + } } } } diff --git a/src/random_ray/flat_source_domain.cpp b/src/random_ray/flat_source_domain.cpp index c2effaa5d4a..62986b570a6 100644 --- a/src/random_ray/flat_source_domain.cpp +++ b/src/random_ray/flat_source_domain.cpp @@ -109,22 +109,29 @@ void FlatSourceDomain::update_single_neutron_source(SourceRegionHandle& srh) // Add scattering + fission source int material = srh.material(); + int temp = srh.temperature_idx(); double density_mult = srh.density_mult(); if (material != MATERIAL_VOID) { double inverse_k_eff = 1.0 / k_eff_; for (int g_out = 0; g_out < negroups_; g_out++) { - double sigma_t = sigma_t_[material * negroups_ + g_out] * density_mult; + double sigma_t = + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g_out] * + density_mult; double scatter_source = 0.0; double fission_source = 0.0; for (int g_in = 0; g_in < negroups_; g_in++) { double scalar_flux = srh.scalar_flux_old(g_in); - double sigma_s = sigma_s_[material * negroups_ * negroups_ + - g_out * negroups_ + g_in] * - density_mult; + double sigma_s = + sigma_s_[((material * ntemperature_ + temp) * negroups_ + g_out) * + negroups_ + + g_in] * + density_mult; double nu_sigma_f = - nu_sigma_f_[material * negroups_ + g_in] * density_mult; - double chi = chi_[material * negroups_ + g_out]; + nu_sigma_f_[(material * ntemperature_ + temp) * negroups_ + g_in] * + density_mult; + double chi = + chi_[(material * ntemperature_ + temp) * negroups_ + g_out]; scatter_source += sigma_s * scalar_flux; if (settings::create_fission_neutrons) { @@ -193,6 +200,7 @@ void FlatSourceDomain::set_flux_to_flux_plus_source( int64_t sr, double volume, int g) { int material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); if (material == MATERIAL_VOID) { source_regions_.scalar_flux_new(sr, g) /= volume; if (settings::run_mode == RunMode::FIXED_SOURCE) { @@ -201,8 +209,9 @@ void FlatSourceDomain::set_flux_to_flux_plus_source( source_regions_.volume_sq(sr); } } else { - double sigma_t = sigma_t_[source_regions_.material(sr) * negroups_ + g] * - source_regions_.density_mult(sr); + double sigma_t = + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + source_regions_.density_mult(sr); source_regions_.scalar_flux_new(sr, g) /= (sigma_t * volume); source_regions_.scalar_flux_new(sr, g) += source_regions_.source(sr, g); } @@ -328,6 +337,7 @@ void FlatSourceDomain::compute_k_eff() } int material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); if (material == MATERIAL_VOID) { continue; } @@ -336,8 +346,9 @@ void FlatSourceDomain::compute_k_eff() double sr_fission_source_new = 0; for (int g = 0; g < negroups_; g++) { - double nu_sigma_f = nu_sigma_f_[material * negroups_ + g] * - source_regions_.density_mult(sr); + double nu_sigma_f = + nu_sigma_f_[(material * ntemperature_ + temp) * negroups_ + g] * + source_regions_.density_mult(sr); sr_fission_source_old += nu_sigma_f * source_regions_.scalar_flux_old(sr, g); sr_fission_source_new += @@ -560,6 +571,7 @@ double FlatSourceDomain::compute_fixed_source_normalization_factor() const #pragma omp parallel for reduction(+ : simulation_external_source_strength) for (int64_t sr = 0; sr < n_source_regions(); sr++) { int material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); double volume = source_regions_.volume(sr) * simulation_volume_; for (int g = 0; g < negroups_; g++) { // For non-void regions, we store the external source pre-divided by @@ -567,8 +579,8 @@ double FlatSourceDomain::compute_fixed_source_normalization_factor() const // to get the total source strength in the expected units. double sigma_t = 1.0; if (material != MATERIAL_VOID) { - sigma_t = - sigma_t_[material * negroups_ + g] * source_regions_.density_mult(sr); + sigma_t = sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + source_regions_.density_mult(sr); } simulation_external_source_strength += source_regions_.external_source(sr, g) * sigma_t * volume; @@ -624,9 +636,9 @@ void FlatSourceDomain::random_ray_tally() // source strength. double volume = source_regions_.volume(sr) * simulation_volume_; - int material = source_regions_.material(sr); + double material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); double density_mult = source_regions_.density_mult(sr); - for (int g = 0; g < negroups_; g++) { double flux = source_regions_.scalar_flux_new(sr, g) * source_normalization_factor; @@ -643,21 +655,27 @@ void FlatSourceDomain::random_ray_tally() case SCORE_TOTAL: if (material != MATERIAL_VOID) { score = - flux * volume * sigma_t_[material * negroups_ + g] * density_mult; + flux * volume * + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + density_mult; } break; case SCORE_FISSION: if (material != MATERIAL_VOID) { score = - flux * volume * sigma_f_[material * negroups_ + g] * density_mult; + flux * volume * + sigma_f_[(material * ntemperature_ + temp) * negroups_ + g] * + density_mult; } break; case SCORE_NU_FISSION: if (material != MATERIAL_VOID) { - score = flux * volume * nu_sigma_f_[material * negroups_ + g] * - density_mult; + score = + flux * volume * + nu_sigma_f_[(material * ntemperature_ + temp) * negroups_ + g] * + density_mult; } break; @@ -666,8 +684,10 @@ void FlatSourceDomain::random_ray_tally() break; case SCORE_KAPPA_FISSION: - score = flux * volume * kappa_fission_[material * negroups_ + g] * - density_mult; + score = + flux * volume * + kappa_fission_[(material * ntemperature_ + temp) * negroups_ + g] * + density_mult; break; default: @@ -925,12 +945,14 @@ void FlatSourceDomain::output_to_vtk() const float total_fission = 0.0; if (fsr >= 0) { int mat = source_regions_.material(fsr); + int temp = source_regions_.temperature_idx(fsr); if (mat != MATERIAL_VOID) { for (int g = 0; g < negroups_; g++) { int64_t source_element = fsr * negroups_ + g; float flux = evaluate_flux_at_point(voxel_positions[i], fsr, g); - double sigma_f = sigma_f_[mat * negroups_ + g] * - source_regions_.density_mult(fsr); + double sigma_f = + sigma_f_[(mat * ntemperature_ + temp) * negroups_ + g] * + source_regions_.density_mult(fsr); total_fission += sigma_f * flux; } } @@ -944,6 +966,7 @@ void FlatSourceDomain::output_to_vtk() const for (int i = 0; i < Nx * Ny * Nz; i++) { int64_t fsr = voxel_indices[i]; int mat = source_regions_.material(fsr); + int temp = source_regions_.temperature_idx(fsr); float total_external = 0.0f; if (fsr >= 0) { for (int g = 0; g < negroups_; g++) { @@ -951,7 +974,7 @@ void FlatSourceDomain::output_to_vtk() const // multiply it back to get the true external source. double sigma_t = 1.0; if (mat != MATERIAL_VOID) { - sigma_t = sigma_t_[mat * negroups_ + g] * + sigma_t = sigma_t_[(mat * ntemperature_ + temp) * negroups_ + g] * source_regions_.density_mult(fsr); } total_external += source_regions_.external_source(fsr, g) * sigma_t; @@ -1131,67 +1154,74 @@ void FlatSourceDomain::flatten_xs() { // Temperature and angle indices, if using multiple temperature // data sets and/or anisotropic data sets. - // TODO: Currently assumes we are only using single temp/single angle data. - const int t = 0; + // TODO: Currently assumes we are only using single angle data. const int a = 0; n_materials_ = data::mg.macro_xs_.size(); + ntemperature_ = 1; + for (int i = 0; i < n_materials_; i++) { + ntemperature_ = + std::max(ntemperature_, data::mg.macro_xs_[i].n_temperature_points()); + } + for (int i = 0; i < n_materials_; i++) { auto& m = data::mg.macro_xs_[i]; - for (int g_out = 0; g_out < negroups_; g_out++) { - if (m.exists_in_model) { - double sigma_t = - m.get_xs(MgxsType::TOTAL, g_out, NULL, NULL, NULL, t, a); - sigma_t_.push_back(sigma_t); - - if (sigma_t < MINIMUM_MACRO_XS) { - Material* mat = model::materials[i].get(); - warning(fmt::format( - "Material \"{}\" (id: {}) has a group {} total cross section " - "({:.3e}) below the minimum threshold " - "({:.3e}). Material will be treated as pure void.", - mat->name(), mat->id(), g_out, sigma_t, MINIMUM_MACRO_XS)); - } + for (int t = 0; t < ntemperature_; t++) { + for (int g_out = 0; g_out < negroups_; g_out++) { + if (m.exists_in_model && t < m.n_temperature_points()) { + double sigma_t = + m.get_xs(MgxsType::TOTAL, g_out, NULL, NULL, NULL, t, a); + sigma_t_.push_back(sigma_t); + + if (sigma_t < MINIMUM_MACRO_XS) { + Material* mat = model::materials[i].get(); + warning(fmt::format( + "Material \"{}\" (id: {}) has a group {} total cross section " + "({:.3e}) below the minimum threshold " + "({:.3e}). Material will be treated as pure void.", + mat->name(), mat->id(), g_out, sigma_t, MINIMUM_MACRO_XS)); + } - double nu_sigma_f = - m.get_xs(MgxsType::NU_FISSION, g_out, NULL, NULL, NULL, t, a); - nu_sigma_f_.push_back(nu_sigma_f); + double nu_sigma_f = + m.get_xs(MgxsType::NU_FISSION, g_out, NULL, NULL, NULL, t, a); + nu_sigma_f_.push_back(nu_sigma_f); - double sigma_f = - m.get_xs(MgxsType::FISSION, g_out, NULL, NULL, NULL, t, a); - sigma_f_.push_back(sigma_f); + double sigma_f = + m.get_xs(MgxsType::FISSION, g_out, NULL, NULL, NULL, t, a); + sigma_f_.push_back(sigma_f); - double chi = - m.get_xs(MgxsType::CHI_PROMPT, g_out, &g_out, NULL, NULL, t, a); - if (!std::isfinite(chi)) { - // MGXS interface may return NaN in some cases, such as when material - // is fissionable but has very small sigma_f. - chi = 0.0; - } - chi_.push_back(chi); - - double kappa_fission = - m.get_xs(MgxsType::KAPPA_FISSION, g_out, NULL, NULL, NULL, t, a); - kappa_fission_.push_back(kappa_fission); - - for (int g_in = 0; g_in < negroups_; g_in++) { - double sigma_s = - m.get_xs(MgxsType::NU_SCATTER, g_in, &g_out, NULL, NULL, t, a); - sigma_s_.push_back(sigma_s); - // For transport corrected XS data, diagonal elements may be negative. - // In this case, set a flag to enable transport stabilization for the - // simulation. - if (g_out == g_in && sigma_s < 0.0) - is_transport_stabilization_needed_ = true; - } - } else { - sigma_t_.push_back(0); - nu_sigma_f_.push_back(0); - sigma_f_.push_back(0); - chi_.push_back(0); - kappa_fission_.push_back(0); - for (int g_in = 0; g_in < negroups_; g_in++) { - sigma_s_.push_back(0); + double chi = + m.get_xs(MgxsType::CHI_PROMPT, g_out, &g_out, NULL, NULL, t, a); + if (!std::isfinite(chi)) { + // MGXS interface may return NaN in some cases, such as when + // material is fissionable but has very small sigma_f. + chi = 0.0; + } + chi_.push_back(chi); + + double kappa_fission = + m.get_xs(MgxsType::KAPPA_FISSION, g_out, NULL, NULL, NULL, t, a); + kappa_fission_.push_back(kappa_fission); + + for (int g_in = 0; g_in < negroups_; g_in++) { + double sigma_s = + m.get_xs(MgxsType::NU_SCATTER, g_in, &g_out, NULL, NULL, t, a); + sigma_s_.push_back(sigma_s); + // For transport corrected XS data, diagonal elements may be + // negative. In this case, set a flag to enable transport + // stabilization for the simulation. + if (g_out == g_in && sigma_s < 0.0) + is_transport_stabilization_needed_ = true; + } + } else { + sigma_t_.push_back(0); + nu_sigma_f_.push_back(0); + sigma_f_.push_back(0); + chi_.push_back(0); + kappa_fission_.push_back(0); + for (int g_in = 0; g_in < negroups_; g_in++) { + sigma_s_.push_back(0); + } } } } @@ -1263,12 +1293,14 @@ void FlatSourceDomain::set_adjoint_sources() #pragma omp parallel for for (int64_t sr = 0; sr < n_source_regions(); sr++) { int material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); if (material == MATERIAL_VOID) { continue; } for (int g = 0; g < negroups_; g++) { double sigma_t = - sigma_t_[material * negroups_ + g] * source_regions_.density_mult(sr); + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + source_regions_.density_mult(sr); source_regions_.external_source(sr, g) /= sigma_t; } } @@ -1278,15 +1310,17 @@ void FlatSourceDomain::transpose_scattering_matrix() { // Transpose the inner two dimensions for each material for (int m = 0; m < n_materials_; ++m) { - int material_offset = m * negroups_ * negroups_; - for (int i = 0; i < negroups_; ++i) { - for (int j = i + 1; j < negroups_; ++j) { - // Calculate indices of the elements to swap - int idx1 = material_offset + i * negroups_ + j; - int idx2 = material_offset + j * negroups_ + i; - - // Swap the elements to transpose the matrix - std::swap(sigma_s_[idx1], sigma_s_[idx2]); + for (int t = 0; t < ntemperature_; t++) { + int material_offset = (m * ntemperature_ + t) * negroups_ * negroups_; + for (int i = 0; i < negroups_; ++i) { + for (int j = i + 1; j < negroups_; ++j) { + // Calculate indices of the elements to swap + int idx1 = material_offset + i * negroups_ + j; + int idx2 = material_offset + j * negroups_ + i; + + // Swap the elements to transpose the matrix + std::swap(sigma_s_[idx1], sigma_s_[idx2]); + } } } } @@ -1506,18 +1540,26 @@ SourceRegionHandle FlatSourceDomain::get_subdivided_source_region_handle( int gs_i_cell = gs.lowest_coord().cell(); Cell& cell = *model::cells[gs_i_cell]; int material = cell.material(gs.cell_instance()); + int temp = 0; // If material total XS is extremely low, just set it to void to avoid // problems with 1/Sigma_t - for (int g = 0; g < negroups_; g++) { - double sigma_t = sigma_t_[material * negroups_ + g]; - if (sigma_t < MINIMUM_MACRO_XS) { - material = MATERIAL_VOID; - break; + if (material != MATERIAL_VOID) { + temp = data::mg.macro_xs_[material].get_temperature_index( + cell.sqrtkT(gs.cell_instance())); + for (int g = 0; g < negroups_; g++) { + double sigma_t = + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g]; + if (sigma_t < MINIMUM_MACRO_XS) { + material = MATERIAL_VOID; + temp = 0; + break; + } } } handle.material() = material; + handle.temperature_idx() = temp; handle.density_mult() = cell.density_mult(gs.cell_instance()); @@ -1550,7 +1592,8 @@ SourceRegionHandle FlatSourceDomain::get_subdivided_source_region_handle( if (material != C_NONE) { for (int g = 0; g < negroups_; g++) { double sigma_t = - sigma_t_[material * negroups_ + g] * handle.density_mult(); + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + handle.density_mult(); handle.external_source(g) /= sigma_t; } } @@ -1625,6 +1668,7 @@ void FlatSourceDomain::apply_transport_stabilization() #pragma omp parallel for for (int64_t sr = 0; sr < n_source_regions(); sr++) { int material = source_regions_.material(sr); + int temp = source_regions_.temperature_idx(sr); double density_mult = source_regions_.density_mult(sr); if (material == MATERIAL_VOID) { continue; @@ -1633,10 +1677,14 @@ void FlatSourceDomain::apply_transport_stabilization() // Only apply stabilization if the diagonal (in-group) scattering XS is // negative double sigma_s = - sigma_s_[material * negroups_ * negroups_ + g * negroups_ + g] * + sigma_s_[((material * ntemperature_ + temp) * negroups_ + g) * + negroups_ + + g] * density_mult; if (sigma_s < 0.0) { - double sigma_t = sigma_t_[material * negroups_ + g] * density_mult; + double sigma_t = + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + density_mult; double phi_new = source_regions_.scalar_flux_new(sr, g); double phi_old = source_regions_.scalar_flux_old(sr, g); diff --git a/src/random_ray/linear_source_domain.cpp b/src/random_ray/linear_source_domain.cpp index 02f4c9e2355..b4701ed1fa9 100644 --- a/src/random_ray/linear_source_domain.cpp +++ b/src/random_ray/linear_source_domain.cpp @@ -43,13 +43,16 @@ void LinearSourceDomain::update_single_neutron_source(SourceRegionHandle& srh) // Add scattering + fission source int material = srh.material(); + int temp = srh.temperature_idx(); double density_mult = srh.density_mult(); if (material != MATERIAL_VOID) { double inverse_k_eff = 1.0 / k_eff_; MomentMatrix invM = srh.mom_matrix().inverse(); for (int g_out = 0; g_out < negroups_; g_out++) { - double sigma_t = sigma_t_[material * negroups_ + g_out] * density_mult; + double sigma_t = + sigma_t_[(material * ntemperature_ + temp) * negroups_ + g_out] * + density_mult; double scatter_flat = 0.0f; double fission_flat = 0.0f; @@ -62,12 +65,16 @@ void LinearSourceDomain::update_single_neutron_source(SourceRegionHandle& srh) MomentArray flux_linear = srh.flux_moments_old(g_in); // Handles for cross sections - double sigma_s = sigma_s_[material * negroups_ * negroups_ + - g_out * negroups_ + g_in] * - density_mult; + double sigma_s = + sigma_s_[((material * ntemperature_ + temp) * negroups_ + g_out) * + negroups_ + + g_in] * + density_mult; double nu_sigma_f = - nu_sigma_f_[material * negroups_ + g_in] * density_mult; - double chi = chi_[material * negroups_ + g_out]; + nu_sigma_f_[(material * ntemperature_ + temp) * negroups_ + g_in] * + density_mult; + double chi = + chi_[(material * ntemperature_ + temp) * negroups_ + g_out]; // Compute source terms for flat and linear components of the flux scatter_flat += sigma_s * flux_flat; diff --git a/src/random_ray/random_ray.cpp b/src/random_ray/random_ray.cpp index 1b61d8c2072..47d901d63a0 100644 --- a/src/random_ray/random_ray.cpp +++ b/src/random_ray/random_ray.cpp @@ -432,11 +432,13 @@ void RandomRay::attenuate_flux_flat_source( // Get material int material = srh.material(); + int temp = srh.temperature_idx(); // MOC incoming flux attenuation + source contribution/attenuation equation for (int g = 0; g < negroups_; g++) { float sigma_t = - domain_->sigma_t_[material * negroups_ + g] * srh.density_mult(); + domain_->sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + srh.density_mult(); float tau = sigma_t * distance; float exponential = cjosey_exponential(tau); // exponential = 1 - exp(-tau) float new_delta_psi = (angular_flux_[g] - srh.source(g)) * exponential; @@ -531,6 +533,7 @@ void RandomRay::attenuate_flux_linear_source( n_event()++; int material = srh.material(); + int temp = srh.temperature_idx(); Position& centroid = srh.centroid(); Position midpoint = r + u() * (distance / 2.0); @@ -560,7 +563,8 @@ void RandomRay::attenuate_flux_linear_source( // Compute tau, the optical thickness of the ray segment float sigma_t = - domain_->sigma_t_[material * negroups_ + g] * srh.density_mult(); + domain_->sigma_t_[(material * ntemperature_ + temp) * negroups_ + g] * + srh.density_mult(); float tau = sigma_t * distance; // If tau is very small, set it to zero to avoid numerical issues. @@ -765,6 +769,7 @@ void RandomRay::attenuate_flux_linear_source_void( void RandomRay::initialize_ray(uint64_t ray_id, FlatSourceDomain* domain) { domain_ = domain; + ntemperature_ = domain->ntemperature_; // Reset particle event counter n_event() = 0; diff --git a/src/random_ray/random_ray_simulation.cpp b/src/random_ray/random_ray_simulation.cpp index 7bab3a9b1b6..7531f3ffe7c 100644 --- a/src/random_ray/random_ray_simulation.cpp +++ b/src/random_ray/random_ray_simulation.cpp @@ -121,10 +121,6 @@ void validate_random_ray_inputs() fatal_error("Anisotropic MGXS detected. Only isotropic XS data sets " "supported in random ray mode."); } - if (material.get_xsdata().size() > 1) { - warning("Non-isothermal MGXS detected. Only isothermal XS data sets " - "supported in random ray mode. Using lowest temperature."); - } for (int g = 0; g < data::mg.num_energy_groups_; g++) { if (material.exists_in_model) { // Temperature and angle indices, if using multiple temperature diff --git a/src/random_ray/source_region.cpp b/src/random_ray/source_region.cpp index 15c65221aa7..78543c5ab53 100644 --- a/src/random_ray/source_region.cpp +++ b/src/random_ray/source_region.cpp @@ -11,11 +11,11 @@ namespace openmc { //============================================================================== SourceRegionHandle::SourceRegionHandle(SourceRegion& sr) : negroups_(sr.scalar_flux_old_.size()), material_(&sr.material_), - density_mult_(&sr.density_mult_), is_small_(&sr.is_small_), - n_hits_(&sr.n_hits_), is_linear_(sr.source_gradients_.size() > 0), - lock_(&sr.lock_), volume_(&sr.volume_), volume_t_(&sr.volume_t_), - volume_sq_(&sr.volume_sq_), volume_sq_t_(&sr.volume_sq_t_), - volume_naive_(&sr.volume_naive_), + temperature_idx_(&sr.temperature_idx_), density_mult_(&sr.density_mult_), + is_small_(&sr.is_small_), n_hits_(&sr.n_hits_), + is_linear_(sr.source_gradients_.size() > 0), lock_(&sr.lock_), + volume_(&sr.volume_), volume_t_(&sr.volume_t_), volume_sq_(&sr.volume_sq_), + volume_sq_t_(&sr.volume_sq_t_), volume_naive_(&sr.volume_naive_), position_recorded_(&sr.position_recorded_), external_source_present_(&sr.external_source_present_), position_(&sr.position_), centroid_(&sr.centroid_), @@ -71,6 +71,7 @@ void SourceRegionContainer::push_back(const SourceRegion& sr) // Scalar fields material_.push_back(sr.material_); + temperature_idx_.push_back(sr.temperature_idx_); density_mult_.push_back(sr.density_mult_); is_small_.push_back(sr.is_small_); n_hits_.push_back(sr.n_hits_); @@ -125,6 +126,7 @@ void SourceRegionContainer::assign( // Clear existing data n_source_regions_ = 0; material_.clear(); + temperature_idx_.clear(); density_mult_.clear(); is_small_.clear(); n_hits_.clear(); @@ -183,6 +185,7 @@ SourceRegionHandle SourceRegionContainer::get_source_region_handle(int64_t sr) SourceRegionHandle handle; handle.negroups_ = negroups(); handle.material_ = &material(sr); + handle.temperature_idx_ = &temperature_idx(sr); handle.density_mult_ = &density_mult(sr); handle.is_small_ = &is_small(sr); handle.n_hits_ = &n_hits(sr); diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/__init__.py b/tests/regression_tests/random_ray_auto_convert_temperature/__init__.py new file mode 100644 index 00000000000..e69de29bb2d diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/inputs_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/inputs_true.dat new file mode 100644 index 00000000000..e88b494ab4c --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/inputs_true.dat @@ -0,0 +1,68 @@ + + + + mgxs.h5 + + + + + + + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 10 + 5 + + + -0.63 -0.63 -1 0.63 0.63 1 + + + true + + + multi-group + nearest + true + 200.0 400.0 + 200.0 + + + + -0.63 -0.63 -1.0 0.63 0.63 1.0 + + + 30.0 + 150.0 + + + + + + linear + + + 2 2 + -0.63 -0.63 + 0.63 0.63 + + + diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/results_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/results_true.dat new file mode 100644 index 00000000000..fee8bf67088 --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/infinite_medium/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +7.499800E-01 1.615317E-02 diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/inputs_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/inputs_true.dat new file mode 100644 index 00000000000..e88b494ab4c --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/inputs_true.dat @@ -0,0 +1,68 @@ + + + + mgxs.h5 + + + + + + + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 10 + 5 + + + -0.63 -0.63 -1 0.63 0.63 1 + + + true + + + multi-group + nearest + true + 200.0 400.0 + 200.0 + + + + -0.63 -0.63 -1.0 0.63 0.63 1.0 + + + 30.0 + 150.0 + + + + + + linear + + + 2 2 + -0.63 -0.63 + 0.63 0.63 + + + diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/results_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/results_true.dat new file mode 100644 index 00000000000..ef0a4b87aef --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/material_wise/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +7.367927E-01 6.850805E-03 diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/inputs_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/inputs_true.dat new file mode 100644 index 00000000000..e88b494ab4c --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/inputs_true.dat @@ -0,0 +1,68 @@ + + + + mgxs.h5 + + + + + + + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 10 + 5 + + + -0.63 -0.63 -1 0.63 0.63 1 + + + true + + + multi-group + nearest + true + 200.0 400.0 + 200.0 + + + + -0.63 -0.63 -1.0 0.63 0.63 1.0 + + + 30.0 + 150.0 + + + + + + linear + + + 2 2 + -0.63 -0.63 + 0.63 0.63 + + + diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/results_true.dat b/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/results_true.dat new file mode 100644 index 00000000000..a702ec851b7 --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/stochastic_slab/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +6.431774E-01 2.076589E-02 diff --git a/tests/regression_tests/random_ray_auto_convert_temperature/test.py b/tests/regression_tests/random_ray_auto_convert_temperature/test.py new file mode 100644 index 00000000000..99c99e6147f --- /dev/null +++ b/tests/regression_tests/random_ray_auto_convert_temperature/test.py @@ -0,0 +1,73 @@ +import os + +import openmc +from openmc.examples import pwr_pin_cell +from openmc import RegularMesh +from openmc.utility_funcs import change_directory +import pytest + +from tests.testing_harness import TolerantPyAPITestHarness + + +class MGXSTestHarness(TolerantPyAPITestHarness): + def _cleanup(self): + super()._cleanup() + f = 'mgxs.h5' + if os.path.exists(f): + os.remove(f) + + +@pytest.mark.parametrize("method", ["material_wise", "stochastic_slab", "infinite_medium"]) +def test_random_ray_auto_convert(method): + with change_directory(method): + openmc.reset_auto_ids() + + # Start with a normal continuous energy model + model = pwr_pin_cell() + + temp_settings = { + 'method' : 'nearest', + 'tolerance' : 200.0, + 'range' : (200.0, 400.0), + 'multipole' : True + } + + # Convert to a multi-group model + model.convert_to_multigroup( + method=method, groups='CASMO-2', nparticles=100, + overwrite_mgxs_library=False, mgxs_path="mgxs.h5", + temperatures=[294.0, 394.0], temperature_settings=temp_settings + ) + + # Convert to a random ray model + model.convert_to_random_ray() + model.settings.temperature = temp_settings + + # Set all material temperatures to room temperature + for mat in model.geometry.get_all_materials().values(): + mat.temperature = 294.0 + + # Set the cell temperature of the fuel such that it moves up to the next + # temperature bin. + for cell in model.geometry.get_all_cells().values(): + if cell.name == "Fuel": + cell.temperature = [395.0] + + # Set the number of particles + model.settings.particles = 100 + + # Overlay a basic 2x2 mesh + n = 2 + mesh = RegularMesh() + mesh.dimension = (n, n) + bbox = model.geometry.bounding_box + mesh.lower_left = (bbox.lower_left[0], bbox.lower_left[1]) + mesh.upper_right = (bbox.upper_right[0], bbox.upper_right[1]) + model.settings.random_ray['source_region_meshes'] = [ + (mesh, [model.geometry.root_universe])] + + # Set the source shape to linear + model.settings.random_ray['source_shape'] = 'linear' + + harness = MGXSTestHarness('statepoint.10.h5', model) + harness.main() diff --git a/tests/regression_tests/random_ray_cell_density/test.py b/tests/regression_tests/random_ray_cell_density/test.py index cb6062cdc18..48ebe0baaa2 100644 --- a/tests/regression_tests/random_ray_cell_density/test.py +++ b/tests/regression_tests/random_ray_cell_density/test.py @@ -22,7 +22,7 @@ def test_random_ray_basic(run_mode): if run_mode == "eigen": openmc.reset_auto_ids() model = random_ray_lattice() - # Double the densities of the lower-left fuel pin -> cell instances [0, 9). + # Double the densities of the lower-left fuel pin -> cell instances [0, 8). for id, cell in model.geometry.get_all_cells().items(): if cell.fill.name == "UO2 fuel": cell.density = [((i < 8) + 1.0) for i in range(24)] diff --git a/tests/regression_tests/random_ray_cell_temperature/__init__.py b/tests/regression_tests/random_ray_cell_temperature/__init__.py new file mode 100644 index 00000000000..e69de29bb2d diff --git a/tests/regression_tests/random_ray_cell_temperature/inputs_true.dat b/tests/regression_tests/random_ray_cell_temperature/inputs_true.dat new file mode 100644 index 00000000000..f451f9e5705 --- /dev/null +++ b/tests/regression_tests/random_ray_cell_temperature/inputs_true.dat @@ -0,0 +1,112 @@ + + + + mgxs.h5 + + + + + + + + + + + + + + + + + + + + + + + + + + + + + 0.126 0.126 + 10 10 + -0.63 -0.63 + +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 +3 3 3 3 3 3 3 3 3 3 + + + 1.26 1.26 + 2 2 + -1.26 -1.26 + +2 2 +2 5 + + + + + + + + + + + + + + + + + + + + + eigenvalue + 100 + 10 + 5 + multi-group + nearest + 200.0 400.0 + 10.0 + + 100.0 + 20.0 + + + -1.26 -1.26 -1 1.26 1.26 1 + + + true + + + + + 2 2 + -1.26 -1.26 + 1.26 1.26 + + + 1 + + + 1e-05 0.0635 10.0 100.0 1000.0 500000.0 1000000.0 20000000.0 + + + 1 2 + flux fission nu-fission + analog + + + diff --git a/tests/regression_tests/random_ray_cell_temperature/results_true.dat b/tests/regression_tests/random_ray_cell_temperature/results_true.dat new file mode 100644 index 00000000000..50476beb614 --- /dev/null +++ b/tests/regression_tests/random_ray_cell_temperature/results_true.dat @@ -0,0 +1,171 @@ +k-combined: +8.721099E-01 6.686066E-03 +tally 1: +1.530055E+00 +4.684582E-01 +2.918032E-01 +1.703772E-02 +7.101906E-01 +1.009209E-01 +7.732582E-01 +1.197247E-01 +5.787213E-02 +6.705930E-04 +1.408492E-01 +3.972179E-03 +4.223488E-01 +3.605283E-02 +6.818063E-03 +9.396342E-06 +1.659380E-02 +5.565812E-05 +5.958551E-01 +7.215230E-02 +9.932163E-03 +2.004773E-05 +2.417290E-02 +1.187504E-04 +1.685106E+00 +5.752729E-01 +9.834826E-03 +1.960236E-05 +2.393629E-02 +1.161151E-04 +4.400457E+00 +3.886563E+00 +3.318560E-03 +2.211175E-06 +8.211544E-03 +1.353859E-05 +2.814971E+00 +1.585202E+00 +1.876190E-02 +7.040316E-05 +5.218528E-02 +5.446713E-04 +1.970406E+00 +7.765020E-01 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +8.643086E-01 +1.494764E-01 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +4.421562E-01 +3.961223E-02 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +6.421077E-01 +8.385255E-02 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +1.717683E+00 +5.974687E-01 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +4.024715E+00 +3.251629E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +2.507885E+00 +1.258277E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +1.237600E+00 +3.066006E-01 +4.632318E-01 +4.294997E-02 +1.127413E+00 +2.544090E-01 +7.085300E-01 +1.005838E-01 +1.066845E-01 +2.280854E-03 +2.596486E-01 +1.351037E-02 +4.080999E-01 +3.359223E-02 +1.334606E-02 +3.590848E-05 +3.248163E-02 +2.126996E-04 +5.587686E-01 +6.339264E-02 +1.868162E-02 +7.084921E-05 +4.546734E-02 +4.196670E-04 +1.647852E+00 +5.503007E-01 +1.940765E-02 +7.635704E-05 +4.723492E-02 +4.523030E-04 +4.672431E+00 +4.381400E+00 +7.228567E-03 +1.048470E-05 +1.788658E-02 +6.419579E-05 +3.054158E+00 +1.866081E+00 +4.202691E-02 +3.534383E-04 +1.168957E-01 +2.734362E-03 +1.313678E+00 +3.454374E-01 +4.964081E-01 +4.934325E-02 +1.208158E+00 +2.922789E-01 +7.298176E-01 +1.067017E-01 +1.106319E-01 +2.452869E-03 +2.692559E-01 +1.452928E-02 +4.129038E-01 +3.441066E-02 +1.357107E-02 +3.713729E-05 +3.302926E-02 +2.199783E-04 +5.677464E-01 +6.543807E-02 +1.908992E-02 +7.390190E-05 +4.646105E-02 +4.377492E-04 +1.651067E+00 +5.522453E-01 +1.956802E-02 +7.754346E-05 +4.762523E-02 +4.593308E-04 +4.583305E+00 +4.217589E+00 +7.135836E-03 +1.022408E-05 +1.765713E-02 +6.260005E-05 +2.988394E+00 +1.786305E+00 +4.141550E-02 +3.431964E-04 +1.151951E-01 +2.655125E-03 diff --git a/tests/regression_tests/random_ray_cell_temperature/test.py b/tests/regression_tests/random_ray_cell_temperature/test.py new file mode 100644 index 00000000000..ef36d9238a0 --- /dev/null +++ b/tests/regression_tests/random_ray_cell_temperature/test.py @@ -0,0 +1,36 @@ +import os + +import openmc +from openmc.examples import random_ray_lattice, random_ray_three_region_cube +from openmc.utility_funcs import change_directory +import pytest + +from tests.testing_harness import TolerantPyAPITestHarness + + +class MGXSTestHarness(TolerantPyAPITestHarness): + def _cleanup(self): + super()._cleanup() + f = 'mgxs.h5' + if os.path.exists(f): + os.remove(f) + + +def test_random_ray_basic(): + openmc.reset_auto_ids() + model = random_ray_lattice(second_temp=True) + # Set the temperature of the lower-left pin to 395 K -> cell instances [0, 8). + # All other pins are set to 295. + for id, cell in model.geometry.get_all_cells().items(): + if cell.fill.name == "UO2 fuel": + cell.temperature = [(100.0 * (i < 8) + 295.0) for i in range(24)] + + model.settings.temperature = { + 'method' : 'nearest', + 'tolerance' : 10.0, + 'range' : (200.0, 400.0) + } + + # Gold file was generated with manually scaled fuel cross sections. + harness = MGXSTestHarness('statepoint.10.h5', model) + harness.main()