From 49efa93ab996ff5b7ce15bfca89c6f3ee9a661b5 Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Fri, 15 Aug 2025 12:04:20 +0200 Subject: [PATCH 1/3] Lint sort_mesh.py --- .../mpas_tools/mesh/creation/sort_mesh.py | 94 +++++++++---------- 1 file changed, 47 insertions(+), 47 deletions(-) diff --git a/conda_package/mpas_tools/mesh/creation/sort_mesh.py b/conda_package/mpas_tools/mesh/creation/sort_mesh.py index 9b7ec3b9f..c413182d9 100644 --- a/conda_package/mpas_tools/mesh/creation/sort_mesh.py +++ b/conda_package/mpas_tools/mesh/creation/sort_mesh.py @@ -1,8 +1,8 @@ +import argparse +import os import numpy as np -import os import xarray -import argparse from scipy.sparse import csr_matrix from scipy.sparse.csgraph import reverse_cuthill_mckee @@ -24,9 +24,9 @@ def sort_mesh(mesh): """ # Authors: Darren Engwirda - ncells = mesh.sizes["nCells"] - nedges = mesh.sizes["nEdges"] - nduals = mesh.sizes["nVertices"] + ncells = mesh.sizes['nCells'] + nedges = mesh.sizes['nEdges'] + nduals = mesh.sizes['nVertices'] cell_fwd = np.arange(0, ncells) + 1 cell_rev = np.arange(0, ncells) + 1 @@ -42,23 +42,20 @@ def sort_mesh(mesh): cell_rev = np.zeros(ncells, dtype=np.int32) cell_rev[cell_fwd - 1] = np.arange(ncells) + 1 - mesh["cellsOnCell"][:] = \ - _sort_rev(mesh["cellsOnCell"], cell_rev) - mesh["cellsOnEdge"][:] = \ - _sort_rev(mesh["cellsOnEdge"], cell_rev) - mesh["cellsOnVertex"][:] = \ - _sort_rev(mesh["cellsOnVertex"], cell_rev) + mesh['cellsOnCell'][:] = _sort_rev(mesh['cellsOnCell'], cell_rev) + mesh['cellsOnEdge'][:] = _sort_rev(mesh['cellsOnEdge'], cell_rev) + mesh['cellsOnVertex'][:] = _sort_rev(mesh['cellsOnVertex'], cell_rev) for var in mesh.keys(): dims = mesh.variables[var].dims - if ("nCells" in dims): + if 'nCells' in dims: mesh[var][:] = _sort_fwd(mesh[var], cell_fwd) - mesh["indexToCellID"][:] = np.arange(ncells) + 1 + mesh['indexToCellID'][:] = np.arange(ncells) + 1 # sort duals via pseudo-linear cell-wise ordering - dual_fwd = np.ravel(mesh["verticesOnCell"].values) + dual_fwd = np.ravel(mesh['verticesOnCell'].values) dual_fwd = dual_fwd[dual_fwd > 0] __, imap = np.unique(dual_fwd, return_index=True) @@ -68,22 +65,20 @@ def sort_mesh(mesh): dual_rev = np.zeros(nduals, dtype=np.int32) dual_rev[dual_fwd - 1] = np.arange(nduals) + 1 - mesh["verticesOnCell"][:] = \ - _sort_rev(mesh["verticesOnCell"], dual_rev) + mesh['verticesOnCell'][:] = _sort_rev(mesh['verticesOnCell'], dual_rev) - mesh["verticesOnEdge"][:] = \ - _sort_rev(mesh["verticesOnEdge"], dual_rev) + mesh['verticesOnEdge'][:] = _sort_rev(mesh['verticesOnEdge'], dual_rev) for var in mesh.keys(): dims = mesh.variables[var].dims - if ("nVertices" in dims): + if 'nVertices' in dims: mesh[var][:] = _sort_fwd(mesh[var], dual_fwd) - mesh["indexToVertexID"][:] = np.arange(nduals) + 1 + mesh['indexToVertexID'][:] = np.arange(nduals) + 1 # sort edges via pseudo-linear cell-wise ordering - edge_fwd = np.ravel(mesh["edgesOnCell"].values) + edge_fwd = np.ravel(mesh['edgesOnCell'].values) edge_fwd = edge_fwd[edge_fwd > 0] __, imap = np.unique(edge_fwd, return_index=True) @@ -93,37 +88,42 @@ def sort_mesh(mesh): edge_rev = np.zeros(nedges, dtype=np.int32) edge_rev[edge_fwd - 1] = np.arange(nedges) + 1 - mesh["edgesOnCell"][:] = \ - _sort_rev(mesh["edgesOnCell"], edge_rev) + mesh['edgesOnCell'][:] = _sort_rev(mesh['edgesOnCell'], edge_rev) - mesh["edgesOnEdge"][:] = \ - _sort_rev(mesh["edgesOnEdge"], edge_rev) + mesh['edgesOnEdge'][:] = _sort_rev(mesh['edgesOnEdge'], edge_rev) - mesh["edgesOnVertex"][:] = \ - _sort_rev(mesh["edgesOnVertex"], edge_rev) + mesh['edgesOnVertex'][:] = _sort_rev(mesh['edgesOnVertex'], edge_rev) for var in mesh.keys(): dims = mesh.variables[var].dims - if ("nEdges" in dims): + if 'nEdges' in dims: mesh[var][:] = _sort_fwd(mesh[var], edge_fwd) - mesh["indexToEdgeID"][:] = np.arange(nedges) + 1 + mesh['indexToEdgeID'][:] = np.arange(nedges) + 1 return mesh def main(): parser = argparse.ArgumentParser( - description=__doc__, - formatter_class=argparse.RawTextHelpFormatter) + description=__doc__, formatter_class=argparse.RawTextHelpFormatter + ) parser.add_argument( - "--mesh-file", dest="mesh_file", type=str, - required=True, help="Path+name to unsorted mesh file.") + '--mesh-file', + dest='mesh_file', + type=str, + required=True, + help='Path+name to unsorted mesh file.', + ) parser.add_argument( - "--sort-file", dest="sort_file", type=str, - required=True, help="Path+name to sorted output file.") + '--sort-file', + dest='sort_file', + type=str, + required=True, + help='Path+name to sorted output file.', + ) args = parser.parse_args() @@ -131,22 +131,23 @@ def main(): sort_mesh(mesh) - with open(os.path.join(os.path.dirname( - args.sort_file), "graph.info"), "w") as fptr: - cellsOnCell = mesh["cellsOnCell"].values + with open( + os.path.join(os.path.dirname(args.sort_file), 'graph.info'), 'w' + ) as fptr: + cellsOnCell = mesh['cellsOnCell'].values - ncells = mesh.sizes["nCells"] + ncells = mesh.sizes['nCells'] nedges = np.count_nonzero(cellsOnCell) // 2 - fptr.write(f"{ncells} {nedges}\n") + fptr.write(f'{ncells} {nedges}\n') for cell in range(ncells): data = cellsOnCell[cell, :] data = data[data > 0] for item in data: - fptr.write(f"{item} ") - fptr.write("\n") + fptr.write(f'{item} ') + fptr.write('\n') - mesh.to_netcdf(args.sort_file, format="NETCDF4") + mesh.to_netcdf(args.sort_file, format='NETCDF4') def _sort_fwd(data, fwd): @@ -212,11 +213,10 @@ def _cell_del2(mesh): ivec = np.array([], dtype=np.int32) jvec = np.array([], dtype=np.int32) - topolOnCell = mesh["nEdgesOnCell"].values - cellsOnCell = mesh["cellsOnCell"].values + topolOnCell = mesh['nEdgesOnCell'].values + cellsOnCell = mesh['cellsOnCell'].values for edge in range(np.max(topolOnCell)): - # cell-to-cell pairs, if edges exist mask = topolOnCell > edge idx_self = np.argwhere(mask).ravel() @@ -237,5 +237,5 @@ def _cell_del2(mesh): return csr_matrix((xvec, (ivec, jvec))) -if (__name__ == "__main__"):\ +if __name__ == '__main__': main() From 85b0413c276388fb012cc9d64d59fd00a48b17cc Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Fri, 15 Aug 2025 12:05:52 +0200 Subject: [PATCH 2/3] Add format argument to `sort_mesh` and use `write_netcdf()` This should hopefully fix an issue with the `Time` dimension. --- conda_package/mpas_tools/mesh/creation/sort_mesh.py | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/conda_package/mpas_tools/mesh/creation/sort_mesh.py b/conda_package/mpas_tools/mesh/creation/sort_mesh.py index c413182d9..3729d3f72 100644 --- a/conda_package/mpas_tools/mesh/creation/sort_mesh.py +++ b/conda_package/mpas_tools/mesh/creation/sort_mesh.py @@ -6,6 +6,8 @@ from scipy.sparse import csr_matrix from scipy.sparse.csgraph import reverse_cuthill_mckee +from mpas_tools.io import write_netcdf + def sort_mesh(mesh): """ @@ -125,6 +127,15 @@ def main(): help='Path+name to sorted output file.', ) + parser.add_argument( + '--format', + dest='format', + type=str, + required=False, + default='NETCDF4', + help='Output format for the sorted mesh file.', + ) + args = parser.parse_args() mesh = xarray.open_dataset(args.mesh_file) @@ -147,7 +158,7 @@ def main(): fptr.write(f'{item} ') fptr.write('\n') - mesh.to_netcdf(args.sort_file, format='NETCDF4') + write_netcdf(mesh, args.sort_file, format=args.format) def _sort_fwd(data, fwd): From e27d1efee1889c4d0dea00d863c1a3a0fe9cce5d Mon Sep 17 00:00:00 2001 From: Xylar Asay-Davis Date: Fri, 15 Aug 2025 12:54:50 +0200 Subject: [PATCH 3/3] Remove encoding unlimited_dims if Time not present When xarray 2025.08.0 opens a file with an empty unlimited Time dimension but no variables with Time, it drops `Time` as a dimension but keeps it in `unlimited_dims`. This is probably not an expected behavior but it's easy enough to address by setting `unlimited_dims = None` if `Time` is not present. --- conda_package/mpas_tools/io.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/conda_package/mpas_tools/io.py b/conda_package/mpas_tools/io.py index f4244b4ab..b6ebec839 100644 --- a/conda_package/mpas_tools/io.py +++ b/conda_package/mpas_tools/io.py @@ -135,6 +135,8 @@ def write_netcdf( # make sure the Time dimension is unlimited because MPAS has trouble # reading Time otherwise ds.encoding['unlimited_dims'] = {'Time'} + else: + ds.encoding['unlimited_dims'] = None # for performance, we have to handle this as a special case convert = format == 'NETCDF3_64BIT_DATA'