Tensor Mesh

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

Bases: discretize.TensorMesh.BaseTensorMesh, discretize.BaseMesh.BaseRectangularMesh, discretize.View.TensorView, discretize.DiffOperators.DiffOperators, discretize.InnerProducts.InnerProducts, discretize.MeshIO.TensorMeshIO

TensorMesh is a mesh class that deals with tensor product meshes.

Any Mesh that has a constant width along the entire axis such that it can defined by a single width vector, called ‘h’.

hx = np.array([1, 1, 1])
hy = np.array([1, 2])
hz = np.array([1, 1, 1, 1])

mesh = Mesh.TensorMesh([hx, hy, hz])

Example of a padded tensor mesh using discretize.utils.meshutils.meshTensor():

import discretize
M = discretize.TensorMesh([
    [(10, 10, -1.3), (10, 40), (10, 10, 1.3)],
    [(10, 10, -1.3), (10, 20)]
])
M.plotGrid()

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

../../_images/api_MeshCode-1.png

For a quick tensor mesh on a (10x12x15) unit cube:

mesh = discretize.TensorMesh([10, 12, 15])

Required Properties:

  • 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 <type ‘float’> with shape (*)) with length between 0 and 3
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘float’> with shape (*)
vol

Construct cell volumes of the 3D model as 1d array.

areaFx

Area of the x-faces

areaFy

Area of the y-faces

areaFz

Area of the z-faces

area

Construct face areas of the 3D model as 1d array.

edgeEx

x-edge lengths

edgeEy

y-edge lengths

edgeEz

z-edge lengths

edge

Construct edge legnths of the 3D model as 1d array.

faceBoundaryInd

Find indices of boundary faces in each direction

cellBoundaryInd

Find indices of boundary faces in each direction

Cylindrical Mesh

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

Bases: discretize.TensorMesh.BaseTensorMesh, discretize.BaseMesh.BaseRectangularMesh, discretize.InnerProducts.InnerProducts, discretize.View.CylView, discretize.DiffOperators.DiffOperators

CylMesh is a mesh class for cylindrical problems

Note

for a cylindrically symmetric mesh use [hx, 1, hz]

cs, nc, npad = 20., 30, 8
hx = utils.meshTensor([(cs,npad+10,-0.7), (cs,nc), (cs,npad,1.3)])
hz = utils.meshTensor([(cs,npad   ,-1.3), (cs,nc), (cs,npad,1.3)])
mesh = Mesh.CylMesh([hx,1,hz], [0.,0,-hz.sum()/2.])

Required Properties:

  • cartesianOrigin (Array): Cartesian origin of the mesh, a list or numpy array of <type ‘float’> with shape (*)
  • 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 <type ‘float’> with shape (*)) with length between 0 and 3
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘float’> with shape (*)
cartesianOrigin

cartesianOrigin (Array) – Cartesian origin of the mesh, a list or numpy array of <type ‘float’> with shape (*)

check_cartesian_origin_shape(change)[source]
isSymmetric

Is the mesh cylindrically symmetric?

Return type:bool
Returns:True if the mesh is cylindrically symmetric, False otherwise
nNx

Number of nodes in the x-direction

Return type:int
Returns:nNx
nNy

Number of nodes in the y-direction

Return type:int
Returns:nNy
nN

Total number of nodes

Return type:int
Returns:nN
vnFx

Number of x-faces in each direction

Return type:numpy.array
Returns:vnFx, (dim, )
vnEy

Number of y-edges in each direction

Return type:numpy.array
Returns:vnEy or None if dim < 2, (dim, )
vnEz

Number of z-edges in each direction

Return type:numpy.array
Returns:vnEz or None if nCy > 1, (dim, )
nEz

Number of z-edges

Return type:int
Returns:nEz
vectorCCx

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

vectorCCy

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

vectorNx

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

vectorNy

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

edgeEx

x-edge lengths - these are the radial edges. Radial edges only exist for a 3D cyl mesh.

Return type:numpy.ndarray
Returns:vector of radial edge lengths
edgeEy

y-edge lengths - these are the azimuthal edges. Azimuthal edges exist for all cylindrical meshes. These are arc-lengths (:math:` heta r`)

Return type:numpy.ndarray
Returns:vector of the azimuthal edges
edgeEz

z-edge lengths - these are the vertical edges. Vertical edges only exist for a 3D cyl mesh.

Return type:numpy.ndarray
Returns:vector of the vertical edges
edge

Edge lengths

Return type:numpy.ndarray
Returns:vector of edge lengths \((r, heta, z)\)
areaFx

Area of the x-faces (radial faces). Radial faces exist on all cylindrical meshes

