pyissm.model.mesh

Functions for building and interacting with an ISSM model mesh.

Functions

`bamg`(md, **kwargs)	Create a triangular mesh using the BAMG (Bidimensional Anisotropic Mesh Generator) algorithm.
`bamg_flowband`(md, x, surf, base, **kwargs)	Create a flowband mesh using BAMG (Bidimensional Anisotropic Mesh Generator).
`compute_hessian`(index, x, y, field, type)	Compute the Hessian matrix from a field.
`compute_metric`(hessian, scale, epsilon, ...)	Calculates anisotropic metric tensors used for adaptive mesh generation based on Hessian matrices.
`depth_average`(md, vector)	Compute depth average of a 3D vector using the trapezoidal rule.
`elements_from_edge`(elements, A, B)	Find elements connected to one edge defined by nodes A and B.
`export_gmsh`()	Export the model mesh to a Gmsh .msh file.
`find_element_types`(md, ice_levelset, ...)	Identify element types (ice, ice-front, ocean, floating ice, grounded ice, grounding line) from level set data.
`find_node_types`(md, ice_levelset, ocean_levelset)	Identify node types (ice, ice-front, ocean, floating ice, grounded ice) from level set data.
`find_segments`()	Build segments model field
`flag_elements`(md[, region, inside])	Flag elements based on their location within the model domain.
`get_element_areas_volumes`(index, x, y[, z])	Computes areas of triangular elements or volumes of pentahedrons.
`get_mesh`(mesh_x, mesh_y, mesh_elements)	Create a triangular mesh from an unstructured model object.
`get_nodal_functions_coeff`(index, x, y)	Compute the coefficients alpha, beta and gamma of 2D triangular elements. For each triangular element, the nodal functions are defined as: N_i(x, y) = alpha_i * x + beta_i * y + gamma_i.
`grid_model_field`(md, model_field, grid_x, grid_y)	Interpolate unstructured model data onto a regular 2D grid.
`make_gridded_domain_mask`(mesh_x, mesh_y, ...)	Generate a binary domain mask on a regular grid based on an unstructured mesh.
`mesh_convert`(md, **kwargs)	Convert mesh to BAMG format for advanced mesh operations.
`model_intersect_3d`()	Intersect a 3D model with a plane defined by points (xs, ys, zs).
`model_merge_3d`()	Merge two 3D models into a single model.
`process_mesh`(md)	Process ISSM model mesh
`round_mesh`(md, radius, resolution[, ...])	Create a structured triangular mesh on a circular domain.
`square_mesh`(md, Lx, Ly, nx, ny)	Create a structured triangular mesh on a rectangular domain.
`triangle`(md, domain_name, resolution[, ...])	Create a triangular mesh for an ISSM model using Triangle mesh generator.
`twod_to_3d`()	Convert 2D mesh to 3D surface mesh.

pyissm.model.mesh.bamg(md, **kwargs)

Create a triangular mesh using the BAMG (Bidimensional Anisotropic Mesh Generator) algorithm.

This function generates high-quality anisotropic triangular meshes for complex geometries using the BAMG mesh generator. It supports various mesh constraints including domain boundaries, holes, subdomains, rifts, and anisotropic metrics.

Parameters:

