diff --git a/SimPEG/Mesh/CylMesh.py b/SimPEG/Mesh/CylMesh.py index 935c8ce0..34f19c1e 100644 --- a/SimPEG/Mesh/CylMesh.py +++ b/SimPEG/Mesh/CylMesh.py @@ -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" diff --git a/SimPEG/Mesh/InnerProducts.py b/SimPEG/Mesh/InnerProducts.py index 427e5ee2..a5257c8d 100644 --- a/SimPEG/Mesh/InnerProducts.py +++ b/SimPEG/Mesh/InnerProducts.py @@ -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! diff --git a/SimPEG/Tests/test_cylMesh.py b/SimPEG/Tests/test_cylMesh.py index c0db7982..3103e1d5 100644 --- a/SimPEG/Tests/test_cylMesh.py +++ b/SimPEG/Tests/test_cylMesh.py @@ -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):