\[A_x = r heta h_z\]
Return type:numpy.ndarray
Returns:area of x-faces
areaFy

Area of y-faces (Azimuthal faces). Azimuthal faces exist only on 3D cylindrical meshes.

\[A_y = h_x h_z\]
Return type:numpy.ndarray
Returns:area of y-faces
areaFz

Area of z-faces.

\[A_z = \frac{ heta}{2} (r_2^2 - r_1^2)z\]
Return type:numpy.ndarray
Returns:area of the z-faces
area

Face areas

For a 3D cyl mesh: [radial, azimuthal, vertical], while a cylindrically symmetric mesh doesn’t have y-Faces, so it returns [radial, vertical]

Return type:numpy.ndarray
Returns:face areas
vol

Volume of each cell

Return type:numpy.ndarray
Returns:cell volumes
gridN

Nodal grid in cylindrical coordinates \((r, heta, z)\). Nodes do not exist in a cylindrically symmetric mesh.

Return type:numpy.ndarray
Returns:grid locations of nodes
gridFx

Grid of x-faces (radial-faces) in cylindrical coordinates \((r, heta, z)\).

Return type:numpy.ndarray
Returns:grid locations of radial faces
gridEy

Grid of y-edges (azimuthal-faces) in cylindrical coordinates \((r, heta, z)\).

Return type:numpy.ndarray
Returns:grid locations of azimuthal faces
gridEz

Grid of z-faces (vertical-faces) in cylindrical coordinates \((r, heta, z)\).

Return type:numpy.ndarray
Returns:grid locations of radial faces
faceDiv

Construct divergence operator (faces to cell-centres).

faceDivx

Construct divergence operator in the x component (faces to cell-centres).

faceDivy

Construct divergence operator in the y component (faces to cell-centres).

faceDivz

Construct divergence operator in the z component (faces to cell-centres).

cellGradx
cellGrad

The cell centered Gradient, takes you to cell faces.

nodalGrad

Construct gradient operator (nodes to edges).

nodalLaplacian

Construct laplacian operator (nodes to edges).

edgeCurl

The edgeCurl (edges to faces)

Return type:scipy.sparse.csr_matrix
Returns:edge curl operator
aveEx2CC

averaging operator of x-edges (radial) to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from x-edges to cell centers
aveEy2CC

averaging operator of y-edges (azimuthal) to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from y-edges to cell centers
aveEz2CC

averaging operator of z-edges to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from z-edges to cell centers
aveE2CC

averaging operator of edges to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from edges to cell centers
aveE2CCV

averaging operator of edges to a cell centered vector

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from edges to cell centered vectors
aveFx2CC

averaging operator of x-faces (radial) to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from x-faces to cell centers
aveFy2CC

averaging operator of y-faces (azimuthal) to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from y-faces to cell centers
aveFz2CC

averaging operator of z-faces (vertical) to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from z-faces to cell centers
aveF2CC

averaging operator of faces to cell centers

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from faces to cell centers
aveF2CCV

averaging operator of x-faces (radial) to cell centered vectors

Return type:scipy.sparse.csr_matrix
Returns:matrix that averages from faces to cell centered vectors
getInterpolationMat(loc, locType='CC', zerosOutside=False)[source]

Produces interpolation matrix

Parameters:
  • loc (numpy.ndarray) – Location of points to interpolate to
  • locType (str) – What to interpolate (see below)
Return type:

scipy.sparse.csr_matrix

Returns:

M, the interpolation matrix

locType can be:

'Ex'    -> x-component of field defined on edges
'Ey'    -> y-component of field defined on edges
'Ez'    -> z-component of field defined on edges
'Fx'    -> x-component of field defined on faces
'Fy'    -> y-component of field defined on faces
'Fz'    -> z-component of field defined on faces
'N'     -> scalar field defined on nodes
'CC'    -> scalar field defined on cell centers
'CCVx'  -> x-component of vector field defined on cell centers
'CCVy'  -> y-component of vector field defined on cell centers
'CCVz'  -> z-component of vector field defined on cell centers
cartesianGrid(locType='CC', theta_shift=None)[source]

Takes a grid location (‘CC’, ‘N’, ‘Ex’, ‘Ey’, ‘Ez’, ‘Fx’, ‘Fy’, ‘Fz’) and returns that grid in cartesian coordinates

Parameters:locType (str) – grid location
Return type:numpy.ndarray
Returns:cartesian coordinates for the cylindrical grid
getInterpolationMatCartMesh(Mrect, locType='CC', locTypeTo=None)[source]

Takes a cartesian mesh and returns a projection to translate onto the cartesian grid.

