pygimli.meshtools

Mesh generation and modification.

Note

Although we discriminate here between grids (structured meshes) and meshes (unstructured), both objects are treated the same internally.

Overview

Functions

appendTetrahedronBoundary(mesh[, xbound, …]) Return new mesh surrounded by tetrahedron boundary box.
appendTriangleBoundary(mesh[, xbound, …]) Add a triangle mesh boundary to a given mesh.
cellDataToBoundaryData(mesh, data) TODO DOCUMENT_ME
cellDataToNodeData(mesh, data[, style]) Convert cell data to node data.
convertHDF5Mesh(h5Mesh[, group, indices, …]) Converts instance of a hdf5 mesh to a GIMLI::Mesh.
createCircle([pos, radius, segments, start, end]) Create simple circle polygon.
createCylinder(radius, height[, nSegments, …]) Create plc of a cylinder.
createGrid([x, y, z]) Create grid style mesh.
createLine(start, end[, segments]) Create simple line polygon.
createMesh(poly[, quality, area, smooth, …]) Create a mesh for a given geometry polygon.
createParaDomain2D(*args, **kwargs) API change here ..
createParaMesh(*args, **kwargs) Create parameter mesh from list of sensor positions.
createParaMesh2DGrid(sensors[, paraDX, …]) Create a grid style mesh for an inversion parameter mesh.
createParaMeshPLC(sensors[, paraDX, …]) Create a PLC mesh for an inversion parameter mesh.
createPolygon(verts[, isClosed, isHole]) Create a polygon from a list of vertices.
createRectangle([start, end, pos, size]) Create rectangle polygon.
createWorld(start, end[, marker, area, …]) Create simple rectangular world.
exportFenicsHDF5Mesh(mesh, exportname, **kwargs) Exports Gimli mesh in HDF5 format suitable for Fenics.
exportHDF5Mesh(mesh, exportname[, group, …]) Writes given GIMLI::Mesh in a hdf5 format file.
exportPLC(poly, fname, **kwargs) General writer to save piece-wise linear complex (PLC) as poly file.
exportSTL(mesh, fileName[, ascii]) Write STL surface mesh and returns a GIMLI::Mesh.
fillEmptyToCellArray(mesh, vals[, slope]) Prolongate empty cell values to complete cell attributes.
interpolate(*args, **kwargs) Interpolation convinience function.
interpolateAlongCurve(curve, t, **kwargs) Interpolate along curve.
merge2Meshes(m1, m2) Merge two meshes into one new mesh and return the combined mesh.
mergeMeshes(meshList[, verbose]) Merge several meshes into one new mesh and return the new mesh.
mergePLC(pols[, tol]) Merge multiply polygons.
nodeDataToBoundaryData(mesh, data) Assuming [NodeCount, dim] data DOCUMENT_ME
nodeDataToCellData(mesh, data) Convert node data to cell data.
quality(mesh[, measure]) Return the quality of a given triangular mesh.
readFenicsHDF5Mesh(fileName[, group, verbose]) Reads FEniCS mesh from file format .h5 and returns a GIMLI::Mesh.
readGmsh(fname[, verbose]) Read Gmsh ASCII file and return instance of GIMLI::Mesh class.
readHDF5Mesh(fileName[, group, indices, …]) Function for loading a mesh from HDF5 file format.
readHydrus2dMesh([fname]) Import mesh from Hydrus2D.
readHydrus3dMesh([fileName]) Import mesh from Hydrus3D.
readPLC(filename[, comment]) Generic PLC reader.
readSTL(fileName[, ascii]) Read STL surface mesh and returns a GIMLI::Mesh.
readTetgen(fname[, comment, verbose, …]) Reads and converts a mesh from the basic Tetgen output.
readTriangle(fname[, verbose]) Read Triangle [She96] mesh.
refineQuad2Tri(mesh[, style]) Refine mesh of quadrangles into a mesh of triangle cells.
syscallTetgen(filename[, quality, area, …]) Create a mesh with Tetgen from file.
tapeMeasureToCoordinates(tape, pos) Interpolate 2D tape measured topography to 2D Cartesian coordinates.
writePLC(*args, **kwargs) Backward compatibility.

Functions

appendTetrahedronBoundary

pygimli.meshtools.appendTetrahedronBoundary(mesh, xbound=100, ybound=100, zbound=100, marker=1, quality=2, isSubSurface=False, verbose=False)[source]

Return new mesh surrounded by tetrahedron boundary box.

Creates a tetrahedral box around a given mesh suitable for geo-simulation (surface boundary at top).

Parameters:
mesh : mesh object

Mesh to which the tetrahedron boundary should be appended.

xbound : float, optional

Horizontal prolongation distance in x-direction.

ybound : float, optional

Horizonal prolongation distance in y-direction.

zbound : float, optional

Vertical prolongation distance.

marker : int, optional

Marker of new cells.

quality : float, optional

Triangle quality.

isSubSurface : boolean, optional

Apply boundary conditions suitable for geo-simulaion and prolongate mesh to the surface if necessary.

verbose : boolean, optional

Be verbose.

Notes

Boundaries of mesh need marker 1.

appendTriangleBoundary

pygimli.meshtools.appendTriangleBoundary(mesh, xbound=10, ybound=10, marker=1, quality=34.0, area=0.0, smooth=False, markerBoundary=1, isSubSurface=False, verbose=False)[source]

Add a triangle mesh boundary to a given mesh.

Returns a new mesh that contains a triangulated box around a given mesh suitable for geo-simulation (surface boundary at top).

Parameters:
mesh : mesh object

Mesh to which the triangle boundary should be appended.

xbound : float, optional

Horizontal prolongation distance. Minimal mesh 0.5 x extension.

ybound : float, optional

Vertical prolongation distance. Minimal mesh 0.5 y extension.

marker : int, optional

Marker of new cells.

markerBoundary : int, optional

Marker of the inner boundary edges between mesh and new boundary.

quality : float, optional

Triangle quality.

area: float, optional

Triangle max size within the boundary.

