wassname/simpeg
1
0
Fork 0
mirror of https://github.com/wassname/simpeg.git synced 2026-07-15 11:26:09 +08:00
Code Issues Packages Projects Releases Wiki Activity

update tensors and tests

Browse Source
This commit is contained in:
rowanc1
2014-02-21 11:26:14 -08:00
parent ba52e16bea
commit 0fb229a3c0
3 changed files with 46 additions and 47 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SimPEG/Mesh/CylMesh.py
+7 -43
View File
@@ -115,9 +115,8 @@ class CylMesh(BaseTensorMesh):
@property
def vectorCCx(self):
"""Cell-centered grid vector (1D) in the x direction."""
if self.nCy == 1:
return np.r_[-self.hx[0]*0.5, self.hx[:-1].cumsum()] + self.hx*0.5 # - self.hx[0]/2
return np.r_[0, self.hx[:-1].cumsum()] + self.hx*0.5
firstEl = -self.hx[0]*0.5 if self.nCy == 1 else 0
return np.r_[firstEl, self.hx[:-1].cumsum()] + self.hx*0.5
@property
def vectorCCy(self):
@@ -128,55 +127,18 @@ class CylMesh(BaseTensorMesh):
def vectorNx(self):
"""Nodal grid vector (1D) in the x direction."""
if self.nCy == 1:
return self.hx.cumsum()# - self.hx[0]/2
return self.hx.cumsum()
return np.r_[0, self.hx].cumsum()
@property
def vectorNy(self):
"""Nodal grid vector (1D) in the y direction."""
if self.nCy == 1:
# There aren't really any nodes, but all the grids need
# somewhere to live, why not zero?!
return np.r_[0]
return np.r_[0, self.hy[:-1].cumsum()] + self.hy[0]*0.5
def getTensor(self, locType):
""" Returns a tensor list.
:param str locType: What tensor (see below)
:rtype: list
:return: list of the tensors that make up the mesh.
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
"""
if locType is 'Fx':
ten = [self.vectorNx , self.vectorCCy, self.vectorCCz]
elif locType is 'Fy':
ten = [self.vectorCCx, self.vectorNy , self.vectorCCz]
elif locType is 'Fz':
ten = [self.vectorCCx, self.vectorCCy, self.vectorNz ]
elif locType is 'Ex':
ten = [self.vectorCCx, self.vectorNy , self.vectorNz ]
elif locType is 'Ey':
ten = [self.vectorNx , self.vectorCCy, self.vectorNz ]
elif locType is 'Ez':
ten = [self.vectorNx , self.vectorNy , self.vectorCCz]
elif locType is 'CC':
ten = [self.vectorCCx, self.vectorCCy, self.vectorCCz]
elif locType is 'N':
ten = [self.vectorNx , self.vectorNy , self.vectorNz ]
return [t for t in ten if t is not None]
@property
def edge(self):
"""Edge lengths"""
@@ -270,6 +232,8 @@ class CylMesh(BaseTensorMesh):
@property
def nodalGrad(self):
"""Construct gradient operator (nodes to edges)."""
if self.nCy == 1:
raise Exception('Nodal grad does not make sense for cylindrically symmetric mesh.')
raise NotImplementedError('nodalGrad not yet implemented')
@property
SimPEG/Mesh/TensorMesh.py
+6
View File
@@ -95,6 +95,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridFx(self):
if self.nFx == 0: return
"""Face staggered grid in the x direction."""
if getattr(self, '_gridFx', None) is None:
self._gridFx = Utils.ndgrid(self.getTensor('Fx'))
@@ -102,6 +103,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridFy(self):
if self.nFy == 0: return
"""Face staggered grid in the y direction."""
if getattr(self, '_gridFy', None) is None and self.dim > 1:
self._gridFy = Utils.ndgrid(self.getTensor('Fy'))
@@ -109,6 +111,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridFz(self):
if self.nFz == 0: return
"""Face staggered grid in the z direction."""
if getattr(self, '_gridFz', None) is None and self.dim > 2:
self._gridFz = Utils.ndgrid(self.getTensor('Fz'))
@@ -116,6 +119,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridEx(self):
if self.nEx == 0: return
"""Edge staggered grid in the x direction."""
if getattr(self, '_gridEx', None) is None:
self._gridEx = Utils.ndgrid(self.getTensor('Ex'))
@@ -123,6 +127,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridEy(self):
if self.nEy == 0: return
"""Edge staggered grid in the y direction."""
if getattr(self, '_gridEy', None) is None and self.dim > 1:
self._gridEy = Utils.ndgrid(self.getTensor('Ey'))
@@ -130,6 +135,7 @@ class BaseTensorMesh(BaseRectangularMesh):
@property
def gridEz(self):
if self.nEz == 0: return
"""Edge staggered grid in the z direction."""
if getattr(self, '_gridEz', None) is None and self.dim > 2:
self._gridEz = Utils.ndgrid(self.getTensor('Ez'))
SimPEG/Tests/test_cylMesh.py
+33 -4
View File
@@ -78,15 +78,44 @@ class TestCyl2DMesh(unittest.TestCase):
def test_gridSizes(self):
self.assertTrue(self.mesh.gridCC.shape == (self.mesh.nC, 3))
# self.assertTrue(self.mesh.gridN.shape == (self.mesh.nN, 3))
self.assertTrue(self.mesh.gridN.shape == (9, 3))
self.assertTrue(self.mesh.gridFx.shape == (self.mesh.nFx, 3))
# self.assertTrue(self.mesh.gridFy.shape == (self.mesh.nFy, 3))
self.assertTrue(self.mesh.gridFy is None)
self.assertTrue(self.mesh.gridFz.shape == (self.mesh.nFz, 3))
# self.assertTrue(self.mesh.gridEx.shape == (self.mesh.nEx, 3))
self.assertTrue(self.mesh.gridEx is None)
self.assertTrue(self.mesh.gridEy.shape == (self.mesh.nEy, 3))
# self.assertTrue(self.mesh.gridEz.shape == (self.mesh.nEz, 3))
self.assertTrue(self.mesh.gridEz is None)
def test_gridCC(self):
x = np.r_[0,1.5,2.25,0,1.5,2.25]
y = np.zeros(6)
z = np.r_[1,1,1,2.5,2.5,2.5]
G = np.c_[x,y,z]
self.assertTrue(np.linalg.norm((G-self.mesh.gridCC).ravel()) == 0)
def test_gridFx(self):
x = np.r_[1,2,2.5,1,2,2.5]
y = np.zeros(6)
z = np.r_[1,1,1,2.5,2.5,2.5]
G = np.c_[x,y,z]
self.assertTrue(np.linalg.norm((G-self.mesh.gridFx).ravel()) == 0)
def test_gridFz(self):
x = np.r_[0,1.5,2.25,0,1.5,2.25,0,1.5,2.25]
y = np.zeros(9)
z = np.r_[0,0,0,2,2,2,3,3,3.]
G = np.c_[x,y,z]
self.assertTrue(np.linalg.norm((G-self.mesh.gridFz).ravel()) == 0)
def test_gridEy(self):
x = np.r_[1,2,2.5,1,2,2.5,1,2,2.5]
y = np.zeros(9)
z = np.r_[0,0,0,2,2,2,3,3,3.]
G = np.c_[x,y,z]
self.assertTrue(np.linalg.norm((G-self.mesh.gridEy).ravel()) == 0)
MESHTYPES = ['uniformCylMesh']