mirror of
https://github.com/wassname/simpeg.git
synced 2026-07-01 16:37:27 +08:00
vol and area tests passing
This commit is contained in:
@@ -157,15 +157,17 @@ class CylMesh(TensorMesh):
|
||||
if self.nCy == 1:
|
||||
self._edge = 2*pi*self.gridN[:,0]
|
||||
else:
|
||||
raise NotImplementedError('edges not implemented for 3D cyl mesh')
|
||||
raise NotImplementedError('edges not yet implemented for 3D cyl mesh')
|
||||
return self._edge
|
||||
|
||||
@property
|
||||
def area(self):
|
||||
"""Face areas"""
|
||||
if getattr(self, '_area', None) is None:
|
||||
areaR = np.kron(self.hz, 2*pi*self.vectorNr)
|
||||
areaZ = np.kron(np.ones_like(self.vectorNz),pi*(self.vectorNr**2 - np.r_[0, self.vectorNr[:-1]]**2))
|
||||
if self.nCy > 1:
|
||||
raise NotImplementedError('area not yet implemented for 3D cyl mesh')
|
||||
areaR = np.kron(self.hz, 2*pi*self.vectorNx)
|
||||
areaZ = np.kron(np.ones_like(self.vectorNz),pi*(self.vectorNx**2 - np.r_[0, self.vectorNx[:-1]]**2))
|
||||
self._area = np.r_[areaR, areaZ]
|
||||
return self._area
|
||||
|
||||
@@ -173,7 +175,9 @@ class CylMesh(TensorMesh):
|
||||
def vol(self):
|
||||
"""Volume of each cell"""
|
||||
if getattr(self, '_vol', None) is None:
|
||||
az = pi*(self.vectorNr**2 - np.r_[0, self.vectorNr[:-1]]**2)
|
||||
if self.nCy > 1:
|
||||
raise NotImplementedError('vol not yet implemented for 3D cyl mesh')
|
||||
az = pi*(self.vectorNx**2 - np.r_[0, self.vectorNx[:-1]]**2)
|
||||
self._vol = np.kron(self.hz,az)
|
||||
return self._vol
|
||||
|
||||
@@ -295,7 +299,7 @@ class CylMesh(TensorMesh):
|
||||
|
||||
loc = np.atleast_2d(loc)
|
||||
|
||||
assert np.all(loc[:,0]<=self.vectorNr.max()) & \
|
||||
assert np.all(loc[:,0]<=self.vectorNx.max()) & \
|
||||
np.all(loc[:,1]>=self.vectorNz.min()) & \
|
||||
np.all(loc[:,1]<=self.vectorNz.max()), \
|
||||
"Points outside of mesh"
|
||||
|
||||
@@ -303,6 +303,8 @@ class InnerProducts(object):
|
||||
def _makeTensor(M, sigma):
|
||||
if sigma is None: # default is ones
|
||||
sigma = np.ones((M.nC, 1))
|
||||
elif type(sigma) is float:
|
||||
sigma = np.ones(self.nC)*sigma
|
||||
|
||||
if M.dim == 2:
|
||||
if sigma.size == M.nC: # Isotropic!
|
||||
|
||||
@@ -15,6 +15,8 @@ class TestCyl2DMesh(unittest.TestCase):
|
||||
|
||||
def test_cylMesh_numbers(self):
|
||||
self.assertTrue(self.mesh.dim == 3)
|
||||
|
||||
self.assertTrue(self.mesh.nC == 6)
|
||||
self.assertTrue(self.mesh.nCx == 3)
|
||||
self.assertTrue(self.mesh.nCy == 1)
|
||||
self.assertTrue(self.mesh.nCz == 2)
|
||||
@@ -42,6 +44,7 @@ class TestCyl2DMesh(unittest.TestCase):
|
||||
self.assertTrue(self.mesh.nEz == 0)
|
||||
self.assertTrue(np.all(self.mesh.vnEz == [3, 0, 2]))
|
||||
self.assertTrue(self.mesh.nE == 9)
|
||||
self.assertTrue(np.all(self.mesh.vnE == [0, 9, 0]))
|
||||
|
||||
def test_vectorsCC(self):
|
||||
v = np.r_[0, 1, 1.75]
|
||||
@@ -60,9 +63,20 @@ class TestCyl2DMesh(unittest.TestCase):
|
||||
self.assertTrue(np.linalg.norm((v-self.mesh.vectorNz)) == 0)
|
||||
|
||||
def test_dimensions(self):
|
||||
v = np.r_[0.5, 1.5, 2, 0.5, 1.5, 2, 0.5, 1.5, 2] * 2 * np.pi
|
||||
self.assertTrue(np.linalg.norm((v-self.mesh.edge)) == 0)
|
||||
edge = np.r_[0.5, 1.5, 2, 0.5, 1.5, 2, 0.5, 1.5, 2] * 2 * np.pi
|
||||
self.assertTrue(np.linalg.norm((edge-self.mesh.edge)) == 0)
|
||||
|
||||
r = np.r_[0, 0.5, 1.5, 2]
|
||||
a = r[1:]*2*np.pi
|
||||
areaX = np.r_[2*a,a]
|
||||
a = (r[1:]**2 - r[:-1]**2)*np.pi
|
||||
areaZ = np.r_[a,a,a]
|
||||
area = np.r_[areaX, areaZ]
|
||||
self.assertTrue(np.linalg.norm((area-self.mesh.area)) == 0)
|
||||
|
||||
a = (r[1:]**2 - r[:-1]**2)*np.pi
|
||||
vol = np.r_[2*a,a]
|
||||
self.assertTrue(np.linalg.norm((vol-self.mesh.vol)) == 0)
|
||||
|
||||
class TestCyl3DMesh(unittest.TestCase):
|
||||
|
||||
|
||||
Reference in New Issue
Block a user