smooth : boolean, optional

Apply mesh smoothing.

isSubSurface : boolean, optional

Apply boundary conditions suitable for geo-simulation and prolongate mesh to the surface if necessary.

verbose : boolean, optional

Be verbose.

Examples

>>> import matplotlib.pyplot as plt
>>> import pygimli as pg
>>> from pygimli.mplviewer import drawMesh, drawModel
>>> from pygimli.meshtools import appendTriangleBoundary
>>> inner = pg.createGrid(range(5), range(5), marker=1)
>>> mesh = appendTriangleBoundary(inner, xbound=3, ybound=6, marker=2)
>>> fig, (ax1, ax2) = plt.subplots(1,2)
>>> p1 = drawMesh(ax1, inner)
>>> p2 = drawModel(ax2, mesh, mesh.cellMarkers(), label='marker')
>>> p3 = drawMesh(ax2, mesh)
>>> txt1 = ax1.set_title("a) Input grid")
>>> txt2 = ax2.set_title("b) With triangle boundary")

(png, pdf)

../../_images/pygimli-meshtools-1.png

Examples using pygimli.meshtools.appendTriangleBoundary

cellDataToBoundaryData

pygimli.meshtools.cellDataToBoundaryData(mesh, data)[source]

TODO DOCUMENT_ME

cellDataToNodeData

pygimli.meshtools.cellDataToNodeData(mesh, data, style='mean')[source]

Convert cell data to node data.

Convert cell data to node data via non-weighted averaging (mean) of common cell data.

Parameters:
mesh : GIMLI::Mesh

2D or 3D GIMLi mesh

data : iterable [float]

Data of len mesh.cellCount(). TODO complex, R3Vector, ndarray

style : str [‘mean’]

Interpolation style. * ‘mean’ : non-weighted averaging TODO harmonic averaging TODO weighted averaging (mean, harmonic) TODO interpolation via cell centered mesh

Examples

>>> import pygimli as pg
>>> grid = pg.createGrid(x=(1,2,3),y=(1,2,3))
>>> celldata = np.array([1, 2, 3, 4])
>>> nodedata = pg.meshtools.cellDataToNodeData(grid, celldata)
>>> print(nodedata.array())
[ 1.   1.5  2.   2.   2.5  3.   3.   3.5  4. ]

convertHDF5Mesh

pygimli.meshtools.convertHDF5Mesh(h5Mesh, group='mesh', indices='cell_indices', pos='coordinates', cells='topology', marker='values', marker_default=0, dimension=3, verbose=True, useFenicsIndices=False)[source]

Converts instance of a hdf5 mesh to a GIMLI::Mesh. For full documentation please see pygimli:meshtools:readHDF5Mesh.

createCircle

pygimli.meshtools.createCircle(pos=None, radius=1, segments=12, start=0, end=6.283185307179586, **kwargs)[source]

Create simple circle polygon.

Create simple circle polygon with given attributes.

Parameters:
pos : [x, y] [[0.0, 0.0]]

Center position

radius : float | [a,b] [1]

radius or halfaxes of the circle

segments : int

Discrete amount of segments for the circle

start : double [0]

Starting angle in radians

end : double [2*pi]

Ending angle in radians

**kwargs:
marker : int [1]

Marker for the resulting triangle cells after mesh generation

area : float [0]

Maximum cell size for resulting triangles after mesh generation

boundaryMarker : int [1]

Marker for the resulting boundary edges

leftDirection : bool [True]

Rotational direction

isHole : bool [False]

The Polygone will become a hole instead of a triangulation

isClosed : bool [True]

Add closing edge between last and first node.

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>> import matplotlib.pyplot as plt
>>> import math
>>> from pygimli.mplviewer import drawMesh
>>> from pygimli.meshtools import polytools as plc
>>> c0 = plc.createCircle(pos=(-5.0, 0.0), radius=2, segments=6)
>>> c1 = plc.createCircle(pos=(0.0, 0.0), segments=5, start=0, end=math.pi)
>>> c2 = plc.createCircle(pos=(5.0, 0.0), segments=3, start=math.pi,
...                       end=1.5*math.pi, isClosed=False)
>>> plc = plc.mergePLC([c0, c1, c2])
>>> fig, ax = plt.subplots()
>>> drawMesh(ax, plc, fillRegion=False)
>>> plt.show()

(png, pdf)

../../_images/pygimli-meshtools-2.png

Examples using pygimli.meshtools.createCircle

createCylinder

pygimli.meshtools.createCylinder(radius, height, nSegments=8, area=0.0, pos=None, marker=1)[source]

Create plc of a cylinder.

Out of core wrapper for dcfemlib::polytools.

Note, there is a bug in the old polytools which ignores the area settings for marker == 0.

Parameters:
radius : float

Radius of the cylinder.

height : float

Height of the cylinder

nSegments : int [8]

Number of segments of the cylinder.

area : float [0.0]

Largest size for the resulting tetrahedrons.

pos : pg.Pos [None]

The center position, default is at the origin.

marker : int [1]

Cell marker the resulting tetrahedrons.

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

createGrid

pygimli.meshtools.createGrid(x=None, y=None, z=None, **kwargs)[source]

Create grid style mesh.

Generate simple grid with defined node positions for each dimension. The resulting grid depends on the amount of given coordinate arguments and consists out of edges (1D - x), quads (2D- x and y), or hexahedrons(3D- x, y, and z).

Parameters:
kwargs :
x : array

x-coordinates for all Nodes (1D, 2D, 3D)

y : array

y-coordinates for all Nodes (2D, 3D)

z : array

z-coordinates for all Nodes (3D)

marker : int = 0

Marker for resulting cells.

worldBoundaryMarker : bool = False

Boundaries are enumerated with world marker, i.e., Top = -1 All remaining = -2. Default marker are left=1, right=2, top=3, bottom=4, front=5, back=6

Examples