md (object) – ISSM model object whose mesh fields will be populated with the generated mesh.
kwargs (dict, optional) – Additional keyword arguments to customize mesh generation. Supported options include: - anisomax : float, maximum anisotropy ratio allowed in the mesh (default= 1e30) - coeff : float, global mesh size coefficient multiplier (default= 1.0) - Crack : int, enable crack processing (0=disabled, 1=enabled) (default= 0) - cutoff : float, cutoff value for metric interpolation (default= 1e-5) - domain : str or list, path to domain outline file or list of domain contours (default= None) - err : float, target interpolation error for mesh adaptation (default= 0.01) - errg : float, target geometric error for mesh adaptation (default= 0.1) - field : array_like, field values for metric computation (default= empty array) - gradation : float, mesh size gradation parameter (default= 1.5) - Hessiantype : int, type of Hessian computation (0=P1, 1=P2) (default= 0) - hmax : float, maximum allowed edge length (default= 1e100) - hmin : float, minimum allowed edge length (default= 1e-100) - hmaxVertices : array_like, maximum edge lengths at specific vertices (default= empty array) - hminVertices : array_like, minimum edge lengths at specific vertices (default= empty array) - holes : str or list, path to holes file or list of hole contours (default= None) - hVertices : array_like, target edge lengths at specific vertices (default= empty array) - KeepVertices : int, keep vertices from previous mesh (0=no, 1=yes) (default= 1) - Markers : array_like, edge markers for boundary identification (default= None) - maxnbv : float, maximum number of vertices allowed (default= 1.0e6) - maxsubdiv : float, maximum number of edge subdivisions (default= 10.0) - metric : array_like, anisotropic metric tensor field (default= empty array) - Metrictype : int, type of metric (0=isotropic, 1=anisotropic) (default= 0) - nbjacobi : int, number of Jacobi smoothing iterations (default= 1) - nbsmooth : int, number of mesh smoothing iterations (default= 3) - NoBoundaryRefinement : int, disable boundary refinement for domain edges (0=allow, 1=disable) (default= 0) - NoBoundaryRefinementAllBoundaries : int, disable boundary refinement for all edges (0=allow, 1=disable) (default= 0) - omega : float, relaxation parameter for smoothing (default= 1.8) - power : float, power for metric computation (default= 1.0) - RequiredVertices : array_like, coordinates of vertices that must be included in the mesh (default= None) - rifts : str, path to rifts file for fracture modeling (default= None) - splitcorners : int, split corners in mesh generation (0=no, 1=yes) (default= 1) - subdomains : str or list, path to subdomains file or list of subdomain contours (default= None) - tol : float, tolerance for geometric operations (default= None) - toltip : float, tolerance for rift tip processing (default= None) - tracks : str or array_like, path to tracks file or track coordinates (default= None) - verbose : int, verbosity level (0=quiet, 1=verbose) (default= 1) - vertical : int, create 2D vertical mesh (0=standard 2D, 1=vertical) (default= 0) - 3dsurface : int, create 3D surface mesh (0=standard 2D, 1=3D surface) (default= 0)

Returns:

md – The input ISSM model object with updated mesh properties including: - mesh.x, mesh.y: Node coordinates - mesh.elements: Element connectivity matrix - mesh.edges: Edge connectivity matrix - mesh.segments: Boundary segment definitions - mesh.segmentmarkers: Boundary segment markers - mesh.numberofvertices: Total number of mesh vertices - mesh.numberofelements: Total number of mesh elements - mesh.numberofedges: Total number of mesh edges - mesh.vertexonboundary: Boolean array indicating boundary vertices - mesh.elementconnectivity: Element-to-element connectivity - private.bamg: BAMG-specific mesh and geometry data

Return type:

object

Raises:

IOError – If specified input files (domain, holes, subdomains, rifts) do not exist.
RuntimeError – If mesh generation fails or if incompatible options are specified.
TypeError – If input arguments are of incorrect type.

Notes

This function is a comprehensive interface to the BAMG mesh generator, supporting complex geometries with multiple constraints. The mesh can be adapted based on metric fields for anisotropic meshing. Special handling is provided for rifts and fractures in ice sheet modeling applications.

Examples

>>> import pyissm
>>> md = pyissm.Model()
>>> # Basic domain meshing
>>> md = pyissm.model.mesh.bamg(md, domain='outline.exp', hmax=1000.0)
>>> # Anisotropic meshing with metric
>>> md = pyissm.model.mesh.bamg(md, domain='outline.exp',
...                             metric=metric_field, err=0.005)
>>> # Mesh with holes and subdomains
>>> md = pyissm.model.mesh.bamg(md, domain='outline.exp',
...                             holes='holes.exp', subdomains='regions.exp')

pyissm.model.mesh.bamg_flowband(md, x, surf, base, **kwargs)

Create a flowband mesh using BAMG (Bidimensional Anisotropic Mesh Generator).

