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

updates to area calculation

Browse Source
This commit is contained in:
Rowan Cockett
2015-02-09 22:12:22 -08:00
parent 9a8f4d60f7
commit 006dcf393d
2 changed files with 40 additions and 6 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SimPEG/Mesh/NewTreeMesh.py
+6 -6
View File
@@ -1,5 +1,5 @@
import numpy as np, scipy.sparse as sp
from SimPEG.Utils import ndgrid, mkvc
from SimPEG.Utils import ndgrid, mkvc, sdiag
NUM, ACTIVE, NX, NY, NZ = range(5)
@@ -265,7 +265,7 @@ class TreeMesh(object):
q = np.sum((n3 - n1)**2,axis=1)**0.5
# Area of an arbitrary quadrilateral (in a plane)
V = 0.25 * (4.0*(p**2)*(q**2) - (a**2 + c**2 - b**2 - d**2))**0.5
V = 0.25 * (4.0*(p**2)*(q**2) - (a**2 + c**2 - b**2 - d**2)**2)**0.5
P = np.argsort(C[activeCells,NUM])
self._vol = V[P]
@@ -391,7 +391,7 @@ class TreeMesh(object):
VOL = self.vol
D = sp.csr_matrix((V,(I,J)), shape=(self.nC, self.nF))
S = self.area
self._faceDiv = Utils.sdiag(1/VOL)*D*Utils.sdiag(S)
self._faceDiv = sdiag(1/VOL)*D*sdiag(S)
return self._faceDiv
def number(self):
@@ -452,10 +452,10 @@ if __name__ == '__main__':
from SimPEG import Mesh, Utils
import matplotlib.pyplot as plt
tM = TreeMesh(np.ones(3),np.ones(2))
tM = TreeMesh([np.ones(3),np.ones(2)])
tM.refineFace(0)
tM.refineFace(9)
# tM.refineFace(0)
# tM.refineFace(9)
# print tM._faces
# print tM._edges[0,:]
SimPEG/Tests/test_NewTreeMesh.py
+34
View File
@@ -69,5 +69,39 @@ class TestQuadTreeMesh(unittest.TestCase):
self.assertTrue(np.linalg.norm((np.r_[ax,ay]-self.M.area)) < TOL)
class SimpleOctreeOperatorTests(unittest.TestCase):
def setUp(self):
h1 = np.random.rand(5)
h2 = np.random.rand(7)
h3 = np.random.rand(3)
# self.tM = TensorMesh([h1,h2,h3])
# self.oM = TreeMesh([h1,h2,h3])
self.tM2 = TensorMesh([h1,h2])
self.oM2 = TreeMesh([h1,h2])
# self.oM2.plotGrid(showIt=True)
def test_faceDiv(self):
# self.assertAlmostEqual((self.tM.faceDiv - self.oM.faceDiv).toarray().sum(), 0)
self.assertAlmostEqual((self.tM2.faceDiv - self.oM2.faceDiv).toarray().sum(), 0)
# def test_nodalGrad(self):
# self.assertAlmostEqual((self.tM.nodalGrad - self.oM.nodalGrad).toarray().sum(), 0)
# self.assertAlmostEqual((self.tM2.nodalGrad - self.oM2.nodalGrad).toarray().sum(), 0)
# def test_edgeCurl(self):
# self.assertAlmostEqual((self.tM.edgeCurl - self.oM.edgeCurl).toarray().sum(), 0)
# # self.assertAlmostEqual((self.tM2.edgeCurl - self.oM2.edgeCurl).toarray().sum(), 0)
# def test_InnerProducts(self):
# self.assertAlmostEqual((self.tM.getFaceInnerProduct() - self.oM.getFaceInnerProduct()).toarray().sum(), 0)
# self.assertAlmostEqual((self.tM2.getFaceInnerProduct() - self.oM2.getFaceInnerProduct()).toarray().sum(), 0)
# self.assertAlmostEqual((self.tM2.getEdgeInnerProduct() - self.oM2.getEdgeInnerProduct()).toarray().sum(), 0)
# self.assertAlmostEqual((self.tM.getEdgeInnerProduct() - self.oM.getEdgeInnerProduct()).toarray().sum(), 0)
if __name__ == '__main__':
unittest.main()