Parameters:
  • Mrect (discretize.BaseMesh.BaseMesh) – the mesh to interpolate on to
  • locType (str) – grid location (‘CC’, ‘N’, ‘Ex’, ‘Ey’, ‘Ez’, ‘Fx’, ‘Fy’, ‘Fz’)
  • locTypeTo (str) – grid location to interpolate to. If None, the same grid type as locType will be assumed

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:

  • 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 <type ‘float’> with shape (*)) with length between 0 and 3
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘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]

Curvilinear Mesh

class discretize.CurvilinearMesh(nodes=None, **kwargs)[source]

Bases: discretize.BaseMesh.BaseRectangularMesh, discretize.DiffOperators.DiffOperators, discretize.InnerProducts.InnerProducts, discretize.View.CurviView

CurvilinearMesh is a mesh class that deals with curvilinear meshes.

Example of a curvilinear mesh:

import discretize
X, Y = discretize.utils.exampleLrmGrid([3,3],'rotate')
M = discretize.CurvilinearMesh([X, Y])
M.plotGrid(showIt=True)

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

../../_images/api_MeshCode-2.png

Required Properties:

  • nodes (a list of Array): List of arrays describing the node locations, a list (each item is a list or numpy array of <type ‘float’>, <type ‘int’> with shape (*, *, *) or (*, *)) with length between 2 and 3
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘float’> with shape (*)
nodes

nodes (a list of Array) – List of arrays describing the node locations, a list (each item is a list or numpy array of <type ‘float’>, <type ‘int’> with shape (*, *, *) or (*, *)) with length between 2 and 3

check_nodes(change)[source]
gridCC

Cell-centered grid

gridN

Nodal grid.

gridFx

Face staggered grid in the x direction.

gridFy

Face staggered grid in the y direction.

gridFz

Face staggered grid in the y direction.

gridEx

Edge staggered grid in the x direction.

gridEy

Edge staggered grid in the y direction.

gridEz

Edge staggered grid in the z direction.

vol

Construct cell volumes of the 3D model as 1d array

area
normals

Face normals – calling this will average the computed normals so that there is one per face. This is especially relevant in 3D, as there are up to 4 different normals for each face that will be different.

To reshape the normals into a matrix and get the y component:

NyX, NyY, NyZ = M.r(M.normals, 'F', 'Fy', 'M')
edge

Edge lengths

tangents

Edge tangents

Base Rectangular Mesh

class discretize.BaseMesh.BaseRectangularMesh(n, x0=None)[source]

Bases: discretize.BaseMesh.BaseMesh

Required Properties:

  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘float’> with shape (*)
nCx

Number of cells in the x direction

Return type:int
Returns:nCx
nCy

Number of cells in the y direction

Return type:int
Returns:nCy or None if dim < 2
nCz

Number of cells in the z direction

Return type:int
Returns:nCz or None if dim < 3
vnC

Total number of cells in each direction

Return type:numpy.array
Returns:[nCx, nCy, nCz]
nNx

Number of nodes in the x-direction

Return type:int
Returns:nNx
nNy

Number of nodes in the y-direction

Return type:int
Returns:nNy or None if dim < 2
nNz

Number of nodes in the z-direction

Return type:int
Returns:nNz or None if dim < 3
vnN

Total number of nodes in each direction

Return type:numpy.array
Returns:[nNx, nNy, nNz]
vnEx

Number of x-edges in each direction

Return type:numpy.array
Returns:vnEx
vnEy

Number of y-edges in each direction

Return type:numpy.array
Returns:vnEy or None if dim < 2
vnEz

Number of z-edges in each direction

Return type:numpy.array
Returns:vnEz or None if dim < 3
vnFx

Number of x-faces in each direction

Return type:numpy.array
Returns:vnFx
vnFy

Number of y-faces in each direction

Return type:numpy.array
Returns:vnFy or None if dim < 2
vnFz

Number of z-faces in each direction

Return type:numpy.array
Returns:vnFz or None if dim < 3
nC

Total number of cells

Return type:int
Returns:nC
nN

Total number of nodes

Return type:int
Returns:nN
nEx

Number of x-edges

Return type:int
Returns:nEx
nEy

Number of y-edges

Return type:int
Returns:nEy
nEz

Number of z-edges

Return type:int
Returns:nEz
nFx

Number of x-faces

Return type:int
Returns:nFx
nFy

Number of y-faces

Return type:int
Returns:nFy
nFz

Number of z-faces

Return type:int
Returns:nFz
r(x, xType='CC', outType='CC', format='V')[source]

r is a quick reshape command that will do the best it can at giving you what you want.