This function generates a triangular mesh for a flowband (vertical 2D slice) using the BAMG mesh generator. The flowband is defined by surface and base profiles along a specified coordinate path, creating a vertical mesh suitable for ice flow modeling in the vertical plane.

Parameters:

md (object) – ISSM model object whose mesh fields will be populated with the generated flowband mesh.
x (array_like) – 1D array of coordinates along the flowband path. These represent the horizontal positions where the surface and base elevations are defined.
surf (array_like) – 1D array of surface elevations corresponding to each x coordinate. Must have the same length as x.
base (array_like) – 1D array of base (bed) elevations corresponding to each x coordinate. Must have the same length as x.
**kwargs (dict, optional) – Additional keyword arguments passed to the bamg function. See bamg() documentation for supported options.

Returns:

md – A new ISSM model object with a vertical 2D mesh populated, including: - mesh.x, mesh.y: Node coordinates in the flowband coordinate system - mesh.elements: Element connectivity matrix for triangular elements - mesh.edges: Edge connectivity matrix - mesh.segments: Boundary segment definitions with markers - mesh.segmentmarkers: Boundary segment markers (1=base, 2=right, 3=surface, 4=left) - mesh.numberofvertices: Total number of mesh vertices - mesh.numberofelements: Total number of mesh elements - mesh.numberofedges: Total number of mesh edges - mesh.vertexonboundary: Boolean array indicating boundary vertices - mesh.vertexonbase: Boolean array indicating vertices on the base boundary - mesh.vertexonsurface: Boolean array indicating vertices on the surface boundary - mesh.elementconnectivity: Element-to-element connectivity

Return type:

object

Raises:

ValueError – If x, surf, and base arrays do not have the same length.
RuntimeError – If BAMG mesh generation fails or if incompatible meshing options are specified.

Notes

This function creates a vertical 2D mesh by: 1. Constructing a closed domain from the surface and base profiles 2. Assigning boundary markers: 1=base, 2=right side, 3=surface, 4=left side 3. Calling the BAMG mesh generator with vertical=1 option (to convert to 2D vertical mesh) 4. Post-processing to identify vertices on base and surface boundaries

The resulting mesh is suitable for flowband modeling where ice flow is assumed to be primarily in the vertical plane defined by the x-coordinate path.

Examples

>>> md = pyissm.Model()
>>> x = np.arange(1, 3001, 100).T
>>> h = np.linspace(1000, 300, np.size(x)).T
>>> b = -917. / 1023. * h
>>> md = pyissm.model.mesh.bamg_flowband(md, x = x, surf = b + h, base = b, hmax = 80.)

pyissm.model.mesh.compute_hessian(index, x, y, field, type)

Compute the Hessian matrix from a field.

Computes the Hessian matrix of a given field and returns the three components Hxx, Hxy, Hyy for each element or each node.

Parameters:

index (ndarray) – Element connectivity matrix defining the triangular mesh elements. Shape: (num_elements, 3) with 1-based indexing.
x (ndarray) – X-coordinates of the mesh nodes. Shape: (num_nodes,).
y (ndarray) – Y-coordinates of the mesh nodes. Shape: (num_nodes,).
field (ndarray) – Field values defined either on nodes or elements. Shape: (num_nodes,) or (num_elements,).
type (str) – Type of output desired. Must be either ‘node’ or ‘element’.

Returns:

Hessian matrix components. Shape depends on type:

If type is ‘element’: (num_elements, 3) with columns [Hxx, Hxy, Hyy] for each element.
If type is ‘node’: (num_nodes, 3) with columns [Hxx, Hxy, Hyy] for each node.

Return type:

ndarray

Raises:

TypeError – If the field is not defined on nodes or elements, or if type is not ‘node’ or ‘element’.

Examples

>>> hessian = compute_hessian(md.mesh.elements, md.mesh.x, md.mesh.y,
...                          md.inversion.vel_obs, 'node')
>>> hessian = compute_hessian(md.mesh.elements, md.mesh.x, md.mesh.y,
...                          md.thermal.temperature, 'element')

Notes

The Hessian computation uses finite element nodal functions and area-weighted averaging for nodal values. For element-based fields, values are first interpolated to nodes before computing gradients and Hessian components.