>>> # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> mesh = pg.meshtools.createGrid(x=[0, 1, 1.5, 2],
...                                y=[-1, -0.5, -0.25, 0], marker=2)
>>> print(mesh)
Mesh: Nodes: 16 Cells: 9 Boundaries: 24
>>> _ = pg.show(mesh)
>>> pg.wait()

(png, pdf)

../../_images/pygimli-meshtools-3.png

createLine

pygimli.meshtools.createLine(start, end, segments=1, **kwargs)[source]

Create simple line polygon.

Create simple line polygon from start to end.

Parameters:
start : [x, y]

start position

end : [x, y]

end position

segments : int

Discrete amount of segments for the line

**kwargs:
boundaryMarker : int [1]

Marker for the resulting boundary edges

leftDirection : bool [True]

Rotational direction

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>>  # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> import pygimli.meshtools as mt
>>>
>>> w = mt.createWorld(start=[0, 0], end=[3, 3])
>>> l1 = mt.createLine(start=[1, 1], end=[1, 2], segments=1,
...                    leftDirection=False)
>>> l2 = mt.createLine(start=[1, 1], end=[2, 1], segments=20,
...                    leftDirection=True)
>>>
>>> ax, _ = pg.show(mt.createMesh([w, l1, l2,]))
>>> ax, _ = pg.show([w, l1, l2,], ax=ax, fillRegion=False)
>>> pg.wait()

(png, pdf)

../../_images/pygimli-meshtools-4.png

createMesh

pygimli.meshtools.createMesh(poly, quality=30, area=0.0, smooth=None, switches=None, verbose=False, **kwargs)[source]

Create a mesh for a given geometry polygon.

The mesh is created by triangle or tetgen if the gimli support for these mesh generators are installed. The geometry needs to contain nodes and boundaries and should be valid in the sense that the boundaries are non intersecting.

If poly is a list of coordinates a simple Delaunay mesh of the convex hull will be created. TODO: Tetgen support need to be implemented

Parameters:
poly: :gimliapi:`GIMLI::Mesh` or list or ndarray
  • 2D or 3D gimli mesh that contains the PLC.
  • 2D mesh needs edges
  • 3D mesh needs a plc and tetgen as system component
  • List of x y pairs [[x0, y0], … ,[xN, yN]]
  • ndarray [x_i, y_i]
  • PLC or list of PLC
quality: float

2D triangle quality sets a minimum angle constraint. Be careful with values above 34 degrees.

area: float

2D maximum triangle size in m*m

smooth: tuple

[smoothing algorithm, number of iterations] 0, no smoothing 1, node center 2, weighted node center

switches: str

Force triangle to use the gives command switches.

Returns:
mesh: :gimliapi:`GIMLI::Mesh`

Examples

>>> # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> import pygimli.meshtools as mt
>>> rect = mt.createRectangle(start=[0, 0], end=[4, 1])
>>> ax, _ = pg.show(mt.createMesh(rect, quality=10))
>>> ax, _ = pg.show(mt.createMesh(rect, quality=33))
>>> ax, _ = pg.show(mt.createMesh(rect, quality=33, area=0.01))
>>> pg.wait()
../../_images/pygimli-meshtools-5_00.png

(png, pdf)

../../_images/pygimli-meshtools-5_01.png

(png, pdf)

../../_images/pygimli-meshtools-5_02.png

(png, pdf)

Examples using pygimli.meshtools.createMesh

createParaDomain2D

pygimli.meshtools.createParaDomain2D(*args, **kwargs)[source]

API change here .. use createParaMeshPLC instead.

createParaMesh

pygimli.meshtools.createParaMesh(*args, **kwargs)[source]

Create parameter mesh from list of sensor positions.

Create parameter mesh from list of sensor positions.

Parameters:
sensors : list of RVector3 objects

Sensor positions. Must be sorted and unique in positive x direction. Depth need to be y-coordinate.

paraDX : float

Relativ distance for refinement nodes between two electrodes (1=none), e.g., 0.5 means 1 additional node between two neighboring electrodes e.g., 0.33 means 2 additional equidistant nodes between two electrodes

paraDepth : float, optional

Maximum depth for parametric domain, 0 (default) means 0.4 * maximum sensor range.

paraBoundary : float, optional

Margin for parameter domain in absolute sensor distances. 2 (default).

paraMaxCellSize: double [0], optional

Maximum size for parametric size in m*m (0-no constraint)

boundaryMaxCellSize: double [0], optional

Maximum cells size in the boundary region in m*m (0-no constraint)

boundary : float, optional

Boundary width to be appended for domain prolongation in absolute para domain width. Values lover 0 force the boundary to be 4 times para domain width.

Returns:
poly: :gimliapi:`GIMLI::Mesh`

createParaMesh2DGrid

pygimli.meshtools.createParaMesh2DGrid(sensors, paraDX=1, paraDZ=1, paraDepth=0, nLayers=11, boundary=-1, paraBoundary=2, **kwargs)[source]

Create a grid style mesh for an inversion parameter mesh.

Create a grid style mesh for an inversion parameter mesh. Return parameter grid for a given list of sensor positions. Uses and forwards arguments to pygimli.meshtools.appendTriangleBoundary.

Parameters:
sensors : list of RVector3 objects or data container with sensorPositions

Sensor positions. Must be sorted in positive x direction

paraDX : float, optional

Horizontal distance between sensors, relative regarding sensor distance. Value must be greater than 0 otherwise 1 is assumed.

paraDZ : float, optional

Vertical distance to the first depth layer, relative regarding sensor distance. Value must be greater than 0 otherwise 1 is assumed.

paraDepth : float, optional

Maximum depth for parametric domain, 0 (default) means 0.4 * maximum sensor range.

nLayers : int, optional [11]

Number of depth layers.

boundary : int, optional [-1]

Boundary width to be appended for domain prolongation in absolute para domain width. Values lower than 0 force the boundary to be 4 times para domain width.

paraBoundary : int, optional [2]

Offset to the parameter domain boundary in absolute sensor spacing.

Returns:
mesh: :gimliapi:`GIMLI::Mesh`

Examples

