Tree Mesh

class discretize.TreeMesh(h, x0=None, **kwargs)[source]

Bases: discretize.tree_ext._TreeMesh, discretize.TensorMesh.BaseTensorMesh, discretize.InnerProducts.InnerProducts, discretize.MeshIO.TreeMeshIO

Required Properties:

  • axis_u (Vector3): Vector orientation of u-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: X
  • axis_v (Vector3): Vector orientation of v-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: Y
  • axis_w (Vector3): Vector orientation of w-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: Z
  • h (a list of Array): h is a list containing the cell widths of the tensor mesh in each dimension., a list (each item is a list or numpy array of <class ‘float’> with shape (*)) with length between 0 and 3
  • reference_system (String): The type of coordinate reference frame. Can take on the values cartesian, cylindrical, or spherical. Abbreviations of these are allowed., a unicode string, Default: cartesian
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <class ‘float’> with shape (*)
vntF
vntE
cellGradStencil
cellGrad

Cell centered Gradient operator built off of the faceDiv operator. Grad = - (Mf)^{-1} * Div * diag (volume)

cellGradx

Cell centered Gradient operator in x-direction (Gradx) Grad = sp.vstack((Gradx, Grady, Gradz))

cellGrady

Cell centered Gradient operator in y-direction (Gradx) Grad = sp.vstack((Gradx, Grady, Gradz))

cellGradz

Cell centered Gradient operator in z-direction (Gradz) Grad = sp.vstack((Gradx, Grady, Gradz))

faceDivx
faceDivy
faceDivz
point2index(locs)[source]
permuteCC
permuteF
permuteE
plotSlice(v, vType='CC', normal='Z', ind=None, grid=True, view='real', ax=None, clim=None, showIt=False, pcolorOpts=None, streamOpts=None, gridOpts=None)[source]
save(*args, **kwargs)[source]
load(*args, **kwargs)[source]
copy(*args, **kwargs)[source]
area
aveCC2F
aveCC2Fx
aveCC2Fy

Construct the averaging operator on cell centers to cell y-faces.

aveCC2Fz

Construct the averaging operator on cell centers to cell z-faces.

aveCCV2F
aveE2CC
aveE2CCV
aveEx2CC
aveEy2CC
aveEz2CC
aveF2CC

Construct the averaging operator on cell faces to cell centers.

aveF2CCV

Construct the averaging operator on cell faces to cell centers.

aveFx2CC
aveFy2CC
aveFz2CC
aveN2CC
aveN2E

Construct the averaging operator on cell nodes to cell edges, keeping each dimension separate.

aveN2Ex
aveN2Ey
aveN2Ez
aveN2F

Construct the averaging operator on cell nodes to cell edges, keeping each dimension separate.

aveN2Fx
aveN2Fy
aveN2Fz
axis_u

axis_u (Vector3): Vector orientation of u-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: X

axis_v

axis_v (Vector3): Vector orientation of v-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: Y

axis_w

axis_w (Vector3): Vector orientation of w-direction. For more details see the docs for the rotation_matrix property., a 3D Vector of <class ‘float’> with shape (3), Default: Z

cellBoundaryInd
check_n_shape(change)
check_x0(change)
deserialize(value, trusted=False, strict=False, assert_valid=False, **kwargs)

Creates HasProperties instance from serialized dictionary

This uses the Property deserializers to deserialize all JSON-compatible dictionary values into their corresponding Property values on a new instance of a HasProperties class. Extra keys in the dictionary that do not correspond to Properties will be ignored.

Parameters:

  • value - Dictionary to deserialize new instance from.
  • trusted - If True (and if the input dictionary has '__class__' keyword and this class is in the registry), the new HasProperties class will come from the dictionary. If False (the default), only the HasProperties class this method is called on will be constructed.
  • strict - Requires '__class__', if present on the input dictionary, to match the deserialized instance’s class. Also disallows unused properties in the input dictionary. Default is False.
  • assert_valid - Require deserialized instance to be valid. Default is False.
  • Any other keyword arguments will be passed through to the Property deserializers.
dim

The dimension of the mesh (1, 2, or 3).

Return type:int
Returns:dim
edge
edgeCurl
equal(other)

Determine if two HasProperties instances are equivalent

Equivalence is determined by checking if all Property values on two instances are equal, using Property.equal.

faceBoundaryInd
faceDiv
fill

How filled is the mesh compared to a TensorMesh? As a fraction: [0, 1].

finalize()
getEdgeInnerProduct(prop=None, invProp=False, invMat=False, doFast=True)
Parameters:
  • prop (numpy.ndarray) – material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
  • invProp (bool) – inverts the material property
  • invMat (bool) – inverts the matrix
  • doFast (bool) – do a faster implementation if available.
Return type:

scipy.sparse.csr_matrix