The Hessian matrix H has components: H = [[Hxx, Hxy], [Hxy, Hyy]]

This function returns the three unique components as [Hxx, Hxy, Hyy].

pyissm.model.mesh.compute_metric(hessian, scale, epsilon, hmin, hmax, pos)

Calculates anisotropic metric tensors used for adaptive mesh generation based on Hessian matrices. The metric tensor controls element size and orientation in the mesh by analyzing eigenvalues and eigenvectors of the Hessian matrix.

Parameters:

hessian (numpy.ndarray) – Array of shape (n, 3) containing Hessian matrix components for each node. Columns represent [H11, H12, H22] where H is the 2x2 Hessian matrix.
scale (float) – Scaling factor for the metric computation.
epsilon (float) – Tolerance parameter used in the metric scaling calculation.
hmin (float) – Minimum allowed element size in the mesh.
hmax (float) – Maximum allowed element size in the mesh.
pos (numpy.ndarray) – Array of indices corresponding to water elements or special boundary conditions that require uniform metric treatment.

Returns:

Array of shape (n, 3) containing the computed metric tensor components [M11, M12, M22] for each node, where M is the 2x2 symmetric metric tensor.

Return type:

numpy.ndarray

Raises:

RuntimeError – If NaN values persist in the metric tensor after all correction attempts.

Notes

The function performs the following key operations: 1. Computes eigenvalues and eigenvectors of the Hessian matrix 2. Applies size constraints using hmin and hmax parameters 3. Handles special cases (zero eigenvalues, water elements) 4. Uses numpy.linalg.eig as a fallback for numerical issues 5. Ensures the resulting metric is free of NaN values The metric tensor M is used in adaptive mesh generation where element sizes are controlled by the relationship: h^T * M * h = 1, where h represents the edge vector in the mesh.

Examples

>>> hessian = compute_hessian(md.mesh.elements, md.mesh.x, md.mesh.y,
...                          md.inversion.vel_obs, 'node')
>>> metric = compute_metric(hessian, 1.0, 0.01, 0.1, 10.0, np.array([]))

pyissm.model.mesh.depth_average(md, vector)

Compute depth average of a 3D vector using the trapezoidal rule.

This function computes the depth-averaged value of a 3D field defined on an extruded mesh and returns the result projected onto the corresponding 2D mesh. The depth averaging is performed using the trapezoidal integration rule.

Parameters:

md (ISSM Model object) – ISSM Model object containing a 3D mesh. Must have valid 3D mesh attributes (e.g., md.mesh.numberofvertices, md.mesh.numberoflayers, md.mesh.z) and geometry information (md.geometry.thickness).
vector (ndarray or array_like) – 3D field values to be depth-averaged. Can be: - (md.mesh.numberofvertices,) for node-based data - (md.mesh.numberofelements,) for element-based data

Returns:

vector_average – Depth-averaged field values projected onto the 2D mesh. Shape depends on input: - (md.mesh.numberofvertices2d,) for node-based input - (md.mesh.numberofelements2d,) for element-based input

Return type:

ndarray

Raises:

TypeError – If md is not provided or does not contain a 3D mesh.
ValueError – If vector size does not match expected dimensions.

Notes

The function uses the trapezoidal rule for integration:

\[\bar{v} = \frac{1}{H} \int_0^H v(z) dz\]

where H is the total thickness and the integral is approximated using layer-by-layer trapezoidal integration.

Examples

Compute depth-averaged velocity:

>>> vel_bar = depth_average(md, md.initialization.vel)

Compute depth-averaged temperature:

>>> temp_bar = depth_average(md, md.thermal.temperature)

pyissm.model.mesh.elements_from_edge(elements, A, B)

Find elements connected to one edge defined by nodes A and B.

Parameters:

elements (array_like) – Array of element connectivity information where each row represents an element and columns represent the nodes that define the element.
A (int) – First node ID defining the edge.
B (int) – Second node ID defining the edge.

Returns:

1D array of element IDs (1-based indexing) that contain the edge defined by nodes A and B.

Return type:

ndarray

Examples

