diff --git a/openmc/__init__.py b/openmc/__init__.py index 8985099ba8..7366e8383e 100644 --- a/openmc/__init__.py +++ b/openmc/__init__.py @@ -30,8 +30,10 @@ from openmc.plotter import * from openmc.search import * from openmc.polynomial import * from . import examples +from openmc.mpi import DummyCommunicator, MPI, comm # Import a few names from the model module from openmc.model import rectangular_prism, hexagonal_prism, Model + __version__ = '0.13.0-dev' diff --git a/openmc/deplete/__init__.py b/openmc/deplete/__init__.py index 00b29e3196..b8c1cdffff 100644 --- a/openmc/deplete/__init__.py +++ b/openmc/deplete/__init__.py @@ -4,17 +4,6 @@ openmc.deplete A depletion front-end tool. """ -import sys -from unittest.mock import Mock - -from .dummy_comm import DummyCommunicator - -try: - from mpi4py import MPI - comm = MPI.COMM_WORLD -except ImportError: - MPI = Mock() - comm = DummyCommunicator() from .nuclide import * from .chain import * diff --git a/openmc/deplete/abc.py b/openmc/deplete/abc.py index ca16926eaf..bf9f4da827 100644 --- a/openmc/deplete/abc.py +++ b/openmc/deplete/abc.py @@ -22,7 +22,7 @@ from uncertainties import ufloat from openmc.data import DataLibrary from openmc.lib import MaterialFilter, Tally from openmc.checkvalue import check_type, check_greater_than -from . import comm +from openmc import comm from .results import Results from .chain import Chain from .results_list import ResultsList diff --git a/openmc/deplete/helpers.py b/openmc/deplete/helpers.py index 6c3badc090..a8cf5212ee 100644 --- a/openmc/deplete/helpers.py +++ b/openmc/deplete/helpers.py @@ -10,7 +10,7 @@ import sys from numpy import dot, zeros, newaxis, asarray -from . import comm +from openmc import comm from openmc.checkvalue import check_type, check_greater_than from openmc.data import JOULE_PER_EV, REACTION_MT from openmc.lib import ( diff --git a/openmc/deplete/operator.py b/openmc/deplete/operator.py index 5b000a3891..e4db13f95d 100644 --- a/openmc/deplete/operator.py +++ b/openmc/deplete/operator.py @@ -20,7 +20,7 @@ from uncertainties import ufloat import openmc from openmc.checkvalue import check_value import openmc.lib -from . import comm +from openmc import comm from .abc import TransportOperator, OperatorResult from .atom_number import AtomNumber from .reaction_rates import ReactionRates @@ -176,6 +176,9 @@ class Operator(TransportOperator): results are to be used. diff_burnable_mats : bool Whether to differentiate burnable materials with multiple instances + cleanup_when_done : bool + Whether to finalize and clear the shared library memory when the + depletion operation is complete. Defaults to clearing the library. """ _fission_helpers = { "average": AveragedFissionYieldHelper, @@ -202,7 +205,7 @@ class Operator(TransportOperator): self.settings = settings self.geometry = geometry self.diff_burnable_mats = diff_burnable_mats - self._cleanup_when_done = True + self.cleanup_when_done = True # Reduce the chain before we create more materials if reduce_chain: @@ -318,6 +321,10 @@ class Operator(TransportOperator): # Reset results in OpenMC openmc.lib.reset() + # Update the number densities regardless of the source rate + self.number.set_density(vec) + self._update_materials() + # If the source rate is zero, return zero reaction rates without running # a transport solve if source_rate == 0.0: @@ -328,11 +335,7 @@ class Operator(TransportOperator): # Prevent OpenMC from complaining about re-creating tallies openmc.reset_auto_ids() - # Update status - self.number.set_density(vec) - - # Update material compositions and tally nuclides - self._update_materials() + # Update tally nuclides data in preparation for transport solve nuclides = self._get_tally_nuclides() self._rate_helper.nuclides = nuclides self._normalization_helper.nuclides = nuclides @@ -556,7 +559,7 @@ class Operator(TransportOperator): def finalize(self): """Finalize a depletion simulation and release resources.""" - if self._cleanup_when_done: + if self.cleanup_when_done: openmc.lib.finalize() def _update_materials(self): diff --git a/openmc/deplete/results.py b/openmc/deplete/results.py index ce8331e825..05d78f6293 100644 --- a/openmc/deplete/results.py +++ b/openmc/deplete/results.py @@ -9,7 +9,7 @@ import copy import h5py import numpy as np -from . import comm, MPI +from openmc import comm, MPI from .reaction_rates import ReactionRates VERSION_RESULTS = (1, 1) diff --git a/openmc/lib/core.py b/openmc/lib/core.py index 6b72dde555..6f08f76b2a 100644 --- a/openmc/lib/core.py +++ b/openmc/lib/core.py @@ -2,6 +2,7 @@ from contextlib import contextmanager from ctypes import (c_bool, c_int, c_int32, c_int64, c_double, c_char_p, c_char, POINTER, Structure, c_void_p, create_string_buffer) import sys +import os import numpy as np from numpy.ctypeslib import as_array @@ -110,9 +111,14 @@ def global_bounding_box(): return llc, urc -def calculate_volumes(): +def calculate_volumes(output=True): """Run stochastic volume calculation""" - _dll.openmc_calculate_volumes() + + if output: + _dll.openmc_calculate_volumes() + else: + with quiet_dll(): + _dll.openmc_calculate_volumes() def current_batch(): @@ -127,13 +133,17 @@ def current_batch(): return c_int.in_dll(_dll, 'current_batch').value -def export_properties(filename=None): +def export_properties(filename=None, output=True): """Export physical properties. + .. versionchanged:: 0.13.0 + Parameters ---------- filename : str or None Filename to export properties to (defaults to "properties.h5") + output: bool, optional + Whether or not to show output. Defaults to showing output See Also -------- @@ -142,7 +152,11 @@ def export_properties(filename=None): """ if filename is not None: filename = c_char_p(filename.encode()) - _dll.openmc_properties_export(filename) + if output: + _dll.openmc_properties_export(filename) + else: + with quiet_dll(): + _dll.openmc_properties_export(filename) def finalize(): @@ -220,15 +234,19 @@ def import_properties(filename): _dll.openmc_properties_import(filename.encode()) -def init(args=None, intracomm=None): +def init(args=None, intracomm=None, output=True): """Initialize OpenMC + .. versionchanged:: 0.13.0 + Parameters ---------- - args : list of str + args : list of str, optional Command-line arguments - intracomm : mpi4py.MPI.Intracomm or None + intracomm : mpi4py.MPI.Intracomm or None, optional MPI intracommunicator + output: bool, optional + Whether or not to show output. Defaults to showing output """ if args is not None: @@ -254,7 +272,11 @@ def init(args=None, intracomm=None): address = MPI._addressof(intracomm) intracomm = c_void_p(address) - _dll.openmc_init(argc, argv, intracomm) + if output: + _dll.openmc_init(argc, argv, intracomm) + else: + with quiet_dll(): + _dll.openmc_init(argc, argv, intracomm) openmc.lib.is_initialized = True @@ -345,9 +367,22 @@ def next_batch(): return status.value -def plot_geometry(): - """Plot geometry""" - _dll.openmc_plot_geometry() +def plot_geometry(output=True): + """Plot geometry + + .. versionchanged:: 0.13.0 + + Parameters + ---------- + output: bool, optional + Whether or not to show output. Defaults to showing output + """ + + if output: + _dll.openmc_plot_geometry() + else: + with quiet_dll(): + _dll.openmc_plot_geometry() def reset(): @@ -360,9 +395,22 @@ def reset_timers(): _dll.openmc_reset_timers() -def run(): - """Run simulation""" - _dll.openmc_run() +def run(output=True): + """Run simulation + + .. versionchanged:: 0.13.0 + + Parameters + ---------- + output: bool, optional + Whether or not to show output. Defaults to showing output + """ + + if output: + _dll.openmc_run() + else: + with quiet_dll(): + _dll.openmc_run() def simulation_init(): @@ -478,3 +526,30 @@ class _FortranObjectWithID(_FortranObject): # assigned. If the array index of the object is out of bounds, an # OutOfBoundsError will be raised here by virtue of referencing self.id self.id + + +@contextmanager +def quiet_dll(): + """This context manager allows us to suppress standard output from DLLs""" + sys.stdout.flush() + # Save the initial file descriptor states + initial_stdout = sys.stdout + initial_stdout_fno = os.dup(sys.stdout.fileno()) + # Get a garbage descriptor so we can throw away output + devnull = os.open(os.devnull, os.O_WRONLY) + + # Get the current stdout stream and make a duplicate of it + new_stdout = os.dup(1) + # Copy the garbage output to the stdout stream + os.dup2(devnull, 1) + os.close(devnull) + # Now point stdout to the re-defined stdout + sys.stdout = os.fdopen(new_stdout, 'w') + + try: + yield + finally: + # Now we just clean up after ourselves and reset the streams + sys.stdout = initial_stdout + sys.stdout.flush() + os.dup2(initial_stdout_fno, 1) diff --git a/openmc/model/model.py b/openmc/model/model.py index e560ef9fb5..c65470e3e3 100644 --- a/openmc/model/model.py +++ b/openmc/model/model.py @@ -191,17 +191,11 @@ class Model: @intracomm.setter def intracomm(self, intracomm): - try: - from mpi4py import MPI - mpi_avail = True - except ImportError: - mpi_avail = False - if intracomm is None or not mpi_avail: - # TODO: move dummy_comm from openmc.deplete to openmc - from openmc.deplete.dummy_comm import DummyCommunicator - self._intracomm = DummyCommunicator() + if intracomm is None: + self._intracomm = openmc.comm else: - check_type('intracomm', intracomm, MPI.Comm) + check_type('intracomm', intracomm, + (openmc.MPI.Comm, openmc.DummyCommunicator)) self._intracomm = intracomm @classmethod @@ -263,25 +257,12 @@ class Model: # where it exists (here in init), but in the future the functionality # should be exposed so that it can be accessed via model.run(...) - # TODO: the output flag is not yet implemented for the openmc.lib.init - # command. This will be the subject of future work. - - args = [] - - if isinstance(threads, Integral) and threads > 0: - args += ['-s', str(threads)] - - if geometry_debug: - args.append('-g') - - if event_based: - args.append('-e') - - if isinstance(restart_file, str): - args += ['-r', restart_file] - - if tracks: - args.append('-t') + args = openmc.process_CLI_arguments( + volume=False, geometry_debug=geometry_debug, + restart_file=restart_file, threads=threads, tracks=tracks, + event_based=event_based) + # Args adds the openmc_exec command in the first entry; remove it + args = args[1:] self.finalize_lib() @@ -289,7 +270,13 @@ class Model: self.export_to_xml() self.intracomm.barrier() - openmc.lib.init(args=args, intracomm=self.intracomm) + if isinstance(self.intracomm, openmc.DummyCommunicator): + # openmc.lib.init does not accept DummyCommunicator, and importing + # the DummyCommunicator class is overkill. Filter it here + intracomm = None + else: + intracomm = self.intracomm + openmc.lib.init(args=args, intracomm=intracomm, output=output) def finalize_lib(self): """Finalize simulation and free memory allocated for the C API @@ -301,7 +288,7 @@ class Model: openmc.lib.finalize() def deplete(self, timesteps, method='cecm', final_step=True, - operator_kwargs=None, integrator_kwargs=None, directory='.'): + operator_kwargs=None, directory='.', **integrator_kwargs): """Deplete model using specified timesteps/power .. versionchanged:: 0.13.0 @@ -320,25 +307,34 @@ class Model: operator_kwargs : dict Keyword arguments passed to the depletion Operator initializer (e.g., :func:`openmc.deplete.Operator`) - integrator_kwargs : dict - Keyword arguments passed to the depletion Operator initializer - (e.g., :func:`openmc.deplete.integrator.cecm`) directory : str, optional Directory to write XML files to. If it doesn't exist already, it will be created. Defaults to the current working directory - **kwargs - Keyword arguments passed to integration function (e.g., - :func:`openmc.deplete.integrator.cecm`) + integrator_kwargs : dict + Remaining keyword arguments passed to the depletion Integrator + initializer (e.g., :func:`openmc.deplete.integrator.cecm`). """ + if operator_kwargs is None: + op_kwargs = {} + elif isinstance(operator_kwargs, dict): + op_kwargs = operator_kwargs + else: + msg = "operator_kwargs must be a dict or None" + raise ValueError(msg) + # Import openmc.deplete here so the Model can be used even if the # shared library is unavailable. import openmc.deplete as dep + # Store whether or not the library was initialized when we started + started_initialized = self.is_initialized + with _change_directory(Path(directory)): depletion_operator = \ - dep.Operator(self.geometry, self.settings, **operator_kwargs) + dep.Operator(self.geometry, self.settings, **op_kwargs) + # Tell depletion_operator.finalize NOT to clear C API memory when # it is done depletion_operator.cleanup_when_done = False @@ -351,13 +347,9 @@ class Model: **integrator_kwargs) # Now perform the depletion + # TODO: add output parameter to integrate integrator.integrate(final_step) - # If we did not perform a transport calculation on the final step, - # then make the code update the C API material inventory - if not final_step: - depletion_operator._update_materials() - # Now make the python Materials match the C API material data for mat_id, mat in self._materials_by_id.items(): if mat.depletable: @@ -370,6 +362,11 @@ class Model: mat.add_nuclide(nuc, density) mat.set_density('atom/b-cm', sum(densities)) + # If we didnt start intialized, we should cleanup after ourselves + if not started_initialized: + depletion_operator.cleanup_when_done = True + depletion_operator.finalize() + def export_to_xml(self, directory='.'): """Export model to XML files. @@ -531,25 +528,16 @@ class Model: if isinstance(particles, Integral) and particles > 0: openmc.lib.settings.particles = particles - # If we dont want output, make the verbosity quiet - if not output: - init_verbosity = openmc.lib.settings.verbosity - if not output: - openmc.lib.settings.verbosity = 1 - # Event-based can be set at init-time or on a case-basis. # Handle the argument here. # TODO This will be dealt with in a future change to the C API # Then run using the C API - openmc.lib.run() + openmc.lib.run(output) # Reset changes for the openmc.run kwargs handling if particles is not None: openmc.lib.settings.particles = init_particles - if not output: - # Then re-set the initial verbosity - openmc.lib.settings.verbosity = init_verbosity else: # Then run via the command line @@ -603,7 +591,7 @@ class Model: to the model. Defaults to applying the volumes. """ - if len(self.settings.volume_calculation) == 0: + if len(self.settings.volume_calculations) == 0: # Then there is no volume calculation specified raise ValueError("The Settings.volume_calculation attribute must" " be specified before executing this method!") @@ -618,18 +606,9 @@ class Model: "with the call to mpi4py" raise ValueError(msg) - # Apply the output settings - if not output: - init_verbosity = openmc.lib.settings.verbosity - if not output: - openmc.lib.settings.verbosity = 1 - # Compute the volumes - openmc.lib.calculate_volumes() + openmc.lib.calculate_volumes(output) - # Reset the output verbosity - if not output: - openmc.lib.settings.verbosity = init_verbosity else: self.export_to_xml() openmc.calculate_volumes(threads=threads, output=output, @@ -638,23 +617,19 @@ class Model: # Now we apply the volumes if apply_volumes: - # Load the results - f_names = \ - [f"volume_{i + 1}.h5" - for i in range(len(self.settings.volume_calculations))] - vol_calcs = [openmc.VolumeCalculation.load_results(f_name) - for f_name in f_names] - # And now we can add them to the model - for vol_calc in vol_calcs: + # Load the results and add them to the model + for i, vol_calc in enumerate(self.settings.volume_calculations): + f_name = f"volume_{i + 1}.h5" + vol_calc.load_results(f_name) # First add them to the Python side self.geometry.add_volume_information(vol_calc) # And now repeat for the C API - if vol_calc.domain == 'material': + if self.is_initialized and vol_calc.domain_type == 'material': # Then we can do this in the C API - for domain_id in vol_calc.domains: - self.update_material_volumes( - [domain_id], vol_calc.volumes[domain_id]) + for domain_id in vol_calc.ids: + openmc.lib.materials[domain_id].volume = \ + vol_calc.volumes[domain_id].n def plot_geometry(self, output=True, cwd='.', openmc_exec='openmc', convert=True, convert_exec='convert'): @@ -689,25 +664,15 @@ class Model: "before executing this method!") with _change_directory(Path(cwd)): - # TODO: openmis_initialized doesnt read plots.xml unless it is in plot mode - # so the following will not work. Commented out for now and - # replacing with non-C API code + # TODO: openmc.is_initialized doesnt read plots.xml unless it is + # in plot mode so the following will not work. Commented out for + # now and replacing with non-C API code self.export_to_xml() openmc.plot_geometry(output=output, openmc_exec=openmc_exec) # if self.is_initialized: - # # Apply the output settings - # if not output: - # init_verbosity = openmc.lib.settings.verbosity - # if not output: - # openmc.lib.settings.verbosity = 1 - # # Compute the volumes - # openmc.lib.plot_geometry() - - # # Reset the output verbosity - # if not output: - # openmc.lib.settings.verbosity = init_verbosity + # openmc.lib.plot_geometry(output) # else: # self.export_to_xml() # openmc.plot_geometry(output=output, openmc_exec=openmc_exec) @@ -721,8 +686,8 @@ class Model: png_file = ppm_file.replace('.ppm', '.png') subprocess.check_call([convert_exec, ppm_file, png_file]) - def _change_py_C_attribs(self, names_or_ids, value, obj_type, attrib_name, - density_units='atom/b-cm'): + def _change_py_lib_attribs(self, names_or_ids, value, obj_type, + attrib_name, density_units='atom/b-cm'): # Method to do the same work whether it is a cell or material and # a temperature or volume check_type('names_or_ids', names_or_ids, Iterable, (Integral, str)) @@ -819,7 +784,7 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, vector, 'cell', 'rotation') + self._change_py_lib_attribs(names_or_ids, vector, 'cell', 'rotation') def translate_cells(self, names_or_ids, vector): """Translate the identified cell(s) by the specified translation vector. @@ -841,7 +806,7 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, vector, 'cell', 'translation') + self._change_py_lib_attribs(names_or_ids, vector, 'cell', 'translation') def update_densities(self, names_or_ids, density, density_units='atom/b-cm'): """Update the density of a given set of materials to a new value @@ -863,7 +828,7 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, density, 'material', 'density', + self._change_py_lib_attribs(names_or_ids, density, 'material', 'density', density_units) def update_cell_temperatures(self, names_or_ids, temperature): @@ -884,7 +849,7 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, temperature, 'cell', + self._change_py_lib_attribs(names_or_ids, temperature, 'cell', 'temperature') def update_material_temperatures(self, names_or_ids, temperature): @@ -905,7 +870,7 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, temperature, 'material', + self._change_py_lib_attribs(names_or_ids, temperature, 'material', 'temperature') def update_material_volumes(self, names_or_ids, volume): @@ -926,4 +891,4 @@ class Model: """ - self._change_py_C_attribs(names_or_ids, volume, 'material', 'volume') + self._change_py_lib_attribs(names_or_ids, volume, 'material', 'volume') diff --git a/openmc/deplete/dummy_comm.py b/openmc/mpi.py similarity index 79% rename from openmc/deplete/dummy_comm.py rename to openmc/mpi.py index 7ac9be6c36..b64a84bc6f 100644 --- a/openmc/deplete/dummy_comm.py +++ b/openmc/mpi.py @@ -1,4 +1,5 @@ import sys +from unittest.mock import Mock class DummyCommunicator: @@ -31,3 +32,10 @@ class DummyCommunicator: def Abort(self, exit_code_or_msg): sys.exit(exit_code_or_msg) + +try: + from mpi4py import MPI + comm = MPI.COMM_WORLD +except ImportError: + MPI = Mock() + comm = DummyCommunicator() diff --git a/tests/unit_tests/conftest.py b/tests/unit_tests/conftest.py index 7ddca6dd73..51d1b19a33 100644 --- a/tests/unit_tests/conftest.py +++ b/tests/unit_tests/conftest.py @@ -1,138 +1,9 @@ -from math import pi import openmc import pytest from tests.regression_tests import config -@pytest.fixture(scope='module') -def pin_model_attributes(): - uo2 = openmc.Material(name='UO2') - uo2.set_density('g/cm3', 10.29769) - uo2.add_element('U', 1., enrichment=2.4) - uo2.add_element('O', 2.) - - zirc = openmc.Material(name='Zirc') - zirc.set_density('g/cm3', 6.55) - zirc.add_element('Zr', 1.) - - borated_water = openmc.Material(name='Borated water') - borated_water.set_density('g/cm3', 0.740582) - borated_water.add_element('B', 4.0e-5) - borated_water.add_element('H', 5.0e-2) - borated_water.add_element('O', 2.4e-2) - borated_water.add_s_alpha_beta('c_H_in_H2O') - - mats = openmc.Materials([uo2, zirc, borated_water]) - - pitch = 1.25984 - fuel_or = openmc.ZCylinder(r=0.39218, name='Fuel OR') - clad_or = openmc.ZCylinder(r=0.45720, name='Clad OR') - box = openmc.model.rectangular_prism(pitch, pitch, - boundary_type='reflective') - - # Define cells - fuel_inf_cell = openmc.Cell(name='inf fuel', fill=uo2) - fuel_inf_univ = openmc.Universe(cells=[fuel_inf_cell]) - fuel = openmc.Cell(name='fuel', fill=fuel_inf_univ, region=-fuel_or) - clad = openmc.Cell(fill=zirc, region=+fuel_or & -clad_or) - water = openmc.Cell(fill=borated_water, region=+clad_or & box) - - # Define overall geometry - geom = openmc.Geometry([fuel, clad, water]) - uo2.volume = pi * fuel_or.r**2 - - settings = openmc.Settings() - settings.batches = 100 - settings.inactive = 10 - settings.particles = 1000 - - # Create an initial uniform spatial source distribution over fissionable zones - bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1] - uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:], only_fissionable=True) - settings.source = openmc.source.Source(space=uniform_dist) - - entropy_mesh = openmc.RegularMesh() - entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50] - entropy_mesh.upper_right = [0.39218, 0.39218, 1.e50] - entropy_mesh.dimension = [10, 10, 1] - settings.entropy_mesh = entropy_mesh - - tals = openmc.Tallies() - tal = openmc.Tally(name='test') - tal.scores = ['flux'] - tals.append(tal) - - plot1 = openmc.Plot() - plot1.origin = (0., 0., 0.) - plot1.width = (pitch, pitch) - plot1.pixels = (300, 300) - plot1.color_by = 'material' - plot1.filename = 'test' - plot2 = openmc.Plot() - plot2.origin = (0., 0., 0.) - plot2.width = (pitch, pitch) - plot2.pixels = (300, 300) - plot2.color_by = 'cell' - plots = openmc.Plots((plot1, plot2)) - - chain = './chain_simple.xml' - fission_q = {'U235': 200e6} - - chain_file_xml = """ - - - - - - - - - - - - - - - - - - - - - 2.53000e-02 - - Gd157 Gd156 I135 Xe135 Xe136 Cs135 - 1.093250e-04 2.087260e-04 2.780820e-02 6.759540e-03 2.392300e-02 4.356330e-05 - - - - - - - 2.53000e-02 - - Gd157 Gd156 I135 Xe135 Xe136 Cs135 - 6.142710e-5 1.483250e-04 0.0292737 0.002566345 0.0219242 4.9097e-6 - - - - - - - 2.53000e-02 - - Gd157 Gd156 I135 Xe135 Xe136 Cs135 - 4.141120e-04 7.605360e-04 0.0135457 0.00026864 0.0024432 3.7100E-07 - - - - -""" - operator_kwargs = {'chain_file': chain, 'fission_q': fission_q} - - return (mats, geom, settings, tals, plots, operator_kwargs, chain_file_xml) - @pytest.fixture(scope='module') def mpi_intracomm(): if config['mpi']: diff --git a/tests/unit_tests/test_deplete_chain.py b/tests/unit_tests/test_deplete_chain.py index fa90b35f2e..0904ff39ae 100644 --- a/tests/unit_tests/test_deplete_chain.py +++ b/tests/unit_tests/test_deplete_chain.py @@ -7,7 +7,8 @@ import os from pathlib import Path import numpy as np -from openmc.deplete import comm, Chain, reaction_rates, nuclide, cram, pool +from openmc import comm +from openmc.deplete import Chain, reaction_rates, nuclide, cram, pool import pytest from tests import cdtemp diff --git a/tests/unit_tests/test_deplete_integrator.py b/tests/unit_tests/test_deplete_integrator.py index 3195c59768..ec424b05b8 100644 --- a/tests/unit_tests/test_deplete_integrator.py +++ b/tests/unit_tests/test_deplete_integrator.py @@ -14,11 +14,11 @@ import numpy as np from uncertainties import ufloat import pytest +from openmc import comm from openmc.deplete import ( - ReactionRates, Results, ResultsList, comm, OperatorResult, - PredictorIntegrator, CECMIntegrator, CF4Integrator, CELIIntegrator, - EPCRK4Integrator, LEQIIntegrator, SICELIIntegrator, SILEQIIntegrator, - cram) + ReactionRates, Results, ResultsList, OperatorResult, PredictorIntegrator, + CECMIntegrator, CF4Integrator, CELIIntegrator, EPCRK4Integrator, + LEQIIntegrator, SICELIIntegrator, SILEQIIntegrator, cram) from tests import dummy_operator diff --git a/tests/unit_tests/test_lib.py b/tests/unit_tests/test_lib.py index 41613832f3..079b4f9eb3 100644 --- a/tests/unit_tests/test_lib.py +++ b/tests/unit_tests/test_lib.py @@ -713,7 +713,6 @@ def test_trigger_set_n_batches(uo2_trigger_model, mpi_intracomm): def test_cell_translation(pincell_model_w_univ, mpi_intracomm): openmc.lib.finalize() openmc.lib.init(intracomm=mpi_intracomm) - openmc.lib.simulation_init() # Cell 1 is filled with a material so it has a translation, but we can't # set it. cell = openmc.lib.cells[1] @@ -727,9 +726,12 @@ def test_cell_translation(pincell_model_w_univ, mpi_intracomm): # This time we *can* set it cell.translation = (1., 0., -1.) assert cell.translation == pytest.approx([1., 0., -1.]) + openmc.lib.finalize() -def test_cell_rotation(pincell_model_w_univ): +def test_cell_rotation(pincell_model_w_univ, mpi_intracomm): + openmc.lib.finalize() + openmc.lib.init(intracomm=mpi_intracomm) # Cell 1 is filled with a material so we cannot rotate it, but we can get # its rotation matrix (which will be the identity matrix) cell = openmc.lib.cells[1] @@ -742,3 +744,4 @@ def test_cell_rotation(pincell_model_w_univ): assert cell.rotation == pytest.approx([0., 0., 0.]) cell.rotation = (180., 0., 0.) assert cell.rotation == pytest.approx([180., 0., 0.]) + openmc.lib.finalize() diff --git a/tests/unit_tests/test_model.py b/tests/unit_tests/test_model.py index 6dc135056b..4ab4c7eb78 100644 --- a/tests/unit_tests/test_model.py +++ b/tests/unit_tests/test_model.py @@ -1,10 +1,110 @@ +from math import pi import numpy as np import pytest from pathlib import Path from shutil import which import openmc import openmc.lib -from openmc.deplete.dummy_comm import DummyCommunicator +from openmc.mpi import DummyCommunicator + + +@pytest.fixture(scope='function') +def pin_model_attributes(): + uo2 = openmc.Material(material_id=1, name='UO2') + uo2.set_density('g/cm3', 10.29769) + uo2.add_element('U', 1., enrichment=2.4) + uo2.add_element('O', 2.) + uo2.depletable = True + + zirc = openmc.Material(material_id=2, name='Zirc') + zirc.set_density('g/cm3', 6.55) + zirc.add_element('Zr', 1.) + zirc.depletable = False + + borated_water = openmc.Material(material_id=3, name='Borated water') + borated_water.set_density('g/cm3', 0.740582) + borated_water.add_element('B', 4.0e-5) + borated_water.add_element('H', 5.0e-2) + borated_water.add_element('O', 2.4e-2) + borated_water.add_s_alpha_beta('c_H_in_H2O') + borated_water.depletable = False + + mats = openmc.Materials([uo2, zirc, borated_water]) + + pitch = 1.25984 + fuel_or = openmc.ZCylinder(r=0.39218, name='Fuel OR') + clad_or = openmc.ZCylinder(r=0.45720, name='Clad OR') + box = openmc.model.rectangular_prism(pitch, pitch, + boundary_type='reflective') + + # Define cells + fuel_inf_cell = openmc.Cell(cell_id=1, name='inf fuel', fill=uo2) + fuel_inf_univ = openmc.Universe(universe_id=1, cells=[fuel_inf_cell]) + fuel = openmc.Cell(cell_id=2, name='fuel', + fill=fuel_inf_univ, region=-fuel_or) + clad = openmc.Cell(cell_id=3, fill=zirc, region=+fuel_or & -clad_or) + water = openmc.Cell(cell_id=4, fill=borated_water, region=+clad_or & box) + + # Define overall geometry + geom = openmc.Geometry([fuel, clad, water]) + uo2.volume = pi * fuel_or.r**2 + + settings = openmc.Settings() + settings.batches = 100 + settings.inactive = 10 + settings.particles = 1000 + + # Create a uniform spatial source distribution over fissionable zones + bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1] + uniform_dist = openmc.stats.Box( + bounds[:3], bounds[3:], only_fissionable=True) + settings.source = openmc.source.Source(space=uniform_dist) + + entropy_mesh = openmc.RegularMesh() + entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50] + entropy_mesh.upper_right = [0.39218, 0.39218, 1.e50] + entropy_mesh.dimension = [10, 10, 1] + settings.entropy_mesh = entropy_mesh + + tals = openmc.Tallies() + tal = openmc.Tally(tally_id=1, name='test') + tal.filters = [openmc.MaterialFilter(bins=[uo2])] + tal.scores = ['flux', 'fission'] + tals.append(tal) + + plot1 = openmc.Plot(plot_id=1) + plot1.origin = (0., 0., 0.) + plot1.width = (pitch, pitch) + plot1.pixels = (300, 300) + plot1.color_by = 'material' + plot1.filename = 'test' + plot2 = openmc.Plot(plot_id=2) + plot2.origin = (0., 0., 0.) + plot2.width = (pitch, pitch) + plot2.pixels = (300, 300) + plot2.color_by = 'cell' + plots = openmc.Plots((plot1, plot2)) + + chain = './test_chain.xml' + + chain_file_xml = """ + + + + + + 2.53000e-02 + + Xe136 + 1.0 + + + + +""" + operator_kwargs = {'chain_file': chain} + + return (mats, geom, settings, tals, plots, operator_kwargs, chain_file_xml) def test_init(run_in_tmpdir, pin_model_attributes, mpi_intracomm): @@ -130,7 +230,7 @@ def test_init_finalize_lib(run_in_tmpdir, pin_model_attributes, mpi_intracomm): # We are going to init and then make sure data is loaded mats, geom, settings, tals, plots, _, _ = pin_model_attributes test_model = openmc.Model(geom, mats, settings, tals, plots, mpi_intracomm) - test_model.init_lib() + test_model.init_lib(output=False) # First check that the API is advertised as initialized assert openmc.lib.is_initialized is True @@ -157,13 +257,13 @@ def test_import_properties(run_in_tmpdir, mpi_intracomm): # Create PWR pin cell model and write XML files openmc.reset_auto_ids() model = openmc.examples.pwr_pin_cell() - model.init_lib() + model.init_lib(output=False) # Change fuel temperature and density and export properties cell = openmc.lib.cells[1] cell.set_temperature(600.0) cell.fill.set_density(5.0, 'g/cm3') - openmc.lib.export_properties() + openmc.lib.export_properties(output=False) # Import properties to existing model model.import_properties("properties.h5") @@ -203,17 +303,20 @@ def test_run(run_in_tmpdir, pin_model_attributes, mpi_intracomm): sp_path = test_model.run(output=False) with openmc.StatePoint(sp_path) as sp: cli_keff = sp.k_combined - cli_flux = sp.get_tally(id=1).get_values()[0, 0, 0] + cli_flux = sp.get_tally(id=1).get_values(scores=['flux'])[0, 0, 0] + cli_fiss = sp.get_tally(id=1).get_values(scores=['fission'])[0, 0, 0] - test_model.init_lib() + test_model.init_lib(output=False) sp_path = test_model.run(output=False) with openmc.StatePoint(sp_path) as sp: - C_keff = sp.k_combined - C_flux = sp.get_tally(id=1).get_values()[0, 0, 0] + lib_keff = sp.k_combined + lib_flux = sp.get_tally(id=1).get_values(scores=['flux'])[0, 0, 0] + lib_fiss = sp.get_tally(id=1).get_values(scores=['fission'])[0, 0, 0] # and lets compare results - assert abs(C_keff - cli_keff) < 1e-13 - assert abs(C_flux - cli_flux) < 1e-13 + assert abs(lib_keff - cli_keff) < 1e-13 + assert abs(lib_flux - cli_flux) < 1e-13 + assert abs(lib_fiss - cli_fiss) < 1e-13 # Now we should make sure that the flags for items which should be handled # by init are properly set @@ -248,8 +351,8 @@ def test_plots(run_in_tmpdir, pin_model_attributes, mpi_intracomm): # We will run the test twice, the first time without C API, the second with for i in range(2): if i == 1: - test_model.init_lib() - test_model.plot_geometry(output=True, convert=convert) + test_model.init_lib(output=False) + test_model.plot_geometry(output=False, convert=convert) # Now look for the files, expect to find test.ppm, plot_2.ppm, and if # convert is True, test.png, plot_2.png @@ -262,11 +365,11 @@ def test_plots(run_in_tmpdir, pin_model_attributes, mpi_intracomm): test_model.finalize_lib() -def test_py_C_attributes(run_in_tmpdir, pin_model_attributes, mpi_intracomm): +def test_py_lib_attributes(run_in_tmpdir, pin_model_attributes, mpi_intracomm): mats, geom, settings, tals, plots, _, _ = pin_model_attributes test_model = openmc.Model(geom, mats, settings, tals, plots, mpi_intracomm) - test_model.init_lib() + test_model.init_lib(output=False) # Now we can call rotate_cells, translate_cells, update_densities, # update_cell_temperatures, and update_material_temperatures and make sure @@ -358,3 +461,106 @@ def test_py_C_attributes(run_in_tmpdir, pin_model_attributes, mpi_intracomm): assert abs(test_model.materials[0].volume - 2.) < 1e-13 test_model.finalize_lib() + + +# def test_deplete(run_in_tmpdir, pin_model_attributes, mpi_intracomm): +# mats, geom, settings, tals, plots, op_kwargs, chain_file_xml = \ +# pin_model_attributes +# with open('test_chain.xml', 'w') as f: +# f.write(chain_file_xml) +# test_model = openmc.Model(geom, mats, settings, tals, plots, mpi_intracomm) + +# initial_mat = mats[0].clone() +# initial_u = initial_mat.get_nuclide_atom_densities()['U235'][1] + +# # Note that the chain file includes only U-235 fission to a stable Xe136 w/ +# # a yield of 100%. Thus all the U235 we lose becomes Xe136 + +# # In this test we first run without pre-initializing the shared library +# # data and then compare. Then we repeat with the C API already initialized +# # and make sure we get the same answer +# test_model.deplete([1e6], 'predictor', final_step=False, +# operator_kwargs=op_kwargs, +# power=1.) +# # Get the new Xe136 and U235 atom densities +# after_xe = mats[0].get_nuclide_atom_densities()['Xe136'][1] +# after_u = mats[0].get_nuclide_atom_densities()['U235'][1] +# assert abs((after_xe + after_u) - initial_u) < 1e-15 +# assert test_model.is_initialized is False + +# # Reset the initial material densities +# mats[0].nuclides.clear() +# densities = initial_mat.get_nuclide_atom_densities() +# tot_density = 0. +# for nuc, density in densities.values(): +# mats[0].add_nuclide(nuc, density) +# tot_density += density +# mats[0].set_density('atom/b-cm', tot_density) + +# # Now we can re-run with the pre-initialized API +# test_model.init_lib(output=False) +# test_model.deplete([1e6], 'predictor', final_step=False, +# operator_kwargs=op_kwargs, +# power=1.) +# # Get the new Xe136 and U235 atom densities +# after_lib_xe = mats[0].get_nuclide_atom_densities()['Xe136'][1] +# after_lib_u = mats[0].get_nuclide_atom_densities()['U235'][1] +# assert abs((after_lib_xe + after_lib_u) - initial_u) < 1e-15 +# assert test_model.is_initialized is True + +# # And end by comparing to the previous case +# assert abs(after_xe - after_lib_xe) < 1e-15 +# assert abs(after_u - after_lib_u) < 1e-15 + +# test_model.finalize_lib() + + +def test_calc_volumes(run_in_tmpdir, pin_model_attributes, mpi_intracomm): + mats, geom, settings, tals, plots, _, _ = pin_model_attributes + + test_model = openmc.Model(geom, mats, settings, tals, plots, mpi_intracomm) + + # With no vol calcs, it should fail + with pytest.raises(ValueError): + test_model.calculate_volumes(output=False) + + # Add a cell and mat volume calc + material_vol_calc = openmc.VolumeCalculation( + [mats[2]], samples=1000, lower_left=(-.63, -.63, -100.), + upper_right=(.63, .63, 100.)) + cell_vol_calc = openmc.VolumeCalculation( + [geom.root_universe.cells[3]], samples=1000, + lower_left=(-.63, -.63, -100.), upper_right=(.63, .63, 100.)) + test_model.settings.volume_calculations = \ + [material_vol_calc, cell_vol_calc] + + # Now lets compute the volumes and check to see if it was applied + # First lets do without using the C-API + # Make sure the volumes are unassigned first + assert mats[2].volume is None + assert geom.root_universe.cells[3].volume is None + test_model.calculate_volumes(output=False, apply_volumes=True) + + # Now let's test that we have volumes assigned; we arent checking the + # value, just that the value was changed + assert mats[2].volume > 0. + assert geom.root_universe.cells[3].volume > 0. + + # Now reset the values + mats[2].volume = None + geom.root_universe.cells[3].volume = None + + # And do again with an initialized library + for file in ['volume_1.h5', 'volume_2.h5']: + file = Path(file) + file.unlink() + test_model.init_lib(output=False) + test_model.calculate_volumes(output=False, apply_volumes=True) + assert mats[2].volume > 0. + assert geom.root_universe.cells[3].volume > 0. + assert openmc.lib.materials[3].volume == mats[2].volume + + test_model.finalize_lib() + + +