Returns:

M, the inner product matrix (nE, nE)

getEdgeInnerProductDeriv(prop, doFast=True, invProp=False, invMat=False)
Parameters:
  • prop (numpy.ndarray) – material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
  • doFast (bool) – do a faster implementation if available.
  • invProp (bool) – inverts the material property
  • invMat (bool) – inverts the matrix
Return type:

scipy.sparse.csr_matrix

Returns:

dMdm, the derivative of the inner product matrix (nE, nC*nA)

getFaceInnerProduct(prop=None, invProp=False, invMat=False, doFast=True)
Parameters:
  • prop (numpy.ndarray) – material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
  • invProp (bool) – inverts the material property
  • invMat (bool) – inverts the matrix
  • doFast (bool) – do a faster implementation if available.
Return type:

scipy.sparse.csr_matrix

Returns:

M, the inner product matrix (nF, nF)

getFaceInnerProductDeriv(prop, doFast=True, invProp=False, invMat=False)
Parameters:
  • prop (numpy.ndarray) – material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
  • doFast (bool) – do a faster implementation if available.
  • invProp (bool) – inverts the material property
  • invMat (bool) – inverts the matrix
Returns:

dMdmu(u), the derivative of the inner product matrix (u)

Given u, dMdmu returns (nF, nC*nA)

Parameters:u (numpy.ndarray) – vector that multiplies dMdmu
Return type:scipy.sparse.csr_matrix
Returns:dMdmu, the derivative of the inner product matrix for a certain u
getInterpolationMat()
getTensor(key)

Returns a tensor list.

Parameters:key (str) – What tensor (see below)
Return type:list
Returns:list of the tensors that make up the mesh.

key can be:

'CC'    -> scalar field defined on cell centers
'N'     -> scalar field defined on nodes
'Fx'    -> x-component of field defined on faces
'Fy'    -> y-component of field defined on faces
'Fz'    -> z-component of field defined on faces
'Ex'    -> x-component of field defined on edges
'Ey'    -> y-component of field defined on edges
'Ez'    -> z-component of field defined on edges
get_boundary_cells()

Returns the indices of boundary cells in a given direction given an active index array.

Optional Input: :param numpy.array active_ind: None or Boolean array of active indexes in the mesh :param str direction: one of (‘zu’, ‘zd’, ‘xu’, ‘xd’, ‘yu’, ‘yd’)

Output: :rtype: numpy.array :return: Array of indices for the boundary cells in a given direction

gridCC

Returns an M by N numpy array with the center locations of all cells in order. M is the number of cells and N=2,3 is the dimension of the mesh.

gridEx
gridEy
gridEz
gridFx
gridFy
gridFz
gridN

Returns an M by N numpy array with the widths of all cells in order. M is the number of nodes and N=2,3 is the dimension of the mesh.

gridhEx
gridhEy
gridhEz
gridhFx
gridhFy
gridhFz
gridhN
h

h (a list of Array): h is a list containing the cell widths of the tensor mesh in each dimension., a list (each item is a list or numpy array of <class ‘float’> with shape (*)) with length between 0 and 3

h_gridded

Returns an (nC, dim) numpy array with the widths of all cells in order

hx

Width of cells in the x direction

hy

Width of cells in the y direction

hz

Width of cells in the z direction

insert_cells()
isInside(pts, locType='N')

Determines if a set of points are inside a mesh.

Parameters:pts (numpy.ndarray) – Location of points to test
Return type:numpy.ndarray
Returns:inside, numpy array of booleans
maxLevel

The maximum level used, which may be less than levels.

max_level
nC
nE
nEx
nEy
nEz
nF
nFx
nFy
nFz
nN
nhE
nhEx
nhEy
nhEz
nhF
nhFx
nhFy
nhFz
nhN
nodalGrad
normals

Face Normals

Return type:numpy.ndarray
Returns:normals, (sum(nF), dim)
ntE
ntEx
ntEy
ntEz
ntF
ntFx
ntFy
ntFz
ntN
number()
plotGrid()
plotImage()
projectEdgeVector(eV)

Given a vector, eV, in cartesian coordinates, this will project it onto the mesh using the tangents

Parameters:eV (numpy.ndarray) – edge vector with shape (nE, dim)
Return type:numpy.ndarray
Returns:projected edge vector, (nE, )
projectFaceVector(fV)

Given a vector, fV, in cartesian coordinates, this will project it onto the mesh using the normals

Parameters:fV (numpy.ndarray) – face vector with shape (nF, dim)
Return type:numpy.ndarray
Returns:projected face vector, (nF, )
readModelUBC(mesh, fileName)

Read UBC OcTree model and get vector :param string fileName: path to the UBC GIF model file to read :rtype: numpy.ndarray :return: OcTree model

readUBC(TreeMesh, meshFile)