>>> edgeelements = elements_from_edge(md.mesh.elements, node1, node2)

pyissm.model.mesh.export_gmsh()

Export the model mesh to a Gmsh .msh file.

Raises:: NotImplementedError – Function is not yet implemented.

pyissm.model.mesh.find_element_types(md, ice_levelset, ocean_levelset)

Identify element types (ice, ice-front, ocean, floating ice, grounded ice, grounding line) from level set data.

This function processes level set fields for ice and ocean and classifies mesh elements into several categories based on their sign.

For 3D meshes, only the surface layer (where vertexonsurface == 1) is processed (see find_node_types()).

Parameters:

md (ISSM Model object) – ISSM Model object containing the mesh and geometry information. Must have attributes: ‘md.mesh.*’ used by process_mesh().
ice_levelset (ndarray) –

1D array of values from the ice level set:
Ice < 0 Ice front = 0 No Ice > 0
ocean_levelset (ndarray) –

1D array of values from the ocean level set:
Ocean < 0 No ocean > 0

Returns:

Dictionary with boolean arrays (same length as number of surface nodes), indicating:

’ice_elements’ : Elements with ice
’ice_front_elements’ : Elements on the ice front
’ocean_elements’ : Elements with ocean
’floating_ice_elements’ : Elements with floating ice
’grounded_ice_elements’ : Elements with grounded ice
’grounding_line_elements’ : Elements on the grounding line

Return type:

dict of str -> ndarray

Warning

If a 3D mesh is detected, only the surface layer is used. A warning is issued.

Example

model_element_types = find_element_types(md, md.mask.ice_levelset, md.mask.ocean_levelset) model_element_types = find_element_types(md, md.results.TransientSolution.MaskIceLevelset[34], md.results.TransientSolution.MaskOceanLevelset[34])

pyissm.model.mesh.find_node_types(md, ice_levelset, ocean_levelset)

Identify node types (ice, ice-front, ocean, floating ice, grounded ice) from level set data.

This function processes level set fields for ice and ocean and classifies mesh nodes into several categories based on their sign.

For 3D meshes, only the surface layer (where vertexonsurface == 1) is processed.

Parameters:

md (ISSM Model object) – ISSM Model object containing the mesh and geometry information. Must have attributes: ‘md.mesh.*’ used by process_mesh().
ice_levelset (ndarray) –

1D array of values from the ice level set:
Ice < 0 Ice front = 0 No Ice > 0
ocean_levelset (ndarray) –

1D array of values from the ocean level set:
Ocean < 0 No ocean > 0

Returns:

Dictionary with boolean arrays (same length as number of surface nodes), indicating:

’ice_nodes’ : Nodes with ice (ice_levelset < 0)
’ice_front_nodes’ : Nodes on the ice front (ice_levelset == 0)
’ocean_nodes’ : Nodes with ocean (ocean_levelset < 0)
’floating_ice_nodes’ : Nodes with floating ice (ice_levelset < 0 & ocean_levelset < 0)
’grounded_ice_nodes’ : Nodes with grounded ice (ice_levelset < 0 & ocean_levelset >= 0)

Return type:

dict of str -> ndarray

Warning

If a 3D mesh is detected, only the surface layer is used. A warning is issued.

Example

model_node_types = find_node_types(md, md.mask.ice_levelset, md.mask.ocean_levelset) model_node_types = find_node_types(md, md.results.TransientSolution.MaskIceLevelset[34], md.results.TransientSolution.MaskOceanLevelset[34])

pyissm.model.mesh.find_segments()

Build segments model field

Raises:: NotImplementedError – Function is not yet implemented.

pyissm.model.mesh.flag_elements(md, region='all', inside=True)

Flag elements based on their location within the model domain.

This function allows users to flag elements in the mesh based on whether they are inside or outside a specified domain. The region can be the entire mesh, no elements, a region specified by a provided *.exp file, or defined by boolean arrays.

Parameters:

md (ISSM Model object) – ISSM Model object containing the mesh
region (str or ndarray, optional) –
Region specification. Options are:
- ’all’ (default): Flag all elements in the mesh.
- ’’: Flag no elements.
- Path to a *.exp file: Flag elements inside or outside the polygon defined in the file.
- ndarray: Boolean array of size (numberofelements,) or (numberofvertices,). If vertices array, elements are flagged when all vertices are flagged.
inside (bool, optional) – If region is a polygon file or array, this parameter specifies whether to flag elements inside (True, default) or outside (False) the region. Default is True.

Returns:

Boolean array of length md.mesh.numberofelements where True indicates flagged elements and False indicates unflagged elements.

Return type:

ndarray of bool

Raises:

RuntimeError – If python wrappers are not installed and a *.exp file is provided.
ValueError – If region array does not match number of elements or vertices.
TypeError – If region is neither a string nor an array.

Examples

Flag all elements in the mesh:

>>> flags = flag_elements(md)
>>> flags = flag_elements(md, region='all')

Flag no elements:

>>> flags = flag_elements(md, region='')

Flag elements inside a polygon:

>>> flags = flag_elements(md, region='path/to/polygon.exp')

Flag elements outside a polygon:

>>> flags = flag_elements(md, region='path/to/polygon.exp', inside=False)

pyissm.model.mesh.get_element_areas_volumes(index, x, y, z=array([], dtype=float64))

Computes areas of triangular elements or volumes of pentahedrons.

Parameters:

index (ndarray) – Element connectivity array. For 2D meshes, should have 3 columns. For 3D meshes, should have 6 columns.
x (ndarray) – 1D array of x-coordinates of mesh nodes.
y (ndarray) – 1D array of y-coordinates of mesh nodes.
z (ndarray, optional) – 1D array of z-coordinates of mesh nodes. If provided, volumes are computed. Default is empty array (areas computed).

Returns:

areas – 1D array of element areas (2D) or volumes (3D).

Return type:

ndarray

Raises:

TypeError – If x, y, and z arrays don’t have the same length. If index contains values above the number of nodes. If index doesn’t have the correct number of columns for the mesh type.

Examples

Compute areas of triangular elements:

>>> areas = get_element_areas(md.mesh.elements, md.mesh.x, md.mesh.y)

Compute volumes of pentahedral elements:

>>> volumes = get_element_areas(md.mesh.elements, md.mesh.x, md.mesh.y, md.mesh.z)

pyissm.model.mesh.get_mesh(mesh_x, mesh_y, mesh_elements)

Create a triangular mesh from an unstructured model object.

Extracts node coordinates and element connectivity from the model and constructs a matplotlib.tri.Triangulation object for downstream operations.

Parameters:

mesh_x (1d array) – x-coordinates of the mesh nodes.
mesh_y (1d array) – y-coordinates of the mesh nodes.
mesh_elements (2d array) – element connectivity.

Returns:

mesh – Triangular mesh object representing the model domain.

Return type:

matplotlib.tri.Triangulation

Notes

If necessary, the function adjusts the element indexing from 1-based to 0-based indexing

required by Triangulation.

See also

grid_data: Interpolates data onto a regular grid, using this mask by default.

pyissm.model.mesh.mesh_convert(md, **kwargs)

Convert mesh to BAMG format for advanced mesh operations.

This function converts an existing mesh to BAMG (Bidimensional Anisotropic Mesh Generator) format, enabling access to BAMG’s advanced mesh manipulation capabilities. The conversion creates internal BAMG data structures while preserving the original mesh geometry and connectivity.

Parameters:

md (object) – ISSM model object containing the mesh to be converted. The mesh should have valid elements, coordinates, and connectivity information.
**kwargs (dict, optional) –
Additional keyword arguments to customize the conversion: - index : array_like, optional

Element connectivity matrix. Defaults to md.mesh.elements.
- xarray_like, optional
  X-coordinates of mesh vertices. Defaults to md.mesh.x.
- yarray_like, optional
  Y-coordinates of mesh vertices. Defaults to md.mesh.y.

Returns:

md – The input ISSM model object with updated mesh properties and BAMG data structures: - mesh.x, mesh.y: Node coordinates - mesh.elements: Element connectivity matrix - mesh.edges: Edge connectivity matrix - mesh.segments: Boundary segment definitions - mesh.segmentmarkers: Boundary segment markers - mesh.numberofvertices: Total number of mesh vertices - mesh.numberofelements: Total number of mesh elements - mesh.numberofedges: Total number of mesh edges - mesh.vertexonboundary: Boolean array indicating boundary vertices - mesh.elementconnectivity: Element-to-element connectivity - private.bamg: BAMG-specific mesh and geometry data structures

Return type:

object

Notes

This function is primarily used to prepare existing meshes for advanced BAMG operations such as mesh adaptation, refinement, or anisotropic meshing. The conversion creates internal BAMG data structures that enable seamless integration with other BAMG-based mesh operations.

The function preserves all original mesh properties while adding BAMG-specific data structures to the model’s private fields. This allows subsequent BAMG operations to work efficiently without data conversion overhead.

Examples

>>> md = pyissm.Model()
>>> md = pyissm.model.mesh.triangle(md, 'domain.exp', 1000.0)
>>> md = pyissm.model.mesh.meshconvert(md)
>>> md = pyissm.model.mesh.meshconvert(md, x=custom_x, y=custom_y)

pyissm.model.mesh.model_intersect_3d()

Intersect a 3D model with a plane defined by points (xs, ys, zs).

Raises:: NotImplementedError – Function is not yet implemented.

pyissm.model.mesh.model_merge_3d()

Merge two 3D models into a single model.

Raises:: NotImplementedError – Function is not yet implemented.

pyissm.model.mesh.process_mesh(md)

Process ISSM model mesh

This function processes key elements of an ISSM model mesh and is used in several core pyISSM functions to ensure consistency

Parameters:

md (ISSM Model object) – ISSM Model object from which the mesh should be extracted/processed.

Returns:

mesh (matplotlib.tri.Triangulation) – Triangular mesh object representing the model domai
mesh_x (1d array) – x-coordinates of the mesh nodes.
mesh_y (1d array) – y-coordinates of the mesh nodes.
mesh_elements (2d array) – element connectivity.
is3d (bool) – ‘True’ if elements2d exists (and the model is 3D), ‘False’ otherwise

Example

mesh, mesh_x, mesh_y, mesh_elements, is3d = process_mesh(md)

pyissm.model.mesh.round_mesh(md, radius, resolution, exp_output_name=None, keep_exp=False)

Create a structured triangular mesh on a circular domain.

This function generates a triangular mesh on a circular domain by first creating a domain outline file and then using the Triangle mesh generator. The mesh is created with approximately uniform resolution around the circle perimeter and moves the closest node to the origin.

Parameters:

md (object) – ISSM model object whose mesh fields will be populated.
radius (float) – Radius of the circular domain in meters.
resolution (float) – Target mesh resolution in meters. This represents the characteristic edge length for mesh elements around the circle perimeter.
exp_output_name (str, optional) – Path for the output domain outline file (.exp format). If None, defaults to ‘round_mesh.exp’. Default is None.
keep_exp (bool, optional) – Whether to keep the temporary domain outline file after mesh creation. If False, the file is automatically deleted. Default is False.

Returns:

md – The input ISSM model object with updated mesh properties including: - mesh.x, mesh.y: Node coordinates - mesh.elements: Element connectivity matrix - mesh.segments: Boundary segment definitions - mesh.segmentmarkers: Boundary segment markers - mesh.numberofvertices: Total number of mesh vertices - mesh.numberofelements: Total number of mesh elements - mesh.vertexonboundary: Boolean array indicating boundary vertices - mesh.vertexconnectivity: Vertex-to-vertex connectivity - mesh.elementconnectivity: Element-to-element connectivity

Return type:

object

Raises:

IOError – If the specified exp_output_name file already exists.
RuntimeError – If the Triangle Python wrappers are not installed.

Warning

UserWarning: If the model mesh is not empty and will be overwritten.

Notes

The function creates a circular mesh by: 1. Generating points uniformly distributed around the circle perimeter 2. Writing these points to a domain outline file (.exp format) 3. Using the Triangle mesh generator to create the triangular mesh 4. Moving the closest node to the origin (0,0) for convenience 5. Optionally removing the temporary domain outline file