For example, you have a face variable, and you want the x component of it reshaped to a 3D matrix.

r can fulfil your dreams:

mesh.r(V, 'F', 'Fx', 'M')
       |   |     |    |
       |   |     |    {
       |   |     |      How: 'M' or ['V'] for a matrix
       |   |     |      (ndgrid style) or a vector (n x dim)
       |   |     |    }
       |   |     {
       |   |       What you want: ['CC'], 'N',
       |   |                       'F', 'Fx', 'Fy', 'Fz',
       |   |                       'E', 'Ex', 'Ey', or 'Ez'
       |   |     }
       |   {
       |     What is it: ['CC'], 'N',
       |                  'F', 'Fx', 'Fy', 'Fz',
       |                  'E', 'Ex', 'Ey', or 'Ez'
       |   }
       {
         The input: as a list or ndarray
       }

For example:

..code:

# Separates each component of the Ex grid into 3 matrices
Xex, Yex, Zex = r(mesh.gridEx, 'Ex', 'Ex', 'M')

# Given an edge vector, return just the x edges as a vector
XedgeVector = r(edgeVector, 'E', 'Ex', 'V')

# Separates each component of the edgeVector into 3 vectors
eX, eY, eZ = r(edgeVector, 'E', 'E', 'V')

Base Tensor Mesh

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

Bases: discretize.BaseMesh.BaseMesh

Required Properties:

  • 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 <type ‘float’> with shape (*)) with length between 0 and 3
  • x0 (Array): origin of the mesh (dim, ), a list or numpy array of <type ‘float’> with shape (*)
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 <type ‘float’> with shape (*)) with length between 0 and 3

hx

Width of cells in the x direction

hy

Width of cells in the y direction

hz

Width of cells 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.

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.

gridCC

Cell-centered grid.

gridN

Nodal grid.

h_gridded

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

gridFx

Face staggered grid in the x direction.

gridFy

Face staggered grid in the y direction.

gridFz

Face staggered grid in the z direction.

gridEx

Edge staggered grid in the x direction.

gridEy

Edge staggered grid in the y direction.

gridEz

Edge staggered grid in the z direction.

getTensor(key)[source]

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
isInside(pts, locType='N')[source]

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
getInterpolationMat(loc, locType='CC', zerosOutside=False)[source]

Produces interpolation matrix

Parameters:
  • loc (numpy.ndarray) – Location of points to interpolate to
  • locType (str) – What to interpolate (see below)
Return type:

scipy.sparse.csr_matrix

Returns:

M, the interpolation matrix

locType can be:

'Ex'    -> x-component of field defined on edges
'Ey'    -> y-component of field defined on edges
'Ez'    -> z-component of field defined on edges
'Fx'    -> x-component of field defined on faces
'Fy'    -> y-component of field defined on faces
'Fz'    -> z-component of field defined on faces
'N'     -> scalar field defined on nodes
'CC'    -> scalar field defined on cell centers
'CCVx'  -> x-component of vector field defined on cell centers
'CCVy'  -> y-component of vector field defined on cell centers
'CCVz'  -> z-component of vector field defined on cell centers

Mesh IO

discretize.MeshIO.load_mesh(filename)[source]

Open a json file and load the mesh into the target class

As long as there are no namespace conflicts, the target __class__ will be stored on the properties.HasProperties registry and may be fetched from there.

Parameters:filename (str) – name of file to read in
class discretize.MeshIO.TensorMeshIO[source]

Bases: object

classmethod readUBC(TensorMesh, fileName, directory='')[source]

Wrapper to Read UBC GIF 2D and 3D tensor mesh and generate same dimension TensorMesh.

Input: :param str fileName: path to the UBC GIF mesh file or just its name if directory is specified :param str directory: directory where the UBC GIF file lives

Output: :rtype: TensorMesh :return: The tensor mesh for the fileName.

classmethod readVTK(TensorMesh, fileName, directory='')[source]

Read VTK Rectilinear (vtr xml file) and return Tensor mesh and model

Input: :param str fileName: path to the vtr model file to read or just its name if directory is specified :param str directory: directory where the UBC GIF file lives

Output: :rtype: tuple :return: (TensorMesh, modelDictionary)

writeVTK(mesh, fileName, models=None, directory='')[source]

Makes and saves a VTK rectilinear file (vtr) for a Tensor mesh and model.

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

readModelUBC(mesh, fileName, directory='')[source]
Read UBC 2D or 3D Tensor mesh model
and generate Tensor mesh model

Input: :param str fileName: path to the UBC GIF mesh file to read or just its name if directory is specified :param str directory: directory where the UBC GIF file lives

