fipy.meshes package

Submodules

fipy.meshes.abstractMesh module

class fipy.meshes.abstractMesh.AbstractMesh(communicator, _RepresentationClass=<class 'fipy.meshes.representations.abstractRepresentation._AbstractRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.abstractTopology._AbstractTopology'>)

Bases: object

A class encapsulating all commonalities among meshes in FiPy.

property VTKCellDataSet

Returns a TVTK DataSet representing the cells of this mesh

property VTKFaceDataSet

Returns a TVTK DataSet representing the face centers of this mesh

__add__(other)

Either translate a Mesh or concatenate two Mesh objects.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

If a vector is added to a Mesh, a translated Mesh is returned

>>> translatedMesh = baseMesh + ((5,), (10,))
>>> print(translatedMesh.cellCenters)
[[  5.5   6.5   5.5   6.5]
 [ 10.5  10.5  11.5  11.5]]

If a Mesh is added to a Mesh, a concatenation of the two Mesh objects is returned

>>> addedMesh = baseMesh + (baseMesh + ((2,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  2.5  3.5  2.5  3.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]

The two Mesh objects need not be properly aligned in order to concatenate them but the resulting mesh may not have the intended connectivity

>>> addedMesh = baseMesh + (baseMesh + ((3,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  3.5  4.5  3.5  4.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]
>>> addedMesh = baseMesh + (baseMesh + ((2,), (2,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  2.5  3.5  2.5  3.5]
 [ 0.5  0.5  1.5  1.5  2.5  2.5  3.5  3.5]]

No provision is made to avoid or consolidate overlapping Mesh objects

>>> addedMesh = baseMesh + (baseMesh + ((1,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  1.5  2.5  1.5  2.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]

Different Mesh classes can be concatenated

>>> from fipy.meshes import Tri2D
>>> triMesh = Tri2D(dx = 1.0, dy = 1.0, nx = 2, ny = 1)
>>> triMesh = triMesh + ((2,), (0,))
>>> triAddedMesh = baseMesh + triMesh
>>> cellCenters = [[0.5, 1.5, 0.5, 1.5, 2.83333333,  3.83333333,
...                 2.5, 3.5, 2.16666667, 3.16666667, 2.5, 3.5],
...                [0.5, 0.5, 1.5, 1.5, 0.5, 0.5, 0.83333333, 0.83333333,
...                 0.5, 0.5, 0.16666667, 0.16666667]]
>>> print(numerix.allclose(triAddedMesh.cellCenters,
...                        cellCenters))
True

again, their faces need not align, but the mesh may not have the desired connectivity

>>> triMesh = Tri2D(dx = 1.0, dy = 2.0, nx = 2, ny = 1)
>>> triMesh = triMesh + ((2,), (0,))
>>> triAddedMesh = baseMesh + triMesh
>>> cellCenters = [[ 0.5, 1.5, 0.5, 1.5, 2.83333333, 3.83333333,
...                  2.5, 3.5, 2.16666667, 3.16666667, 2.5, 3.5],
...                [ 0.5, 0.5, 1.5, 1.5, 1., 1.,
...                  1.66666667, 1.66666667, 1., 1., 0.33333333, 0.33333333]]
>>> print(numerix.allclose(triAddedMesh.cellCenters,
...                        cellCenters))
True

Mesh concatenation is not limited to 2D meshes

>>> from fipy.meshes import Grid3D
>>> threeDBaseMesh = Grid3D(dx = 1.0, dy = 1.0, dz = 1.0,
...                         nx = 2, ny = 2, nz = 2)
>>> threeDSecondMesh = Grid3D(dx = 1.0, dy = 1.0, dz = 1.0,
...                           nx = 1, ny = 1, nz = 1)
>>> threeDAddedMesh = threeDBaseMesh + (threeDSecondMesh + ((2,), (0,), (0,)))
>>> print(threeDAddedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  0.5  1.5  0.5  1.5  2.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5  0.5]
 [ 0.5  0.5  0.5  0.5  1.5  1.5  1.5  1.5  0.5]]

but the different Mesh objects must, of course, have the same dimensionality.

>>> InvalidMesh = threeDBaseMesh + baseMesh 
Traceback (most recent call last):
...
MeshAdditionError: Dimensions do not match
__dict__ = mappingproxy({'__module__': 'fipy.meshes.abstractMesh', '__doc__': '\n A class encapsulating all commonalities among meshes in FiPy.\n ', '__init__': <function AbstractMesh.__init__>, '_setTopology': <function AbstractMesh._setTopology>, '_setGeometry': <function AbstractMesh._setGeometry>, '_setScale': <function AbstractMesh._setScale>, 'scale': <property object>, '_calcScaleArea': <function AbstractMesh._calcScaleArea>, '_calcScaleVolume': <function AbstractMesh._calcScaleVolume>, '_getPointToCellDistances': <function AbstractMesh._getPointToCellDistances>, 'getNearestCell': <function AbstractMesh.getNearestCell>, '_getCellFaceIDsInternal': <function AbstractMesh._getCellFaceIDsInternal>, '_setCellFaceIDsInternal': <function AbstractMesh._setCellFaceIDsInternal>, 'cellFaceIDs': <property object>, 'interiorFaces': <property object>, '_setExteriorFaces': <function AbstractMesh._setExteriorFaces>, 'exteriorFaces': <property object>, '_isOrthogonal': <property object>, 'faceCenters': <property object>, 'cellToFaceDistanceVectors': <property object>, 'cellDistanceVectors': <property object>, 'cellVolumes': <property object>, 'cellCenters': <property object>, 'x': <property object>, 'y': <property object>, 'z': <property object>, 'extents': <property object>, 'scaledFaceAreas': <property object>, 'scaledCellVolumes': <property object>, 'scaledFaceToCellDistances': <property object>, 'scaledCellDistances': <property object>, 'scaledCellToCellDistances': <property object>, '_connectFaces': <function AbstractMesh._connectFaces>, '_concatenableMesh': <property object>, '_translate': <function AbstractMesh._translate>, '_getAddedMeshValues': <function AbstractMesh._getAddedMeshValues>, 'interiorFaceIDs': <property object>, 'interiorFaceCellIDs': <property object>, '_numberOfFacesPerCell': <property object>, '_maxFacesPerCell': <property object>, '_numberOfVertices': <property object>, '_globalNonOverlappingCellIDs': <property object>, '_globalOverlappingCellIDs': <property object>, '_localNonOverlappingCellIDs': <property object>, '_localOverlappingCellIDs': <property object>, '_globalNonOverlappingFaceIDs': <property object>, '_globalOverlappingFaceIDs': <property object>, '_localNonOverlappingFaceIDs': <property object>, '_localOverlappingFaceIDs': <property object>, '_vertexCellIDs': <property object>, '_vertexFaceIDs': <property object>, 'facesLeft': <property object>, 'facesRight': <property object>, 'facesBottom': <property object>, 'facesDown': <property object>, 'facesTop': <property object>, 'facesUp': <property object>, 'facesBack': <property object>, 'facesFront': <property object>, '_cellVertexIDs': <property object>, '_calcOrderedCellVertexIDs': <function AbstractMesh._calcOrderedCellVertexIDs>, '_orderedCellVertexIDs': <property object>, '_cellDistanceNormals': <property object>, '_cellAreaProjections': <property object>, '_concatenatedClass': <property object>, '__add__': <function AbstractMesh.__add__>, '__radd__': <function AbstractMesh.__add__>, '__mul__': <function AbstractMesh.__mul__>, '__rmul__': <function AbstractMesh.__mul__>, '__sub__': <function AbstractMesh.__sub__>, '__truediv__': <function AbstractMesh.__truediv__>, '__div__': <function AbstractMesh.__truediv__>, '__getstate__': <function AbstractMesh.__getstate__>, '__setstate__': <function AbstractMesh.__setstate__>, '__repr__': <function AbstractMesh.__repr__>, '_VTKCellType': <property object>, 'VTKCellDataSet': <property object>, 'VTKFaceDataSet': <property object>, '_toVTK3D': <function AbstractMesh._toVTK3D>, '_makePeriodic': <function AbstractMesh._makePeriodic>, '_facesPerCell': <property object>, '_cellFaceVertices': <property object>, '_unsortedNodesPerFace': <property object>, '_nodesPerFace': <property object>, '_faceOrder': <property object>, '_cellTopology': <property object>, 'aspect2D': <property object>, '_minmin': <function AbstractMesh._minmin>, '_maxmax': <function AbstractMesh._maxmax>, '__dict__': <attribute '__dict__' of 'AbstractMesh' objects>, '__weakref__': <attribute '__weakref__' of 'AbstractMesh' objects>, '__annotations__': {}})
__div__(other)

Tests. >>> from fipy import * >>> print((Grid1D(nx=1) / 2.).cellCenters) [[ 0.25]] >>> AbstractMesh(communicator=None) / 2. Traceback (most recent call last): … NotImplementedError

__getstate__()
__init__(communicator, _RepresentationClass=<class 'fipy.meshes.representations.abstractRepresentation._AbstractRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.abstractTopology._AbstractTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.abstractMesh'
__mul__(other)
__radd__(other)

Either translate a Mesh or concatenate two Mesh objects.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

If a vector is added to a Mesh, a translated Mesh is returned

>>> translatedMesh = baseMesh + ((5,), (10,))
>>> print(translatedMesh.cellCenters)
[[  5.5   6.5   5.5   6.5]
 [ 10.5  10.5  11.5  11.5]]

If a Mesh is added to a Mesh, a concatenation of the two Mesh objects is returned

>>> addedMesh = baseMesh + (baseMesh + ((2,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  2.5  3.5  2.5  3.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]

The two Mesh objects need not be properly aligned in order to concatenate them but the resulting mesh may not have the intended connectivity

>>> addedMesh = baseMesh + (baseMesh + ((3,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  3.5  4.5  3.5  4.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]
>>> addedMesh = baseMesh + (baseMesh + ((2,), (2,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  2.5  3.5  2.5  3.5]
 [ 0.5  0.5  1.5  1.5  2.5  2.5  3.5  3.5]]

No provision is made to avoid or consolidate overlapping Mesh objects

>>> addedMesh = baseMesh + (baseMesh + ((1,), (0,)))
>>> print(addedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  1.5  2.5  1.5  2.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]]

Different Mesh classes can be concatenated

>>> from fipy.meshes import Tri2D
>>> triMesh = Tri2D(dx = 1.0, dy = 1.0, nx = 2, ny = 1)
>>> triMesh = triMesh + ((2,), (0,))
>>> triAddedMesh = baseMesh + triMesh
>>> cellCenters = [[0.5, 1.5, 0.5, 1.5, 2.83333333,  3.83333333,
...                 2.5, 3.5, 2.16666667, 3.16666667, 2.5, 3.5],
...                [0.5, 0.5, 1.5, 1.5, 0.5, 0.5, 0.83333333, 0.83333333,
...                 0.5, 0.5, 0.16666667, 0.16666667]]
>>> print(numerix.allclose(triAddedMesh.cellCenters,
...                        cellCenters))
True

again, their faces need not align, but the mesh may not have the desired connectivity

>>> triMesh = Tri2D(dx = 1.0, dy = 2.0, nx = 2, ny = 1)
>>> triMesh = triMesh + ((2,), (0,))
>>> triAddedMesh = baseMesh + triMesh
>>> cellCenters = [[ 0.5, 1.5, 0.5, 1.5, 2.83333333, 3.83333333,
...                  2.5, 3.5, 2.16666667, 3.16666667, 2.5, 3.5],
...                [ 0.5, 0.5, 1.5, 1.5, 1., 1.,
...                  1.66666667, 1.66666667, 1., 1., 0.33333333, 0.33333333]]
>>> print(numerix.allclose(triAddedMesh.cellCenters,
...                        cellCenters))
True

Mesh concatenation is not limited to 2D meshes

>>> from fipy.meshes import Grid3D
>>> threeDBaseMesh = Grid3D(dx = 1.0, dy = 1.0, dz = 1.0,
...                         nx = 2, ny = 2, nz = 2)
>>> threeDSecondMesh = Grid3D(dx = 1.0, dy = 1.0, dz = 1.0,
...                           nx = 1, ny = 1, nz = 1)
>>> threeDAddedMesh = threeDBaseMesh + (threeDSecondMesh + ((2,), (0,), (0,)))
>>> print(threeDAddedMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5  0.5  1.5  0.5  1.5  2.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5  0.5]
 [ 0.5  0.5  0.5  0.5  1.5  1.5  1.5  1.5  0.5]]

but the different Mesh objects must, of course, have the same dimensionality.

>>> InvalidMesh = threeDBaseMesh + baseMesh 
Traceback (most recent call last):
...
MeshAdditionError: Dimensions do not match
__repr__()

Return repr(self).

__rmul__(other)
__setstate__(state)
__sub__(other)

Tests. >>> from fipy import * >>> m = Grid1D() >>> print((m - ((1,))).cellCenters) [[-0.5]] >>> ((1,)) - m Traceback (most recent call last): … TypeError: unsupported operand type(s) for -: ‘tuple’ and ‘UniformGrid1D’

__truediv__(other)

Tests. >>> from fipy import * >>> print((Grid1D(nx=1) / 2.).cellCenters) [[ 0.25]] >>> AbstractMesh(communicator=None) / 2. Traceback (most recent call last): … NotImplementedError

__weakref__

list of weak references to the object (if defined)

property aspect2D

The physical y vs x aspect ratio of a 2D mesh

property cellCenters
property cellDistanceVectors
property cellFaceIDs

Topology properties

property cellToFaceDistanceVectors
property cellVolumes
property extents
property exteriorFaces
property faceCenters
property facesBack

Return list of faces on back boundary of 3D Mesh with the z-axis running from front to back.

>>> from fipy import Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((6, 7, 8, 9, 10, 11),
...                        numerix.nonzero(mesh.facesBack)[0])) 
True
>>> ignore = mesh.facesBack.value 
property facesBottom

Return list of faces on bottom boundary of 2D or 3D Mesh with the y-axis running from bottom to top.

>>> from fipy import Grid2D, Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((12, 13, 14),
...                        numerix.nonzero(mesh.facesBottom)[0])) 
True
>>> ignore = mesh.facesBottom.value 
>>> x, y, z = mesh.faceCenters
>>> print(numerix.allequal((12, 13),
...                        numerix.nonzero(mesh.facesBottom & (x < 1))[0])) 
True
>>> ignore = mesh.facesBottom.value 
property facesDown

Return list of faces on bottom boundary of 2D or 3D Mesh with the y-axis running from bottom to top.

>>> from fipy import Grid2D, Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((12, 13, 14),
...                        numerix.nonzero(mesh.facesBottom)[0])) 
True
>>> ignore = mesh.facesBottom.value 
>>> x, y, z = mesh.faceCenters
>>> print(numerix.allequal((12, 13),
...                        numerix.nonzero(mesh.facesBottom & (x < 1))[0])) 
True
>>> ignore = mesh.facesBottom.value 
property facesFront

Return list of faces on front boundary of 3D Mesh with the z-axis running from front to back.

>>> from fipy import Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((0, 1, 2, 3, 4, 5),
...                        numerix.nonzero(mesh.facesFront)[0])) 
True
>>> ignore = mesh.facesFront.value 
property facesLeft

Return face on left boundary of Mesh as list with the x-axis running from left to right.

>>> from fipy import Grid2D, Grid3D
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((21, 25),
...                        numerix.nonzero(mesh.facesLeft)[0])) 
True
>>> ignore = mesh.facesLeft.value 
>>> mesh = Grid2D(nx = 3, ny = 2, dx = 0.5, dy = 2.)
>>> print(numerix.allequal((9, 13),
...                        numerix.nonzero(mesh.facesLeft)[0])) 
True
>>> ignore = mesh.facesLeft.value 
property facesRight

Return list of faces on right boundary of Mesh with the x-axis running from left to right.

>>> from fipy import Grid2D, Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((24, 28),
...                        numerix.nonzero(mesh.facesRight)[0])) 
True
>>> ignore = mesh.facesRight.value 
>>> mesh = Grid2D(nx = 3, ny = 2, dx = 0.5, dy = 2.)
>>> print(numerix.allequal((12, 16),
...                        numerix.nonzero(mesh.facesRight)[0])) 
True
>>> ignore = mesh.facesRight.value 
property facesTop

Return list of faces on top boundary of 2D or 3D Mesh with the y-axis running from bottom to top.

>>> from fipy import Grid2D, Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((18, 19, 20),
...                        numerix.nonzero(mesh.facesTop)[0])) 
True
>>> ignore = mesh.facesTop.value 
>>> mesh = Grid2D(nx = 3, ny = 2, dx = 0.5, dy = 2.)
>>> print(numerix.allequal((6, 7, 8),
...                        numerix.nonzero(mesh.facesTop)[0])) 
True
>>> ignore = mesh.facesTop.value 
property facesUp

Return list of faces on top boundary of 2D or 3D Mesh with the y-axis running from bottom to top.

>>> from fipy import Grid2D, Grid3D, numerix
>>> mesh = Grid3D(nx = 3, ny = 2, nz = 1, dx = 0.5, dy = 2., dz = 4.)
>>> print(numerix.allequal((18, 19, 20),
...                        numerix.nonzero(mesh.facesTop)[0])) 
True
>>> ignore = mesh.facesTop.value 
>>> mesh = Grid2D(nx = 3, ny = 2, dx = 0.5, dy = 2.)
>>> print(numerix.allequal((6, 7, 8),
...                        numerix.nonzero(mesh.facesTop)[0])) 
True
>>> ignore = mesh.facesTop.value 
getNearestCell(point)
property interiorFaceCellIDs
property interiorFaceIDs
property interiorFaces
property scale
property scaledCellDistances
property scaledCellToCellDistances
property scaledCellVolumes
property scaledFaceAreas
property scaledFaceToCellDistances
property x

Equivalent to using cellCenters[0].

>>> from fipy import *
>>> print(Grid1D(nx=2).x)
[ 0.5  1.5]
property y

Equivalent to using cellCenters[1].

>>> from fipy import *
>>> print(Grid2D(nx=2, ny=2).y)
[ 0.5  0.5  1.5  1.5]
>>> print(Grid1D(nx=2).y)
Traceback (most recent call last):
  ...
AttributeError: 1D meshes do not have a "y" attribute.
property z

Equivalent to using cellCenters[2].

>>> from fipy import *
>>> print(Grid3D(nx=2, ny=2, nz=2).z)
[ 0.5  0.5  0.5  0.5  1.5  1.5  1.5  1.5]
>>> print(Grid2D(nx=2, ny=2).z)
Traceback (most recent call last):
  ...
AttributeError: 1D and 2D meshes do not have a "z" attribute.

fipy.meshes.cylindricalGrid1D module

fipy.meshes.cylindricalGrid2D module

fipy.meshes.cylindricalNonUniformGrid1D module

1D Mesh

class fipy.meshes.cylindricalNonUniformGrid1D.CylindricalNonUniformGrid1D(dx=1.0, nx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.nonUniformGrid1D.NonUniformGrid1D

Creates a 1D cylindrical grid mesh.

>>> mesh = CylindricalNonUniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
>>> mesh = CylindricalNonUniformGrid1D(dx = (1, 2, 3))
>>> print(mesh.cellCenters)
[[ 0.5  2.   4.5]]
>>> print(numerix.allclose(mesh.cellVolumes, (0.5, 4., 13.5))) 
True
>>> mesh = CylindricalNonUniformGrid1D(nx = 2, dx = (1, 2, 3))
Traceback (most recent call last):
...
IndexError: nx != len(dx)
>>> mesh = CylindricalNonUniformGrid1D(nx=2, dx=(1., 2.)) + ((1.,),)
>>> print(mesh.cellCenters)
[[ 1.5  3. ]]
>>> print(numerix.allclose(mesh.cellVolumes, (1.5, 6))) 
True
__init__(dx=1.0, nx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.cylindricalNonUniformGrid1D'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape

fipy.meshes.cylindricalNonUniformGrid2D module

2D rectangular Mesh

class fipy.meshes.cylindricalNonUniformGrid2D.CylindricalNonUniformGrid2D(dx=1.0, dy=1.0, nx=None, ny=None, origin=((0.0), (0.0)), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.nonUniformGrid2D.NonUniformGrid2D

Creates a 2D cylindrical grid mesh with horizontal faces numbered first and then vertical faces.

__init__(dx=1.0, dy=1.0, nx=None, ny=None, origin=((0.0), (0.0)), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.cylindricalNonUniformGrid2D'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape

fipy.meshes.cylindricalUniformGrid1D module

1D Mesh

class fipy.meshes.cylindricalUniformGrid1D.CylindricalUniformGrid1D(dx=1.0, nx=1, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.uniformGrid1D.UniformGrid1D

Creates a 1D cylindrical grid mesh.

>>> mesh = CylindricalUniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
__init__(dx=1.0, nx=1, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.cylindricalUniformGrid1D'
property cellVolumes

fipy.meshes.cylindricalUniformGrid2D module

2D cylindrical rectangular Mesh with constant spacing in x and constant spacing in y

class fipy.meshes.cylindricalUniformGrid2D.CylindricalUniformGrid2D(dx=1.0, dy=1.0, nx=1, ny=1, origin=((0), (0)), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.uniformGrid2D.UniformGrid2D

Creates a 2D cylindrical grid in the radial and axial directions, appropriate for axial symmetry.

__init__(dx=1.0, dy=1.0, nx=1, ny=1, origin=((0), (0)), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.cylindricalUniformGrid2D'
property cellVolumes

fipy.meshes.factoryMeshes module

fipy.meshes.factoryMeshes.CylindricalGrid1D(dr=None, nr=None, Lr=None, dx=1.0, nx=None, Lx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between CylindricalUniformGrid1D and CylindricalNonUniformGrid1D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.factoryMeshes.CylindricalGrid2D(dr=None, dz=None, nr=None, nz=None, Lr=None, Lz=None, dx=1.0, dy=1.0, nx=None, ny=None, Lx=None, Ly=None, origin=((0), (0)), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between CylindricalUniformGrid2D and CylindricalNonUniformGrid2D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • dz (float) – grid spacing in the vertical direction. Alternative: dy.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • nz (int) – Number of cells in the vertical direction. Alternative: ny.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • Lz (float) – Domain length in the vertical direction. Alternative: Ly.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.factoryMeshes.Grid1D(dx=1.0, nx=None, Lx=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid1D and NonUniformGrid1D. If Lx is specified the length of the domain is always Lx regardless of dx, unless dx is a list of spacings, in which case Lx will be the sum of dx and nx will be the count of dx.

Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • nx (int) – Number of cells in the horizontal direction

  • Lx (float) – Domain length in the horizontal direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.factoryMeshes.Grid2D(dx=1.0, dy=1.0, nx=None, ny=None, Lx=None, Ly=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid2D and NonUniformGrid2D. If L{x,y} is specified, the length of the domain is always L{x,y} regardless of d{x,y}, unless d{x,y} is a list of spacings, in which case L{x,y} will be the sum of d{x,y} and n{x,y} will be the count of d{x,y}.

>>> print(Grid2D(Lx=3., nx=2).dx)
1.5
Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • dy (float) – Grid spacing in the vertical direction

  • nx (int) – Number of cells in the horizontal direction

  • ny (int) – Number of cells in the vertical direction

  • Lx (float) – Domain length in the horizontal direction

  • Ly (float) – Domain length in the vertical direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.factoryMeshes.Grid3D(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, Lx=None, Ly=None, Lz=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid3D and NonUniformGrid3D. If L{x,y,z} is specified, the length of the domain is always L{x,y,z} regardless of d{x,y,z}, unless d{x,y,z} is a list of spacings, in which case L{x,y,z} will be the sum of d{x,y,z} and n{x,y,z} will be the count of d{x,y,z}.

Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • dy (float) – Grid spacing in the vertical direction

  • dz (float) – Grid spacing in the depth direction

  • nx (int) – Number of cells in the horizontal direction

  • ny (int) – Number of cells in the vertical direction

  • nz (int) – Number of cells in the depth direction

  • Lx (float) – Domain length in the horizontal direction

  • Ly (float) – Domain length in the vertical direction

  • Lz (float) – Domain length in the depth direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.factoryMeshes.SphericalGrid1D(dr=None, nr=None, Lr=None, dx=1.0, nx=None, Lx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between SphericalUniformGrid1D and SphericalNonUniformGrid1D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.gmshMesh module

class fipy.meshes.gmshMesh.Gmsh2D(arg, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.mesh2D.Mesh2D

Construct a 2D Mesh using Gmsh

>>> radius = 5.
>>> side = 4.
>>> squaredCircle = Gmsh2D('''
... // A mesh consisting of a square inside a circle inside a circle
...
... // define the basic dimensions of the mesh
...
... cellSize = 1;
... radius = %(radius)g;
... side = %(side)g;
...
... // define the compass points of the inner circle
...
... Point(1) = {0, 0, 0, cellSize};
... Point(2) = {-radius, 0, 0, cellSize};
... Point(3) = {0, radius, 0, cellSize};
... Point(4) = {radius, 0, 0, cellSize};
... Point(5) = {0, -radius, 0, cellSize};
...
... // define the compass points of the outer circle
...
... Point(6) = {-2*radius, 0, 0, cellSize};
... Point(7) = {0, 2*radius, 0, cellSize};
... Point(8) = {2*radius, 0, 0, cellSize};
... Point(9) = {0, -2*radius, 0, cellSize};
...
... // define the corners of the square
...
... Point(10) = {side/2, side/2, 0, cellSize/2};
... Point(11) = {-side/2, side/2, 0, cellSize/2};
... Point(12) = {-side/2, -side/2, 0, cellSize/2};
... Point(13) = {side/2, -side/2, 0, cellSize/2};
...
... // define the inner circle
...
... Circle(1) = {2, 1, 3};
... Circle(2) = {3, 1, 4};
... Circle(3) = {4, 1, 5};
... Circle(4) = {5, 1, 2};
...
... // define the outer circle
...
... Circle(5) = {6, 1, 7};
... Circle(6) = {7, 1, 8};
... Circle(7) = {8, 1, 9};
... Circle(8) = {9, 1, 6};
...
... // define the square
...
... Line(9) = {10, 13};
... Line(10) = {13, 12};
... Line(11) = {12, 11};
... Line(12) = {11, 10};
...
... // define the three boundaries
...
... Line Loop(1) = {1, 2, 3, 4};
... Line Loop(2) = {5, 6, 7, 8};
... Line Loop(3) = {9, 10, 11, 12};
...
... // define the three domains
...
... Plane Surface(1) = {2, 1};
... Plane Surface(2) = {1, 3};
... Plane Surface(3) = {3};
...
... // label the three domains
...
... // attention: if you use any "Physical" labels, you *must* label
... // all elements that correspond to FiPy Cells (Physical Surface in 2D
... // and Physical Volume in 3D) or Gmsh will not include them and FiPy
... // will not be able to include them in the Mesh.
...
... // note: if you do not use any labels, all Cells will be included.
...
... Physical Surface("Outer") = {1};
... Physical Surface("Middle") = {2};
... Physical Surface("Inner") = {3};
...
... // label the "north-west" part of the exterior boundary
...
... // note: you only need to label the Face elements
... // (Physical Line in 2D and Physical Surface in 3D) that correspond
... // to boundaries you are interested in. FiPy does not need them to
... // construct the Mesh.
...
... Physical Line("NW") = {5};
... ''' % locals()) 

It can be easier to specify certain domains and boundaries within Gmsh than it is to define the same domains and boundaries with FiPy expressions.

Here we compare obtaining the same Cells and Faces using FiPy’s parametric descriptions and Gmsh’s labels.

>>> x, y = squaredCircle.cellCenters 
>>> middle = ((x**2 + y**2 <= radius**2)
...           & ~((x > -side/2) & (x < side/2)
...               & (y > -side/2) & (y < side/2))) 
>>> print((middle == squaredCircle.physicalCells["Middle"]).all()) 
True
>>> X, Y = squaredCircle.faceCenters 
>>> NW = ((X**2 + Y**2 > (1.99*radius)**2)
...       & (X**2 + Y**2 < (2.01*radius)**2)
...       & (X <= 0) & (Y >= 0)) 
>>> print((NW == squaredCircle.physicalFaces["NW"]).all()) 
True

It is possible to direct Gmsh to give the mesh different densities in different locations

>>> geo = '''
... // A mesh consisting of a square
...
... // define the corners of the square
...
... Point(1) = {1, 1, 0, 1};
... Point(2) = {0, 1, 0, 1};
... Point(3) = {0, 0, 0, 1};
... Point(4) = {1, 0, 0, 1};
...
... // define the square
...
... Line(1) = {1, 2};
... Line(2) = {2, 3};
... Line(3) = {3, 4};
... Line(4) = {4, 1};
...
... // define the boundary
...
... Line Loop(1) = {1, 2, 3, 4};
...
... // define the domain
...
... Plane Surface(1) = {1};
... '''
>>> from fipy import CellVariable, numerix
>>> error = []
>>> bkg = None
>>> from builtins import range
>>> for refine in range(4):
...     square = Gmsh2D(geo, background=bkg) 
...     x, y = square.cellCenters 
...     bkg = CellVariable(mesh=square, value=abs(x / 4) + 0.01) 
...     error.append(((2 * numerix.sqrt(square.cellVolumes) / bkg - 1)**2).cellVolumeAverage) 

Check that the mesh is (semi)monotonically approaching the desired density (the first step may increase, depending on the number of partitions)

>>> print(numerix.greater(error[:-1], error[1:]).all()) 
True

and that the final density is close enough to the desired density

>>> print(error[-1] < 0.02) 
True

The initial mesh doesn’t have to be from Gmsh

>>> from fipy import Tri2D
>>> trisquare = Tri2D(nx=1, ny=1)
>>> x, y = trisquare.cellCenters
>>> bkg = CellVariable(mesh=trisquare, value=abs(x / 4) + 0.01)
>>> std1 = (numerix.sqrt(2 * trisquare.cellVolumes) / bkg).std()
>>> square = Gmsh2D(geo, background=bkg) 
>>> x, y = square.cellCenters 
>>> bkg = CellVariable(mesh=square, value=abs(x / 4) + 0.01) 
>>> std2 = (numerix.sqrt(2 * square.cellVolumes) / bkg).std() 
>>> print(std1 > std2) 
True
Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
__setstate__(state)
class fipy.meshes.gmshMesh.Gmsh2DIn3DSpace(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.gmshMesh.Gmsh2D

Create a topologically 2D Mesh in 3D coordinates using Gmsh

Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
class fipy.meshes.gmshMesh.Gmsh3D(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.mesh.Mesh

Create a 3D Mesh using Gmsh

Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
__setstate__(state)
class fipy.meshes.gmshMesh.GmshGrid2D(dx=1.0, dy=1.0, nx=1, ny=None, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1)

Bases: fipy.meshes.gmshMesh.Gmsh2D

Should serve as a drop-in replacement for Grid2D

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=1, ny=None, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
class fipy.meshes.gmshMesh.GmshGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=None, nz=None, communicator=SerialPETScCommWrapper(), overlap=1)

Bases: fipy.meshes.gmshMesh.Gmsh3D

Should serve as a drop-in replacement for Grid3D

__annotations__ = {}
__init__(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=None, nz=None, communicator=SerialPETScCommWrapper(), overlap=1)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
fipy.meshes.gmshMesh.openMSHFile(name, dimensions=None, coordDimensions=None, communicator=SerialPETScCommWrapper(), overlap=1, mode='r', background=None)

Open a Gmsh MSH file

Parameters
  • filename (str) – Gmsh output file

  • dimensions (int) – Dimension of mesh

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • mode (str) – Beginning with r for reading and w for writing. The file will be created if it doesn’t exist when opened for writing; it will be truncated when opened for writing. Add a b to the mode for binary files.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

fipy.meshes.gmshMesh.openPOSFile(name, communicator=SerialPETScCommWrapper(), mode='w')

Open a Gmsh POS post-processing file

fipy.meshes.grid1D module

fipy.meshes.grid2D module

fipy.meshes.grid3D module

fipy.meshes.mesh module

class fipy.meshes.mesh.Mesh(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._MeshTopology'>)

Bases: fipy.meshes.abstractMesh.AbstractMesh

Generic mesh class using numerix to do the calculations

Meshes contain cells, faces, and vertices.

This is built for a non-mixed element mesh.

__init__(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._MeshTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.mesh'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape
__rmul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape
exception fipy.meshes.mesh.MeshAdditionError

Bases: Exception

__module__ = 'fipy.meshes.mesh'
__weakref__

list of weak references to the object (if defined)

fipy.meshes.mesh1D module

Generic mesh class using numerix to do the calculations

Meshes contain cells, faces, and vertices.

This is built for a non-mixed element mesh.

class fipy.meshes.mesh1D.Mesh1D(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh1DTopology'>)

Bases: fipy.meshes.mesh.Mesh

__init__(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh1DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.mesh1D'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape

fipy.meshes.mesh2D module

Generic mesh class using numerix to do the calculations

Meshes contain cells, faces, and vertices.

This is built for a non-mixed element mesh.

class fipy.meshes.mesh2D.Mesh2D(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Bases: fipy.meshes.mesh.Mesh

__init__(vertexCoords, faceVertexIDs, cellFaceIDs, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.meshRepresentation._MeshRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.mesh2D'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape
extrude(extrudeFunc=<function Mesh2D.<lambda>>, layers=1)

This function returns a new 3D mesh. The 2D mesh is extruded using the extrudeFunc and the number of layers.

>>> from fipy.meshes.nonUniformGrid2D import NonUniformGrid2D
>>> print(NonUniformGrid2D(nx=2, ny=2).extrude(layers=2).cellCenters)
[[ 0.5  1.5  0.5  1.5  0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5  0.5  0.5  1.5  1.5]
 [ 0.5  0.5  0.5  0.5  1.5  1.5  1.5  1.5]]
>>> from fipy.meshes.tri2D import Tri2D
>>> print(Tri2D().extrude(layers=2).cellCenters.allclose([[ 0.83333333, 0.5,        0.16666667, 0.5,       0.83333333, 0.5,
...                                                      0.16666667, 0.5       ],
...                                                       [ 0.5,        0.83333333, 0.5,        0.16666667, 0.5,        0.83333333,
...                                                      0.5,        0.16666667],
...                                                      [ 0.5,        0.5,        0.5,        0.5,        1.5,        1.5,        1.5,
...                                                      1.5       ]]))
True
Parameters
  • extrudeFunc (function) – Takes the vertex coordinates and returns the displaced values

  • layers (int) – Number of layers in the extruded mesh (number of times extrudeFunc will be called)

fipy.meshes.nonUniformGrid1D module

1D Mesh

class fipy.meshes.nonUniformGrid1D.NonUniformGrid1D(dx=1.0, nx=None, overlap=2, communicator=SerialPETScCommWrapper(), _BuilderClass=<class 'fipy.meshes.builders.grid1DBuilder._NonuniformGrid1DBuilder'>, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid1DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid1DTopology'>)

Bases: fipy.meshes.mesh1D.Mesh1D

Creates a 1D grid mesh.

>>> mesh = NonUniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
>>> mesh = NonUniformGrid1D(dx = (1, 2, 3))
>>> print(mesh.cellCenters)
[[ 0.5  2.   4.5]]
>>> mesh = NonUniformGrid1D(nx = 2, dx = (1, 2, 3))
Traceback (most recent call last):
...
IndexError: nx != len(dx)
__init__(dx=1.0, nx=None, overlap=2, communicator=SerialPETScCommWrapper(), _BuilderClass=<class 'fipy.meshes.builders.grid1DBuilder._NonuniformGrid1DBuilder'>, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid1DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid1DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.nonUniformGrid1D'

fipy.meshes.nonUniformGrid2D module

2D rectangular Mesh

class fipy.meshes.nonUniformGrid2D.NonUniformGrid2D(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid2DTopology'>)

Bases: fipy.meshes.mesh2D.Mesh2D

Creates a 2D grid mesh with horizontal faces numbered first and then vertical faces.

__init__(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid2DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.nonUniformGrid2D'

fipy.meshes.nonUniformGrid3D module

class fipy.meshes.nonUniformGrid3D.NonUniformGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid3DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid3DTopology'>)

Bases: fipy.meshes.mesh.Mesh

3D rectangular-prism Mesh

X axis runs from left to right. Y axis runs from bottom to top. Z axis runs from front to back.

Numbering System:

Vertices: Numbered in the usual way. X coordinate changes most quickly, then Y, then Z.

Cells: Same numbering system as vertices.

Faces: XY faces numbered first, then XZ faces, then YZ faces. Within each subcategory, it is numbered in the usual way.

__init__(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid3DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid3DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.nonUniformGrid3D'

fipy.meshes.periodicGrid1D module

Periodic 1D Mesh

class fipy.meshes.periodicGrid1D.PeriodicGrid1D(dx=1.0, nx=None, overlap=2, *args, **kwargs)

Bases: fipy.meshes.nonUniformGrid1D.NonUniformGrid1D

Creates a Periodic grid mesh.

>>> mesh = PeriodicGrid1D(dx = (1, 2, 3))
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [3])) 
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-999),
...                        [[2, 0, 1, 2],
...                         [0, 1, 2, -999]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [ 2., 1.5, 2.5, 1.5])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[ 2.,   1.5,  2.5],
...                         [ 1.5,  2.5,  2. ]])) 
True
>>> print(numerix.allclose(mesh.faceNormals,
...                        [[ 1.,  1.,  1.,  1.]])) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[1, 2, 2],
...                        [0, 1, 0]])) 
True
__annotations__ = {}
__init__(dx=1.0, nx=None, overlap=2, *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.periodicGrid1D'
property cellCenters

Defined outside of a geometry class since we need the CellVariable version of cellCenters; that is, the cellCenters defined in fipy.meshes.mesh and not in any geometry (since a CellVariable requires a reference to a mesh).

fipy.meshes.periodicGrid2D module

2D periodic rectangular Mesh

class fipy.meshes.periodicGrid2D.PeriodicGrid2D(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

Creates a periodic 2D grid mesh with horizontal faces numbered first and then vertical faces. Vertices and cells are numbered in the usual way.

>>> from fipy import numerix
>>> mesh = PeriodicGrid2D(dx = 1., dy = 0.5, nx = 2, ny = 2)
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [ 4,  5,  8, 11]))  
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-1),
...                        [[2, 3, 0, 1, 2, 3, 1, 0, 1, 3, 2, 3],
...                         [0, 1, 2, 3, -1, -1, 0, 1, -1, 2, 3, -1]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [ 0.5, 0.5, 0.5, 0.5, 0.25, 0.25, 1., 1., 0.5, 1., 1., 0.5])) 
True
>>> print(numerix.allclose(mesh.cellFaceIDs,
...                        [[0, 1, 2, 3],
...                         [7, 6, 10, 9],
...                         [2, 3, 0, 1],
...                         [6, 7, 9, 10]])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[ 0.5, 0.5, 0.5, 0.5],
...                         [ 1., 1., 1., 1. ],
...                         [ 0.5, 0.5, 0.5, 0.5],
...                         [ 1., 1., 1., 1. ]])) 
True
>>> normals = [[0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1],
...            [1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0]]
>>> print(numerix.allclose(mesh.faceNormals, normals)) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[4, 5, 7, 8],
...                         [3, 4, 6, 7],
...                         [1, 2, 4, 5],
...                         [0, 1, 3, 4]])) 
True
__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'
class fipy.meshes.periodicGrid2D.PeriodicGrid2DLeftRight(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'
class fipy.meshes.periodicGrid2D.PeriodicGrid2DTopBottom(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'

fipy.meshes.periodicGrid3D module

3D periodic rectangular Mesh

class fipy.meshes.periodicGrid3D.PeriodicGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

Creates a periodic 3D grid mesh with horizontal faces numbered first and then vertical faces. Vertices and cells are numbered in the usual way.

>>> from fipy import numerix
>>> mesh = PeriodicGrid3D(dx=1., dy=0.5, dz=2., nx=2, ny=2, nz=1)
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [4, 5, 6, 7, 12, 13, 16, 19]))  
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-1),
...                        [[0, 1, 2, 3, 0, 1, 2, 3, 2, 3,
...                          0, 1, 2, 3, 1, 0, 1, 3, 2, 3],
...                         [0, 1, 2, 3, -1, -1, -1, -1, 0, 1,
...                          2, 3, -1, -1, 0, 1, -1, 2, 3, -1]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [2., 2., 2., 2., 1., 1., 1., 1., 0.5, 0.5,
...                         0.5, 0.5, 0.25, 0.25, 1., 1., 0.5, 1., 1., 0.5])) 
True
>>> print(numerix.allclose(mesh.cellFaceIDs,
...                        [[14, 15, 17, 18],
...                         [15, 14, 18, 17],
...                         [8, 9, 10, 11],
...                         [10, 11, 8, 9],
...                         [0, 1, 2, 3],
...                         [0, 1, 2, 3]])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[1., 1., 1., 1.],
...                         [1., 1., 1., 1.],
...                         [0.5, 0.5, 0.5, 0.5],
...                         [0.5, 0.5, 0.5, 0.5],
...                         [2., 2., 2., 2.],
...                         [2., 2., 2., 2.]])) 
True
>>> normals = [[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1],
...            [0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0],
...            [1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]
>>> print(numerix.allclose(mesh.faceNormals, normals)) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[13, 14, 16, 17],
...                         [12, 13, 15, 16],
...                         [10, 11, 13, 14],
...                         [9, 10, 12, 13],
...                         [4, 5, 7, 8],
...                         [3, 4, 6, 7],
...                         [1, 2, 4, 5],
...                         [0, 1, 3, 4]])) 
True
__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DLeftRight(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DLeftRightFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DLeftRightTopBottom(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DTopBottom(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.periodicGrid3D.PeriodicGrid3DTopBottomFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'

fipy.meshes.skewedGrid2D module

class fipy.meshes.skewedGrid2D.SkewedGrid2D(dx=1.0, dy=1.0, nx=None, ny=1, rand=0, *args, **kwargs)

Bases: fipy.meshes.mesh2D.Mesh2D

Creates a 2D grid mesh with horizontal faces numbered first and then vertical faces. The points are skewed by a random amount (between rand and -rand) in the X and Y directions.

Note

This Mesh only operates in serial

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=None, ny=1, rand=0, *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.skewedGrid2D'
property physicalShape

Return physical dimensions of Grid2D.

property shape

fipy.meshes.sphericalNonUniformGrid1D module

1D Mesh

class fipy.meshes.sphericalNonUniformGrid1D.SphericalNonUniformGrid1D(dx=1.0, nx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.nonUniformGrid1D.NonUniformGrid1D

Creates a 1D spherical grid mesh.

>>> mesh = SphericalNonUniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
>>> mesh = SphericalNonUniformGrid1D(dx = (1, 2, 3))
>>> print(mesh.cellCenters)
[[ 0.5  2.   4.5]]
>>> print(numerix.allclose(mesh.cellVolumes, (0.5, 13., 94.5))) 
True
>>> mesh = SphericalNonUniformGrid1D(nx = 2, dx = (1, 2, 3))
Traceback (most recent call last):
...
IndexError: nx != len(dx)
>>> mesh = SphericalNonUniformGrid1D(nx=2, dx=(1., 2.)) + ((1.,),)
>>> print(mesh.cellCenters)
[[ 1.5  3. ]]
>>> print(numerix.allclose(mesh.cellVolumes, (3.5, 28))) 
True
__init__(dx=1.0, nx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.sphericalNonUniformGrid1D'
__mul__(factor)

Dilate a Mesh by factor.

>>> from fipy.meshes import Grid2D
>>> baseMesh = Grid2D(dx = 1.0, dy = 1.0, nx = 2, ny = 2)
>>> print(baseMesh.cellCenters)
[[ 0.5  1.5  0.5  1.5]
 [ 0.5  0.5  1.5  1.5]]

The factor can be a scalar

>>> dilatedMesh = baseMesh * 3
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.5  1.5  4.5  4.5]]

or a vector

>>> dilatedMesh = baseMesh * ((3,), (2,))
>>> print(dilatedMesh.cellCenters)
[[ 1.5  4.5  1.5  4.5]
 [ 1.   1.   3.   3. ]]

but the vector must have the same dimensionality as the Mesh

>>> dilatedMesh = baseMesh * ((3,), (2,), (1,)) 
Traceback (most recent call last):
...
ValueError: shape mismatch: objects cannot be broadcast to a single shape

fipy.meshes.sphericalUniformGrid1D module

1D Mesh

class fipy.meshes.sphericalUniformGrid1D.SphericalUniformGrid1D(dx=1.0, nx=1, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.uniformGrid1D.UniformGrid1D

Creates a 1D spherical grid mesh.

>>> mesh = SphericalUniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
__init__(dx=1.0, nx=1, origin=(0), overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.sphericalUniformGrid1D'
property cellVolumes

fipy.meshes.test module

Test implementation of the mesh

fipy.meshes.tri2D module

class fipy.meshes.tri2D.Tri2D(dx=1.0, dy=1.0, nx=1, ny=1, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Bases: fipy.meshes.mesh2D.Mesh2D

This class creates a mesh made out of triangles. It does this by starting with a standard Cartesian mesh (Grid2D) and dividing each cell in that mesh (hereafter referred to as a “box”) into four equal parts with the dividing lines being the diagonals.

Creates a 2D triangular mesh with horizontal faces numbered first then vertical faces, then diagonal faces. Vertices are numbered starting with the vertices at the corners of boxes and then the vertices at the centers of boxes. Cells on the right of boxes are numbered first, then cells on the top of boxes, then cells on the left of boxes, then cells on the bottom of boxes. Within each of the “sub-categories” in the above, the vertices, cells and faces are numbered in the usual way.

Parameters
  • dx (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • dy (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • nx (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

  • ny (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=1, ny=1, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Creates a 2D triangular mesh with horizontal faces numbered first then vertical faces, then diagonal faces. Vertices are numbered starting with the vertices at the corners of boxes and then the vertices at the centers of boxes. Cells on the right of boxes are numbered first, then cells on the top of boxes, then cells on the left of boxes, then cells on the bottom of boxes. Within each of the “sub-categories” in the above, the vertices, cells and faces are numbered in the usual way.

Parameters
  • dx (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • dy (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • nx (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

  • ny (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

__module__ = 'fipy.meshes.tri2D'
property physicalShape

Return physical dimensions of Grid2D.

property shape

fipy.meshes.uniformGrid module

class fipy.meshes.uniformGrid.UniformGrid(communicator, _RepresentationClass=<class 'fipy.meshes.representations.abstractRepresentation._AbstractRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.abstractTopology._AbstractTopology'>)

Bases: fipy.meshes.abstractMesh.AbstractMesh

Wrapped scaled geometry properties

__module__ = 'fipy.meshes.uniformGrid'

fipy.meshes.uniformGrid1D module

1D Mesh

class fipy.meshes.uniformGrid1D.UniformGrid1D(dx=1.0, nx=1, origin=(0, ), overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid1DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid1DTopology'>)

Bases: fipy.meshes.uniformGrid.UniformGrid

Creates a 1D grid mesh.

>>> mesh = UniformGrid1D(nx = 3)
>>> print(mesh.cellCenters)
[[ 0.5  1.5  2.5]]
__init__(dx=1.0, nx=1, origin=(0, ), overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid1DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid1DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.uniformGrid1D'
__mul__(factor)
property exteriorFaces

Geometry set and calc

property faceCellIDs
property faceNormals
property vertexCoords

fipy.meshes.uniformGrid2D module

2D rectangular Mesh with constant spacing in x and constant spacing in y

class fipy.meshes.uniformGrid2D.UniformGrid2D(dx=1.0, dy=1.0, nx=1, ny=1, origin=((0, ), (0, )), overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid2DTopology'>)

Bases: fipy.meshes.uniformGrid.UniformGrid

Creates a 2D grid mesh with horizontal faces numbered first and then vertical faces.

__init__(dx=1.0, dy=1.0, nx=1, ny=1, origin=((0, ), (0, )), overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid2DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.uniformGrid2D'
__mul__(factor)
property faceCellIDs
property faceNormals
property faceVertexIDs
property vertexCoords

fipy.meshes.uniformGrid3D module

class fipy.meshes.uniformGrid3D.UniformGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=1, nz=1, origin=[[0], [0], [0]], overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid3DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid3DTopology'>)

Bases: fipy.meshes.uniformGrid.UniformGrid

3D rectangular-prism Mesh with uniform grid spacing in each dimension.

X axis runs from left to right. Y axis runs from bottom to top. Z axis runs from front to back.

Numbering System:

Vertices: Numbered in the usual way. X coordinate changes most quickly, then Y, then Z.

* arrays are arranged Z, Y, X because in numerix, the final index is the one that changes the most quickly *

Cells: Same numbering system as vertices.

Faces: XY faces numbered first, then XZ faces, then YZ faces. Within each subcategory, it is numbered in the usual way.

__init__(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=1, nz=1, origin=[[0], [0], [0]], overlap=2, communicator=SerialPETScCommWrapper(), _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid3DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.gridTopology._Grid3DTopology'>)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.uniformGrid3D'
__mul__(factor)
property faceCellIDs
property faceNormals
property faceVertexIDs
property vertexCoords

Module contents

fipy.meshes.CylindricalGrid1D(dr=None, nr=None, Lr=None, dx=1.0, nx=None, Lx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between CylindricalUniformGrid1D and CylindricalNonUniformGrid1D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.CylindricalGrid2D(dr=None, dz=None, nr=None, nz=None, Lr=None, Lz=None, dx=1.0, dy=1.0, nx=None, ny=None, Lx=None, Ly=None, origin=((0), (0)), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between CylindricalUniformGrid2D and CylindricalNonUniformGrid2D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • dz (float) – grid spacing in the vertical direction. Alternative: dy.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • nz (int) – Number of cells in the vertical direction. Alternative: ny.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • Lz (float) – Domain length in the vertical direction. Alternative: Ly.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

class fipy.meshes.Gmsh2D(arg, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.mesh2D.Mesh2D

Construct a 2D Mesh using Gmsh

>>> radius = 5.
>>> side = 4.
>>> squaredCircle = Gmsh2D('''
... // A mesh consisting of a square inside a circle inside a circle
...
... // define the basic dimensions of the mesh
...
... cellSize = 1;
... radius = %(radius)g;
... side = %(side)g;
...
... // define the compass points of the inner circle
...
... Point(1) = {0, 0, 0, cellSize};
... Point(2) = {-radius, 0, 0, cellSize};
... Point(3) = {0, radius, 0, cellSize};
... Point(4) = {radius, 0, 0, cellSize};
... Point(5) = {0, -radius, 0, cellSize};
...
... // define the compass points of the outer circle
...
... Point(6) = {-2*radius, 0, 0, cellSize};
... Point(7) = {0, 2*radius, 0, cellSize};
... Point(8) = {2*radius, 0, 0, cellSize};
... Point(9) = {0, -2*radius, 0, cellSize};
...
... // define the corners of the square
...
... Point(10) = {side/2, side/2, 0, cellSize/2};
... Point(11) = {-side/2, side/2, 0, cellSize/2};
... Point(12) = {-side/2, -side/2, 0, cellSize/2};
... Point(13) = {side/2, -side/2, 0, cellSize/2};
...
... // define the inner circle
...
... Circle(1) = {2, 1, 3};
... Circle(2) = {3, 1, 4};
... Circle(3) = {4, 1, 5};
... Circle(4) = {5, 1, 2};
...
... // define the outer circle
...
... Circle(5) = {6, 1, 7};
... Circle(6) = {7, 1, 8};
... Circle(7) = {8, 1, 9};
... Circle(8) = {9, 1, 6};
...
... // define the square
...
... Line(9) = {10, 13};
... Line(10) = {13, 12};
... Line(11) = {12, 11};
... Line(12) = {11, 10};
...
... // define the three boundaries
...
... Line Loop(1) = {1, 2, 3, 4};
... Line Loop(2) = {5, 6, 7, 8};
... Line Loop(3) = {9, 10, 11, 12};
...
... // define the three domains
...
... Plane Surface(1) = {2, 1};
... Plane Surface(2) = {1, 3};
... Plane Surface(3) = {3};
...
... // label the three domains
...
... // attention: if you use any "Physical" labels, you *must* label
... // all elements that correspond to FiPy Cells (Physical Surface in 2D
... // and Physical Volume in 3D) or Gmsh will not include them and FiPy
... // will not be able to include them in the Mesh.
...
... // note: if you do not use any labels, all Cells will be included.
...
... Physical Surface("Outer") = {1};
... Physical Surface("Middle") = {2};
... Physical Surface("Inner") = {3};
...
... // label the "north-west" part of the exterior boundary
...
... // note: you only need to label the Face elements
... // (Physical Line in 2D and Physical Surface in 3D) that correspond
... // to boundaries you are interested in. FiPy does not need them to
... // construct the Mesh.
...
... Physical Line("NW") = {5};
... ''' % locals()) 

It can be easier to specify certain domains and boundaries within Gmsh than it is to define the same domains and boundaries with FiPy expressions.

Here we compare obtaining the same Cells and Faces using FiPy’s parametric descriptions and Gmsh’s labels.

>>> x, y = squaredCircle.cellCenters 
>>> middle = ((x**2 + y**2 <= radius**2)
...           & ~((x > -side/2) & (x < side/2)
...               & (y > -side/2) & (y < side/2))) 
>>> print((middle == squaredCircle.physicalCells["Middle"]).all()) 
True
>>> X, Y = squaredCircle.faceCenters 
>>> NW = ((X**2 + Y**2 > (1.99*radius)**2)
...       & (X**2 + Y**2 < (2.01*radius)**2)
...       & (X <= 0) & (Y >= 0)) 
>>> print((NW == squaredCircle.physicalFaces["NW"]).all()) 
True

It is possible to direct Gmsh to give the mesh different densities in different locations

>>> geo = '''
... // A mesh consisting of a square
...
... // define the corners of the square
...
... Point(1) = {1, 1, 0, 1};
... Point(2) = {0, 1, 0, 1};
... Point(3) = {0, 0, 0, 1};
... Point(4) = {1, 0, 0, 1};
...
... // define the square
...
... Line(1) = {1, 2};
... Line(2) = {2, 3};
... Line(3) = {3, 4};
... Line(4) = {4, 1};
...
... // define the boundary
...
... Line Loop(1) = {1, 2, 3, 4};
...
... // define the domain
...
... Plane Surface(1) = {1};
... '''
>>> from fipy import CellVariable, numerix
>>> error = []
>>> bkg = None
>>> from builtins import range
>>> for refine in range(4):
...     square = Gmsh2D(geo, background=bkg) 
...     x, y = square.cellCenters 
...     bkg = CellVariable(mesh=square, value=abs(x / 4) + 0.01) 
...     error.append(((2 * numerix.sqrt(square.cellVolumes) / bkg - 1)**2).cellVolumeAverage) 

Check that the mesh is (semi)monotonically approaching the desired density (the first step may increase, depending on the number of partitions)

>>> print(numerix.greater(error[:-1], error[1:]).all()) 
True

and that the final density is close enough to the desired density

>>> print(error[-1] < 0.02) 
True

The initial mesh doesn’t have to be from Gmsh

>>> from fipy import Tri2D
>>> trisquare = Tri2D(nx=1, ny=1)
>>> x, y = trisquare.cellCenters
>>> bkg = CellVariable(mesh=trisquare, value=abs(x / 4) + 0.01)
>>> std1 = (numerix.sqrt(2 * trisquare.cellVolumes) / bkg).std()
>>> square = Gmsh2D(geo, background=bkg) 
>>> x, y = square.cellCenters 
>>> bkg = CellVariable(mesh=square, value=abs(x / 4) + 0.01) 
>>> std2 = (numerix.sqrt(2 * square.cellVolumes) / bkg).std() 
>>> print(std1 > std2) 
True
Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
__setstate__(state)
class fipy.meshes.Gmsh2DIn3DSpace(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.gmshMesh.Gmsh2D

Create a topologically 2D Mesh in 3D coordinates using Gmsh

Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
class fipy.meshes.Gmsh3D(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Bases: fipy.meshes.mesh.Mesh

Create a 3D Mesh using Gmsh

Parameters
  • arg (str) – (i) the path to an MSH file, (ii) a path to a Gmsh geometry (.geo) file, or (iii) a Gmsh geometry script

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

__annotations__ = {}
__init__(arg, communicator=SerialPETScCommWrapper(), overlap=1, background=None)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
__setstate__(state)
class fipy.meshes.GmshGrid2D(dx=1.0, dy=1.0, nx=1, ny=None, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1)

Bases: fipy.meshes.gmshMesh.Gmsh2D

Should serve as a drop-in replacement for Grid2D

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=1, ny=None, coordDimensions=2, communicator=SerialPETScCommWrapper(), overlap=1)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
class fipy.meshes.GmshGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=None, nz=None, communicator=SerialPETScCommWrapper(), overlap=1)

Bases: fipy.meshes.gmshMesh.Gmsh3D

Should serve as a drop-in replacement for Grid3D

__annotations__ = {}
__init__(dx=1.0, dy=1.0, dz=1.0, nx=1, ny=None, nz=None, communicator=SerialPETScCommWrapper(), overlap=1)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.gmshMesh'
fipy.meshes.Grid1D(dx=1.0, nx=None, Lx=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid1D and NonUniformGrid1D. If Lx is specified the length of the domain is always Lx regardless of dx, unless dx is a list of spacings, in which case Lx will be the sum of dx and nx will be the count of dx.

Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • nx (int) – Number of cells in the horizontal direction

  • Lx (float) – Domain length in the horizontal direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.Grid2D(dx=1.0, dy=1.0, nx=None, ny=None, Lx=None, Ly=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid2D and NonUniformGrid2D. If L{x,y} is specified, the length of the domain is always L{x,y} regardless of d{x,y}, unless d{x,y} is a list of spacings, in which case L{x,y} will be the sum of d{x,y} and n{x,y} will be the count of d{x,y}.

>>> print(Grid2D(Lx=3., nx=2).dx)
1.5
Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • dy (float) – Grid spacing in the vertical direction

  • nx (int) – Number of cells in the horizontal direction

  • ny (int) – Number of cells in the vertical direction

  • Lx (float) – Domain length in the horizontal direction

  • Ly (float) – Domain length in the vertical direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

fipy.meshes.Grid3D(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, Lx=None, Ly=None, Lz=None, overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between UniformGrid3D and NonUniformGrid3D. If L{x,y,z} is specified, the length of the domain is always L{x,y,z} regardless of d{x,y,z}, unless d{x,y,z} is a list of spacings, in which case L{x,y,z} will be the sum of d{x,y,z} and n{x,y,z} will be the count of d{x,y,z}.

Parameters
  • dx (float) – Grid spacing in the horizontal direction

  • dy (float) – Grid spacing in the vertical direction

  • dz (float) – Grid spacing in the depth direction

  • nx (int) – Number of cells in the horizontal direction

  • ny (int) – Number of cells in the vertical direction

  • nz (int) – Number of cells in the depth direction

  • Lx (float) – Domain length in the horizontal direction

  • Ly (float) – Domain length in the vertical direction

  • Lz (float) – Domain length in the depth direction

  • overlap (int) – Number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

class fipy.meshes.PeriodicGrid1D(dx=1.0, nx=None, overlap=2, *args, **kwargs)

Bases: fipy.meshes.nonUniformGrid1D.NonUniformGrid1D

Creates a Periodic grid mesh.

>>> mesh = PeriodicGrid1D(dx = (1, 2, 3))
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [3])) 
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-999),
...                        [[2, 0, 1, 2],
...                         [0, 1, 2, -999]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [ 2., 1.5, 2.5, 1.5])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[ 2.,   1.5,  2.5],
...                         [ 1.5,  2.5,  2. ]])) 
True
>>> print(numerix.allclose(mesh.faceNormals,
...                        [[ 1.,  1.,  1.,  1.]])) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[1, 2, 2],
...                        [0, 1, 0]])) 
True
__annotations__ = {}
__init__(dx=1.0, nx=None, overlap=2, *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.periodicGrid1D'
property cellCenters

Defined outside of a geometry class since we need the CellVariable version of cellCenters; that is, the cellCenters defined in fipy.meshes.mesh and not in any geometry (since a CellVariable requires a reference to a mesh).

class fipy.meshes.PeriodicGrid2D(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

Creates a periodic 2D grid mesh with horizontal faces numbered first and then vertical faces. Vertices and cells are numbered in the usual way.

>>> from fipy import numerix
>>> mesh = PeriodicGrid2D(dx = 1., dy = 0.5, nx = 2, ny = 2)
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [ 4,  5,  8, 11]))  
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-1),
...                        [[2, 3, 0, 1, 2, 3, 1, 0, 1, 3, 2, 3],
...                         [0, 1, 2, 3, -1, -1, 0, 1, -1, 2, 3, -1]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [ 0.5, 0.5, 0.5, 0.5, 0.25, 0.25, 1., 1., 0.5, 1., 1., 0.5])) 
True
>>> print(numerix.allclose(mesh.cellFaceIDs,
...                        [[0, 1, 2, 3],
...                         [7, 6, 10, 9],
...                         [2, 3, 0, 1],
...                         [6, 7, 9, 10]])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[ 0.5, 0.5, 0.5, 0.5],
...                         [ 1., 1., 1., 1. ],
...                         [ 0.5, 0.5, 0.5, 0.5],
...                         [ 1., 1., 1., 1. ]])) 
True
>>> normals = [[0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1],
...            [1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0]]
>>> print(numerix.allclose(mesh.faceNormals, normals)) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[4, 5, 7, 8],
...                         [3, 4, 6, 7],
...                         [1, 2, 4, 5],
...                         [0, 1, 3, 4]])) 
True
__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'
class fipy.meshes.PeriodicGrid2DLeftRight(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'
class fipy.meshes.PeriodicGrid2DTopBottom(dx=1.0, dy=1.0, nx=None, ny=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid2D._BasePeriodicGrid2D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid2D'
class fipy.meshes.PeriodicGrid3D(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

Creates a periodic 3D grid mesh with horizontal faces numbered first and then vertical faces. Vertices and cells are numbered in the usual way.

>>> from fipy import numerix
>>> mesh = PeriodicGrid3D(dx=1., dy=0.5, dz=2., nx=2, ny=2, nz=1)
>>> print(numerix.allclose(numerix.nonzero(mesh.exteriorFaces)[0],
...                        [4, 5, 6, 7, 12, 13, 16, 19]))  
True
>>> print(numerix.allclose(mesh.faceCellIDs.filled(-1),
...                        [[0, 1, 2, 3, 0, 1, 2, 3, 2, 3,
...                          0, 1, 2, 3, 1, 0, 1, 3, 2, 3],
...                         [0, 1, 2, 3, -1, -1, -1, -1, 0, 1,
...                          2, 3, -1, -1, 0, 1, -1, 2, 3, -1]])) 
True
>>> print(numerix.allclose(mesh._cellDistances,
...                        [2., 2., 2., 2., 1., 1., 1., 1., 0.5, 0.5,
...                         0.5, 0.5, 0.25, 0.25, 1., 1., 0.5, 1., 1., 0.5])) 
True
>>> print(numerix.allclose(mesh.cellFaceIDs,
...                        [[14, 15, 17, 18],
...                         [15, 14, 18, 17],
...                         [8, 9, 10, 11],
...                         [10, 11, 8, 9],
...                         [0, 1, 2, 3],
...                         [0, 1, 2, 3]])) 
True
>>> print(numerix.allclose(mesh._cellToCellDistances,
...                        [[1., 1., 1., 1.],
...                         [1., 1., 1., 1.],
...                         [0.5, 0.5, 0.5, 0.5],
...                         [0.5, 0.5, 0.5, 0.5],
...                         [2., 2., 2., 2.],
...                         [2., 2., 2., 2.]])) 
True
>>> normals = [[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1],
...            [0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0],
...            [1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]
>>> print(numerix.allclose(mesh.faceNormals, normals)) 
True
>>> print(numerix.allclose(mesh._cellVertexIDs,
...                        [[13, 14, 16, 17],
...                         [12, 13, 15, 16],
...                         [10, 11, 13, 14],
...                         [9, 10, 12, 13],
...                         [4, 5, 7, 8],
...                         [3, 4, 6, 7],
...                         [1, 2, 4, 5],
...                         [0, 1, 3, 4]])) 
True
__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DLeftRight(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DLeftRightFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DLeftRightTopBottom(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DTopBottom(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.PeriodicGrid3DTopBottomFrontBack(dx=1.0, dy=1.0, dz=1.0, nx=None, ny=None, nz=None, overlap=2, communicator=SerialPETScCommWrapper(), *args, **kwargs)

Bases: fipy.meshes.periodicGrid3D._BasePeriodicGrid3D

__annotations__ = {}
__module__ = 'fipy.meshes.periodicGrid3D'
class fipy.meshes.SkewedGrid2D(dx=1.0, dy=1.0, nx=None, ny=1, rand=0, *args, **kwargs)

Bases: fipy.meshes.mesh2D.Mesh2D

Creates a 2D grid mesh with horizontal faces numbered first and then vertical faces. The points are skewed by a random amount (between rand and -rand) in the X and Y directions.

Note

This Mesh only operates in serial

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=None, ny=1, rand=0, *args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

__module__ = 'fipy.meshes.skewedGrid2D'
property physicalShape

Return physical dimensions of Grid2D.

property shape
fipy.meshes.SphericalGrid1D(dr=None, nr=None, Lr=None, dx=1.0, nx=None, Lx=None, origin=(0), overlap=2, communicator=SerialPETScCommWrapper())

Factory function to select between SphericalUniformGrid1D and SphericalNonUniformGrid1D. If Lr is specified the length of the domain is always Lr regardless of dr, unless dr is a list of spacings, in which case Lr will be the sum of dr.

Parameters
  • dr (float) – Grid spacing in the radial direction. Alternative: dx.

  • nr (int) – Number of cells in the radial direction. Alternative: nx.

  • Lr (float) – Domain length in the radial direction. Alternative: Lx.

  • overlap (int) – the number of overlapping cells for parallel simulations. Generally 2 is adequate. Higher order equations or discretizations require more.

  • communicator (CommWrapper) – Generally, fipy.tools.serialComm or fipy.tools.parallelComm. Select ~fipy.tools.serialComm to create a serial mesh when running in parallel; mostly used for test purposes.

class fipy.meshes.Tri2D(dx=1.0, dy=1.0, nx=1, ny=1, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Bases: fipy.meshes.mesh2D.Mesh2D

This class creates a mesh made out of triangles. It does this by starting with a standard Cartesian mesh (Grid2D) and dividing each cell in that mesh (hereafter referred to as a “box”) into four equal parts with the dividing lines being the diagonals.

Creates a 2D triangular mesh with horizontal faces numbered first then vertical faces, then diagonal faces. Vertices are numbered starting with the vertices at the corners of boxes and then the vertices at the centers of boxes. Cells on the right of boxes are numbered first, then cells on the top of boxes, then cells on the left of boxes, then cells on the bottom of boxes. Within each of the “sub-categories” in the above, the vertices, cells and faces are numbered in the usual way.

Parameters
  • dx (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • dy (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • nx (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

  • ny (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

__annotations__ = {}
__init__(dx=1.0, dy=1.0, nx=1, ny=1, _RepresentationClass=<class 'fipy.meshes.representations.gridRepresentation._Grid2DRepresentation'>, _TopologyClass=<class 'fipy.meshes.topologies.meshTopology._Mesh2DTopology'>)

Creates a 2D triangular mesh with horizontal faces numbered first then vertical faces, then diagonal faces. Vertices are numbered starting with the vertices at the corners of boxes and then the vertices at the centers of boxes. Cells on the right of boxes are numbered first, then cells on the top of boxes, then cells on the left of boxes, then cells on the bottom of boxes. Within each of the “sub-categories” in the above, the vertices, cells and faces are numbered in the usual way.

Parameters
  • dx (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • dy (float) – The X and Y dimensions of each “box”. If dx <> dy, the line segments connecting the cell centers will not be orthogonal to the faces.

  • nx (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

  • ny (int) – The number of boxes in the X direction and the Y direction. The total number of boxes will be equal to nx * ny, and the total number of cells will be equal to 4 * nx * ny.

__module__ = 'fipy.meshes.tri2D'
property physicalShape

Return physical dimensions of Grid2D.

property shape
fipy.meshes.openMSHFile(name, dimensions=None, coordDimensions=None, communicator=SerialPETScCommWrapper(), overlap=1, mode='r', background=None)

Open a Gmsh MSH file

Parameters
  • filename (str) – Gmsh output file

  • dimensions (int) – Dimension of mesh

  • coordDimensions (int) – Dimension of shapes

  • overlap (int) – The number of overlapping cells for parallel simulations. Generally 1 is adequate. Higher order equations or discretizations require more. If overlap is greater than one, communication reverts to serial, as Gmsh only provides one layer of ghost cells.

  • mode (str) – Beginning with r for reading and w for writing. The file will be created if it doesn’t exist when opened for writing; it will be truncated when opened for writing. Add a b to the mode for binary files.

  • background (CellVariable) – Specifies the desired characteristic lengths of the mesh cells

fipy.meshes.openPOSFile(name, communicator=SerialPETScCommWrapper(), mode='w')

Open a Gmsh POS post-processing file

Last updated on Jan 14, 2021. Created using Sphinx 3.4.3.