The number of points on the circle perimeter is calculated based on the target resolution to ensure approximately uniform spacing.

Examples

>>> import pyissm
>>> md = pyissm.Model()
>>> md = pyissm.model.mesh.round_mesh(md, radius=5000.0, resolution=500.0)
>>> md = pyissm.model.mesh.round_mesh(md, radius=1000.0, resolution=100.0,
...                                   exp_output_name='circle.exp', keep_exp=True)

pyissm.model.mesh.square_mesh(md, Lx, Ly, nx, ny)

Create a structured triangular mesh on a rectangular domain.

This function generates a structured triangular mesh on a rectangular domain and populates the mesh fields of an ISSM model object. The mesh consists of triangular elements arranged in a regular grid pattern.

Parameters:

md (object) – ISSM model object whose mesh fields will be populated.
Lx (float) – Length of the domain in the x-direction.
Ly (float) – Length of the domain in the y-direction.
nx (int) – Number of nodes in the x-direction.
ny (int) – Number of nodes in the y-direction.

Raises:

RuntimeError – If Python wrappers are not installed.

Warning

UserWarning: If the model mesh is not empty and will be overwritten.

Notes

The function creates a structured triangular mesh by: 1. Generating node coordinates on a regular grid 2. Creating triangular elements by splitting each grid cell into two triangles 3. Defining boundary segments around the domain perimeter 4. Computing connectivity arrays for mesh topology The resulting mesh will have (nx-1)*(ny-1)*2 triangular elements and nx*ny nodes.

Examples

>>> import pyissm
>>> md = pyissm.Model()
>>> md = pyissm.model.mesh.square_mesh(md, Lx=100.0, Ly=50.0, nx=11, ny=6)

pyissm.model.mesh.triangle(md, domain_name, resolution, rift_name=None)

Create a triangular mesh for an ISSM model using Triangle mesh generator.

This function generates a triangular mesh based on a domain outline and optional rift constraints. It uses the Triangle mesh generator to create high-quality Delaunay triangulations with specified resolution constraints.

Parameters:

md (object) – ISSM model object containing the mesh structure to be populated.
domain_name (str) – Path to the file containing the domain outline geometry.
resolution (float) – Target mesh resolution in meters. This represents the characteristic edge length for mesh elements. The actual mesh area constraint is calculated as resolution squared.
rift_name (str, optional) – Path to the file containing rift constraint geometry. If provided, these constraints will be incorporated into the mesh generation. Default is None.

Returns:

md – The input ISSM model object with updated mesh properties including: - mesh.x, mesh.y: Node coordinates - mesh.elements: Element connectivity matrix - mesh.segments: Boundary segment definitions - mesh.segmentmarkers: Boundary segment markers - mesh.numberofvertices: Total number of mesh vertices - mesh.numberofelements: Total number of mesh elements - mesh.vertexonboundary: Boolean array indicating boundary vertices - mesh.vertexconnectivity: Vertex-to-vertex connectivity - mesh.elementconnectivity: Element-to-element connectivity

Return type:

object

Raises:

IOError – If the domain outline file or rift file (when specified) does not exist.
RuntimeError – If the Triangle Python wrappers are not installed and mesh creation fails.

Warning

Issues a warning if the existing mesh is not empty and will be overwritten.
Issues a warning if orphaned nodes are found and removed from the mesh.

Notes

This function requires the Triangle Python wrappers to be installed for mesh generation. If wrappers are not available, the function raises a RuntimeError. The function automatically handles orphaned nodes (nodes not belonging to any element) by removing them and updating the connectivity accordingly.

Examples

>>> import pyissm
>>> md = pyissm.Model()
>>> md = pyissm.model.mesh.triangle(md, 'domain.exp', 1000.0)
>>> md = pyissm.model.mesh.triangle(md, 'domain.exp', 500.0, 'rifts.exp')

pyissm.model.mesh.twod_to_3d()

Convert 2D mesh to 3D surface mesh.

Raises:: NotImplementedError – Function is not yet implemented.