Output: :rtype: numpy.ndarray :return: model with TensorMesh ordered

writeModelUBC(mesh, fileName, model, directory='')[source]

Writes a model associated with a TensorMesh to a UBC-GIF format model file.

Input: :param str fileName: File to write to or just its name if directory is specified :param str directory: directory where the UBC GIF file lives :param numpy.ndarray model: The model

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

Writes a TensorMesh to a UBC-GIF format mesh file.

Input: :param str fileName: File to write to :param str directory: directory where to save model :param dict models: A dictionary of the models :param str comment_lines: comment lines preceded with ‘!’ to add

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

writeVTK(fileName, models=None)[source]

Function to write a VTU file from a TreeMesh and model.

Mesh Viewing

class discretize.View.TensorView[source]

Bases: object

Provides viewing functions for TensorMesh

This class is inherited by TensorMesh

plotImage(v)[source]

Plots scalar fields on the given mesh.

Input:

Parameters:v (numpy.array) – vector

Optional Inputs:

Parameters:
  • vType (str) – type of vector (‘CC’, ‘N’, ‘F’, ‘Fx’, ‘Fy’, ‘Fz’, ‘E’, ‘Ex’, ‘Ey’, ‘Ez’)
  • ax (matplotlib.axes.Axes) – axis to plot to
  • showIt (bool) – call plt.show()

3D Inputs:

Parameters:
  • numbering (bool) – show numbering of slices, 3D only
  • annotationColor (str) – color of annotation, e.g. ‘w’, ‘k’, ‘b’
import discretize
import numpy as np
M = discretize.TensorMesh([20, 20])
v = np.sin(M.gridCC[:, 0]*2*np.pi)*np.sin(M.gridCC[:, 1]*2*np.pi)
M.plotImage(v, showIt=True)

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

../../_images/api_MeshCode-3.png
import discretize
import numpy as np
M = discretize.TensorMesh([20, 20, 20])
v = np.sin(M.gridCC[:, 0]*2*np.pi)*np.sin(M.gridCC[:, 1]*2*np.pi)*np.sin(M.gridCC[:, 2]*2*np.pi)
M.plotImage(v, annotationColor='k', showIt=True)

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

../../_images/api_MeshCode-4.png
plotSlice(v, vType='CC', normal='Z', ind=None, grid=False, view='real', ax=None, clim=None, showIt=False, pcolorOpts=None, streamOpts=None, gridOpts=None, range_x=None, range_y=None, sample_grid=None, stream_threshold=None)[source]

Plots a slice of a 3D mesh.

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

../../_images/api_MeshCode-5.png
plotGrid(ax=None, nodes=False, faces=False, centers=False, edges=False, lines=True, showIt=False)[source]

Plot the nodal, cell-centered and staggered grids for 1,2 and 3 dimensions.

Parameters:
  • nodes (bool) – plot nodes
  • faces (bool) – plot faces
  • centers (bool) – plot centers
  • edges (bool) – plot edges
  • lines (bool) – plot lines connecting nodes
  • showIt (bool) – call plt.show()
import discretize
import numpy as np
h1 = np.linspace(.1, .5, 3)
h2 = np.linspace(.1, .5, 5)
mesh = discretize.TensorMesh([h1, h2])
mesh.plotGrid(nodes=True, faces=True, centers=True, lines=True, showIt=True)

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

../../_images/api_MeshCode-6.png
import discretize
import numpy as np
h1 = np.linspace(.1, .5, 3)
h2 = np.linspace(.1, .5, 5)
h3 = np.linspace(.1, .5, 3)
mesh = discretize.TensorMesh([h1, h2, h3])
mesh.plotGrid(nodes=True, faces=True, centers=True, lines=True, showIt=True)

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

../../_images/api_MeshCode-7.png
class discretize.View.CylView[source]

Bases: object

plotGrid(*args, **kwargs)[source]
plotImage(*args, **kwargs)[source]
class discretize.View.CurviView[source]

Bases: object

Provides viewing functions for CurvilinearMesh

This class is inherited by CurvilinearMesh

plotGrid(ax=None, nodes=False, faces=False, centers=False, edges=False, lines=True, showIt=False)[source]

Plot the nodal, cell-centered and staggered grids for 1, 2 and 3 dimensions.

import discretize
X, Y = discretize.utils.exampleLrmGrid([3, 3], 'rotate')
M = discretize.CurvilinearMesh([X, Y])
M.plotGrid(showIt=True)

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

../../_images/api_MeshCode-8.png
plotImage(I, ax=None, showIt=False, grid=False, clim=None)[source]