mirror of
https://github.com/wassname/simpeg.git
synced 2026-07-01 08:54:09 +08:00
Delete old implementations
This commit is contained in:
Rowan Cockett
File diff suppressed because it is too large
Load Diff
File diff suppressed because it is too large
Load Diff
+1628
-1102
File diff suppressed because it is too large
Load Diff
@@ -1,367 +0,0 @@
|
||||
from SimPEG.Mesh import TensorMesh
|
||||
from SimPEG.Mesh.NewTreeMesh import TreeMesh, TreeCell
|
||||
import numpy as np
|
||||
import unittest
|
||||
import matplotlib.pyplot as plt
|
||||
|
||||
TOL = 1e-10
|
||||
|
||||
class TestQuadTreeMesh(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
M = TreeMesh([np.ones(x) for x in [3,2]])
|
||||
M.refineFace(0)
|
||||
self.M = M
|
||||
M.number()
|
||||
# M.plotGrid(showIt=True)
|
||||
|
||||
def test_numbering(self):
|
||||
M = TreeMesh([2,2,2])
|
||||
M.number()
|
||||
M.refineCell(0)
|
||||
M.refineCell(3)
|
||||
assert M.isNumbered is False
|
||||
|
||||
def test_MeshSizes(self):
|
||||
self.assertTrue(self.M.nC==9)
|
||||
self.assertTrue(self.M.nF==25)
|
||||
self.assertTrue(self.M.nFx==12)
|
||||
self.assertTrue(self.M.nFy==13)
|
||||
self.assertTrue(self.M.nE==25)
|
||||
self.assertTrue(self.M.nEx==13)
|
||||
self.assertTrue(self.M.nEy==12)
|
||||
|
||||
def test_gridCC(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.5,1.5,2.5]
|
||||
y = np.r_[0.25,0.25,0.5,0.5,0.75,0.75,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridCC).flatten()) == 0)
|
||||
|
||||
def test_gridN(self):
|
||||
x = np.r_[0,0.5,1,2,3,0,0.5,1,0,0.5,1,2,3,0,1,2,3]
|
||||
y = np.r_[0,0,0,0,0,.5,.5,.5,1,1,1,1,1,2,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridN).flatten()) == 0)
|
||||
|
||||
def test_gridFx(self):
|
||||
x = np.r_[0.0,0.5,1.0,2.0,3.0,0.0,0.5,1.0,0.0,1.0,2.0,3.0]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.5,0.75,0.75,0.75,1.5,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridFx).flatten()) == 0)
|
||||
|
||||
def test_gridFy(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.25,0.75,1.5,2.5,0.5,1.5,2.5]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,1,1,1,1,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridFy).flatten()) == 0)
|
||||
|
||||
def test_gridEx(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.25,0.75,1.5,2.5,0.5,1.5,2.5]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,1,1,1,1,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridEx).flatten()) == 0)
|
||||
|
||||
def test_gridEy(self):
|
||||
x = np.r_[0.0,0.5,1.0,2.0,3.0,0.0,0.5,1.0,0.0,1.0,2.0,3.0]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.5,0.75,0.75,0.75,1.5,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridEy).flatten()) == 0)
|
||||
|
||||
def test_vol(self):
|
||||
v = np.r_[0.25,0.25,1,1,0.25,0.25,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((v-self.M.vol)) < TOL)
|
||||
|
||||
def test_edge(self):
|
||||
ex = np.r_[0.5,0.5,1,1,0.5,0.5,0.5,0.5,1,1,1,1,1]
|
||||
ey = np.r_[0.5,0.5,0.5,1,1,0.5,0.5,0.5,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.r_[ex,ey]-self.M.edge)) < TOL)
|
||||
|
||||
def test_area(self):
|
||||
ax = np.r_[0.5,0.5,0.5,1,1,0.5,0.5,0.5,1,1,1,1]
|
||||
ay = np.r_[0.5,0.5,1,1,0.5,0.5,0.5,0.5,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.r_[ax,ay]-self.M.area)) < TOL)
|
||||
|
||||
class TestOcTreeConnectivity(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
self.oM = TreeMesh([1,1,1])
|
||||
self.oM.refine(lambda c: 1)
|
||||
|
||||
def test_setup(self):
|
||||
C = TreeCell(self.oM, 0)
|
||||
children = C.children
|
||||
assert not C.isleaf
|
||||
assert len(children.index) == 8
|
||||
# assert not TreeCell(self.oM, 0).isleaf
|
||||
c0, c1, c2, c3, c4, c5, c6, c7 = [TreeCell(self.oM, i) for i in range(1,9)]
|
||||
|
||||
|
||||
# .----------------.----------------.
|
||||
# /| /| /|
|
||||
# / | / | / |
|
||||
# / | c6 / | c7 / |
|
||||
# / | / | / |
|
||||
# .----------------.----+-----------. |
|
||||
# /| . ---------/|----.----------/|----.
|
||||
# fZp / | /| / | /| / | /|
|
||||
# | / | / | c4 / | / | c5 / | X |
|
||||
# 6 ------eX3------ 7 / | / | / | / | / | / |
|
||||
# /| | / | . -------------- .----------------. |/ |
|
||||
# /eZ2 . / eZ3 | . ---+------|----.----+------|----. |
|
||||
# eY2 | fYp eY3 | | /| .______|___/|____.______|___/|____.
|
||||
# / | / fXp| | / | / c2 | / | / c3 | / | /
|
||||
# 4 ------eX2----- 5 | | / | / | / | / | / | /
|
||||
# |fXm 2 -----eX1--|---- 3 z . ---+---------- . ---+---------- . | /
|
||||
# eZ0 / | eY1 ^ y | |/ | |/ | |/
|
||||
# | eY0 . fYm eZ1 / | / | . ----------|----.-----------|----.
|
||||
# | / | | / | / | / c0 | / c1 | /
|
||||
# 0 ------eX0------1 o----> x | / | / | /
|
||||
# | | / | / | /
|
||||
# fZm . -------------- . -------------- .
|
||||
#
|
||||
#
|
||||
# fX fY fZ
|
||||
# 2___________3 2___________3 2___________3
|
||||
# | e1 | | e1 | | e1 |
|
||||
# | | | | | |
|
||||
# e2 | x | e3 z e2 | x | e3 z e2 | x | e3 y
|
||||
# | | ^ | | ^ | | ^
|
||||
# |___________| |___> y |___________| |___> x |___________| |___> x
|
||||
# 0 e0 1 0 e0 1 0 e0 1
|
||||
#
|
||||
|
||||
|
||||
# there are two faces for each edge
|
||||
for ii, c in enumerate([c0, c1, c2, c3, c4, c5, c6, c7]):
|
||||
assert c.fZm.e0.index == c.fYm.e0.index, "Cell %d: fZm.e0 and fYm.e0"%ii
|
||||
assert c.fZm.e1.index == c.fYp.e0.index, "Cell %d: fZm.e1 and fYp.e0"%ii
|
||||
assert c.fZp.e0.index == c.fYm.e1.index, "Cell %d: fZp.e0 and fYm.e1"%ii
|
||||
assert c.fZp.e1.index == c.fYp.e1.index, "Cell %d: fZp.e1 and fYp.e1"%ii
|
||||
assert c.fZm.e2.index == c.fXm.e0.index, "Cell %d: fZm.e2 and fXm.e0"%ii
|
||||
assert c.fZm.e3.index == c.fXp.e0.index, "Cell %d: fZm.e3 and fXp.e0"%ii
|
||||
assert c.fZp.e2.index == c.fXm.e1.index, "Cell %d: fZp.e2 and fXm.e1"%ii
|
||||
assert c.fZp.e3.index == c.fXp.e1.index, "Cell %d: fZp.e3 and fXp.e1"%ii
|
||||
assert c.fYm.e2.index == c.fXm.e2.index, "Cell %d: fYm.e2 and fXm.e2"%ii
|
||||
assert c.fYm.e3.index == c.fXp.e2.index, "Cell %d: fYm.e3 and fXp.e2"%ii
|
||||
assert c.fYp.e2.index == c.fXm.e3.index, "Cell %d: fYp.e2 and fXm.e3"%ii
|
||||
assert c.fYp.e3.index == c.fXp.e3.index, "Cell %d: fYp.e3 and fXp.e3"%ii
|
||||
|
||||
assert c0.eZ1.index == c1.eZ0.index
|
||||
assert c0.eZ3.index == c1.eZ2.index
|
||||
assert c2.eZ1.index == c3.eZ0.index
|
||||
assert c2.eZ3.index == c3.eZ2.index
|
||||
|
||||
assert c4.eZ1.index == c5.eZ0.index
|
||||
assert c4.eZ3.index == c5.eZ2.index
|
||||
assert c6.eZ1.index == c7.eZ0.index
|
||||
assert c6.eZ3.index == c7.eZ2.index
|
||||
|
||||
assert c0.n7.index == c7.n0.index
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
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.tM.getEdgeInnerProduct() - self.oM.getEdgeInnerProduct()).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)
|
||||
|
||||
|
||||
class SimpleOctreeOperatorTestsRefined(unittest.TestCase):
|
||||
def setUp(self):
|
||||
self.tM = TreeMesh([1,1,1])
|
||||
self.tM.refine(lambda c:2)
|
||||
self.M = TensorMesh([4,4,4])
|
||||
|
||||
self.tM2 = TreeMesh([1,1])
|
||||
self.tM2.refine(lambda c:2)
|
||||
self.M2 = TensorMesh([4,4])
|
||||
|
||||
def test_grids(self):
|
||||
self.assertAlmostEqual((self.tM2.gridN - self.M2.gridN).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM2.gridCC - self.M2.gridCC).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM2.gridFx - self.M2.gridFx).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM2.gridFy - self.M2.gridFy).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM2.gridEx - self.M2.gridEx).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM2.gridEy - self.M2.gridEy).sum(), 0)
|
||||
|
||||
self.assertAlmostEqual((self.tM.gridN - self.M.gridN).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridCC - self.M.gridCC).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridFx - self.M.gridFx).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridFy - self.M.gridFy).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridFz - self.M.gridFz).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridEx - self.M.gridEx).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridEy - self.M.gridEy).sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.gridEz - self.M.gridEz).sum(), 0)
|
||||
|
||||
def test_InnerProducts(self):
|
||||
self.assertAlmostEqual((self.tM.getFaceInnerProduct() - self.M.getFaceInnerProduct()).toarray().sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.getEdgeInnerProduct() - self.M.getEdgeInnerProduct()).toarray().sum(), 0)
|
||||
|
||||
# self.assertAlmostEqual((self.tM2.getFaceInnerProduct() - self.M2.getFaceInnerProduct()).toarray().sum(), 0)
|
||||
# self.assertAlmostEqual((self.tM2.getEdgeInnerProduct() - self.M2.getEdgeInnerProduct()).toarray().sum(), 0)
|
||||
|
||||
def test_faceDiv(self):
|
||||
self.assertAlmostEqual((self.tM2.faceDiv - self.M2.faceDiv).toarray().sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.faceDiv - self.M.faceDiv).toarray().sum(), 0)
|
||||
|
||||
def test_nodalGrad(self):
|
||||
self.assertAlmostEqual((self.tM2.nodalGrad - self.M2.nodalGrad).toarray().sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.nodalGrad - self.M.nodalGrad).toarray().sum(), 0)
|
||||
|
||||
def test_edgeCurl(self):
|
||||
self.assertAlmostEqual((self.tM.edgeCurl - self.M.edgeCurl).toarray().sum(), 0)
|
||||
# self.assertAlmostEqual((self.tM2.edgeCurl - self.M2.edgeCurl).toarray().sum(), 0)
|
||||
|
||||
def test_InnerProducts(self):
|
||||
self.assertAlmostEqual((self.tM.getFaceInnerProduct() - self.M.getFaceInnerProduct()).toarray().sum(), 0)
|
||||
self.assertAlmostEqual((self.tM.getEdgeInnerProduct() - self.M.getEdgeInnerProduct()).toarray().sum(), 0)
|
||||
# self.assertAlmostEqual((self.tM2.getFaceInnerProduct() - self.M2.getFaceInnerProduct()).toarray().sum(), 0)
|
||||
# self.assertAlmostEqual((self.tM2.getEdgeInnerProduct() - self.M2.getEdgeInnerProduct()).toarray().sum(), 0)
|
||||
|
||||
|
||||
class TestOcTreeObjects(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
self.M = TreeMesh([2,1,1])
|
||||
self.M.number()
|
||||
|
||||
self.Mr = TreeMesh([2,1,1])
|
||||
self.Mr.refineCell(0)
|
||||
self.Mr.number()
|
||||
|
||||
def test_counts(self):
|
||||
self.assertTrue(self.M.nC == 2)
|
||||
self.assertTrue(self.M.nFx == 3)
|
||||
self.assertTrue(self.M.nFy == 4)
|
||||
self.assertTrue(self.M.nFz == 4)
|
||||
self.assertTrue(self.M.nF == 11)
|
||||
self.assertTrue(self.M.nEx == 8)
|
||||
self.assertTrue(self.M.nEy == 6)
|
||||
self.assertTrue(self.M.nEz == 6)
|
||||
self.assertTrue(self.M.nE == 20)
|
||||
self.assertTrue(self.M.nN == 12)
|
||||
|
||||
self.assertTrue(self.Mr.nC == 9)
|
||||
self.assertTrue(self.Mr.nFx == 13)
|
||||
self.assertTrue(self.Mr.nFy == 14)
|
||||
self.assertTrue(self.Mr.nFz == 14)
|
||||
self.assertTrue(self.Mr.nF == 41)
|
||||
|
||||
self.assertTrue(self.Mr.nN == 31)
|
||||
self.assertTrue(self.Mr.nEx == 22)
|
||||
self.assertTrue(self.Mr.nEy == 20)
|
||||
self.assertTrue(self.Mr.nEz == 20)
|
||||
|
||||
|
||||
def test_gridCC(self):
|
||||
x = np.r_[0.25,0.75]
|
||||
y = np.r_[0.5,0.5]
|
||||
z = np.r_[0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridCC).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375]
|
||||
y = np.r_[0.25,0.25,0.5,0.75,0.75,0.25,0.25,0.75,0.75]
|
||||
z = np.r_[0.25,0.25,0.5,0.25,0.25,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridCC).flatten()) == 0)
|
||||
|
||||
def test_gridN(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1,0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0,0,0,1,1,1,0,0,0,1,1,1.]
|
||||
z = np.r_[0,0,0,0,0,0,1,1,1,1,1,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridN).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,1]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,0.5,1,1,1,1,0,0,0,0.5,0.5,0.5,1,1,1,0,0,0,0,0.5,0.5,0.5,1,1,1,1]
|
||||
z = np.r_[0,0,0,0,0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridN).flatten()) == 0)
|
||||
|
||||
def test_gridFx(self):
|
||||
x = np.r_[0.0,0.5,1.0]
|
||||
y = np.r_[0.5,0.5,0.5]
|
||||
z = np.r_[0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFx).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.0,0.25,0.5,1.0,0.0,0.25,0.5,0.0,0.25,0.5,0.0,0.25,0.5]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.25,0.25,0.25,0.75,0.75,0.75]
|
||||
z = np.r_[0.25,0.25,0.25,0.5,0.25,0.25,0.25,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFx).flatten()) == 0)
|
||||
|
||||
def test_gridFy(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0,0,1.,1.]
|
||||
z = np.r_[0.5,0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFy).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375]
|
||||
y = np.r_[0,0,0,0.5,0.5,1,1,1,0,0,0.5,0.5,1,1]
|
||||
z = np.r_[0.25,0.25,0.5,0.25,0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFy).flatten()) == 0)
|
||||
|
||||
def test_gridFz(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0.5,0.5,0.5,0.5]
|
||||
z = np.r_[0,0,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFz).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375,0.125,0.375,0.75,0.125,0.375]
|
||||
y = np.r_[0.25,0.25,0.5,0.75,0.75,0.25,0.25,0.75,0.75,0.25,0.25,0.5,0.75,0.75]
|
||||
z = np.r_[0,0,0,0,0,0.5,0.5,0.5,0.5,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFz).flatten()) == 0)
|
||||
|
||||
|
||||
def test_gridEx(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75,0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0,0,1.,1.,0,0,1.,1.]
|
||||
z = np.r_[0,0,0,0,1.,1.,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEx).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75]
|
||||
y = np.r_[0,0,0,0.5,0.5,1,1,1,0,0,0.5,0.5,1,1,0,0,0,0.5,0.5,1,1,1]
|
||||
z = np.r_[0,0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEx).flatten()) == 0)
|
||||
|
||||
def test_gridEy(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0.5,0.5,0.5,0.5,0.5,0.5]
|
||||
z = np.r_[0,0,0,1.,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEy).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.25,0.25,0.25,0.75,0.75,0.75,0.25,0.25,0.25,0.5,0.75,0.75,0.75]
|
||||
z = np.r_[0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEy).flatten()) == 0)
|
||||
|
||||
def test_gridEz(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0,0,0,1.,1.,1.]
|
||||
z = np.r_[0.5,0.5,0.5,0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEz).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0 ,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0 ,0.25,0.5,0 ,0.25,0.5]
|
||||
y = np.r_[0,0 ,0 ,0,0.5,0.5 ,0.5,1,1 ,1 ,1,0,0 ,0 ,0.5,0.5 ,0.5,1 ,1 ,1 ]
|
||||
z = np.r_[0.25,0.25,0.25,0.5,0.25,0.25,0.25,0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEz).flatten()) == 0)
|
||||
|
||||
if __name__ == '__main__':
|
||||
unittest.main()
|
||||
+135
-481
@@ -1,504 +1,158 @@
|
||||
from SimPEG.Mesh import TensorMesh
|
||||
from SimPEG.Mesh.TreeMesh import TreeMesh, TreeFace, TreeCell
|
||||
from SimPEG import Mesh
|
||||
import numpy as np
|
||||
import unittest
|
||||
import matplotlib.pyplot as plt
|
||||
import unittest
|
||||
|
||||
TOL = 1e-10
|
||||
|
||||
class TestOcTreeObjects(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
self.M = TreeMesh([2,1,1])
|
||||
self.M.number()
|
||||
|
||||
self.Mr = TreeMesh([2,1,1])
|
||||
self.Mr.children[0,0,0].refine()
|
||||
self.Mr.number()
|
||||
|
||||
def q(s):
|
||||
if s[0] == 'M':
|
||||
m = self.M
|
||||
s = s[1:]
|
||||
else:
|
||||
m = self.Mr
|
||||
c = m.sortedCells[int(s[1])]
|
||||
if len(s) == 2: return c
|
||||
if s[2] == 'f' and len(s) == 5: return c.faceDict[s[2:]]
|
||||
if s[2] == 'f': return getattr(c.faceDict[s[2:5]], 'edg' +s[5:])
|
||||
if s[2] == 'e': return getattr(c,s[2:])
|
||||
if s[2] == 'n': return getattr(c,'node'+s[3:])
|
||||
|
||||
self.q = q
|
||||
class TestSimpleQuadTree(unittest.TestCase):
|
||||
|
||||
def test_counts(self):
|
||||
self.assertTrue(self.M.nC == 2)
|
||||
self.assertTrue(self.M.nFx == 3)
|
||||
self.assertTrue(self.M.nFy == 4)
|
||||
self.assertTrue(self.M.nFz == 4)
|
||||
self.assertTrue(self.M.nF == 11)
|
||||
self.assertTrue(self.M.nEx == 8)
|
||||
self.assertTrue(self.M.nEy == 6)
|
||||
self.assertTrue(self.M.nEz == 6)
|
||||
self.assertTrue(self.M.nE == 20)
|
||||
self.assertTrue(self.M.nN == 12)
|
||||
|
||||
self.assertTrue(self.Mr.nC == 9)
|
||||
self.assertTrue(self.Mr.nFx == 13)
|
||||
self.assertTrue(self.Mr.nFy == 14)
|
||||
self.assertTrue(self.Mr.nFz == 14)
|
||||
self.assertTrue(self.Mr.nF == 41)
|
||||
|
||||
|
||||
for cell in self.Mr.sortedCells:
|
||||
for e in cell.edgeDict:
|
||||
self.assertTrue(cell.edgeDict[e].edgeType==e[1].lower())
|
||||
|
||||
self.assertTrue(self.Mr.nN == 31)
|
||||
self.assertTrue(self.Mr.nEx == 22)
|
||||
self.assertTrue(self.Mr.nEy == 20)
|
||||
self.assertTrue(self.Mr.nEz == 20)
|
||||
|
||||
def test_sizes(self):
|
||||
q = self.q
|
||||
|
||||
for key in ['Mc0','Mc1']:
|
||||
self.assertTrue(q(key).vol == 0.5)
|
||||
self.assertTrue(q(key+'fXm').area == 1.)
|
||||
self.assertTrue(q(key+'fXp').area == 1.)
|
||||
self.assertTrue(q(key+'fYm').area == 0.5)
|
||||
self.assertTrue(q(key+'fYp').area == 0.5)
|
||||
self.assertTrue(q(key+'fZm').area == 0.5)
|
||||
self.assertTrue(q(key+'fZp').area == 0.5)
|
||||
|
||||
def test_pointersM(self):
|
||||
q = self.q
|
||||
|
||||
self.assertTrue(q('Mc0fXp') is q('Mc1fXm'))
|
||||
self.assertTrue(q('Mc0fXpe0') is q('Mc1fXme0'))
|
||||
self.assertTrue(q('Mc0fXpe1') is q('Mc1fXme1'))
|
||||
self.assertTrue(q('Mc0fXpe2') is q('Mc1fXme2'))
|
||||
self.assertTrue(q('Mc0fXpe3') is q('Mc1fXme3'))
|
||||
self.assertTrue(q('Mc0fYp') is not q('c1fYm'))
|
||||
self.assertTrue(q('Mc0fXm') is not q('c1fXm'))
|
||||
|
||||
# Test connectivity of shared edges
|
||||
self.assertTrue(q('Mc0fZpe3') is not q('c1fZpe0'))
|
||||
self.assertTrue(q('Mc0fZpe3') is not q('c1fZpe1'))
|
||||
self.assertTrue(q('Mc0fZpe3') is q('Mc1fZpe2'))
|
||||
self.assertTrue(q('Mc0fZpe3') is not q('c1fZpe3'))
|
||||
|
||||
self.assertTrue(q('Mc0fZme3') is not q('c1fZme0'))
|
||||
self.assertTrue(q('Mc0fZme3') is not q('c1fZme1'))
|
||||
self.assertTrue(q('Mc0fZme3') is q('Mc1fZme2'))
|
||||
self.assertTrue(q('Mc0fZme3') is not q('c1fZme3'))
|
||||
|
||||
self.assertTrue(q('Mc0fYpe3') is not q('c1fYpe0'))
|
||||
self.assertTrue(q('Mc0fYpe3') is not q('c1fYpe1'))
|
||||
self.assertTrue(q('Mc0fYpe3') is q('Mc1fYpe2'))
|
||||
self.assertTrue(q('Mc0fYpe3') is not q('c1fYpe3'))
|
||||
|
||||
self.assertTrue(q('Mc0fYme3') is not q('c1fYme0'))
|
||||
self.assertTrue(q('Mc0fYme3') is not q('c1fYme1'))
|
||||
self.assertTrue(q('Mc0fYme3') is q('Mc1fYme2'))
|
||||
self.assertTrue(q('Mc0fYme3') is not q('c1fYme3'))
|
||||
|
||||
self.assertTrue(q('Mc0fZme3') is q('Mc1fXme0'))
|
||||
self.assertTrue(q('Mc0fZpe3') is q('Mc1fXme1'))
|
||||
self.assertTrue(q('Mc0fYme3') is q('Mc1fXme2'))
|
||||
self.assertTrue(q('Mc0fYpe3') is q('Mc1fXme3'))
|
||||
|
||||
self.assertTrue(q('Mc0fZme3') is q('Mc0fXpe0'))
|
||||
self.assertTrue(q('Mc0fZpe3') is q('Mc0fXpe1'))
|
||||
self.assertTrue(q('Mc0fYme3') is q('Mc0fXpe2'))
|
||||
self.assertTrue(q('Mc0fYpe3') is q('Mc0fXpe3'))
|
||||
|
||||
self.assertTrue(q('Mc1fZme2') is q('Mc1fXme0'))
|
||||
self.assertTrue(q('Mc1fZpe2') is q('Mc1fXme1'))
|
||||
self.assertTrue(q('Mc1fYme2') is q('Mc1fXme2'))
|
||||
self.assertTrue(q('Mc1fYpe2') is q('Mc1fXme3'))
|
||||
|
||||
self.assertTrue(q('Mc1fZme2') is q('Mc0fXpe0'))
|
||||
self.assertTrue(q('Mc1fZpe2') is q('Mc0fXpe1'))
|
||||
self.assertTrue(q('Mc1fYme2') is q('Mc0fXpe2'))
|
||||
self.assertTrue(q('Mc1fYpe2') is q('Mc0fXpe3'))
|
||||
|
||||
|
||||
def test_nodePointers(self):
|
||||
q = self.q
|
||||
c0 = self.Mr.sortedCells[0]
|
||||
c0n0 = c0.node0
|
||||
self.assertTrue(c0n0 is q('c0n0'))
|
||||
self.assertTrue(np.all(q('c0n0').center == np.r_[0,0,0.]))
|
||||
self.assertTrue(q('c0n0').num == 0)
|
||||
self.assertTrue(q('c0n1').num == 1)
|
||||
self.assertTrue(q('c0n2').num == 4)
|
||||
self.assertTrue(q('c0n3').num == 5)
|
||||
self.assertTrue(q('c0n4').num == 11)
|
||||
self.assertTrue(q('c0n5').num == 12)
|
||||
self.assertTrue(q('c0n6').num == 14)
|
||||
self.assertTrue(q('c0n7').num == 15)
|
||||
|
||||
def test_pointersMr(self):
|
||||
q = self.q
|
||||
|
||||
c0 = self.Mr.sortedCells[0]
|
||||
c0fXm = c0.fXm
|
||||
c0eX0 = c0.eX0
|
||||
c0fYme0 = c0.fYm.edge0
|
||||
self.assertTrue(c0 is q('c0'))
|
||||
self.assertTrue(c0fXm is q('c0fXm'))
|
||||
self.assertTrue(c0eX0 is q('c0eX0'))
|
||||
self.assertTrue(c0fYme0 is q('c0fYme0'))
|
||||
|
||||
self.assertTrue(q('c0').depth == 1)
|
||||
self.assertTrue(q('c1').depth == 1)
|
||||
self.assertTrue(q('c2').depth == 0)
|
||||
|
||||
# Make sure we know where the center of the cells are.
|
||||
self.assertTrue(np.all(q('c0').center == np.r_[0.125,0.25,0.25]))
|
||||
self.assertTrue(np.all(q('c1').center == np.r_[0.375,0.25,0.25]))
|
||||
self.assertTrue(np.all(q('c2').center == np.r_[0.75,0.5,0.5]))
|
||||
self.assertTrue(np.all(q('c3').center == np.r_[0.125,0.75,0.25]))
|
||||
self.assertTrue(np.all(q('c4').center == np.r_[0.375,0.75,0.25]))
|
||||
self.assertTrue(np.all(q('c5').center == np.r_[0.125,0.25,0.75]))
|
||||
self.assertTrue(np.all(q('c6').center == np.r_[0.375,0.25,0.75]))
|
||||
self.assertTrue(np.all(q('c7').center == np.r_[0.125,0.75,0.75]))
|
||||
self.assertTrue(np.all(q('c8').center == np.r_[0.375,0.75,0.75]))
|
||||
|
||||
# Test X face connectivity and locations and stuff...
|
||||
self.assertTrue(np.all(q('c0fXm').center == np.r_[0,0.25,0.25]))
|
||||
self.assertTrue(np.all(q('c0fXp').center == np.r_[0.25,0.25,0.25]))
|
||||
self.assertTrue(q('c0fXp') is q('c1fXm'))
|
||||
self.assertTrue(np.all(q('c1fXp').center == np.r_[0.5,0.25,0.25]))
|
||||
self.assertTrue(np.all(q('c2fXm').center == np.r_[0.5,0.5,0.5]))
|
||||
self.assertTrue(q('c2fXm').branchdepth == 1)
|
||||
self.assertTrue(q('c2fXm').children[0,0] is q('c1fXp'))
|
||||
self.assertTrue(np.all(q('c3fXm').center == np.r_[0,0.75,0.25]))
|
||||
self.assertTrue(np.all(q('c3fXp').center == np.r_[0.25,0.75,0.25]))
|
||||
self.assertTrue(q('c4fXm') is q('c3fXp'))
|
||||
self.assertTrue(q('c2fXm').children[1,0] is q('c4fXp'))
|
||||
|
||||
#Test some internal stuff (edges held by cell should be same as inside)
|
||||
for key in ['Mc0', 'Mc1'] + ['c%d'%i for i in range(9)]:
|
||||
self.assertTrue(q(key+'eX0') is q(key+'fZme0'))
|
||||
self.assertTrue(q(key+'eX1') is q(key+'fZme1'))
|
||||
self.assertTrue(q(key+'eX2') is q(key+'fZpe0'))
|
||||
self.assertTrue(q(key+'eX3') is q(key+'fZpe1'))
|
||||
|
||||
self.assertTrue(q(key+'eX0') is q(key+'fYme0'))
|
||||
self.assertTrue(q(key+'eX1') is q(key+'fYpe0'))
|
||||
self.assertTrue(q(key+'eX2') is q(key+'fYme1'))
|
||||
self.assertTrue(q(key+'eX3') is q(key+'fYpe1'))
|
||||
|
||||
self.assertTrue(q(key+'eY0') is q(key+'fXme0'))
|
||||
self.assertTrue(q(key+'eY1') is q(key+'fXpe0'))
|
||||
self.assertTrue(q(key+'eY2') is q(key+'fXme1'))
|
||||
self.assertTrue(q(key+'eY3') is q(key+'fXpe1'))
|
||||
|
||||
self.assertTrue(q(key+'eY0') is q(key+'fZme2'))
|
||||
self.assertTrue(q(key+'eY1') is q(key+'fZme3'))
|
||||
self.assertTrue(q(key+'eY2') is q(key+'fZpe2'))
|
||||
self.assertTrue(q(key+'eY3') is q(key+'fZpe3'))
|
||||
|
||||
self.assertTrue(q(key+'eZ0') is q(key+'fXme2'))
|
||||
self.assertTrue(q(key+'eZ1') is q(key+'fXpe2'))
|
||||
self.assertTrue(q(key+'eZ2') is q(key+'fXme3'))
|
||||
self.assertTrue(q(key+'eZ3') is q(key+'fXpe3'))
|
||||
|
||||
self.assertTrue(q(key+'eZ0') is q(key+'fYme2'))
|
||||
self.assertTrue(q(key+'eZ1') is q(key+'fYme3'))
|
||||
self.assertTrue(q(key+'eZ2') is q(key+'fYpe2'))
|
||||
self.assertTrue(q(key+'eZ3') is q(key+'fYpe3'))
|
||||
|
||||
#Test some edge stuff
|
||||
self.assertTrue(np.all(q('c0eX0').center == np.r_[0.125,0,0]))
|
||||
self.assertTrue(np.all(q('c0eX1').center == np.r_[0.125,0.5,0]))
|
||||
self.assertTrue(np.all(q('c0eX2').center == np.r_[0.125,0,0.5]))
|
||||
self.assertTrue(np.all(q('c0eX3').center == np.r_[0.125,0.5,0.5]))
|
||||
|
||||
self.assertTrue(np.all(q('c5eX0').center == np.r_[0.125,0,0.5]))
|
||||
self.assertTrue(np.all(q('c5eX1').center == np.r_[0.125,0.5,0.5]))
|
||||
self.assertTrue(q('c5eX0') is q('c0eX2'))
|
||||
self.assertTrue(q('c5eX1') is q('c0eX3'))
|
||||
|
||||
self.assertTrue(np.all(q('c0eY0').center == np.r_[0,0.25,0]))
|
||||
self.assertTrue(np.all(q('c0eY1').center == np.r_[0.25,0.25,0]))
|
||||
self.assertTrue(np.all(q('c0eY2').center == np.r_[0,0.25,0.5]))
|
||||
self.assertTrue(np.all(q('c0eY3').center == np.r_[0.25,0.25,0.5]))
|
||||
|
||||
self.assertTrue(np.all(q('c1eY0').center == np.r_[0.25,0.25,0]))
|
||||
self.assertTrue(np.all(q('c1eY2').center == np.r_[0.25,0.25,0.5]))
|
||||
self.assertTrue(q('c1eY0') is q('c0eY1'))
|
||||
self.assertTrue(q('c1eY2') is q('c0eY3'))
|
||||
|
||||
|
||||
self.assertTrue(np.all(q('c0eZ0').center == np.r_[0,0,0.25]))
|
||||
self.assertTrue(np.all(q('c0eZ1').center == np.r_[0.25,0,0.25]))
|
||||
self.assertTrue(np.all(q('c0eZ2').center == np.r_[0,0.5,0.25]))
|
||||
self.assertTrue(np.all(q('c0eZ3').center == np.r_[0.25,0.5,0.25]))
|
||||
|
||||
self.assertTrue(np.all(q('c3eZ0').center == np.r_[0,0.5,0.25]))
|
||||
self.assertTrue(np.all(q('c3eZ1').center == np.r_[0.25,0.5,0.25]))
|
||||
self.assertTrue(q('c3eZ0') is q('c0eZ2'))
|
||||
self.assertTrue(q('c3eZ1') is q('c0eZ3'))
|
||||
|
||||
|
||||
self.assertTrue(q('c0fXp') is q('c1fXm'))
|
||||
self.assertTrue(q('c0fYp') is not q('c1fYm'))
|
||||
self.assertTrue(q('c0fXm') is not q('c1fXm'))
|
||||
|
||||
self.assertTrue(q('c1fXp') is q('c2fXm').children[0,0])
|
||||
|
||||
self.assertTrue(q('c1fYp') is q('c4fYm'))
|
||||
self.assertTrue(q('c1fZp') is q('c6fZm'))
|
||||
|
||||
self.assertTrue(q('c6fXp') is q('c2fXm').children[0,1])
|
||||
|
||||
self.assertTrue(q('c4fXp') is q('c2fXm').children[1,0])
|
||||
|
||||
|
||||
def test_gridCC(self):
|
||||
x = np.r_[0.25,0.75]
|
||||
y = np.r_[0.5,0.5]
|
||||
z = np.r_[0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridCC).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375]
|
||||
y = np.r_[0.25,0.25,0.5,0.75,0.75,0.25,0.25,0.75,0.75]
|
||||
z = np.r_[0.25,0.25,0.5,0.25,0.25,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridCC).flatten()) == 0)
|
||||
|
||||
def test_gridN(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1,0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0,0,0,1,1,1,0,0,0,1,1,1.]
|
||||
z = np.r_[0,0,0,0,0,0,1,1,1,1,1,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridN).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,1]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,0.5,1,1,1,1,0,0,0,0.5,0.5,0.5,1,1,1,0,0,0,0,0.5,0.5,0.5,1,1,1,1]
|
||||
z = np.r_[0,0,0,0,0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridN).flatten()) == 0)
|
||||
|
||||
def test_gridFx(self):
|
||||
x = np.r_[0.0,0.5,1.0]
|
||||
y = np.r_[0.5,0.5,0.5]
|
||||
z = np.r_[0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFx).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.0,0.25,0.5,1.0,0.0,0.25,0.5,0.0,0.25,0.5,0.0,0.25,0.5]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.25,0.25,0.25,0.75,0.75,0.75]
|
||||
z = np.r_[0.25,0.25,0.25,0.5,0.25,0.25,0.25,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFx).flatten()) == 0)
|
||||
|
||||
def test_gridFy(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0,0,1.,1.]
|
||||
z = np.r_[0.5,0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFy).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375]
|
||||
y = np.r_[0,0,0,0.5,0.5,1,1,1,0,0,0.5,0.5,1,1]
|
||||
z = np.r_[0.25,0.25,0.5,0.25,0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFy).flatten()) == 0)
|
||||
|
||||
def test_gridFz(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0.5,0.5,0.5,0.5]
|
||||
z = np.r_[0,0,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridFz).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375,0.125,0.375,0.75,0.125,0.375]
|
||||
y = np.r_[0.25,0.25,0.5,0.75,0.75,0.25,0.25,0.75,0.75,0.25,0.25,0.5,0.75,0.75]
|
||||
z = np.r_[0,0,0,0,0,0.5,0.5,0.5,0.5,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridFz).flatten()) == 0)
|
||||
|
||||
|
||||
def test_gridEx(self):
|
||||
x = np.r_[0.25,0.75,0.25,0.75,0.25,0.75,0.25,0.75]
|
||||
y = np.r_[0,0,1.,1.,0,0,1.,1.]
|
||||
z = np.r_[0,0,0,0,1.,1.,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEx).flatten()) == 0)
|
||||
|
||||
x = np.r_[0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.125,0.375,0.125,0.375,0.75,0.125,0.375,0.125,0.375,0.75]
|
||||
y = np.r_[0,0,0,0.5,0.5,1,1,1,0,0,0.5,0.5,1,1,0,0,0,0.5,0.5,1,1,1]
|
||||
z = np.r_[0,0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEx).flatten()) == 0)
|
||||
|
||||
def test_gridEy(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0.5,0.5,0.5,0.5,0.5,0.5]
|
||||
z = np.r_[0,0,0,1.,1.,1.]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEy).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.25,0.25,0.25,0.75,0.75,0.75,0.25,0.25,0.25,0.5,0.75,0.75,0.75]
|
||||
z = np.r_[0,0,0,0,0,0,0,0.5,0.5,0.5,0.5,0.5,0.5,1,1,1,1,1,1,1]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEy).flatten()) == 0)
|
||||
|
||||
def test_gridEz(self):
|
||||
x = np.r_[0,0.5,1,0,0.5,1]
|
||||
y = np.r_[0,0,0,1.,1.,1.]
|
||||
z = np.r_[0.5,0.5,0.5,0.5,0.5,0.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.M.gridEz).flatten()) == 0)
|
||||
|
||||
x = np.r_[0,0.25,0.5,1,0 ,0.25,0.5,0,0.25,0.5,1,0,0.25,0.5,0 ,0.25,0.5,0 ,0.25,0.5]
|
||||
y = np.r_[0,0 ,0 ,0,0.5,0.5 ,0.5,1,1 ,1 ,1,0,0 ,0 ,0.5,0.5 ,0.5,1 ,1 ,1 ]
|
||||
z = np.r_[0.25,0.25,0.25,0.5,0.25,0.25,0.25,0.25,0.25,0.25,0.5,0.75,0.75,0.75,0.75,0.75,0.75,0.75,0.75,0.75]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y,z]-self.Mr.gridEz).flatten()) == 0)
|
||||
|
||||
|
||||
class TestQuadTreeObjects(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
self.M = TreeMesh([2,1])
|
||||
self.Mr = TreeMesh([2,1])
|
||||
self.Mr.children[0,0].refine()
|
||||
self.Mr.number()
|
||||
# self.Mr.plotGrid(showIt=True)
|
||||
|
||||
def test_pointersM(self):
|
||||
c0 = self.M.children[0,0]
|
||||
c0fXm = c0.fXm
|
||||
c0fXp = c0.fXp
|
||||
c0fYm = c0.fYm
|
||||
c0fYp = c0.fYp
|
||||
|
||||
c1 = self.M.children[1,0]
|
||||
c1fXm = c1.fXm
|
||||
c1fXp = c1.fXp
|
||||
c1fYm = c1.fYm
|
||||
c1fYp = c1.fYp
|
||||
|
||||
self.assertTrue(c0fXp is c1fXm)
|
||||
self.assertTrue(c0fYp is not c1fYm)
|
||||
self.assertTrue(c0fXm is not c1fXm)
|
||||
|
||||
self.assertTrue(c0fXm.area == 1)
|
||||
self.assertTrue(c0fYm.area == 0.5)
|
||||
|
||||
self.assertTrue(c0.node1 is c1.node0)
|
||||
self.assertTrue(c0.node3 is c1.node2)
|
||||
self.assertTrue(self.M.nN == 6)
|
||||
|
||||
|
||||
def test_pointersMr(self):
|
||||
c0 = self.Mr.sortedCells[0]
|
||||
c0fXm = c0.fXm
|
||||
c0fXp = c0.fXp
|
||||
c0fYm = c0.fYm
|
||||
c0fYp = c0.fYp
|
||||
|
||||
c1 = self.Mr.sortedCells[1]
|
||||
c1fXm = c1.fXm
|
||||
c1fXp = c1.fXp
|
||||
c1fYm = c1.fYm
|
||||
c1fYp = c1.fYp
|
||||
|
||||
c2 = self.Mr.sortedCells[2]
|
||||
c2fXm = c2.fXm
|
||||
c2fXp = c2.fXp
|
||||
c2fYm = c2.fYm
|
||||
c2fYp = c2.fYp
|
||||
|
||||
c4 = self.Mr.sortedCells[4]
|
||||
c4fXm = c4.fXm
|
||||
c4fXp = c4.fXp
|
||||
c4fYm = c4.fYm
|
||||
c4fYp = c4.fYp
|
||||
|
||||
self.assertTrue(c0fXp is c1fXm)
|
||||
self.assertTrue(c1fXp.node0 is c2fXm.node0)
|
||||
self.assertTrue(c1fXp.node0 is c2fXm.node0)
|
||||
self.assertTrue(c4fYm is c1fYp)
|
||||
self.assertTrue(c4fXp.node1 is c2fXm.node1)
|
||||
self.assertTrue(c4fXp.node0 is c1fYp.node1)
|
||||
self.assertTrue(c0fXp.node1 is c4fYm.node0)
|
||||
|
||||
self.assertTrue(self.Mr.nN == 11)
|
||||
|
||||
self.assertTrue(np.all(c1fXp.node0.x0 == np.r_[0.5,0]))
|
||||
self.assertTrue(np.all(c1fYp.node0.x0 == np.r_[0.25,0.5]))
|
||||
|
||||
|
||||
class TestQuadTreeMesh(unittest.TestCase):
|
||||
|
||||
def setUp(self):
|
||||
M = TreeMesh([np.ones(x) for x in [3,2]])
|
||||
for ii in range(1):
|
||||
M.children[ii,ii].refine()
|
||||
self.M = M
|
||||
nc = 8
|
||||
h1 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h2 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h = [hi/np.sum(hi) for hi in [h1, h2]] # normalize
|
||||
M = Mesh.TreeMesh(h)
|
||||
M._refineCell([0,0,0])
|
||||
M._refineCell([0,0,1])
|
||||
M.number()
|
||||
# M.plotGrid(showIt=True)
|
||||
assert M.nhFx == 2
|
||||
assert M.nFx == 9
|
||||
|
||||
def test_MeshSizes(self):
|
||||
self.assertTrue(self.M.nC==9)
|
||||
self.assertTrue(self.M.nF==25)
|
||||
self.assertTrue(self.M.nFx==12)
|
||||
self.assertTrue(self.M.nFy==13)
|
||||
self.assertTrue(self.M.nE==25)
|
||||
self.assertTrue(self.M.nEx==13)
|
||||
self.assertTrue(self.M.nEy==12)
|
||||
assert np.allclose(M.vol.sum(), 1.0)
|
||||
|
||||
def test_gridCC(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.5,1.5,2.5]
|
||||
y = np.r_[0.25,0.25,0.5,0.5,0.75,0.75,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridCC).flatten()) == 0)
|
||||
assert np.allclose(np.r_[M._areaFxFull, M._areaFyFull], M._deflationMatrix('F') * M.area)
|
||||
|
||||
def test_gridN(self):
|
||||
x = np.r_[0,0.5,1,2,3,0,0.5,1,0,0.5,1,2,3,0,1,2,3]
|
||||
y = np.r_[0,0,0,0,0,.5,.5,.5,1,1,1,1,1,2,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridN).flatten()) == 0)
|
||||
|
||||
def test_gridFx(self):
|
||||
x = np.r_[0.0,0.5,1.0,2.0,3.0,0.0,0.5,1.0,0.0,1.0,2.0,3.0]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.5,0.75,0.75,0.75,1.5,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridFx).flatten()) == 0)
|
||||
|
||||
def test_gridFy(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.25,0.75,1.5,2.5,0.5,1.5,2.5]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,1,1,1,1,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridFy).flatten()) == 0)
|
||||
|
||||
def test_gridEx(self):
|
||||
x = np.r_[0.25,0.75,1.5,2.5,0.25,0.75,0.25,0.75,1.5,2.5,0.5,1.5,2.5]
|
||||
y = np.r_[0,0,0,0,0.5,0.5,1,1,1,1,2,2,2]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridEx).flatten()) == 0)
|
||||
|
||||
def test_gridEy(self):
|
||||
x = np.r_[0.0,0.5,1.0,2.0,3.0,0.0,0.5,1.0,0.0,1.0,2.0,3.0]
|
||||
y = np.r_[0.25,0.25,0.25,0.5,0.5,0.75,0.75,0.75,1.5,1.5,1.5,1.5]
|
||||
self.assertTrue(np.linalg.norm((np.c_[x,y]-self.M.gridEy).flatten()) == 0)
|
||||
|
||||
|
||||
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])
|
||||
|
||||
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)
|
||||
hx, hy = np.r_[1.,2,3,4], np.r_[5.,6,7,8]
|
||||
T = Mesh.TreeMesh([hx, hy], levels=2)
|
||||
T.refine(lambda xc:2)
|
||||
# T.plotGrid(showIt=True)
|
||||
M = Mesh.TensorMesh([hx, hy])
|
||||
assert M.nC == T.nC
|
||||
assert M.nF == T.nF
|
||||
assert M.nFx == T.nFx
|
||||
assert M.nFy == T.nFy
|
||||
assert M.nE == T.nE
|
||||
assert M.nEx == T.nEx
|
||||
assert M.nEy == T.nEy
|
||||
assert np.allclose(M.area, T.permuteF*T.area)
|
||||
assert np.allclose(M.edge, T.permuteE*T.edge)
|
||||
assert np.allclose(M.vol, T.permuteCC*T.vol)
|
||||
|
||||
# plt.subplot(211).spy(M.faceDiv)
|
||||
# plt.subplot(212).spy(T.permuteCC*T.faceDiv*T.permuteF.T)
|
||||
# plt.show()
|
||||
|
||||
assert (M.faceDiv - T.permuteCC*T.faceDiv*T.permuteF.T).nnz == 0
|
||||
|
||||
|
||||
class TestOcTree(unittest.TestCase):
|
||||
|
||||
def test_counts(self):
|
||||
nc = 8
|
||||
h1 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h2 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h3 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h = [hi/np.sum(hi) for hi in [h1, h2, h3]] # normalize
|
||||
M = Mesh.TreeMesh(h, levels=3)
|
||||
M._refineCell([0,0,0,0])
|
||||
M._refineCell([0,0,0,1])
|
||||
M.number()
|
||||
# M.plotGrid(showIt=True)
|
||||
# assert M.nhFx == 2
|
||||
# assert M.nFx == 9
|
||||
|
||||
assert np.allclose(M.vol.sum(), 1.0)
|
||||
|
||||
# assert np.allclose(M._areaFxFull, (M._deflationMatrix('F') * M.area)[:M.ntFx])
|
||||
# assert np.allclose(M._areaFyFull, (M._deflationMatrix('F') * M.area)[M.ntFx:(M.ntFx+M.ntFy)])
|
||||
# assert np.allclose(M._areaFzFull, (M._deflationMatrix('F') * M.area)[(M.ntFx+M.ntFy):])
|
||||
|
||||
# assert np.allclose(M._edgeExFull, (M._deflationMatrix('E') * M.edge)[:M.ntEx])
|
||||
# assert np.allclose(M._edgeEyFull, (M._deflationMatrix('E') * M.edge)[M.ntEx:(M.ntEx+M.ntEy)])
|
||||
# assert np.allclose(M._edgeEzFull, (M._deflationMatrix('E') * M.edge)[(M.ntEx+M.ntEy):])
|
||||
|
||||
def test_faceDiv(self):
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
M = Mesh.TreeMesh([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
assert M.nC == Mr.nC
|
||||
assert M.nF == Mr.nF
|
||||
assert M.nFx == Mr.nFx
|
||||
assert M.nFy == Mr.nFy
|
||||
assert M.nE == Mr.nE
|
||||
assert M.nEx == Mr.nEx
|
||||
assert M.nEy == Mr.nEy
|
||||
assert np.allclose(Mr.area, M.permuteF*M.area)
|
||||
assert np.allclose(Mr.edge, M.permuteE*M.edge)
|
||||
assert np.allclose(Mr.vol, M.permuteCC*M.vol)
|
||||
|
||||
# plt.subplot(211).spy(Mr.faceDiv)
|
||||
# plt.subplot(212).spy(M.permuteCC*M.faceDiv*M.permuteF.T)
|
||||
# plt.show()
|
||||
|
||||
assert (Mr.faceDiv - M.permuteCC*M.faceDiv*M.permuteF.T).nnz == 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)
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
M = Mesh.TreeMesh([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
# plt.subplot(211).spy(Mr.faceDiv)
|
||||
# plt.subplot(212).spy(M.permuteCC.T*M.faceDiv*M.permuteF)
|
||||
# plt.show()
|
||||
|
||||
assert (Mr.edgeCurl - M.permuteF*M.edgeCurl*M.permuteE.T).nnz == 0
|
||||
|
||||
def test_faceInnerProduct(self):
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
# hx, hy, hz = [[(1,4)], [(1,4)], [(1,4)]]
|
||||
|
||||
M = Mesh.TreeMesh([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
# plt.subplot(211).spy(Mr.getFaceInnerProduct())
|
||||
# plt.subplot(212).spy(M.getFaceInnerProduct())
|
||||
# plt.show()
|
||||
|
||||
# print M.nC, M.nF, M.getFaceInnerProduct().shape, M.permuteF.shape
|
||||
|
||||
assert np.allclose(Mr.getFaceInnerProduct().todense(), (M.permuteF * M.getFaceInnerProduct() * M.permuteF.T).todense())
|
||||
assert np.allclose(Mr.getEdgeInnerProduct().todense(), (M.permuteE * M.getEdgeInnerProduct() * M.permuteE.T).todense())
|
||||
|
||||
def test_VectorIdenties(self):
|
||||
hx, hy, hz = [[(1,4)], [(1,4)], [(1,4)]]
|
||||
|
||||
M = Mesh.TreeMesh([hx, hy, hz], levels=2)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
assert (M.faceDiv * M.edgeCurl).nnz == 0
|
||||
assert (Mr.faceDiv * Mr.edgeCurl).nnz == 0
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
|
||||
M = Mesh.TreeMesh([hx, hy, hz], levels=2)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
assert np.max(np.abs((M.faceDiv * M.edgeCurl).todense().flatten())) < TOL
|
||||
assert np.max(np.abs((Mr.faceDiv * Mr.edgeCurl).todense().flatten())) < TOL
|
||||
|
||||
|
||||
|
||||
if __name__ == '__main__':
|
||||
|
||||
@@ -1,160 +0,0 @@
|
||||
from SimPEG import Mesh
|
||||
from SimPEG.Mesh.PointerTree import Tree
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
import unittest
|
||||
|
||||
TOL = 1e-10
|
||||
|
||||
|
||||
|
||||
class TestSimpleQuadTree(unittest.TestCase):
|
||||
|
||||
def test_counts(self):
|
||||
nc = 8
|
||||
h1 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h2 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h = [hi/np.sum(hi) for hi in [h1, h2]] # normalize
|
||||
M = Tree(h)
|
||||
M._refineCell([0,0,0])
|
||||
M._refineCell([0,0,1])
|
||||
M.number()
|
||||
# M.plotGrid(showIt=True)
|
||||
assert M.nhFx == 2
|
||||
assert M.nFx == 9
|
||||
|
||||
assert np.allclose(M.vol.sum(), 1.0)
|
||||
|
||||
assert np.allclose(np.r_[M._areaFxFull, M._areaFyFull], M._deflationMatrix('F') * M.area)
|
||||
|
||||
|
||||
def test_faceDiv(self):
|
||||
|
||||
hx, hy = np.r_[1.,2,3,4], np.r_[5.,6,7,8]
|
||||
T = Tree([hx, hy], levels=2)
|
||||
T.refine(lambda xc:2)
|
||||
# T.plotGrid(showIt=True)
|
||||
M = Mesh.TensorMesh([hx, hy])
|
||||
assert M.nC == T.nC
|
||||
assert M.nF == T.nF
|
||||
assert M.nFx == T.nFx
|
||||
assert M.nFy == T.nFy
|
||||
assert M.nE == T.nE
|
||||
assert M.nEx == T.nEx
|
||||
assert M.nEy == T.nEy
|
||||
assert np.allclose(M.area, T.permuteF*T.area)
|
||||
assert np.allclose(M.edge, T.permuteE*T.edge)
|
||||
assert np.allclose(M.vol, T.permuteCC*T.vol)
|
||||
|
||||
# plt.subplot(211).spy(M.faceDiv)
|
||||
# plt.subplot(212).spy(T.permuteCC*T.faceDiv*T.permuteF.T)
|
||||
# plt.show()
|
||||
|
||||
assert (M.faceDiv - T.permuteCC*T.faceDiv*T.permuteF.T).nnz == 0
|
||||
|
||||
|
||||
class TestOcTree(unittest.TestCase):
|
||||
|
||||
def test_counts(self):
|
||||
nc = 8
|
||||
h1 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h2 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h3 = np.random.rand(nc)*nc*0.5 + nc*0.5
|
||||
h = [hi/np.sum(hi) for hi in [h1, h2, h3]] # normalize
|
||||
M = Tree(h, levels=3)
|
||||
M._refineCell([0,0,0,0])
|
||||
M._refineCell([0,0,0,1])
|
||||
M.number()
|
||||
# M.plotGrid(showIt=True)
|
||||
# assert M.nhFx == 2
|
||||
# assert M.nFx == 9
|
||||
|
||||
assert np.allclose(M.vol.sum(), 1.0)
|
||||
|
||||
# assert np.allclose(M._areaFxFull, (M._deflationMatrix('F') * M.area)[:M.ntFx])
|
||||
# assert np.allclose(M._areaFyFull, (M._deflationMatrix('F') * M.area)[M.ntFx:(M.ntFx+M.ntFy)])
|
||||
# assert np.allclose(M._areaFzFull, (M._deflationMatrix('F') * M.area)[(M.ntFx+M.ntFy):])
|
||||
|
||||
# assert np.allclose(M._edgeExFull, (M._deflationMatrix('E') * M.edge)[:M.ntEx])
|
||||
# assert np.allclose(M._edgeEyFull, (M._deflationMatrix('E') * M.edge)[M.ntEx:(M.ntEx+M.ntEy)])
|
||||
# assert np.allclose(M._edgeEzFull, (M._deflationMatrix('E') * M.edge)[(M.ntEx+M.ntEy):])
|
||||
|
||||
def test_faceDiv(self):
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
M = Tree([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
assert M.nC == Mr.nC
|
||||
assert M.nF == Mr.nF
|
||||
assert M.nFx == Mr.nFx
|
||||
assert M.nFy == Mr.nFy
|
||||
assert M.nE == Mr.nE
|
||||
assert M.nEx == Mr.nEx
|
||||
assert M.nEy == Mr.nEy
|
||||
assert np.allclose(Mr.area, M.permuteF*M.area)
|
||||
assert np.allclose(Mr.edge, M.permuteE*M.edge)
|
||||
assert np.allclose(Mr.vol, M.permuteCC*M.vol)
|
||||
|
||||
# plt.subplot(211).spy(Mr.faceDiv)
|
||||
# plt.subplot(212).spy(M.permuteCC*M.faceDiv*M.permuteF.T)
|
||||
# plt.show()
|
||||
|
||||
assert (Mr.faceDiv - M.permuteCC*M.faceDiv*M.permuteF.T).nnz == 0
|
||||
|
||||
|
||||
def test_edgeCurl(self):
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
M = Tree([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
# plt.subplot(211).spy(Mr.faceDiv)
|
||||
# plt.subplot(212).spy(M.permuteCC.T*M.faceDiv*M.permuteF)
|
||||
# plt.show()
|
||||
|
||||
assert (Mr.edgeCurl - M.permuteF*M.edgeCurl*M.permuteE.T).nnz == 0
|
||||
|
||||
def test_faceInnerProduct(self):
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
# hx, hy, hz = [[(1,4)], [(1,4)], [(1,4)]]
|
||||
|
||||
M = Tree([hx, hy, hz], levels=2)
|
||||
M.refine(lambda xc:2)
|
||||
# M.plotGrid(showIt=True)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
# plt.subplot(211).spy(Mr.getFaceInnerProduct())
|
||||
# plt.subplot(212).spy(M.getFaceInnerProduct())
|
||||
# plt.show()
|
||||
|
||||
# print M.nC, M.nF, M.getFaceInnerProduct().shape, M.permuteF.shape
|
||||
|
||||
assert np.allclose(Mr.getFaceInnerProduct().todense(), (M.permuteF * M.getFaceInnerProduct() * M.permuteF.T).todense())
|
||||
assert np.allclose(Mr.getEdgeInnerProduct().todense(), (M.permuteE * M.getEdgeInnerProduct() * M.permuteE.T).todense())
|
||||
|
||||
def test_VectorIdenties(self):
|
||||
hx, hy, hz = [[(1,4)], [(1,4)], [(1,4)]]
|
||||
|
||||
M = Tree([hx, hy, hz], levels=2)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
assert (M.faceDiv * M.edgeCurl).nnz == 0
|
||||
assert (Mr.faceDiv * Mr.edgeCurl).nnz == 0
|
||||
|
||||
hx, hy, hz = np.r_[1.,2,3,4], np.r_[5.,6,7,8], np.r_[9.,10,11,12]
|
||||
|
||||
M = Tree([hx, hy, hz], levels=2)
|
||||
Mr = Mesh.TensorMesh([hx, hy, hz])
|
||||
|
||||
assert np.max(np.abs((M.faceDiv * M.edgeCurl).todense().flatten())) < TOL
|
||||
assert np.max(np.abs((Mr.faceDiv * Mr.edgeCurl).todense().flatten())) < TOL
|
||||
|
||||
|
||||
|
||||
if __name__ == '__main__':
|
||||
unittest.main()
|
||||
Reference in New Issue
Block a user