>>> import pygimli as pg
>>> import matplotlib.pyplot as plt
>>>
>>> from pygimli.meshtools import createParaMesh2DGrid
>>> mesh = createParaMesh2DGrid(sensors=pg.RVector(range(10)),
...                             boundary=1, paraDX=1,
...                             paraDZ=1, paraDepth=5)
>>> ax, _ = pg.show(mesh, markers=True, showMesh=True)

(png, pdf)

../../_images/pygimli-meshtools-6.png

createParaMeshPLC

pygimli.meshtools.createParaMeshPLC(sensors, paraDX=1, paraDepth=0, paraBoundary=2, paraMaxCellSize=0.0, boundary=-1, boundaryMaxCellSize=0, **kwargs)[source]

Create a PLC mesh for an inversion parameter mesh.

Create a PLC mesh for an inversion parameter mesh with for a given list of sensor positions. Sensor positions are assumed to lie on the surface and must be sorted and unique.

You can create a parameter mesh without sensors by setting [xmin, xmax] as sensors.

The PLC is a GIMLI::Mesh and contains nodes, edges and 2 region markers, one for the parameters domain (marker=2) and a larger boundary around the outside (marker=1)

TODO:

  • closed domains (boundary == 0)
  • additional topopoints
  • spline interpolations between sensorpoints or addpoints
  • subsurface sensors (partly .. see example)
Parameters:
sensors : [RVector3] | DataContainer with sensorPositions() | [xmin, xmax]

Sensor positions. Must be sorted and unique in positive x direction. Depth need to be y-coordinate.

paraDX : float [1]

Relativ distance for refinement nodes between two sensors (1=none), e.g., 0.5 means 1 additional node between two neighboring sensors e.g., 0.33 means 2 additional equidistant nodes between two sensors

paraDepth : float, optional

Maximum depth for parametric domain, 0 (default) means 0.4 * maximum sensor range.

paraBoundary : float, optional

Margin for parameter domain in absolute sensor distances. 2 (default).

paraMaxCellSize: double, optional

Maximum size for parametric size in m*m

boundaryMaxCellSize: double, optional

Maximum cells size in the boundary region in m*m

boundary : float, optional

Boundary width to be appended for domain prolongation in absolute para domain width. Values lover 0 force the boundary to be 4 times para domain width.

Returns:
poly: :gimliapi:`GIMLI::Mesh`

piecewise linear complex (PLC) containing nodes and edges

Examples

>>>  # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> import pygimli.meshtools as plc
>>> # Create the simplest paramesh PLC with a para box of 10 m without
>>> # sensors
>>> p = plc.createParaMeshPLC([0,10])
>>> # you can add subsurface sensors now with
>>> for z in range(1,4):
...     n = p.createNode((5,-z), -99)
>>> ax,_ = pg.show(p)

(png, pdf)

../../_images/pygimli-meshtools-7.png

createPolygon

pygimli.meshtools.createPolygon(verts, isClosed=False, isHole=False, **kwargs)[source]

Create a polygon from a list of vertices.

All vertices needs to be unique and duplicate vertices will be ignored. If you want the polygon be a closed region you can set the ‘isCloses’ flag. Closed region can be attributed by assigning a region marker. The automatic region marker is set in the center of all vertices.

Parameters:
verts : []
  • List of x y pairs [[x0, y0], … ,[xN, yN]]
**kwargs:
  • boundaryMarker : int [1]
    Marker for the resulting boundary edges
  • leftDirection : bool [True]
    Rotational direction
  • marker : int [None]
    Marker for the resulting triangle cells after mesh generation.
  • area : float [0]
    Maximum cell size for resulting triangles after mesh generation
  • isHole : bool [False]
    The Polygone will become a hole instead of a triangulation
isClosed : bool [True]

Add closing edge between last and first node.

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>>  # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> import pygimli.meshtools as mt
>>> p1 = mt.createPolygon([[0.0, 0.0], [1.0, 0.0], [1.0, 1.0]],
...                       isClosed=True, marker=3, area=0.1)
>>> p2 = mt.createPolygon([[0.3, 0.15], [0.85, 0.15], [0.85, 0.7]],
...                       isClosed=True, marker=3, area=0.1, isHole=True)
>>> ax, _ = pg.show(mt.mergePLC([p1,p2]))

(png, pdf)

../../_images/pygimli-meshtools-8.png

Examples using pygimli.meshtools.createPolygon

createRectangle

pygimli.meshtools.createRectangle(start=None, end=None, pos=None, size=None, **kwargs)[source]

Create rectangle polygon.

Create rectangle with start position and a given size. Give either start and end OR pos and size.

Parameters:
start : [x, y]

Left upper corner. Default [-0.5, 0.5]

end : [x, y]

Right lower corner. Default [0.5, -0.5]

pos : [x, y]

Center position. The rectangle will be moved.

size : [x, y]

width and height. The rectangle will be scaled.

**kwargs:
marker : int [1]

Marker for the resulting triangle cells after mesh generation

area : float [0]

Maximum cell size for resulting triangles after mesh generation

boundaryMarker : int [1]

Marker for the resulting boundary edges

leftDirection : bool [True]

TODO Rotational direction

isHole : bool [False]

The Polygone will become a hole instead of a triangulation

isClosed : bool [True]

Add closing edge between last and first node.

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>> import matplotlib.pyplot as plt
>>> import pygimli as pg
>>> from pygimli.meshtools import createRectangle
>>> rectangle = createRectangle(start=[4, -4], end=[6, -6],
...                             marker=4, area=0.1)
>>> _ = pg.show(rectangle)

(png, pdf)

../../_images/pygimli-meshtools-9.png

Examples using pygimli.meshtools.createRectangle

createWorld

pygimli.meshtools.createWorld(start, end, marker=1, area=0.0, layers=None, worldMarker=True)[source]

Create simple rectangular world.

Create simple rectangular world with appropriate boundary conditions. Surface boundary is set do pg.MARKER_BOUND_HOMOGEN_NEUMANN, i.e, -1 and inner subsurface is set to pg.MARKER_BOUND_MIXED, i.e., -2 or Numbered in ascending order in left direction starting upper left if worldMarker is set to false.