Read UBC 3D OcTree mesh file Input: :param str meshFile: path to the UBC GIF OcTree mesh file to read :rtype: discretize.TreeMesh :return: The octree mesh

reference_is_rotated

True if the axes are rotated from the traditional <X,Y,Z> system with vectors of \((1,0,0)\), \((0,1,0)\), and \((0,0,1)\)

reference_system

reference_system (String): The type of coordinate reference frame. Can take on the values cartesian, cylindrical, or spherical. Abbreviations of these are allowed., a unicode string, Default: cartesian

refine()
rotation_matrix

Builds a rotation matrix to transform coordinates from their coordinate system into a conventional cartesian system. This is built off of the three axis_u, axis_v, and axis_w properties; these mapping coordinates use the letters U, V, and W (the three letters preceding X, Y, and Z in the alphabet) to define the projection of the X, Y, and Z durections. These UVW vectors describe the placement and transformation of the mesh’s coordinate sytem assuming at most 3 directions.

Why would you want to use these UVW mapping vectors the this rotation_matrix property? They allow us to define the realationship between local and global coordinate systems and provide a tool for switching between the two while still maintaing the connectivity of the mesh’s cells. For a visual example of this, please see the figure in the docs for the vtkInterface.

serialize(include_class=True, save_dynamic=False, **kwargs)

Serializes a HasProperties instance to dictionary

This uses the Property serializers to serialize all Property values to a JSON-compatible dictionary. Properties that are undefined are not included. If the HasProperties instance contains a reference to itself, a properties.SelfReferenceError will be raised.

Parameters:

  • include_class - If True (the default), the name of the class will also be saved to the serialized dictionary under key '__class__'
  • save_dynamic - If True, dynamic properties are written to the serialized dict (default: False).
  • Any other keyword arguments will be passed through to the Property serializers.
tangents

Edge Tangents

Return type:numpy.ndarray
Returns:normals, (sum(nE), dim)
toVTK(mesh, models=None)

Convert this mesh object to it’s proper VTK data object with the given model dictionary as the cell data of that dataset.

Input:

Parameters:models (dict(numpy.ndarray)) – Name(‘s) and array(‘s). Match number of cells
validate()

Call all registered class validator methods

These are all methods decorated with @properties.validator. Validator methods are expected to raise a ValidationError if they fail.

vectorCCx

Cell-centered grid vector (1D) in the x direction.

vectorCCy

Cell-centered grid vector (1D) in the y direction.

vectorCCz

Cell-centered grid vector (1D) in the z direction.

vectorNx

Nodal grid vector (1D) in the x direction.

vectorNy

Nodal grid vector (1D) in the y direction.

vectorNz

Nodal grid vector (1D) in the z direction.

vnE

Total number of edges in each direction

Return type:numpy.ndarray
Returns:[nEx, nEy, nEz], (dim, )
import discretize
import numpy as np
M = discretize.TensorMesh([np.ones(n) for n in [2,3]])
M.plotGrid(edges=True, showIt=True)

(Source code, png, hires.png, pdf)

../_images/mesh_tree-1.png
vnF

Total number of faces in each direction

Return type:numpy.ndarray
Returns:[nFx, nFy, nFz], (dim, )
import discretize
import numpy as np
M = discretize.TensorMesh([np.ones(n) for n in [2,3]])
M.plotGrid(faces=True, showIt=True)

(Source code, png, hires.png, pdf)

../_images/mesh_tree-2.png
vol
writeUBC(mesh, fileName, models=None)

Write UBC ocTree mesh and model files from a octree mesh and model. :param string fileName: File to write to :param dict models: Models in a dict, where each key is the filename

writeVTK(mesh, fileName, models=None, directory='')

Makes and saves a VTK object from this mesh and given models

Input:

Parameters:
  • fileName (str) – path to the output vtk file or just its name if directory is specified
  • models (dict) – dictionary of numpy.array - Name(‘s) and array(‘s). Match number of cells
  • directory (str) – directory where the UBC GIF file lives
x0

Mesh IO

class discretize.MeshIO.TreeMeshIO[source]

Bases: object

classmethod readUBC(TreeMesh, meshFile)[source]

Read UBC 3D OcTree mesh file Input: :param str meshFile: path to the UBC GIF OcTree mesh file to read :rtype: discretize.TreeMesh :return: The octree mesh

readModelUBC(mesh, fileName)[source]

Read UBC OcTree model and get vector :param string fileName: path to the UBC GIF model file to read :rtype: numpy.ndarray :return: OcTree model

writeUBC(mesh, fileName, models=None)[source]

Write UBC ocTree mesh and model files from a octree mesh and model. :param string fileName: File to write to :param dict models: Models in a dict, where each key is the filename

_TreeMesh

class discretize.tree_ext._TreeMesh