Skip to content

Commit ede21b4

Browse files
committed
Some reverting and small fix in filter.py
1 parent 75dcb51 commit ede21b4

3 files changed

Lines changed: 44 additions & 23 deletions

File tree

openmc/filter.py

Lines changed: 16 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -37,7 +37,6 @@
3737
)
3838

3939

40-
4140
class FilterMeta(ABCMeta):
4241
"""Metaclass for filters that ensures class names are appropriate."""
4342

@@ -1113,12 +1112,20 @@ def from_volumes(cls, mesh: openmc.MeshBase, volumes: openmc.MeshMaterialVolumes
11131112
A new MeshMaterialFilter instance
11141113
11151114
"""
1116-
# Get flat arrays of material IDs and element indices
1117-
mat_ids = volumes._materials[volumes._materials > -1]
1118-
elems, _ = np.where(volumes._materials > -1)
1119-
1120-
# Stack them into a 2D array of (element, material) pairs
1121-
bins = np.column_stack((elems, mat_ids))
1115+
# Build bins sorted by volume ascending (mat ID ascending as tiebreaker)
1116+
# so the bin ordering is deterministic regardless of the ray-traversal
1117+
# order stored in the raw _materials array.
1118+
bins = []
1119+
for i in range(volumes.num_elements):
1120+
entries = volumes.by_element(i)
1121+
# sort by (volume asc, mat_id asc), excluding void (None)
1122+
entries = sorted(
1123+
((mat_id, vol)
1124+
for mat_id, vol in entries if mat_id is not None),
1125+
key=lambda t: (t[1], t[0])
1126+
)
1127+
for mat_id, _ in entries:
1128+
bins.append((i, mat_id))
11221129
return cls(mesh, bins)
11231130

11241131
def __hash__(self):
@@ -1861,7 +1868,7 @@ def __init__(self, particles, energies=None, filter_id=None):
18611868
def __repr__(self):
18621869
string = type(self).__name__ + '\n'
18631870
string += '{: <16}=\t{}\n'.format('\tParticles',
1864-
[str(p) for p in self.particles])
1871+
[str(p) for p in self.particles])
18651872
if self.energies is not None:
18661873
string += '{: <16}=\t{}\n'.format('\tEnergies', self.energies)
18671874
string += '{: <16}=\t{}\n'.format('\tID', self.id)
@@ -2171,7 +2178,7 @@ def paths(self):
21712178
if self._paths is None:
21722179
if not hasattr(self, '_geometry'):
21732180
raise ValueError(
2174-
"Model must be exported before the 'paths' attribute is" \
2181+
"Model must be exported before the 'paths' attribute is"
21752182
"available for a DistribcellFilter.")
21762183

21772184
# Determine paths for cell instances

openmc/mesh.py

Lines changed: 24 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -693,7 +693,7 @@ def num_mesh_cells(self):
693693
def write_data_to_vtk(self,
694694
filename: PathLike,
695695
datasets: dict | None = None,
696-
volume_normalization: bool = False,
696+
volume_normalization: bool | None = None,
697697
curvilinear: bool = False):
698698
"""Creates a VTK object of the mesh
699699
@@ -709,9 +709,10 @@ def write_data_to_vtk(self,
709709
with structured indexing in "C" ordering. See the "expand_dims" flag
710710
of :meth:`~openmc.Tally.get_reshaped_data` on reshaping tally data when using
711711
:class:`~openmc.MeshFilter`'s.
712-
volume_normalization : bool, optional
712+
volume_normalization : bool or None, optional
713713
Whether or not to normalize the data by the volume of the mesh
714-
elements.
714+
elements. When None (default), the format-appropriate default is
715+
used: True for legacy ASCII .vtk files, False for .vtkhdf files.
715716
curvilinear : bool
716717
Whether or not to write curvilinear elements. Only applies to
717718
``SphericalMesh`` and ``CylindricalMesh``.
@@ -748,10 +749,14 @@ def write_data_to_vtk(self,
748749
if write_impl is None:
749750
raise NotImplementedError(
750751
f"VTKHDF output not implemented for {type(self).__name__}")
751-
# write_impl is a bound method – do NOT pass self again
752-
write_impl(filename, datasets, volume_normalization)
752+
# vtkhdf default is False; explicit value is respected
753+
norm = volume_normalization if volume_normalization is not None else False
754+
write_impl(filename, datasets, norm)
753755
return None
754756

757+
# legacy ASCII .vtk default is True
758+
norm = volume_normalization if volume_normalization is not None else True
759+
755760
# vtk is an optional dependency only needed for the legacy ASCII path
756761
import vtk
757762
from vtk.util import numpy_support as nps
@@ -769,20 +774,26 @@ def write_data_to_vtk(self,
769774
# maintain a list of the datasets as added to the VTK arrays to
770775
# ensure they persist in memory until the file is written
771776
datasets_out = []
777+
# Regular/Rectilinear meshes store data in C (ijk) order which
778+
# matches VTK's expected ordering directly — no transpose needed.
779+
# Curvilinear meshes (Cylindrical, Spherical) require a transpose
780+
# to convert from C ordering to the Fortran (kji) ordering that VTK
781+
# expects.
782+
# TODO: update to "C" ordering throughout
783+
needs_transpose = not isinstance(
784+
self, (RegularMesh, RectilinearMesh))
772785
for label, dataset in datasets.items():
773786
dataset = self._reshape_vtk_dataset(dataset)
774787
self._check_vtk_dataset(label, dataset)
775-
# If the array data is 3D, assume is in C ordering and transpose
776-
# before flattening to match the ordering expected by the VTK
777-
# array based on the way mesh indices are ordered in the Python
778-
# API
779-
# TODO: update to "C" ordering throughout
780788
if dataset.ndim == 3:
781-
dataset = dataset.ravel()
789+
dataset = dataset.T.ravel() if needs_transpose else dataset.ravel()
782790
datasets_out.append(dataset)
783791

784-
if volume_normalization:
785-
dataset /= self.volumes.ravel()
792+
if norm:
793+
if needs_transpose:
794+
dataset /= self.volumes.T.ravel()
795+
else:
796+
dataset /= self.volumes.ravel()
786797

787798
dataset_array = vtk.vtkDoubleArray()
788799
dataset_array.SetName(label)

tests/unit_tests/test_mesh.py

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1217,7 +1217,10 @@ def test_write_ascii_vtk_unchanged(run_in_tmpdir):
12171217

12181218
# This should work without requiring vtk module changes
12191219
vtkIOLegacy = pytest.importorskip("vtkmodules.vtkIOLegacy")
1220-
mesh.write_data_to_vtk(datasets={"data": ref_data}, filename=filename)
1220+
mesh.write_data_to_vtk(datasets={"data": ref_data},
1221+
filename=filename,
1222+
volume_normalization=False
1223+
)
12211224

12221225
assert Path(filename).exists()
12231226

0 commit comments

Comments
 (0)