Parameters:
start : [x, y]

Upper/Left Corner

end : [x, y]

Lower/Right Corner

marker : int

Marker for the resulting triangle cells after mesh generation.

area : float | list

Maximum cell size for resulting triangles after mesh generation. If area is a float set it global, if area is a list set it per layer.

layers : [float]

Add some layers to the world.

worldMarker : [bool]

Specify kind of preset boundary marker [-1, -2] or [1, 2, 3, 4 ..]

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>> from pygimli.meshtools import createWorld
>>> from pygimli.mplviewer import drawMesh
>>> import matplotlib.pyplot as plt
>>> world = createWorld(start=[-5, 0], end=[5, -5], layers=[-1,-2,-3])
>>>
>>> fig, ax = plt.subplots()
>>> drawMesh(ax, world)
>>> plt.show()

(png, pdf)

../../_images/pygimli-meshtools-10.png

Examples using pygimli.meshtools.createWorld

exportFenicsHDF5Mesh

pygimli.meshtools.exportFenicsHDF5Mesh(mesh, exportname, **kwargs)[source]

Exports Gimli mesh in HDF5 format suitable for Fenics.

Equivalent to calling the function pygimli.meshtools.exportHDF5Mesh(mesh, exportname, group=['mesh', 'domains'], indices='cell_indices', pos='coordinates', cells='topology', marker='values').

