mirror of
https://github.com/wassname/simpeg.git
synced 2026-08-04 13:13:40 +08:00
@@ -111,12 +111,12 @@ class DiffOperators(object):
|
||||
"""
|
||||
Function that sets the boundary conditions for cell-centred derivative operators.
|
||||
|
||||
Examples:
|
||||
Examples::
|
||||
|
||||
BC = 'neumann' # Neumann in all directions
|
||||
BC = ['neumann', 'dirichlet', 'neumann'] # 3D, Dirichlet in y Neumann else
|
||||
BC = [['neumann', 'dirichlet'], 'dirichlet', 'dirichlet'] # 3D, Neumann in x on bottom of domain,
|
||||
# Dirichlet else
|
||||
BC = 'neumann' # Neumann in all directions
|
||||
BC = ['neumann', 'dirichlet', 'neumann'] # 3D, Dirichlet in y Neumann else
|
||||
BC = [['neumann', 'dirichlet'], 'dirichlet', 'dirichlet'] # 3D, Neumann in x on bottom of domain,
|
||||
# Dirichlet else
|
||||
|
||||
"""
|
||||
if(type(BC) is str):
|
||||
@@ -266,7 +266,7 @@ class DiffOperators(object):
|
||||
nodalAve = property(**nodalAve())
|
||||
|
||||
def nodalVectorAve():
|
||||
doc = "Construct the averaging operator on cell nodes to cell centers, keeping each dimension seperate."
|
||||
doc = "Construct the averaging operator on cell nodes to cell centers, keeping each dimension separate."
|
||||
|
||||
def fget(self):
|
||||
if(self._nodalVectorAve is None):
|
||||
|
||||
+223
-4
@@ -6,17 +6,97 @@ import numpy as np
|
||||
class InnerProducts(object):
|
||||
"""
|
||||
Class creates the inner product matrices that you need!
|
||||
|
||||
InnerProducts is a base class providing inner product matrices for meshes and cannot run on its own. Inherit to your favorite Mesh class.
|
||||
|
||||
|
||||
**Example problem for DC resistivity**
|
||||
|
||||
.. math::
|
||||
|
||||
\sigma^{-1}\mathbf{J} = \\nabla \phi
|
||||
|
||||
We can define in weak form by integrating with a general face function F:
|
||||
|
||||
.. math::
|
||||
|
||||
\int_{\\text{cell}}{\sigma^{-1}\mathbf{J} \cdot \mathbf{F}} = \int_{\\text{cell}}{\\nabla \phi \cdot \mathbf{F}}
|
||||
|
||||
\int_{\\text{cell}}{\sigma^{-1}\mathbf{J} \cdot \mathbf{F}} = \int_{\\text{cell}}{(\\nabla \cdot \mathbf{F}) \phi } + \int_{\partial \\text{cell}}{ \phi \mathbf{F} \cdot \mathbf{n}}
|
||||
|
||||
We can then discretize for every cell:
|
||||
|
||||
.. math::
|
||||
|
||||
v_{\\text{cell}} \sigma^{-1} (\mathbf{J}_x \mathbf{F}_x +\mathbf{J}_y \mathbf{F}_y + \mathbf{J}_z \mathbf{F}_z ) = -\phi^{\\top} v_{\\text{cell}} (\mathbf{D}_{\\text{cell}} \mathbf{F}) + \\text{BC}
|
||||
|
||||
We can represent this in vector form (again this is for every cell), and will generalize for the case of anisotropic (tensor) sigma.
|
||||
|
||||
.. math::
|
||||
|
||||
\mathbf{F}_c^{\\top} (\sqrt{v_{\\text{cell}}} \Sigma^{-1} \sqrt{v_{\\text{cell}}}) \mathbf{J}_c = -\phi^{\\top} v_{\\text{cell}}( v_\\text{cell}^{-1} \mathbf{D}_{\\text{cell}} \mathbf{A} \mathbf{F}) + \\text{BC}
|
||||
|
||||
We multiply by volume on each side of the tensor conductivity to keep symmetry in the system. Here J_c is the Cartesian J (on the faces) and must be calculated differently depending on the mesh:
|
||||
|
||||
.. math::
|
||||
\mathbf{J}_c = \mathbf{Q}_{(i)}\mathbf{J}_\\text{TENSOR} = \mathbf{N}_{(i)}^{-1}\mathbf{Q}_{(i)}\mathbf{J}_\\text{LOM}
|
||||
|
||||
Here the i index refers to where we choose to approximate this integral.
|
||||
We will approximate this relation at every node of the cell, there are 8 in 3D, using a projection matrix Q_i to pick the appropriate fluxes.
|
||||
We will then average to the cell center. For the TENSOR mesh, this looks like:
|
||||
|
||||
.. math::
|
||||
|
||||
\mathbf{F}^{\\top}
|
||||
{1\over 8}
|
||||
\left(\sum_{i=1}^8
|
||||
\mathbf{Q}_{(i)}^{-\\top} \sqrt{v_{\\text{cell}}} \Sigma^{-1} \sqrt{v_{\\text{cell}}} \mathbf{Q}_{(i)}
|
||||
\\right)
|
||||
\mathbf{J}
|
||||
=
|
||||
-\mathbf{F}^{\\top} \mathbf{A} \mathbf{D}_{\\text{cell}}^{\\top} \phi + \\text{BC}
|
||||
|
||||
\mathbf{M}(\Sigma^{-1}) \mathbf{J}
|
||||
=
|
||||
-\mathbf{A} \mathbf{D}_{\\text{cell}}^{\\top} \phi + \\text{BC}
|
||||
|
||||
\mathbf{M}(\Sigma^{-1}) = {1\over 8}
|
||||
\left(\sum_{i=1}^8
|
||||
\mathbf{Q}_{(i)}^{-\\top} \sqrt{v_{\\text{cell}}} \Sigma^{-1} \sqrt{v_{\\text{cell}}} \mathbf{Q}_{(i)}
|
||||
\\right)
|
||||
|
||||
The M is returned if mu is set equal to \Sigma^{-1}.
|
||||
|
||||
If requested (returnP=True) the projection matricies are returned as well (ordered by nodes).
|
||||
Here each P (3*nC, sum(nF)) is a combination of the projection, volume, and any normalization to Cartesian coordinates:
|
||||
|
||||
.. math::
|
||||
\mathbf{P}_{(i)} = \sqrt{ {1\over 8} v_{\\text{cell}}} \overbrace{\mathbf{N}_{(i)}^{-1}}^{\\text{LOM only}} \mathbf{Q}_{(i)}
|
||||
|
||||
Note that this is completed for each cell in the mesh at the same time.
|
||||
"""
|
||||
def __init__(self):
|
||||
raise Exception('InnerProducts is a base class providing inner product matrices for meshes and cannot run on its own. Inherit to your favorite Mesh class.')
|
||||
|
||||
def getFaceInnerProduct(self, mu=None, returnP=False):
|
||||
"""Wrapper function,
|
||||
|
||||
:py:func:`SimPEG.InnerProducts.getEdgeInnerProduct`
|
||||
|
||||
:py:func:`SimPEG.InnerProducts.getEdgeInnerProduct2D`
|
||||
"""
|
||||
if self.dim == 2:
|
||||
return getFaceInnerProduct2D(self, mu, returnP)
|
||||
elif self.dim == 3:
|
||||
return getFaceInnerProduct(self, mu, returnP)
|
||||
|
||||
def getEdgeInnerProduct(self, sigma=None, returnP=False):
|
||||
"""Wrapper function,
|
||||
|
||||
:py:func:`SimPEG.InnerProducts.getFaceInnerProduct`
|
||||
|
||||
:py:func:`SimPEG.InnerProducts.getFaceInnerProduct2D`
|
||||
"""
|
||||
if self.dim == 2:
|
||||
return getEdgeInnerProduct2D(self, sigma, returnP)
|
||||
elif self.dim == 3:
|
||||
@@ -43,6 +123,38 @@ class InnerProducts(object):
|
||||
|
||||
|
||||
def getFaceInnerProduct(mesh, mu=None, returnP=False):
|
||||
"""
|
||||
:param numpy.array mu: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
|
||||
:param bool returnP: returns the projection matrices
|
||||
:rtype: scipy.csr_matrix
|
||||
:return: M, the inner product matrix (sum(nF), sum(nF))
|
||||
|
||||
Depending on the number of columns (either 1, 3, or 6) of mu, the material property is interpreted as follows:
|
||||
|
||||
.. math::
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & 0 & 0 \\\\ 0 & \mu_{1} & 0 \\\\ 0 & 0 & \mu_{1} \end{matrix}\\right]
|
||||
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & 0 & 0 \\\\ 0 & \mu_{2} & 0 \\\\ 0 & 0 & \mu_{3} \end{matrix}\\right]
|
||||
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & \mu_{4} & \mu_{5} \\\\ \mu_{4} & \mu_{2} & \mu_{6} \\\\ \mu_{5} & \mu_{6} & \mu_{3} \end{matrix}\\right]
|
||||
|
||||
\mathbf{M}(\\vec{\mu}) = {1\over 8}
|
||||
\left(\sum_{i=1}^8
|
||||
\mathbf{J}_c^{-\\top} \sqrt{v_{\\text{cell}}} \\vec{\mu} \sqrt{v_{\\text{cell}}} \mathbf{J}_c
|
||||
\\right)
|
||||
|
||||
If requested (returnP=True) the projection matricies are returned as well (ordered by nodes)::
|
||||
|
||||
P = [P000, P001, P010, P011, P100, P101, P110, P111]
|
||||
|
||||
Here each P (3*nC, sum(nF)) is a combination of the projection, volume, and any normalization to Cartesian coordinates:
|
||||
|
||||
.. math::
|
||||
\mathbf{P}_{(i)} = \sqrt{ {1\over 8} v_{\\text{cell}}} \overbrace{\mathbf{N}_{(i)}^{-1}}^{\\text{LOM only}} \mathbf{Q}_{(i)}
|
||||
|
||||
Note that this is completed for each cell in the mesh at the same time.
|
||||
|
||||
"""
|
||||
|
||||
if mu is None: # default is ones
|
||||
mu = np.ones((mesh.nC, 1))
|
||||
@@ -82,10 +194,10 @@ def getFaceInnerProduct(mesh, mu=None, returnP=False):
|
||||
# 100 | i+1,j ,k | i+1, j, k | i, j , k | i, j, k
|
||||
# 010 | i ,j+1,k | i , j, k | i, j+1, k | i, j, k
|
||||
# 110 | i+1,j+1,k | i+1, j, k | i, j+1, k | i, j, k
|
||||
# 001 | i ,j ,k | i , j, k | i, j , k | i, j, k+1
|
||||
# 101 | i+1,j ,k | i+1, j, k | i, j , k | i, j, k+1
|
||||
# 011 | i ,j+1,k | i , j, k | i, j+1, k | i, j, k+1
|
||||
# 111 | i+1,j+1,k | i+1, j, k | i, j+1, k | i, j, k+1
|
||||
# 001 | i ,j ,k+1 | i , j, k | i, j , k | i, j, k+1
|
||||
# 101 | i+1,j ,k+1 | i+1, j, k | i, j , k | i, j, k+1
|
||||
# 011 | i ,j+1,k+1 | i , j, k | i, j+1, k | i, j, k+1
|
||||
# 111 | i+1,j+1,k+1 | i+1, j, k | i, j+1, k | i, j, k+1
|
||||
|
||||
# Square root of cell volume multiplied by 1/8
|
||||
v = np.sqrt(0.125*mesh.vol)
|
||||
@@ -120,6 +232,42 @@ def getFaceInnerProduct(mesh, mu=None, returnP=False):
|
||||
|
||||
|
||||
def getFaceInnerProduct2D(mesh, mu=None, returnP=False):
|
||||
"""
|
||||
:param numpy.array mu: material property (tensor properties are possible) at each cell center (nC, (1, 2, or 3))
|
||||
:param bool returnP: returns the projection matrices
|
||||
:rtype: scipy.csr_matrix
|
||||
:return: M, the inner product matrix (sum(nF), sum(nF))
|
||||
|
||||
Depending on the number of columns (either 1, 2, or 3) of mu, the material property is interpreted as follows:
|
||||
|
||||
.. math::
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & 0 \\\\ 0 & \mu_{1} \end{matrix}\\right]
|
||||
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & 0 \\\\ 0 & \mu_{2} \end{matrix}\\right]
|
||||
|
||||
\\vec{\mu} = \left[\\begin{matrix} \mu_{1} & \mu_{3} \\\\ \mu_{3} & \mu_{2} \end{matrix}\\right]
|
||||
|
||||
|
||||
.. math::
|
||||
|
||||
\mathbf{M}(\\vec{\mu}) = {1\over 4}
|
||||
\left(\sum_{i=1}^4
|
||||
\mathbf{J}_c^{-\\top} \sqrt{v_{\\text{cell}}} \\vec{\mu} \sqrt{v_{\\text{cell}}} \mathbf{J}_c
|
||||
\\right)
|
||||
|
||||
|
||||
If requested (returnP=True) the projection matricies are returned as well (ordered by nodes)::
|
||||
|
||||
P = [P00, P10, P01, P11]
|
||||
|
||||
Here each P (2*nC, sum(nF)) is a combination of the projection, volume, and any normalization to Cartesian coordinates:
|
||||
|
||||
.. math::
|
||||
\mathbf{P}_{(i)} = \sqrt{ {1\over 4} v_{\\text{cell}}} \overbrace{\mathbf{N}_{(i)}^{-1}}^{\\text{LOM only}} \mathbf{Q}_{(i)}
|
||||
|
||||
Note that this is completed for each cell in the mesh at the same time.
|
||||
|
||||
"""
|
||||
|
||||
if mu is None: # default is ones
|
||||
mu = np.ones((mesh.nC, 1))
|
||||
@@ -185,6 +333,41 @@ def getFaceInnerProduct2D(mesh, mu=None, returnP=False):
|
||||
|
||||
|
||||
def getEdgeInnerProduct(mesh, sigma=None, returnP=False):
|
||||
"""
|
||||
:param numpy.array sigma: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
|
||||
:param bool returnP: returns the projection matrices
|
||||
:rtype: scipy.csr_matrix
|
||||
:return: M, the inner product matrix (sum(nE), sum(nE))
|
||||
|
||||
|
||||
Depending on the number of columns (either 1, 3, or 6) of sigma, the material property is interpreted as follows:
|
||||
|
||||
.. math::
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & 0 & 0 \\\\ 0 & \sigma_{1} & 0 \\\\ 0 & 0 & \sigma_{1} \end{matrix}\\right]
|
||||
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & 0 & 0 \\\\ 0 & \sigma_{2} & 0 \\\\ 0 & 0 & \sigma_{3} \end{matrix}\\right]
|
||||
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & \sigma_{4} & \sigma_{5} \\\\ \sigma_{4} & \sigma_{2} & \sigma_{6} \\\\ \sigma_{5} & \sigma_{6} & \sigma_{3} \end{matrix}\\right]
|
||||
|
||||
What is returned:
|
||||
|
||||
.. math::
|
||||
\mathbf{M}(\Sigma) = {1\over 8}
|
||||
\left(\sum_{i=1}^8
|
||||
\mathbf{J}_c^{-\\top} \sqrt{v_{\\text{cell}}} \Sigma \sqrt{v_{\\text{cell}}} \mathbf{J}_c
|
||||
\\right)
|
||||
|
||||
If requested (returnP=True) the projection matricies are returned as well (ordered by nodes)::
|
||||
|
||||
P = [P000, P001, P010, P011, P100, P101, P110, P111]
|
||||
|
||||
Here each P (3*nC, sum(nE)) is a combination of the projection, volume, and any normalization to Cartesian coordinates:
|
||||
|
||||
.. math::
|
||||
\mathbf{P}_{(i)} = \sqrt{ {1\over 8} v_{\\text{cell}}} \overbrace{\mathbf{N}_{(i)}^{-1}}^{\\text{LOM only}} \mathbf{Q}_{(i)}
|
||||
|
||||
Note that this is completed for each cell in the mesh at the same time.
|
||||
"""
|
||||
|
||||
if sigma is None: # default is ones
|
||||
sigma = np.ones((mesh.nC, 1))
|
||||
@@ -262,6 +445,42 @@ def getEdgeInnerProduct(mesh, sigma=None, returnP=False):
|
||||
|
||||
|
||||
def getEdgeInnerProduct2D(mesh, sigma=None, returnP=False):
|
||||
"""
|
||||
:param numpy.array sigma: material property (tensor properties are possible) at each cell center (nC, (1, 2, or 3))
|
||||
:param bool returnP: returns the projection matrices
|
||||
:rtype: scipy.csr_matrix
|
||||
:return: M, the inner product matrix (sum(nE), sum(nE))
|
||||
|
||||
Depending on the number of columns (either 1, 2, or 3) of sigma, the material property is interpreted as follows:
|
||||
|
||||
.. math::
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & 0 \\\\ 0 & \sigma_{1} \end{matrix}\\right]
|
||||
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & 0 \\\\ 0 & \sigma_{2} \end{matrix}\\right]
|
||||
|
||||
\Sigma = \left[\\begin{matrix} \sigma_{1} & \sigma_{3} \\\\ \sigma_{3} & \sigma_{2} \end{matrix}\\right]
|
||||
|
||||
|
||||
.. math::
|
||||
|
||||
\mathbf{M}(\Sigma) = {1\over 4}
|
||||
\left(\sum_{i=1}^4
|
||||
\mathbf{J}_c^{-\\top} \sqrt{v_{\\text{cell}}} \Sigma \sqrt{v_{\\text{cell}}} \mathbf{J}_c
|
||||
\\right)
|
||||
|
||||
|
||||
If requested (returnP=True) the projection matricies are returned as well (ordered by nodes)::
|
||||
|
||||
P = [P00, P10, P01, P11]
|
||||
|
||||
Here each P (2*nC, sum(nE)) is a combination of the projection, volume, and any normalization to Cartesian coordinates:
|
||||
|
||||
.. math::
|
||||
\mathbf{P}_{(i)} = \sqrt{ {1\over 4} v_{\\text{cell}}} \overbrace{\mathbf{N}_{(i)}^{-1}}^{\\text{LOM only}} \mathbf{Q}_{(i)}
|
||||
|
||||
Note that this is completed for each cell in the mesh at the same time.
|
||||
|
||||
"""
|
||||
|
||||
if sigma is None: # default is ones
|
||||
sigma = np.ones((mesh.nC, 1))
|
||||
|
||||
@@ -16,6 +16,10 @@ class LogicallyOrthogonalMesh(BaseMesh, DiffOperators, InnerProducts, LomView):
|
||||
"""
|
||||
LogicallyOrthogonalMesh is a mesh class that deals with logically orthogonal meshes.
|
||||
|
||||
Example of a logically orthogonal mesh:
|
||||
|
||||
.. plot:: examples/mesh/plot_LogicallyOrthogonalMesh.py
|
||||
|
||||
"""
|
||||
_meshType = 'LOM'
|
||||
|
||||
@@ -254,17 +258,16 @@ class LogicallyOrthogonalMesh(BaseMesh, DiffOperators, InnerProducts, LomView):
|
||||
area = property(**area())
|
||||
|
||||
def normals():
|
||||
doc = """
|
||||
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.
|
||||
doc = """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:
|
||||
To reshape the normals into a matrix and get the y component:
|
||||
|
||||
NyX, NyY, NyZ = M.r(M.normals, 'F', 'Fy', 'M')
|
||||
"""
|
||||
NyX, NyY, NyZ = M.r(M.normals, 'F', 'Fy', 'M')
|
||||
"""
|
||||
|
||||
def fget(self):
|
||||
if(self._normals is None):
|
||||
|
||||
+7
-3
@@ -7,15 +7,19 @@ from utils import mkvc
|
||||
|
||||
class LomView(object):
|
||||
"""
|
||||
Provides viewing functions for TensorMesh
|
||||
Provides viewing functions for LogicallyOrthogonalMesh
|
||||
|
||||
This class is inherited by LogicallyOrthogonalMesh
|
||||
|
||||
This class is inherited by TensorMesh
|
||||
"""
|
||||
def __init__(self):
|
||||
pass
|
||||
|
||||
def plotGrid(self, length=0.05):
|
||||
"""Plot the nodal, cell-centered and staggered grids for 1,2 and 3 dimensions."""
|
||||
"""Plot the nodal, cell-centered and staggered grids for 1,2 and 3 dimensions.
|
||||
|
||||
.. plot:: examples/mesh/plot_LogicallyOrthogonalMesh.py
|
||||
"""
|
||||
NN = self.r(self.gridN, 'N', 'N', 'M')
|
||||
if self.dim == 2:
|
||||
fig = plt.figure(2)
|
||||
|
||||
@@ -21,8 +21,9 @@ class TensorMesh(BaseMesh, TensorView, DiffOperators, InnerProducts):
|
||||
|
||||
mesh = TensorMesh([hx, hy, hz])
|
||||
|
||||
.. math::
|
||||
x^2 = 5
|
||||
Example of a padded tensor mesh:
|
||||
|
||||
.. plot:: examples/mesh/plot_TensorMesh.py
|
||||
|
||||
"""
|
||||
_meshType = 'TENSOR'
|
||||
|
||||
+54
-13
@@ -2,6 +2,7 @@ import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
import matplotlib
|
||||
from mpl_toolkits.mplot3d import Axes3D
|
||||
from utils import mkvc
|
||||
|
||||
|
||||
class TensorView(object):
|
||||
@@ -176,9 +177,16 @@ class TensorView(object):
|
||||
if showIt: plt.show()
|
||||
return ph
|
||||
|
||||
def plotGrid(self, showIt=False):
|
||||
def plotGrid(self, nodes=False, faces=False, centers=False, edges=False, lines=True, showIt=False):
|
||||
"""Plot the nodal, cell-centered and staggered grids for 1,2 and 3 dimensions.
|
||||
|
||||
:param bool nodes: plot nodes
|
||||
:param bool faces: plot faces
|
||||
:param bool centers: plot centers
|
||||
:param bool edges: plot edges
|
||||
:param bool lines: plot lines connecting nodes
|
||||
:param bool showIt: call plt.show()
|
||||
|
||||
.. plot:: examples/mesh/plot_grid_2D.py
|
||||
:include-source:
|
||||
|
||||
@@ -209,10 +217,23 @@ class TensorView(object):
|
||||
xs2 = self.gridFy
|
||||
|
||||
ax.hold(True)
|
||||
ax.plot(xn[:, 0], xn[:, 1], 'bs')
|
||||
ax.plot(xc[:, 0], xc[:, 1], 'ro')
|
||||
ax.plot(xs1[:, 0], xs1[:, 1], 'g>')
|
||||
ax.plot(xs2[:, 0], xs2[:, 1], 'g^')
|
||||
if nodes: ax.plot(xn[:, 0], xn[:, 1], 'bs')
|
||||
if centers: ax.plot(xc[:, 0], xc[:, 1], 'ro')
|
||||
if faces:
|
||||
ax.plot(xs1[:, 0], xs1[:, 1], 'g>')
|
||||
ax.plot(xs2[:, 0], xs2[:, 1], 'g^')
|
||||
|
||||
# Plot the grid lines
|
||||
if lines:
|
||||
NN = self.r(self.gridN, 'N', 'N', 'M')
|
||||
X1 = np.c_[mkvc(NN[0][0, :]), mkvc(NN[0][self.nCx, :]), mkvc(NN[0][0, :])*np.nan].flatten()
|
||||
Y1 = np.c_[mkvc(NN[1][0, :]), mkvc(NN[1][self.nCx, :]), mkvc(NN[1][0, :])*np.nan].flatten()
|
||||
X2 = np.c_[mkvc(NN[0][:, 0]), mkvc(NN[0][:, self.nCy]), mkvc(NN[0][:, 0])*np.nan].flatten()
|
||||
Y2 = np.c_[mkvc(NN[1][:, 0]), mkvc(NN[1][:, self.nCy]), mkvc(NN[1][:, 0])*np.nan].flatten()
|
||||
X = np.r_[X1, X2]
|
||||
Y = np.r_[Y1, Y2]
|
||||
plt.plot(X, Y)
|
||||
|
||||
ax.grid(True)
|
||||
ax.hold(False)
|
||||
ax.set_xlabel('x1')
|
||||
@@ -233,14 +254,34 @@ class TensorView(object):
|
||||
xes3 = self.gridEz
|
||||
|
||||
ax.hold(True)
|
||||
ax.plot(xn[:, 0], xn[:, 1], 'bs', zs=xn[:, 2])
|
||||
ax.plot(xc[:, 0], xc[:, 1], 'ro', zs=xc[:, 2])
|
||||
ax.plot(xfs1[:, 0], xfs1[:, 1], 'g>', zs=xfs1[:, 2])
|
||||
ax.plot(xfs2[:, 0], xfs2[:, 1], 'g<', zs=xfs2[:, 2])
|
||||
ax.plot(xfs3[:, 0], xfs3[:, 1], 'g^', zs=xfs3[:, 2])
|
||||
ax.plot(xes1[:, 0], xes1[:, 1], 'k>', zs=xes1[:, 2])
|
||||
ax.plot(xes2[:, 0], xes2[:, 1], 'k<', zs=xes2[:, 2])
|
||||
ax.plot(xes3[:, 0], xes3[:, 1], 'k^', zs=xes3[:, 2])
|
||||
if nodes: ax.plot(xn[:, 0], xn[:, 1], 'bs', zs=xn[:, 2])
|
||||
if centers: ax.plot(xc[:, 0], xc[:, 1], 'ro', zs=xc[:, 2])
|
||||
if faces:
|
||||
ax.plot(xfs1[:, 0], xfs1[:, 1], 'g>', zs=xfs1[:, 2])
|
||||
ax.plot(xfs2[:, 0], xfs2[:, 1], 'g<', zs=xfs2[:, 2])
|
||||
ax.plot(xfs3[:, 0], xfs3[:, 1], 'g^', zs=xfs3[:, 2])
|
||||
if edges:
|
||||
ax.plot(xes1[:, 0], xes1[:, 1], 'k>', zs=xes1[:, 2])
|
||||
ax.plot(xes2[:, 0], xes2[:, 1], 'k<', zs=xes2[:, 2])
|
||||
ax.plot(xes3[:, 0], xes3[:, 1], 'k^', zs=xes3[:, 2])
|
||||
|
||||
# Plot the grid lines
|
||||
if lines:
|
||||
NN = self.r(self.gridN, 'N', 'N', 'M')
|
||||
X1 = np.c_[mkvc(NN[0][0, :, :]), mkvc(NN[0][self.nCx, :, :]), mkvc(NN[0][0, :, :])*np.nan].flatten()
|
||||
Y1 = np.c_[mkvc(NN[1][0, :, :]), mkvc(NN[1][self.nCx, :, :]), mkvc(NN[1][0, :, :])*np.nan].flatten()
|
||||
Z1 = np.c_[mkvc(NN[2][0, :, :]), mkvc(NN[2][self.nCx, :, :]), mkvc(NN[2][0, :, :])*np.nan].flatten()
|
||||
X2 = np.c_[mkvc(NN[0][:, 0, :]), mkvc(NN[0][:, self.nCy, :]), mkvc(NN[0][:, 0, :])*np.nan].flatten()
|
||||
Y2 = np.c_[mkvc(NN[1][:, 0, :]), mkvc(NN[1][:, self.nCy, :]), mkvc(NN[1][:, 0, :])*np.nan].flatten()
|
||||
Z2 = np.c_[mkvc(NN[2][:, 0, :]), mkvc(NN[2][:, self.nCy, :]), mkvc(NN[2][:, 0, :])*np.nan].flatten()
|
||||
X3 = np.c_[mkvc(NN[0][:, :, 0]), mkvc(NN[0][:, :, self.nCz]), mkvc(NN[0][:, :, 0])*np.nan].flatten()
|
||||
Y3 = np.c_[mkvc(NN[1][:, :, 0]), mkvc(NN[1][:, :, self.nCz]), mkvc(NN[1][:, :, 0])*np.nan].flatten()
|
||||
Z3 = np.c_[mkvc(NN[2][:, :, 0]), mkvc(NN[2][:, :, self.nCz]), mkvc(NN[2][:, :, 0])*np.nan].flatten()
|
||||
X = np.r_[X1, X2, X3]
|
||||
Y = np.r_[Y1, Y2, Y3]
|
||||
Z = np.r_[Z1, Z2, Z3]
|
||||
plt.plot(X, Y, 'b-', zs=Z)
|
||||
|
||||
ax.grid(True)
|
||||
ax.hold(False)
|
||||
ax.set_xlabel('x1')
|
||||
|
||||
@@ -1,8 +1,8 @@
|
||||
.. _api_LOMView:
|
||||
|
||||
LOM View
|
||||
***********
|
||||
********
|
||||
|
||||
.. automodule:: SimPEG.LOMView
|
||||
.. automodule:: SimPEG.LomView
|
||||
:members:
|
||||
:undoc-members:
|
||||
|
||||
@@ -0,0 +1,6 @@
|
||||
from SimPEG import LogicallyOrthogonalMesh, utils
|
||||
import matplotlib.pyplot as plt
|
||||
X, Y = utils.exampleLomGird([3,3],'rotate')
|
||||
M = LogicallyOrthogonalMesh([X, Y])
|
||||
M.plotGrid()
|
||||
plt.show()
|
||||
@@ -0,0 +1,21 @@
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
from SimPEG import TensorMesh
|
||||
|
||||
pad = 7
|
||||
padfactor = 1.4
|
||||
xpad = (np.ones(pad)*padfactor)**np.arange(pad)
|
||||
ypad = (np.ones(pad)*padfactor)**np.arange(pad)
|
||||
|
||||
core = 15
|
||||
xcore = np.ones(core)
|
||||
ycore = np.ones(core)
|
||||
|
||||
h1 = np.r_[xpad[::-1],xcore,xpad]
|
||||
h2 = np.r_[ypad[::-1],ycore,ypad]
|
||||
|
||||
mesh = TensorMesh([h1, h2])
|
||||
mesh.plotGrid()
|
||||
plt.axis('tight')
|
||||
|
||||
plt.show()
|
||||
@@ -5,7 +5,7 @@ from SimPEG import TensorMesh
|
||||
h1 = np.linspace(.1,.5,3)
|
||||
h2 = np.linspace(.1,.5,5)
|
||||
mesh = TensorMesh([h1, h2])
|
||||
mesh.plotGrid()
|
||||
mesh.plotGrid(nodes=True, faces=True, centers=True, lines=True)
|
||||
|
||||
plt.show()
|
||||
|
||||
|
||||
@@ -6,7 +6,7 @@ h1 = np.linspace(.1,.5,3)
|
||||
h2 = np.linspace(.1,.5,5)
|
||||
h3 = np.linspace(.1,.5,3)
|
||||
mesh = TensorMesh([h1,h2,h3])
|
||||
mesh.plotGrid()
|
||||
mesh.plotGrid(nodes=True, faces=True, centers=True, lines=True)
|
||||
|
||||
plt.show()
|
||||
|
||||
|
||||
+1
-1
@@ -26,7 +26,7 @@ Meshing & Operators
|
||||
api_TensorMesh
|
||||
api_TensorView
|
||||
api_LogicallyOrthogonalMesh
|
||||
api_LomView
|
||||
api_LOMView
|
||||
api_DiffOperators
|
||||
api_InnerProducts
|
||||
|
||||
|
||||
Reference in New Issue
Block a user