diff --git a/docs/source/io_formats/geometry.rst b/docs/source/io_formats/geometry.rst index cef3bb79c80..e2f53fc635f 100644 --- a/docs/source/io_formats/geometry.rst +++ b/docs/source/io_formats/geometry.rst @@ -406,24 +406,58 @@ Each ```` element can have the following attributes or sub-eleme *Default*: None - :material_overrides: - This element contains information on material overrides to be applied to the - DAGMC universe. It has the following attributes and sub-elements: - - :cell: - Material override information for a single cell. It contains the following - attributes and sub-elements: - - :id: - The cell ID in the DAGMC geometry for which the material override will - apply. - - :materials: - A list of material IDs that will apply to instances of the cell. If the - list contains only one ID, it will replace the original material - assignment of all instances of the DAGMC cell. If the list contains more - than one material, each material ID of the list will be assigned to the - various instances of the DAGMC cell. + :cell: + Zero or more ```` sub-elements may appear to override properties of + individual DAGMC volumes. Each ```` element supports the following + attributes and sub-elements: + + :id: + The integer cell ID in the DAGMC geometry to override. Required. + + :name: + An optional string label for the cell. + + *Default*: None + + :material: + The material ID to assign to this cell. Use ``void`` for vacuum. Multiple + space-separated IDs may be given to specify a distribmat (distributed + material) assignment. May also be specified as a ```` + sub-element. Required. + + :temperature: + Temperature(s) in Kelvin to assign to the cell. Must be ≥ 0. Multiple + space-separated values may be given. May also be specified as a + ```` sub-element. + + *Default*: None + + :density: + Density in g/cm³ to assign to the cell. Must be > 0. Requires a non-void + material fill. Multiple space-separated values may be given. May also be + specified as a ```` sub-element. + + *Default*: None + + :volume: + Volume of the cell in cm³. + + .. note:: DAGMC can compute cell volumes exactly from the triangulated + mesh surfaces. Specifying a manual volume risks inconsistency + with that capability. + + *Default*: None + + The following standard ```` attributes are **not** supported inside + ```` and will raise an error if present: ``region``, + ``fill``, ``universe``, ``translation``, ``rotation``. + + .. deprecated:: + The ```` sub-element (containing ```` + children with ````) is deprecated. A deprecation warning is + emitted and the overrides are converted to the ```` format at parse + time. It is an error to specify both ```` and + ```` sub-elements on the same ````. *Default*: None diff --git a/include/openmc/cell.h b/include/openmc/cell.h index 26c34ebd41d..9beb56c2d44 100644 --- a/include/openmc/cell.h +++ b/include/openmc/cell.h @@ -145,6 +145,30 @@ class Region { bool simple_; //!< Does the region contain only intersections? }; +//============================================================================== +// XML parsing helpers for nodes +//============================================================================== + +//! Parse material IDs from a XML node. +//! \param node XML node containing a "material" attribute or child element +//! \param cell_id Cell ID used in error messages +//! \return Vector of material IDs (MATERIAL_VOID for "void") +vector parse_cell_material_xml(pugi::xml_node node, int32_t cell_id); + +//! Parse temperatures (in Kelvin) from a XML node. +//! Validates that all values are non-negative and the list is non-empty. +//! \param node XML node containing a "temperature" attribute or child element +//! \param cell_id Cell ID used in error messages +//! \return Vector of temperatures in Kelvin +vector parse_cell_temperature_xml(pugi::xml_node node, int32_t cell_id); + +//! Parse densities (in g/cm³) from a XML node. +//! Validates that all values are positive and the list is non-empty. +//! \param node XML node containing a "density" attribute or child element +//! \param cell_id Cell ID used in error messages +//! \return Vector of densities in g/cm³ +vector parse_cell_density_xml(pugi::xml_node node, int32_t cell_id); + //============================================================================== class Cell { diff --git a/include/openmc/dagmc.h b/include/openmc/dagmc.h index 0e27402a159..858636c0bc7 100644 --- a/include/openmc/dagmc.h +++ b/include/openmc/dagmc.h @@ -94,6 +94,10 @@ class DAGCell : public Cell { class DAGUniverse : public Universe { public: + using MaterialOverrides = std::unordered_map>; + using TemperatureOverrides = std::unordered_map>; + using DensityOverrides = std::unordered_map>; + explicit DAGUniverse(pugi::xml_node node); //! Create a new DAGMC universe @@ -112,6 +116,9 @@ class DAGUniverse : public Universe { //! Initialize the DAGMC accel. data structures, indices, material //! assignments, etc. void initialize(); + void initialize(const MaterialOverrides& material_overrides, + const TemperatureOverrides& temperature_overrides, + const DensityOverrides& density_overrides = {}); //! Reads UWUW materials and returns an ID map void read_uwuw_materials(); @@ -146,7 +153,8 @@ class DAGUniverse : public Universe { //! Assign a material overriding normal assignement to a cell //! \param[in] c The OpenMC cell to which the material is assigned - void override_assign_material(std::unique_ptr& c) const; + void override_assign_material(std::unique_ptr& c, + const MaterialOverrides& material_overrides) const; //! Return the index into the model cells vector for a given DAGMC volume //! handle in the universe @@ -187,7 +195,9 @@ class DAGUniverse : public Universe { void set_id(); //!< Deduce the universe id from model::universes void init_dagmc(); //!< Create and initialise DAGMC pointer void init_metadata(); //!< Create and initialise dagmcMetaData pointer - void init_geometry(); //!< Create cells and surfaces from DAGMC entities + void init_geometry(const MaterialOverrides& material_overrides, + const TemperatureOverrides& temperature_overrides, + const DensityOverrides& density_overrides); std::string filename_; //!< Name of the DAGMC file used to create this universe @@ -201,11 +211,6 @@ class DAGUniverse : public Universe { //!< generate new material IDs for the universe bool has_graveyard_; //!< Indicates if the DAGMC geometry has a "graveyard" //!< volume - std::unordered_map> - material_overrides_; //!< Map of material overrides - //!< keys correspond to the DAGMCCell id - //!< values are a list of material ids used - //!< for the override }; //============================================================================== diff --git a/openmc/dagmc.py b/openmc/dagmc.py index fd4258225f3..fd93a331c24 100644 --- a/openmc/dagmc.py +++ b/openmc/dagmc.py @@ -36,12 +36,6 @@ class DAGMCUniverse(openmc.UniverseBase): auto_mat_ids : bool Set IDs automatically on initialization (True) or report overlaps in ID space between OpenMC and UWUW materials (False) - material_overrides : dict, optional - A dictionary of material overrides. The keys are material name strings - and the values are Iterables of openmc.Material objects. If a material - name is found in the DAGMC file, the material will be replaced with the - openmc.Material object in the value. - Attributes ---------- id : int @@ -78,15 +72,6 @@ class DAGMCUniverse(openmc.UniverseBase): The number of surfaces in the model. .. versionadded:: 0.13.2 - material_overrides : dict - A dictionary of material overrides. Keys are cell IDs; values are - iterables of :class:`openmc.Material` objects. The material assignment - of each DAGMC cell ID key will be replaced with the - :class:`~openmc.Material` object in the value. If the value contains - multiple :class:`~openmc.Material` objects, each Material in the list - will be assigned to the corresponding instance of the cell. - - .. versionadded:: 0.15.1 """ def __init__(self, @@ -94,16 +79,12 @@ def __init__(self, universe_id=None, name='', auto_geom_ids=False, - auto_mat_ids=False, - material_overrides=None): + auto_mat_ids=False): super().__init__(universe_id, name) # Initialize class attributes self.filename = filename self.auto_geom_ids = auto_geom_ids self.auto_mat_ids = auto_mat_ids - self._material_overrides = {} - if material_overrides is not None: - self.material_overrides = material_overrides def __repr__(self): string = super().__repr__() @@ -130,47 +111,17 @@ def filename(self, val: cv.PathLike): @property def material_overrides(self): - return self._material_overrides + raise AttributeError( + "DAGMCUniverse.material_overrides has been removed. Use " + "DAGMCCell objects added via add_cell() to manage per-cell " + "material assignments.") @material_overrides.setter def material_overrides(self, val): - cv.check_type('material overrides', val, Mapping) - for key, value in val.items(): - self.add_material_override(key, value) - - def replace_material_assignment(self, material_name: str, material: openmc.Material): - """Replace a material assignment within the DAGMC universe. - - Replace the material assignment of all cells filled with a material in - the DAGMC universe. The universe must be synchronized in an initialized - Model (see :meth:`~openmc.DAGMCUniverse.sync_dagmc_cells`) before - calling this method. - - .. versionadded:: 0.15.1 - - Parameters - ---------- - material_name : str - Material name to replace - material : openmc.Material - Material to replace the material_name with - - """ - if material_name not in self.material_names: - raise ValueError( - f"No material with name '{material_name}' found in the DAGMC universe") - - if not self.cells: - raise RuntimeError("This DAGMC universe has not been synchronized " - "on an initialized Model.") - - for cell in self.cells.values(): - if cell.fill is None: - continue - if isinstance(cell.fill, openmc.Iterable): - cell.fill = list(map(lambda x: material if x.name == material_name else x, cell.fill)) - else: - cell.fill = material if cell.fill.name == material_name else cell.fill + raise AttributeError( + "DAGMCUniverse.material_overrides has been removed. Use " + "DAGMCCell objects added via add_cell() to manage per-cell " + "material assignments.") def add_material_override(self, key, overrides=None): """Add a material override to the universe. @@ -201,7 +152,10 @@ def add_material_override(self, key, overrides=None): if key not in self.cells: raise ValueError(f"Cell ID '{key}' not found in DAGMC universe") - self._material_overrides[key] = overrides + if len(overrides) == 1: + self.cells[key].fill = overrides[0] + else: + self.cells[key].fill = list(overrides) @property def auto_geom_ids(self): @@ -290,12 +244,6 @@ def create_xml_subelement(self, xml_element, memo=None): memo.add(self) - # Ensure that the material overrides are up-to-date - for cell in self.cells.values(): - if cell.fill is None: - continue - self.add_material_override(cell, cell.fill) - # Set xml element values dagmc_element = ET.Element('dagmc_universe') dagmc_element.set('id', str(self.id)) @@ -307,17 +255,10 @@ def create_xml_subelement(self, xml_element, memo=None): if self.auto_mat_ids: dagmc_element.set('auto_mat_ids', 'true') dagmc_element.set('filename', str(self.filename)) - if self._material_overrides: - mat_element = ET.Element('material_overrides') - for key in self._material_overrides: - cell_overrides = ET.Element('cell_override') - cell_overrides.set("id", str(key)) - material_element = ET.Element('material_ids') - material_element.text = ' '.join( - str(t.id) for t in self._material_overrides[key]) - cell_overrides.append(material_element) - mat_element.append(cell_overrides) - dagmc_element.append(mat_element) + if self.cells: + for cell in self.cells.values(): + cell_element = cell.create_xml_subelement(xml_element, memo) + dagmc_element.append(cell_element) xml_element.append(dagmc_element) def bounding_region( @@ -442,7 +383,7 @@ def from_hdf5(cls, group): return out @classmethod - def from_xml_element(cls, elem, mats = None): + def from_xml_element(cls, elem, mats=None): """Generate DAGMC universe from XML element Parameters @@ -471,21 +412,56 @@ def from_xml_element(cls, elem, mats = None): out.auto_geom_ids = bool(get_text(elem, "auto_geom_ids")) out.auto_mat_ids = bool(get_text(elem, "auto_mat_ids")) - el_mat_override = elem.find('material_overrides') - if el_mat_override is not None: - if mats is None: - raise ValueError("Material overrides found in DAGMC universe " - "but no materials were provided to populate " - "the mapping.") - out._material_overrides = {} - for elem in el_mat_override.findall('cell_override'): - cell_id = int(get_text(elem, 'id')) - mat_ids = get_elem_list(elem, "material_ids", str) or [] - mat_objs = [mats[mat_id] for mat_id in mat_ids] - out._material_overrides[cell_id] = mat_objs + has_overrides = elem.find('material_overrides') is not None + has_cells = elem.find('cell') is not None + + if has_overrides and has_cells: + raise ValueError( + "DAGMCUniverse cannot specify both and " + " sub-elements. Use elements only.") + + if has_overrides: + warnings.warn( + "DAGMCUniverse is deprecated and will be " + "removed in a future version. Use nested elements " + "instead.", DeprecationWarning, stacklevel=2) + out._parse_legacy_material_overrides(elem, mats) + elif has_cells: + out._parse_cell_overrides(elem, mats) return out + def _parse_legacy_material_overrides(self, elem, mats): + """Parse the deprecated XML format and populate + the universe with equivalent DAGMCCell objects.""" + if mats is None: + raise ValueError( + "DAGMC material overrides found but no materials were " + "provided to populate the mapping.") + mo_elem = elem.find('material_overrides') + for co_elem in mo_elem.findall('cell_override'): + cell_id = int(get_text(co_elem, 'id')) + mat_ids = co_elem.find('material_ids').text.split() + fill_objs = [mats[mid] for mid in mat_ids] + fill = fill_objs[0] if len(fill_objs) == 1 else fill_objs + if cell_id in self.cells: + raise ValueError( + f"Duplicate DAGMC cell override specified for cell {cell_id}.") + self.add_cell(DAGMCCell(cell_id=cell_id, fill=fill)) + + def _parse_cell_overrides(self, elem, mats): + if mats is None: + raise ValueError("DAGMC cell overrides found in DAGMC universe but " + "no materials were provided to populate the " + "mapping.") + + for cell_elem in elem.findall('cell'): + cell_id = int(get_text(cell_elem, 'id')) + if cell_id in self.cells: + raise ValueError( + f"Duplicate DAGMC cell override specified for cell {cell_id}.") + DAGMCCell.from_xml_element(cell_elem, mats, self) + def _partial_deepcopy(self): """Clone all of the openmc.DAGMCUniverse object's attributes except for its cells, as they are copied within the clone function. This should @@ -565,7 +541,13 @@ def sync_dagmc_cells(self, mats: Iterable[openmc.Material]): fill = [mats_per_id[mat.id] for mat in dag_cell.fill if mat] else: fill = mats_per_id[dag_cell.fill.id] if dag_cell.fill else None - self.add_cell(openmc.DAGMCCell(cell_id=dag_cell_id, fill=fill)) + name = dag_cell.name + if dag_cell_id in self._cells: + self._cells[dag_cell_id].name = name + self._cells[dag_cell_id].fill = fill + else: + self.add_cell( + openmc.DAGMCCell(cell_id=dag_cell_id, name=name, fill=fill)) @add_plot_params def plot(self, *args, **kwargs): @@ -594,6 +576,14 @@ class DAGMCCell(openmc.Cell): DAG_parent_universe : int The parent universe of the cell. + Notes + ----- + DAGMC geometries are composed of triangulated surfaces, which means cell + volumes can in principle be computed exactly (e.g. via mesh-based + integration). Manually specifying :attr:`volume` overrides any such + calculation and may introduce inconsistencies if the value does not + accurately reflect the true geometric volume. + """ def __init__(self, cell_id=None, name='', fill=None): super().__init__(cell_id, name, fill, None) @@ -625,8 +615,66 @@ def plot(self, *args, **kwargs): raise TypeError("plot is not available for DAGMC cells.") def create_xml_subelement(self, xml_element, memo=None): - raise TypeError("create_xml_subelement is not available for DAGMC cells.") + if self.fill_type not in ('void', 'material', 'distribmat'): + raise TypeError("DAGMC cell overrides currently only support " + "material fills.") + if self.temperature is not None and self.fill_type not in ( + 'material', 'distribmat' + ): + raise TypeError("DAGMC cell temperature overrides require a " + "material fill.") + if self.density is not None and self.fill_type not in ( + 'material', 'distribmat' + ): + raise TypeError("DAGMC cell density overrides require a " + "material fill.") + if any(getattr(self, attr) is not None for attr in ( + 'translation', 'rotation' + )): + raise TypeError("DAGMC cell overrides do not support translation " + "or rotation.") + return super().create_xml_subelement(xml_element, memo) @classmethod - def from_xml_element(cls, elem, surfaces, materials, get_universe): - raise TypeError("from_xml_element is not available for DAGMC cells.") + def from_xml_element(cls, elem, mats, universe): + """Generate a DAGMCCell from an XML override element. + + Parameters + ---------- + elem : lxml.etree._Element + `` element containing a DAGMC cell property override + mats : dict + Dictionary mapping material ID strings to + :class:`openmc.Material` instances + universe : DAGMCUniverse + Universe to add the parsed cell to. + + Returns + ------- + DAGMCCell + DAGMCCell instance + """ + if not isinstance(universe, DAGMCUniverse): + raise TypeError( + f"universe must be a DAGMCUniverse instance, " + f"got {type(universe).__name__}.") + + cell_id = int(get_text(elem, 'id')) + + # Validate attributes that are unsupported for DAGMC cell overrides + for tag in ('region', 'fill', 'universe'): + if get_text(elem, tag) is not None: + raise ValueError( + f"DAGMC cell {cell_id} override cannot specify '{tag}'.") + for tag in ('translation', 'rotation'): + if get_text(elem, tag) is not None: + raise ValueError( + f"DAGMC cell {cell_id} override does not support " + f"'{tag}'.") + if get_elem_list(elem, 'material', str) is None: + raise ValueError( + f"DAGMC cell {cell_id} must specify a material override.") + + return super().from_xml_element( + elem, surfaces={}, materials=mats, + get_universe=lambda _: universe) diff --git a/openmc/model/model.py b/openmc/model/model.py index 884e3ff6f95..1b66628379f 100644 --- a/openmc/model/model.py +++ b/openmc/model/model.py @@ -427,6 +427,8 @@ def sync_dagmc_universes(self): This method iterates over all DAGMC universes in the geometry and synchronizes their cells with the current material assignments. Requires that the model has been initialized via :meth:`Model.init_lib`. + Synchronized DAGMC cells can then be edited and exported as nested + `` overrides inside each `` element. .. versionadded:: 0.15.1 diff --git a/src/cell.cpp b/src/cell.cpp index ebe28c3d2ce..06bc29a6629 100644 --- a/src/cell.cpp +++ b/src/cell.cpp @@ -340,6 +340,63 @@ void Cell::to_hdf5(hid_t cell_group) const close_group(group); } +//============================================================================== +// XML parsing helpers for nodes +//============================================================================== + +vector parse_cell_material_xml(pugi::xml_node node, int32_t cell_id) +{ + vector mats { + get_node_array(node, "material", true)}; + if (mats.empty()) { + fatal_error(fmt::format( + "An empty material element was specified for cell {}", cell_id)); + } + vector material; + material.reserve(mats.size()); + for (const auto& mat : mats) { + if (mat == "void") { + material.push_back(MATERIAL_VOID); + } else { + material.push_back(std::stoi(mat)); + } + } + return material; +} + +vector parse_cell_temperature_xml(pugi::xml_node node, int32_t cell_id) +{ + auto temperatures = get_node_array(node, "temperature"); + if (temperatures.empty()) { + fatal_error(fmt::format( + "An empty temperature element was specified for cell {}", cell_id)); + } + for (auto T : temperatures) { + if (T < 0) { + fatal_error(fmt::format( + "Cell {} was specified with a negative temperature", cell_id)); + } + } + return temperatures; +} + +vector parse_cell_density_xml(pugi::xml_node node, int32_t cell_id) +{ + auto densities = get_node_array(node, "density"); + if (densities.empty()) { + fatal_error(fmt::format( + "An empty density element was specified for cell {}", cell_id)); + } + for (auto rho : densities) { + if (rho <= 0) { + fatal_error(fmt::format( + "Cell {} was specified with a density less than or equal to zero", + cell_id)); + } + } + return densities; +} + //============================================================================== // CSGCell implementation //============================================================================== @@ -390,26 +447,12 @@ CSGCell::CSGCell(pugi::xml_node cell_node) // universe), more than one material (distribmats), and some materials may // be "void". if (material_present) { - vector mats { - get_node_array(cell_node, "material", true)}; - if (mats.size() > 0) { - material_.reserve(mats.size()); - for (std::string mat : mats) { - if (mat.compare("void") == 0) { - material_.push_back(MATERIAL_VOID); - } else { - material_.push_back(std::stoi(mat)); - } - } - } else { - fatal_error(fmt::format( - "An empty material element was specified for cell {}", id_)); - } + material_ = parse_cell_material_xml(cell_node, id_); } // Read the temperature element which may be distributed like materials. if (check_for_node(cell_node, "temperature")) { - sqrtkT_ = get_node_array(cell_node, "temperature"); + sqrtkT_ = parse_cell_temperature_xml(cell_node, id_); sqrtkT_.shrink_to_fit(); // Make sure this is a material-filled cell. @@ -420,14 +463,6 @@ CSGCell::CSGCell(pugi::xml_node cell_node) id_)); } - // Make sure all temperatures are non-negative. - for (auto T : sqrtkT_) { - if (T < 0) { - fatal_error(fmt::format( - "Cell {} was specified with a negative temperature", id_)); - } - } - // Convert to sqrt(k*T). for (auto& T : sqrtkT_) { T = std::sqrt(K_BOLTZMANN * T); @@ -440,7 +475,7 @@ CSGCell::CSGCell(pugi::xml_node cell_node) // Note: calculating the actual density multiplier is deferred until materials // are finalized. density_mult_ contains the true density in the meantime. if (check_for_node(cell_node, "density")) { - density_mult_ = get_node_array(cell_node, "density"); + density_mult_ = parse_cell_density_xml(cell_node, id_); density_mult_.shrink_to_fit(); // Make sure this is a material-filled cell. @@ -461,15 +496,6 @@ CSGCell::CSGCell(pugi::xml_node cell_node) id_)); } } - - // Make sure all densities are non-negative and greater than zero. - for (auto rho : density_mult_) { - if (rho <= 0) { - fatal_error(fmt::format( - "Cell {} was specified with a density less than or equal to zero", - id_)); - } - } } // Read the region specification. diff --git a/src/dagmc.cpp b/src/dagmc.cpp index 571182fa6b2..a5cc71b5230 100644 --- a/src/dagmc.cpp +++ b/src/dagmc.cpp @@ -50,6 +50,10 @@ namespace openmc { DAGUniverse::DAGUniverse(pugi::xml_node node) { + MaterialOverrides material_overrides; + TemperatureOverrides temperature_overrides; + DensityOverrides density_overrides; + if (check_for_node(node, "id")) { id_ = std::stoi(get_node_value(node, "id")); } else { @@ -76,24 +80,77 @@ DAGUniverse::DAGUniverse(pugi::xml_node node) adjust_material_ids_ = get_node_value_bool(node, "auto_mat_ids"); } - // get material assignment overloading - if (check_for_node(node, "material_overrides")) { - auto mat_node = node.child("material_overrides"); - // loop over all subelements (each subelement corresponds to a material) - for (pugi::xml_node cell_node : mat_node.children("cell_override")) { - // Store assignment reference name - int32_t ref_assignment = std::stoi(get_node_value(cell_node, "id")); + // Get material assignment overrides from nested DAGMC cell elements. + if (node.child("cell")) { + for (pugi::xml_node cell_node : node.children("cell")) { + if (!check_for_node(cell_node, "id")) { + fatal_error( + "Must specify id for each DAGMC cell override in ."); + } - // Get mat name for each assignement instances - vector instance_mats = - get_node_array(cell_node, "material_ids"); + int32_t cell_id = std::stoi(get_node_value(cell_node, "id")); - // Store mat name for each instances - material_overrides_.emplace(ref_assignment, instance_mats); + if (check_for_node(cell_node, "region")) { + fatal_error(fmt::format( + "DAGMC cell {} override cannot specify a region.", cell_id)); + } + if (check_for_node(cell_node, "fill")) { + fatal_error(fmt::format( + "DAGMC cell {} override currently only supports material fills.", + cell_id)); + } + if (check_for_node(cell_node, "universe")) { + fatal_error(fmt::format( + "DAGMC cell {} override cannot specify a universe.", cell_id)); + } + if (check_for_node(cell_node, "translation") || + check_for_node(cell_node, "rotation")) { + fatal_error(fmt::format( + "DAGMC cell {} override does not support translation or rotation.", + cell_id)); + } + if (!check_for_node(cell_node, "material")) { + fatal_error(fmt::format( + "DAGMC cell {} override must specify material.", cell_id)); + } + + auto inserted = material_overrides.emplace( + cell_id, parse_cell_material_xml(cell_node, cell_id)); + if (!inserted.second) { + fatal_error(fmt::format( + "Duplicate DAGMC cell override specified for cell {}", cell_id)); + } + + if (check_for_node(cell_node, "temperature")) { + temperature_overrides.emplace( + cell_id, parse_cell_temperature_xml(cell_node, cell_id)); + } + + if (check_for_node(cell_node, "density")) { + density_overrides.emplace( + cell_id, parse_cell_density_xml(cell_node, cell_id)); + } + } + } else if (check_for_node(node, "material_overrides")) { + if (node.child("cell")) { + fatal_error("DAGMCUniverse cannot specify both and " + " sub-elements. Use elements only."); + } + warning("DAGMCUniverse is deprecated. Use nested " + " elements under instead."); + for (pugi::xml_node co : + node.child("material_overrides").children("cell_override")) { + int32_t cell_id = std::stoi(get_node_value(co, "id")); + std::istringstream iss(co.child("material_ids").text().get()); + vector mats; + for (std::string s; iss >> s;) { + mats.push_back(s == "void" ? MATERIAL_VOID : std::stoi(s)); + } + material_overrides.emplace(cell_id, mats); } } - initialize(); + initialize(material_overrides, temperature_overrides, density_overrides); } DAGUniverse::DAGUniverse( @@ -110,9 +167,12 @@ DAGUniverse::DAGUniverse(std::shared_ptr dagmc_ptr, : dagmc_instance_(dagmc_ptr), filename_(filename), adjust_geometry_ids_(auto_geom_ids), adjust_material_ids_(auto_mat_ids) { + MaterialOverrides material_overrides; + TemperatureOverrides temperature_overrides; + DensityOverrides density_overrides; set_id(); init_metadata(); - init_geometry(); + init_geometry(material_overrides, temperature_overrides, density_overrides); } void DAGUniverse::set_id() @@ -130,6 +190,15 @@ void DAGUniverse::set_id() } void DAGUniverse::initialize() +{ + MaterialOverrides material_overrides; + TemperatureOverrides temperature_overrides; + initialize(material_overrides, temperature_overrides); +} + +void DAGUniverse::initialize(const MaterialOverrides& material_overrides, + const TemperatureOverrides& temperature_overrides, + const DensityOverrides& density_overrides) { #ifdef OPENMC_UWUW_ENABLED // read uwuw materials from the .h5m file if present @@ -140,7 +209,7 @@ void DAGUniverse::initialize() init_metadata(); - init_geometry(); + init_geometry(material_overrides, temperature_overrides, density_overrides); } void DAGUniverse::init_dagmc() @@ -176,7 +245,9 @@ void DAGUniverse::init_metadata() MB_CHK_ERR_CONT(rval); } -void DAGUniverse::init_geometry() +void DAGUniverse::init_geometry(const MaterialOverrides& material_overrides, + const TemperatureOverrides& temperature_overrides, + const DensityOverrides& density_overrides) { moab::ErrorCode rval; @@ -202,6 +273,9 @@ void DAGUniverse::init_geometry() : dagmc_instance_->id_by_index(3, c->dag_index()); c->universe_ = this->id_; c->fill_ = C_NONE; // no fill, single universe + if (dagmc_instance_->is_implicit_complement(vol_handle)) { + c->name_ = "implicit complement"; + } auto in_map = model::cell_map.find(c->id_); if (in_map == model::cell_map.end()) { @@ -230,16 +304,68 @@ void DAGUniverse::init_geometry() if (mat_str == "graveyard") { graveyard = vol_handle; } - // material void checks - if (mat_str == "void" || mat_str == "vacuum" || mat_str == "graveyard") { + if (material_overrides.count(c->id_)) { + override_assign_material(c, material_overrides); + } else if (mat_str == "void" || mat_str == "vacuum" || + mat_str == "graveyard") { c->material_.push_back(MATERIAL_VOID); + } else if (uses_uwuw()) { + uwuw_assign_material(vol_handle, c); } else { - if (material_overrides_.count(c->id_)) { - override_assign_material(c); - } else if (uses_uwuw()) { - uwuw_assign_material(vol_handle, c); - } else { - legacy_assign_material(mat_str, c); + legacy_assign_material(mat_str, c); + } + + if (temperature_overrides.count(c->id_)) { + if (c->material_.empty() || c->material_[0] == MATERIAL_VOID) { + fatal_error(fmt::format("DAGMC cell {} was specified with a " + "temperature but no non-void material.", + c->id_)); + } + + c->sqrtkT_.clear(); + const auto& temp_overrides = temperature_overrides.at(c->id_); + c->sqrtkT_.reserve(temp_overrides.size()); + for (auto T : temp_overrides) { + c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * T)); + } + + if (settings::verbosity >= 10) { + std::stringstream override_values; + for (size_t i = 0; i < temp_overrides.size(); ++i) { + if (i > 0) { + override_values << " "; + } + override_values << temp_overrides[i]; + } + auto msg = fmt::format("Overriding DAGMC cell {} property " + "'temperature [K]' with value(s): {}", + c->id_, override_values.str()); + write_message(msg, 10); + } + } + + if (density_overrides.count(c->id_)) { + if (c->material_.empty() || c->material_[0] == MATERIAL_VOID) { + fatal_error(fmt::format("DAGMC cell {} was specified with a density " + "but no non-void material.", + c->id_)); + } + // density_mult_ holds the true density until materials are finalized, + // at which point it is converted to a proper multiplier (same as CSG). + c->density_mult_ = density_overrides.at(c->id_); + + if (settings::verbosity >= 10) { + const auto& dens = density_overrides.at(c->id_); + std::stringstream override_values; + for (size_t i = 0; i < dens.size(); ++i) { + if (i > 0) + override_values << " "; + override_values << dens[i]; + } + write_message(fmt::format("Overriding DAGMC cell {} property " + "'density [g/cm³]' with value(s): {}", + c->id_, override_values.str()), + 10); } } @@ -252,18 +378,21 @@ void DAGUniverse::init_geometry() continue; } - // assign cell temperature - const auto& mat = model::materials[model::material_map.at(c->material_[0])]; - if (dagmc_instance_->has_prop(vol_handle, "temp")) { - rval = dagmc_instance_->prop_value(vol_handle, "temp", temp_value); - MB_CHK_ERR_CONT(rval); - double temp = std::stod(temp_value); - c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * temp)); - } else if (mat->temperature() > 0.0) { - c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * mat->temperature())); - } else { - c->sqrtkT_.push_back( - std::sqrt(K_BOLTZMANN * settings::temperature_default)); + // assign cell temperature if not explicitly overridden + if (c->sqrtkT_.empty()) { + const auto& mat = + model::materials[model::material_map.at(c->material_[0])]; + if (dagmc_instance_->has_prop(vol_handle, "temp")) { + rval = dagmc_instance_->prop_value(vol_handle, "temp", temp_value); + MB_CHK_ERR_CONT(rval); + double temp = std::stod(temp_value); + c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * temp)); + } else if (mat->temperature() > 0.0) { + c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * mat->temperature())); + } else { + c->sqrtkT_.push_back( + std::sqrt(K_BOLTZMANN * settings::temperature_default)); + } } model::cells.emplace_back(std::move(c)); @@ -630,7 +759,8 @@ void DAGUniverse::uwuw_assign_material( #endif // OPENMC_UWUW_ENABLED } -void DAGUniverse::override_assign_material(std::unique_ptr& c) const +void DAGUniverse::override_assign_material(std::unique_ptr& c, + const MaterialOverrides& material_overrides) const { // if Cell ID matches an override key, use it to override the material // assignment else if UWUW is used, get the material assignment from the DAGMC @@ -638,17 +768,30 @@ void DAGUniverse::override_assign_material(std::unique_ptr& c) const // Notify User that an override is being applied on a DAGMCCell write_message(fmt::format("Applying override for DAGMCCell {}", c->id_), 8); + const auto& mat_overrides = material_overrides.at(c->id_); if (settings::verbosity >= 10) { - auto msg = fmt::format("Assigning DAGMC cell {} material(s) based on " - "override information (see input XML).", - c->id_); + std::stringstream override_values; + for (size_t i = 0; i < mat_overrides.size(); ++i) { + if (i > 0) { + override_values << " "; + } + if (mat_overrides[i] == MATERIAL_VOID) { + override_values << "void"; + } else { + override_values << mat_overrides[i]; + } + } + auto msg = fmt::format("Overriding DAGMC cell {} property 'material' " + "with value(s): {}", + c->id_, override_values.str()); write_message(msg, 10); } // Override the material assignment for each cell instance using the legacy // assignement - for (auto mat_id : material_overrides_.at(c->id_)) { - if (model::material_map.find(mat_id) == model::material_map.end()) { + for (auto mat_id : mat_overrides) { + if (mat_id != MATERIAL_VOID && + model::material_map.find(mat_id) == model::material_map.end()) { fatal_error(fmt::format( "Material with ID '{}' not found for DAGMC cell {}", mat_id, c->id_)); } diff --git a/tests/unit_tests/dagmc/test_model.py b/tests/unit_tests/dagmc/test_model.py index 0917f5b237c..5889fa5420f 100644 --- a/tests/unit_tests/dagmc/test_model.py +++ b/tests/unit_tests/dagmc/test_model.py @@ -70,28 +70,20 @@ def model(request): openmc.reset_auto_ids() -def test_dagmc_replace_material_assignment(model): - mats = {} - - mats["foo"] = openmc.Material(name="foo") - mats["foo"].add_nuclide("H1", 2.0) - mats["foo"].add_element("O", 1.0) - mats["foo"].set_density("g/cm3", 1.0) - mats["foo"].add_s_alpha_beta("c_H_in_H2O") - +def test_dagmc_sync_cell_names(model): + dag_univ = None for univ in model.geometry.get_all_universes().values(): - if not isinstance(univ, openmc.DAGMCUniverse): + if isinstance(univ, openmc.DAGMCUniverse): + dag_univ = univ break - cells_with_41 = [] - for cell in univ.cells.values(): - if cell.fill is None: - continue - if cell.fill.name == "41": - cells_with_41.append(cell.id) - univ.replace_material_assignment("41", mats["foo"]) - for cell_id in cells_with_41: - assert univ.cells[cell_id] == mats["foo"] + assert dag_univ is not None + + for cell_id, cell in dag_univ.cells.items(): + assert cell.name == openmc.lib.cells[cell_id].name + + assert any(cell.name == "implicit complement" + for cell in dag_univ.cells.values()) def test_dagmc_add_material_override_with_id(model): @@ -114,7 +106,7 @@ def test_dagmc_add_material_override_with_id(model): cells_with_41.append(cell.id) univ.add_material_override(cell.id, mats["foo"]) for cell_id in cells_with_41: - assert univ.cells[cell_id] == mats["foo"] + assert univ.cells[cell_id].fill == mats["foo"] def test_dagmc_add_material_override_with_cell(model): @@ -137,7 +129,7 @@ def test_dagmc_add_material_override_with_cell(model): cells_with_41.append(cell.id) univ.add_material_override(cell, mats["foo"]) for cell_id in cells_with_41: - assert univ.cells[cell_id] == mats["foo"] + assert univ.cells[cell_id].fill == mats["foo"] def test_model_differentiate_depletable_with_dagmc(model, run_in_tmpdir): @@ -174,57 +166,20 @@ def test_model_differentiate_with_dagmc(model): assert len(model.materials) == 4*2 + 4 -def test_bad_override_cell_id(model): - for univ in model.geometry.get_all_universes().values(): - if isinstance(univ, openmc.DAGMCUniverse): - break - with pytest.raises(ValueError, match="Cell ID '1' not found in DAGMC universe"): - univ.material_overrides = {1: model.materials[0]} - - -def test_bad_override_type(model): - not_a_dag_cell = openmc.Cell() - for univ in model.geometry.get_all_universes().values(): - if isinstance(univ, openmc.DAGMCUniverse): - break - with pytest.raises(ValueError, match="Unrecognized key type. Must be an integer or openmc.DAGMCCell object"): - univ.material_overrides = {not_a_dag_cell: model.materials[0]} - - -def test_bad_replacement_mat_name(model): - for univ in model.geometry.get_all_universes().values(): - if isinstance(univ, openmc.DAGMCUniverse): - break - with pytest.raises(ValueError, match="No material with name 'not_a_mat' found in the DAGMC universe"): - univ.replace_material_assignment("not_a_mat", model.materials[0]) - - def test_dagmc_xml(model): - # Set the environment - mats = {} - mats["no-void fuel"] = openmc.Material(1, name="no-void fuel") - mats["no-void fuel"].add_nuclide("U235", 0.03) - mats["no-void fuel"].add_nuclide("U238", 0.97) - mats["no-void fuel"].add_nuclide("O16", 2.0) - mats["no-void fuel"].set_density("g/cm3", 10.0) - - mats[5] = openmc.Material(name="41") - mats[5].add_nuclide("H1", 2.0) - mats[5].add_element("O", 1.0) - mats[5].set_density("g/cm3", 1.0) - mats[5].add_s_alpha_beta("c_H_in_H2O") + override_mat = openmc.Material(name="41") + override_mat.add_nuclide("H1", 2.0) + override_mat.add_element("O", 1.0) + override_mat.set_density("g/cm3", 1.0) + override_mat.add_s_alpha_beta("c_H_in_H2O") + model.materials.append(override_mat) for univ in model.geometry.get_all_universes().values(): if isinstance(univ, openmc.DAGMCUniverse): dag_univ = univ break - for k, v in mats.items(): - if isinstance(k, int): - dag_univ.add_material_override(k, v) - model.materials.append(v) - elif isinstance(k, str): - dag_univ.replace_material_assignment(k, v) + dag_univ.add_material_override(5, override_mat) # Tesing the XML subelement generation root = ET.Element('dagmc_universe') @@ -236,12 +191,24 @@ def test_dagmc_xml(model): assert dagmc_ele.get('filename') == str(dag_univ.filename) assert dagmc_ele.get('auto_geom_ids') == str(dag_univ.auto_geom_ids).lower() - override_eles = dagmc_ele.find('material_overrides').findall('cell_override') - assert len(override_eles) == 4 - - for i, override_ele in enumerate(override_eles): - cell_id = override_ele.get('id') - assert dag_univ.material_overrides[int(cell_id)][0].id == int(override_ele.find('material_ids').text) + assert dagmc_ele.find('material_overrides') is None + + override_elements = dagmc_ele.findall('cell') + assert len(override_elements) == len(dag_univ.cells) + xml_cells = {int(elem.get('id')): elem for elem in override_elements} + for cell_id, cell in dag_univ.cells.items(): + assert cell_id in xml_cells + xml_cell = xml_cells[cell_id] + if cell.fill_type == 'void': + assert xml_cell.get('material') == 'void' + elif cell.fill_type == 'material': + assert xml_cell.get('material') == str(cell.fill.id) + elif cell.fill_type == 'distribmat': + mat_list = xml_cell.find('material').text.split() + expected = ["void" if m is None else str(m.id) for m in cell.fill] + assert mat_list == expected + else: + pytest.fail(f"Unexpected DAGMC cell fill type: {cell.fill_type}") model.export_to_model_xml() @@ -252,7 +219,147 @@ def test_dagmc_xml(model): xml_dagmc_univ = univ break - assert xml_dagmc_univ._material_overrides.keys() == dag_univ._material_overrides.keys() + assert xml_dagmc_univ.cells.keys() == dag_univ.cells.keys() + + for cell_id, cell in dag_univ.cells.items(): + xml_cell = xml_dagmc_univ.cells[cell_id] + assert xml_cell.fill_type == cell.fill_type + if cell.fill_type == 'void': + assert xml_cell.fill is None + elif cell.fill_type == 'material': + assert xml_cell.fill.id == cell.fill.id + elif cell.fill_type == 'distribmat': + xml_ids = [m.id if m is not None else None for m in xml_cell.fill] + model_ids = [m.id if m is not None else None for m in cell.fill] + assert xml_ids == model_ids + else: + pytest.fail(f"Unexpected DAGMC cell fill type: {cell.fill_type}") + + +def test_dagmc_xml_reject_fill_override(): + mats = {'1': openmc.Material(1), 'void': None} + elem = ET.fromstring( + '' + '' + '' + ) + with pytest.raises(ValueError, match="cannot specify 'fill'"): + openmc.DAGMCUniverse.from_xml_element(elem, mats) + + +def test_dagmc_xml_reject_region_override(): + mats = {'1': openmc.Material(1), 'void': None} + elem = ET.fromstring( + '' + '' + '' + ) + with pytest.raises(ValueError, match="cannot specify 'region'"): + openmc.DAGMCUniverse.from_xml_element(elem, mats) + + +def _legacy_xml(cell_overrides): + """Helper to build a with old-format .""" + inner = ''.join( + f'{mids}' + for cid, mids in cell_overrides.items() + ) + return ET.fromstring( + f'' + f'{inner}' + f'' + ) + + +def test_dagmc_xml_legacy_single_material_compat(): + mat = openmc.Material(1) + mats = {'1': mat, 'void': None} + elem = _legacy_xml({3: '1'}) + with pytest.warns(DeprecationWarning, match="deprecated"): + univ = openmc.DAGMCUniverse.from_xml_element(elem, mats) + assert 3 in univ.cells + assert univ.cells[3].fill is mat + + +def test_dagmc_xml_legacy_distribmat_compat(): + mat1, mat2 = openmc.Material(2), openmc.Material(3) + mats = {'2': mat1, '3': mat2, 'void': None} + elem = _legacy_xml({5: '2 3'}) + with pytest.warns(DeprecationWarning): + univ = openmc.DAGMCUniverse.from_xml_element(elem, mats) + assert univ.cells[5].fill_type == 'distribmat' + assert list(univ.cells[5].fill) == [mat1, mat2] + + +def test_dagmc_xml_legacy_void_compat(): + mats = {'void': None} + elem = _legacy_xml({7: 'void'}) + with pytest.warns(DeprecationWarning): + univ = openmc.DAGMCUniverse.from_xml_element(elem, mats) + assert univ.cells[7].fill_type == 'void' + + +def test_dagmc_xml_legacy_both_raises(): + mat = openmc.Material(1) + mats = {'1': mat, 'void': None} + elem = ET.fromstring( + '' + '' + '1' + '' + '' + '' + ) + with pytest.raises(ValueError, match="both"): + openmc.DAGMCUniverse.from_xml_element(elem, mats) + + +def test_dagmc_xml_legacy_deprecation_warning(): + mats = {'1': openmc.Material(1), 'void': None} + elem = _legacy_xml({3: '1'}) + with pytest.warns(DeprecationWarning): + openmc.DAGMCUniverse.from_xml_element(elem, mats) + + +def test_dagmc_xml_legacy_roundtrip(): + """Old-format XML loads correctly and re-exports using the new format.""" + mat = openmc.Material(1) + mats = {'1': mat, 'void': None} + elem = _legacy_xml({3: '1'}) + with pytest.warns(DeprecationWarning): + univ = openmc.DAGMCUniverse.from_xml_element(elem, mats) + + root = ET.Element('geometry') + univ.create_xml_subelement(root) + dagmc_elem = root.find('dagmc_universe') + + assert dagmc_elem.find('material_overrides') is None + cell_elems = dagmc_elem.findall('cell') + assert len(cell_elems) == 1 + assert int(cell_elems[0].get('id')) == 3 + assert cell_elems[0].get('material') == '1' + + +def test_dagmc_xml_temperature_roundtrip(): + mat = openmc.Material(1) + mats = {'1': mat, 'void': None} + + elem = ET.fromstring( + '' + '' + '' + ) + + dag_univ = openmc.DAGMCUniverse.from_xml_element(elem, mats) + assert dag_univ.cells[7].fill.id == 1 + assert dag_univ.cells[7].temperature == pytest.approx(825.0) + + root = ET.Element('geometry') + dag_univ.create_xml_subelement(root) + dagmc_elem = root.find('dagmc_universe') + xml_cell = dagmc_elem.find('cell') + assert xml_cell.get('temperature') == '825.0' - for xml_mats, model_mats in zip(xml_dagmc_univ._material_overrides.values(), dag_univ._material_overrides.values()): - assert all([xml_mat.id == orig_mat.id for xml_mat, orig_mat in zip(xml_mats, model_mats)]) + dag_univ_roundtrip = openmc.DAGMCUniverse.from_xml_element(dagmc_elem, mats) + assert dag_univ_roundtrip.cells[7].fill.id == 1 + assert dag_univ_roundtrip.cells[7].temperature == pytest.approx(825.0) diff --git a/tests/unit_tests/test_cell.py b/tests/unit_tests/test_cell.py index 2baa59f1bfc..8d39c071b7b 100644 --- a/tests/unit_tests/test_cell.py +++ b/tests/unit_tests/test_cell.py @@ -374,6 +374,49 @@ def test_rotation_from_xml(rotation): np.testing.assert_allclose(new_cell.rotation, cell.rotation) +def test_dagmccell_from_xml_element(): + """DAGMCCell.from_xml_element parses material, temperature, density, + and volume; rejects unsupported attributes.""" + mat = openmc.Material(1) + mat.add_nuclide('U235', 1.0) + mats = {'1': mat} + # In practice, from_xml_element is always called from DAGMCUniverse + # during XML parsing. A placeholder universe is used here so the test + # can exercise the method directly without a real DAGMC model file. + placeholder_univ = openmc.DAGMCUniverse('model.h5m') + + # material + temperature + density round-trip + xml = '' + cell = openmc.DAGMCCell.from_xml_element(ET.fromstring(xml), mats, placeholder_univ) + assert cell.id == 5 + assert cell.name == 'fuel' + assert cell.fill is mat + assert cell.density == 10.5 + assert cell.temperature == 900.0 + + # volume round-trip + xml = '' + cell = openmc.DAGMCCell.from_xml_element(ET.fromstring(xml), mats, placeholder_univ) + assert cell.volume == 42.0 + + # forbidden: region, fill, universe + for tag, val in [('region', '-1'), ('fill', '2'), ('universe', '0')]: + xml = f'' + with pytest.raises(ValueError, match=tag): + openmc.DAGMCCell.from_xml_element(ET.fromstring(xml), mats, placeholder_univ) + + # forbidden: translation, rotation + for tag, val in [('translation', '1 0 0'), ('rotation', '0 0 90')]: + xml = f'' + with pytest.raises(ValueError, match=tag): + openmc.DAGMCCell.from_xml_element(ET.fromstring(xml), mats, placeholder_univ) + + # missing material raises + xml = '' + with pytest.raises(ValueError, match='material'): + openmc.DAGMCCell.from_xml_element(ET.fromstring(xml), mats, placeholder_univ) + + def test_plot(run_in_tmpdir): zcyl = openmc.ZCylinder() c = openmc.Cell(region=-zcyl)