Parameters:
mesh: :gimliapi:GIMLI::Mesh`

Mesh to be saved.

exportname: string

Name under which the mesh is saved.

exportHDF5Mesh

pygimli.meshtools.exportHDF5Mesh(mesh, exportname, group='mesh', indices='cell_indices', pos='coordinates', cells='topology', marker='values')[source]

Writes given GIMLI::Mesh in a hdf5 format file.

3D tetrahedron meshes only! Boundary markers are ignored.

Keywords are explained in pygimli.meshtools.readHDFS

exportPLC

pygimli.meshtools.exportPLC(poly, fname, **kwargs)[source]

General writer to save piece-wise linear complex (PLC) as poly file.

Choose from poly.dimension() and forward appropriate to GIMLI::Mesh::exportAsTetgenPolyFile and pygimli.meshtools.writeTrianglePoly

Parameters:
poly : GIMLI::Mesh

The polygon to be written.

fname : string

Filename of the file to write (\*.n, \*.e).

Examples

>>> import pygimli as pg
>>> import tempfile, os
>>> fname = tempfile.mktemp() # Create temporary string for filename.
>>> world2d = pg.meshtools.createWorld(start=[-10, 0], end=[20, -10])
>>> pg.meshtools.exportPLC(world2d, fname)
>>> read2d = pg.meshtools.readPLC(fname)
>>> print(read2d)
Mesh: Nodes: 4 Cells: 0 Boundaries: 4
>>> world3d = pg.createGrid([0, 1], [0, 1], [-1, 0])
>>> pg.meshtools.exportPLC(world3d, fname)
>>> os.remove(fname)

exportSTL

pygimli.meshtools.exportSTL(mesh, fileName, ascii=True)[source]

Write STL surface mesh and returns a GIMLI::Mesh.

Export a three dimensional boundary GIMLI::Mesh into a STL surface mesh. Boundaries with different marker will be separated into different STL solids.

TODO:
  • ASCII=False, write binary STL
  • QuadrangleFace Boundaries
  • p2 Boundaries
Parameters:
mesh : GIMLI::Mesh

Mesh to be exported. Only Boundaries of type TriangleFace will be exported.

fileName : str

name of the .stl file containing the STL surface mesh

ascii : bool [True]

STL Ascii format

fillEmptyToCellArray

pygimli.meshtools.fillEmptyToCellArray(mesh, vals, slope=True)[source]

Prolongate empty cell values to complete cell attributes.

It is possible that you have zero values that need to be filled with appropriate attributes. This function tries to fill the empty values successive prolongation of the non zeros.

Parameters:
mesh : GIMLI::Mesh

For each cell of mesh a value will be returned.

vals : array

Array of size cellCount().

Returns:
atts : array

Array of length mesh.cellCount()

interpolate

pygimli.meshtools.interpolate(*args, **kwargs)[source]

Interpolation convinience function.

Convenience function to interpolate different kind of data. Currently supported interpolation schemes are:

  • Interpolate mesh based data from one mesh to another

(syntactic sugar for the core based interpolate (see below))

Parameters:
args: GIMLI::Mesh, GIMLI::Mesh, iterable
outData = interpolate(outMesh, inMesh, vals) Interpolate values based on inMesh to outMesh. Values can be of length inMesh.cellCount() interpolated to outMesh.cellCenters() or inMesh.nodeCount() which are interpolated tp outMesh.positions().
Returns:
Interpolated values.
  • Mesh based values to arbitrary points, based on finite element interpolation (from gimli core).

    Parameters:
    args: GIMLI::Mesh, …

    Arguments forwarded to GIMLI::interpolate()

    kwargs:

    Arguments forwarded to GIMLI::interpolate()

    interpolate(srcMesh, destMesh)

    All data from inMesh are interpolated to outMesh

    Returns:

    Interpolated values

  • Interpolate along curve. Forwarded to pygimli.meshtools.interpolateAlongCurve

    Parameters:

    args: curve, t

    kwargs:

    Arguments forwarded to pygimli.meshtools.interpolateAlongCurve

  • 1D point set \(u(x)\) for ascending \(x\). Find interpolation function \(I = u(x)\) and returns \(u_{\text{i}} = I(x_{\text{i}})\) (interpolation methods are [linear via matplotlib, cubic spline via scipy, fit with harmonic functions’ via pygimli]) Note, for ‘linear’ and ‘spline’ the interpolate contains all original coordinates while ‘harmonic’ returns an approximate best fit. The amount of harmonic coefficients can be specfied with the ‘nc’ keyword.

    Parameters:
    args: xi, x, u
    • \(x_{\text{i}}\) - target sample points
    • \(x\) - function sample points
    • \(u\) - function values
    kwargs:
    • method : string
      Specify interpolation method ‘linear, ‘spline’, ‘harmonic’
    • nc : int
      Number of harmonic coefficients for the ‘harmonic’ method.
    Returns:
    ui: array of length xi

    \(u_{\text{i}} = I(x_{\text{i}})\), with \(I = u(x)\)

To use the core functions GIMLI::interpolate() start with a mesh instance as first argument or use the appropriate keyword arguments.

TODO

  • 2D parametric to points (method=[‘linear, ‘spline’, ‘harmonic’])
  • 2D/3D point cloud to points/grids (‘Delauney’, ‘linear, ‘spline’, ‘harmonic’)
  • Mesh to points based on nearest neighbour values (pg.core)

Examples

>>> # no need to import matplotlib. pygimli's show does
>>> import numpy as np
>>> import pygimli as pg
>>> fig, ax = pg.plt.subplots(1, 1, figsize=(10, 5))
>>> u = np.array([1.0, 12.0, 3.0, -4.0, 5.0, 6.0, -1.0])
>>> xu = np.array(range(len(u)))
>>> xi = np.linspace(xu[0], xu[-1], 1000)
>>> _= ax.plot(xu, u, 'o')
>>> _= ax.plot(xi, pg.interpolate(xi, xu, u, method='linear'),
...         color='blue', label='linear')
>>> _= ax.plot(xi, pg.interpolate(xi, xu, u, method='spline'),
...            color='red', label='spline')
>>> _= ax.plot(xi, pg.interpolate(xi, xu, u, method='harmonic'),
...         color='green', label='harmonic')
>>> _= ax.legend()

(png, pdf)

../../_images/pygimli-meshtools-11.png

interpolateAlongCurve

pygimli.meshtools.interpolateAlongCurve(curve, t, **kwargs)[source]

Interpolate along curve.

Return curve coordinates for a piecewise linear curve \(C(t) = {x_i,y_i,z_i}\) at positions \(t\). Curve and \(t\) values are expected to be sorted along distance from the origin of the curve.

Parameters:
curve : [[x,z]] | [[x,y,z]] | [GIMLI::RVector3] | GIMLI::R3Vector

Discrete curve for 2D \(x,z\) curve=[[x,z]], 3D \(x,y,z\)

t : 1D iterable

Query positions along the curve in absolute distance

kwargs :

If kwargs are given an additional curve smoothing is applied using pygimli.meshtools.interpolate. The kwargs will be delegated.

Returns:
p : np.array

Curve positions at query points \(t\). Dimension of p match the size of curve the coordinates.

Examples

>>> # no need to import matplotlib. pygimli's show does
>>> import numpy as np
>>> import pygimli as pg
>>> import pygimli.meshtools as mt
>>> fig, axs = pg.plt.subplots(2,2)
>>> topo = np.array([[-2., 0.], [-1., 0.], [0.5, 0.], [3., 2.], [4., 2.], [6., 1.], [10., 1.], [12., 1.]])
>>> t = np.arange(15.0)
>>> p = mt.interpolateAlongCurve(topo, t)
>>> _= axs[0,0].plot(topo[:,0], topo[:,1], '-x', mew=2)
>>> _= axs[0,1].plot(p[:,0], p[:,1], 'o', color='red') 
>>>
>>> p = mt.interpolateAlongCurve(topo, t, method='spline')
>>> _= axs[1,0].plot(p[:,0], p[:,1], '-o', color='black') 
>>>
>>> p = mt.interpolateAlongCurve(topo, t, method='harmonic', nc=3)
>>> _= axs[1,1].plot(p[:,0], p[:,1], '-o', color='green') 
>>>
>>> pg.plt.show()
>>> pg.wait()

(png, pdf)

../../_images/pygimli-meshtools-12.png

merge2Meshes

pygimli.meshtools.merge2Meshes(m1, m2)[source]

Merge two meshes into one new mesh and return the combined mesh.

Merge two meshes into a new mesh and return the combined mesh. Note, there is a duplication check for all nodes which should reuse existing node but NO cells or boundaries.

Parameters:
m1: :gimliapi:`GIMLI::Mesh`

First mesh.

m2: :gimliapi:`GIMLI::Mesh`

Second mesh.

Returns:
mesh: :gimliapi:`GIMLI::Mesh`

Resulting mesh.

Examples using pygimli.meshtools.merge2Meshes

mergeMeshes

pygimli.meshtools.mergeMeshes(meshList, verbose=False)[source]

Merge several meshes into one new mesh and return the new mesh.

Merge several meshes into one new mesh and return the new mesh.

Parameters:
meshList : [GIMLI::Mesh, …] | [str, …]

List of at least two meshes or (filenames to meshes) to be merged.

verbose : bool

Give some output

See also

merge2Meshes

mergePLC

pygimli.meshtools.mergePLC(pols, tol=0.001)[source]

Merge multiply polygons.

Merge multiply polygons into a single polygon. Common nodes and common edges will be checked and removed. When a node touches and edge the edge will be split.

TODO:
  • Crossing or Node/Edge intersections will NOT be

recognized yet.

Parameters:
pols: [:gimliapi:`GIMLI::Mesh`]

List of polygons that need to be merged

tol : double

Tolerance to check for duplicated nodes. [1e-3]

Returns:
poly : GIMLI::Mesh

The resulting polygon is a GIMLI::Mesh.

Examples

>>> from pygimli.meshtools import polytools as plc
>>> from pygimli.meshtools import createMesh
>>> from pygimli.mplviewer import drawMesh
>>> import matplotlib.pyplot as plt
>>> world = plc.createWorld(start=[-10, 0], end=[10, -10], marker=1)
>>> c1 = plc.createCircle([-1, -4], radius=1.5, area=0.1,
...                       marker=2, segments=4)
>>> c2 = plc.createCircle([-6, -5], radius=[1.5, 3.5], isHole=1)
>>> r1 = plc.createRectangle(pos=[3, -5], size=[2, 2], marker=3)
>>> r2 = plc.createRectangle(start=[4, -4], end=[6, -6],
...                          marker=4, area=0.1)
>>> plc = plc.mergePLC([world, c1, c2, r1, r2])
>>> fig, ax = plt.subplots()
>>> drawMesh(ax, plc)
>>> drawMesh(ax, createMesh(plc))
>>> plt.show()

(png, pdf)

../../_images/pygimli-meshtools-13.png

Examples using pygimli.meshtools.mergePLC

nodeDataToBoundaryData

pygimli.meshtools.nodeDataToBoundaryData(mesh, data)[source]

Assuming [NodeCount, dim] data DOCUMENT_ME

nodeDataToCellData

pygimli.meshtools.nodeDataToCellData(mesh, data)[source]

Convert node data to cell data.

Convert node data to cell data via interpolation to cell centers.

Parameters:
mesh : GIMLI::Mesh

2D or 3D GIMLi mesh

data : iterable [float]

Data of len mesh.nodeCount(). TODO complex, R3Vector, ndarray

quality

pygimli.meshtools.quality(mesh, measure='eta')[source]

Return the quality of a given triangular mesh.

Parameters:
mesh : mesh object

Mesh for which the quality is calculated.

measure : quality measure, str

Can be either “eta”, “nsr”, or “minimumAngle”.

See also

eta, nsr, minimumAngle

Examples

>>> # no need to import matplotlib
>>> import pygimli as pg
>>> from pygimli.meshtools import polytools as plc
>>> from pygimli.meshtools import quality
>>> # Create Mesh
>>> world = plc.createWorld(start=[-10, 0], end=[10, -10],
...                         marker=1, worldMarker=False)
>>> c1 = plc.createCircle(pos=[0.0, -5.0], radius=3.0, area=.3)
>>> mesh = pg.meshtools.createMesh([world, c1], quality=21.3)
>>> # Compute and show quality
>>> q = quality(mesh, measure="nsr")
>>> ax, _ = pg.show(mesh, q, cMap="RdYlGn", showMesh=True, cMin=0.5, cMax=1.0,
...                 label="Normalized shape ratio")

(png, pdf)

../../_images/pygimli-meshtools-14.png

Examples using pygimli.meshtools.quality

readFenicsHDF5Mesh

pygimli.meshtools.readFenicsHDF5Mesh(fileName, group='mesh', verbose=True)[source]

Reads FEniCS mesh from file format .h5 and returns a GIMLI::Mesh.

readGmsh

pygimli.meshtools.readGmsh(fname, verbose=False)[source]

Read Gmsh ASCII file and return instance of GIMLI::Mesh class.

Parameters:
fname : string

Filename of the file to read (\*.msh). The file must conform to the MSH ASCII file version 2 format

verbose : boolean, optional

Be verbose during import.

Notes

Physical groups specified in Gmsh are interpreted as follows:

  • Points with the physical number 99 are interpreted as sensors.
  • Physical Lines and Surfaces define boundaries in 2D and 3D, respectively.
    • Physical Number 1: Homogeneous Neumann condition
    • Physical Number 2: Mixed boundary condition
    • Physical Number 3: Homogeneous Dirichlet condition
    • Physical Number 4: Dirichlet condition
  • Physical Surfaces and Volumes define regions in 2D and 3D, respectively.
    • Physical Number 1: No inversion region
    • Physical Number >= 2: Inversion region

Examples

>>> import tempfile, os
>>> from pygimli.meshtools import readGmsh
>>> gmsh = '''
... $MeshFormat
... 2.2 0 8
... $EndMeshFormat
... $Nodes
... 3
... 1 0 0 0
... 2 0 1 0
... 3 1 1 0
... $EndNodes
... $Elements
... 7
... 1 15 2 0 1 1
... 2 15 2 0 2 2
... 3 15 2 0 3 3
... 4 1 2 0 1 2 3
... 5 1 2 0 2 3 1
... 6 1 2 0 3 1 2
... 7 2 2 0 5 1 2 3
... $EndElements
... '''
>>> fname = tempfile.mktemp()
>>> with open(fname, "w") as f:
...     f.writelines(gmsh)
>>> mesh = readGmsh(fname)
>>> print(mesh)
Mesh: Nodes: 3 Cells: 1 Boundaries: 3
>>> os.remove(fname)

Examples using pygimli.meshtools.readGmsh

readHDF5Mesh

pygimli.meshtools.readHDF5Mesh(fileName, group='mesh', indices='cell_indices', pos='coordinates', cells='topology', marker='values', marker_default=0, dimension=3, verbose=True, useFenicsIndices=False)[source]

Function for loading a mesh from HDF5 file format.

Returns an instance of GIMLI::Mesh class. Default values for keywords are suited for FEniCS syntax .h5 meshes.

Requirements: h5py module

TODO:
  • Fenics hdf5 meshs doesn’t have boundary markers.
Parameters:
fileName: string

Name of the mesh that has to be transformed into pyGIMLi format.

group: string [‘domains’]

hdf group that contains the mesh informations (see other keyword arguments). Default is ‘domains’ for FEniCS compatibility.

indices: string [‘cell_indices’]

Key for the part of the hdf file containing the indices of the cells.

pos: string [‘coordinates’]

Key for the part of the hdf file containing the nodepositions.

cells: string [‘topology’]

Key for the part of the hdf file containing the array which defies the cells. Usually of shape (cellCount, 3) for 2D meshs or (cellCount, 4) for 3D tetrahedra meshes. For each cell the indices of the corresponding node indices is given.

marker: string [‘values’]

If marker is part of the hdf data container, the corresponding array is used as identifier for the cell markers. If not found, the cell markers will be set to marker_default.

marker_default: int or array [0]

Default marker if no markers are found in the hdf file. If array, size has to match the cellCount of the mesh.

dimension: int [3]

Dimension of the in/outpu mesh, no own check for dimensions yet. Fixed on 3 for now.

Returns:
mesh:

GIMLI::Mesh

readHydrus2dMesh

pygimli.meshtools.readHydrus2dMesh(fname='MESHTRIA.TXT')[source]

Import mesh from Hydrus2D.

Parameters:
fname : str, optional

Filename of Hydrus output file.

See also

readHydrus3dMesh
Similar routine for three-dimensional meshes.

References

readHydrus3dMesh

pygimli.meshtools.readHydrus3dMesh(fileName='MESHTRIA.TXT')[source]

Import mesh from Hydrus3D.

Parameters:
fname : str, optional

Filename of Hydrus output file.

See also

readHydrus2dMesh
Similar routine for two-dimensional meshes.

References

readPLC

pygimli.meshtools.readPLC(filename, comment='#')[source]

Generic PLC reader.

Read 2D Triangle or 3D Tetgen PLC files.

Parameters:
filename: string

Filename *.poly

comment: string (‘#’)

String containing all characters that define a comment line. Identified lines will be ignored during import.

Returns:
poly :

GIMLI::Mesh

readSTL

pygimli.meshtools.readSTL(fileName, ascii=True)[source]

Read STL surface mesh and returns a GIMLI::Mesh.

Read STL surface mesh and returns a GIMLI::Mesh of triangle boundary faces. Multiple solids are supported with increasing boundary marker.

TODO: ASCII=False, read binary STL

Parameters:
fileName : str

name of the .stl file containing the STL surface mesh

ascii : bool [True]

STL Ascii format

readTetgen

pygimli.meshtools.readTetgen(fname, comment='#', verbose=True, default_cell_marker=0, load_faces=True, quadratic=False)[source]

Reads and converts a mesh from the basic Tetgen output.

Read Tetgen [Si04] ASCII files and return instance of GIMLI::Mesh class. See: http://tetgen.org/

Parameters:
fname: str

Base name of the tetgen output, without ending. All additional files (.ele and .face respectively) has to be the same basename.

comment: str (‘#’)

String consisting of all symbols that indicates a comment in the input files. Standard for tetgen files is the ‘#’.

verbose: boolean (True)

Enables console output during the import process.

default_cell_marker: int (0)

Tetgen files can contain cell markers, but doesn’t have to. If no marker are found, the given interger is used.

load_faces:

Optional decision weather the faces of the Tetgen output are loaded or not. Note that without the -f in during the tetgen call, the faces in the .face file will only contain the faces of the original input poly file and not all faces. If only part of the faces are imported a createNeighbourInfos call of the mesh will fail.

quadratic: boolean (False)

Returns a P2 refined mesh when True (to be removed, as soon as direct import of quadratic meshs is possible).

Returns:
mesh: :gimliapi:`GIMLI::Mesh`

readTriangle

pygimli.meshtools.readTriangle(fname, verbose=False)[source]

Read Triangle [She96] mesh.

Read Triangle [She96] ASCII mesh files and return an instance of GIMLI::Mesh class. See: ://www.cs.cmu.edu/~quake/triangle.html

Parameters:
fname : string

Filename of the file to read (\*.n, \*.e)

verbose : boolean, optional

Be verbose during import.

refineQuad2Tri

pygimli.meshtools.refineQuad2Tri(mesh, style=1)[source]

Refine mesh of quadrangles into a mesh of triangle cells.

TODO mixed meshes
Parameters:
mesh : GIMLI::Mesh

Mesh containing quadrangle cells.

style: int [1]
  • 1 bisect each quadrangle into 2 triangles
  • 2 bisect each quadrangle into 4 triangles
Returns:
ret : GIMLI::Mesh

Mesh containing triangle cells.

Examples

>>> # no need to import matplotlib. pygimli's show does
>>> import pygimli as pg
>>> import pygimli.meshtools as mt
>>> quads = pg.createGrid(range(10), range(10))
>>> ax, _ = pg.show(quads)
>>> ax, _ = pg.show(mt.refineQuad2Tri(quads, style=1))
>>> ax, _ = pg.show(mt.refineQuad2Tri(quads, style=2))
>>> pg.wait()
../../_images/pygimli-meshtools-15_00.png

(png, pdf)

../../_images/pygimli-meshtools-15_01.png

(png, pdf)

../../_images/pygimli-meshtools-15_02.png

(png, pdf)

syscallTetgen

pygimli.meshtools.syscallTetgen(filename, quality=1.2, area=0, preserveBoundary=False, verbose=False)[source]

Create a mesh with Tetgen from file.

Create a Tetgen [Si04] mesh from a PLC.

Forwards to system call tetgen, which must be known to your system.

Parameters:
filename: str
quality: float [1.2]

Refines mesh (to improve mesh quality). [1.1 … ]

area: float [0.0]

Maximum cell size (m³)

preserveBoundary: bool [False]

Preserve PLC boundary mesh

verbose: bool [False]

be verbose

Returns:
mesh : GIMLI::Mesh

tapeMeasureToCoordinates

pygimli.meshtools.tapeMeasureToCoordinates(tape, pos)[source]

Interpolate 2D tape measured topography to 2D Cartesian coordinates.

Tape and pos value are expected to be sorted along distance to the origin.

DEPRECATED will be removed use pygimli.meshtools.interpolateAlongCurve instead

TODO optional smooth curve with harmfit TODO parametric TODO parametric + Topo: 3d

Parameters:
tape : [[x,z]] | [RVector3] | R3Vector

List of tape measured topography points with measured distance (x) from origin and height (z)

pos : iterable

Array of query positions along the tape measured profile t[0 ..

Returns:
res : ndarray(N, 2)

Same as pos but with interpolated height values. The Distance between pos points and res (along curve) points remains.

writePLC

pygimli.meshtools.writePLC(*args, **kwargs)[source]

Backward compatibility. Please use pygimli.meshtools.exportPLC.



2018 - GIMLi Development Team
Created using Bootstrap, Sphinx and pyGIMLi 1.0.9+87.g68804698 on Dec 16, 2018.