mirror of
https://github.com/wassname/simpeg.git
synced 2026-06-27 19:48:52 +08:00
More tests on more meshes, IPDeriv returns a function.
This commit is contained in:
rowanc1
@@ -114,8 +114,14 @@ class InnerProducts(object):
|
||||
:param TensorType tensorType: type of the tensor: TensorType(mesh, sigma)
|
||||
:param list P: list of projection matrices
|
||||
:param bool doFast: do a faster implementation if available.
|
||||
:rtype: function
|
||||
:return: dMdmu(u), the derivative of the inner product matrix (u)
|
||||
|
||||
Given u, dMdmu returns (nF, nC*nA)
|
||||
|
||||
:param np.ndarray u: vector that multiplies dMdmu
|
||||
:rtype: scipy.csr_matrix
|
||||
:return: dMdm, the derivative of the inner product matrix (nF, nC*nA)
|
||||
:return: dMdmu, the derivative of the inner product matrix for a certain u
|
||||
"""
|
||||
assert isinstance(tensorType, TensorType), 'tensorType must be an instance of TensorType.'
|
||||
fast = None
|
||||
@@ -131,7 +137,7 @@ class InnerProducts(object):
|
||||
|
||||
def innerProductDeriv(v):
|
||||
return self._getInnerProductDeriv(tensorType, P, self.nF, v)
|
||||
return DerivOperator(innerProductDeriv)
|
||||
return innerProductDeriv
|
||||
|
||||
|
||||
def getEdgeInnerProductDeriv(self, tensorType, P=None, doFast=True):
|
||||
@@ -155,7 +161,7 @@ class InnerProducts(object):
|
||||
P = self._getInnerProductProjectionMatrices('E', tensorType=tensorType)
|
||||
def innerProductDeriv(v):
|
||||
return self._getInnerProductDeriv(tensorType, P, self.nE, v)
|
||||
return DerivOperator(innerProductDeriv)
|
||||
return innerProductDeriv
|
||||
|
||||
def _getInnerProductDeriv(self, tensorType, P, n, v):
|
||||
"""
|
||||
|
||||
@@ -319,21 +319,21 @@ class BaseTensorMesh(BaseRectangularMesh):
|
||||
# return self.dim * Av.T * V * ones
|
||||
def scalarInnerProductDeriv(v):
|
||||
return Utils.sdiag(v) * self.dim * Av.T * V * ones
|
||||
return Utils.DerivOperator(scalarInnerProductDeriv)
|
||||
return scalarInnerProductDeriv
|
||||
|
||||
if tensorType == 1:
|
||||
Av = getattr(self, 'ave'+projType+'2CC')
|
||||
V = Utils.sdiag(self.vol)
|
||||
def isotropicInnerProductDeriv(v):
|
||||
return Utils.sdiag(v) * self.dim * Av.T * V
|
||||
return Utils.DerivOperator(isotropicInnerProductDeriv)
|
||||
return isotropicInnerProductDeriv
|
||||
|
||||
if tensorType == 2: # anisotropic
|
||||
Av = getattr(self, 'ave'+projType+'2CCV')
|
||||
V = sp.kron(sp.identity(self.dim), Utils.sdiag(self.vol))
|
||||
def anisotropicInnerProductDeriv(v):
|
||||
return Utils.sdiag(v) * Av.T * V
|
||||
return Utils.DerivOperator(anisotropicInnerProductDeriv)
|
||||
return anisotropicInnerProductDeriv
|
||||
|
||||
|
||||
|
||||
|
||||
@@ -6,102 +6,240 @@ from TestUtils import checkDerivative
|
||||
|
||||
class TestInnerProductsDerivs(unittest.TestCase):
|
||||
|
||||
def doTestFace(self, h, rep, fast):
|
||||
mesh = Mesh.TensorMesh(h)
|
||||
def doTestFace(self, h, rep, fast, meshType):
|
||||
if meshType == 'LRM':
|
||||
hRect = Utils.exampleLrmGrid(h,'rotate')
|
||||
mesh = Mesh.LogicallyRectMesh(hRect)
|
||||
elif meshType == 'Tree':
|
||||
mesh = Mesh.TreeMesh(h)
|
||||
elif meshType == 'Tensor':
|
||||
mesh = Mesh.TensorMesh(h)
|
||||
v = np.random.rand(mesh.nF)
|
||||
sig = np.random.rand(1) if rep is 0 else np.random.rand(mesh.nC*rep)
|
||||
Md = mesh.getFaceInnerProductDeriv(Utils.TensorType(mesh, sig), doFast=fast)
|
||||
def fun(sig):
|
||||
M = mesh.getFaceInnerProduct(sig)
|
||||
return M*v, Md*v
|
||||
return M*v, Md(v)
|
||||
print meshType, 'Face', h, rep, fast
|
||||
return checkDerivative(fun, sig, num=5, plotIt=False)
|
||||
|
||||
def doTestEdge(self, h, rep, fast):
|
||||
mesh = Mesh.TensorMesh(h)
|
||||
def doTestEdge(self, h, rep, fast, meshType):
|
||||
if meshType == 'LRM':
|
||||
hRect = Utils.exampleLrmGrid(h,'rotate')
|
||||
mesh = Mesh.LogicallyRectMesh(hRect)
|
||||
elif meshType == 'Tree':
|
||||
mesh = Mesh.TreeMesh(h)
|
||||
elif meshType == 'Tensor':
|
||||
mesh = Mesh.TensorMesh(h)
|
||||
v = np.random.rand(mesh.nE)
|
||||
sig = np.random.rand(1) if rep is 0 else np.random.rand(mesh.nC*rep)
|
||||
Md = mesh.getEdgeInnerProductDeriv(Utils.TensorType(mesh, sig), doFast=fast)
|
||||
def fun(sig):
|
||||
M = mesh.getEdgeInnerProduct(sig)
|
||||
return M*v, Md*v
|
||||
return M*v, Md(v)
|
||||
print meshType, 'Edge', h, rep, fast
|
||||
return checkDerivative(fun, sig, num=5, plotIt=False)
|
||||
|
||||
def test_FaceIP_1D_float(self):
|
||||
self.assertTrue(self.doTestFace([10],0, False))
|
||||
self.assertTrue(self.doTestFace([10],0, False, 'Tensor'))
|
||||
def test_FaceIP_2D_float(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, False))
|
||||
self.assertTrue(self.doTestFace([10, 4],0, False, 'Tensor'))
|
||||
def test_FaceIP_3D_float(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, False))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, False, 'Tensor'))
|
||||
def test_FaceIP_1D_isotropic(self):
|
||||
self.assertTrue(self.doTestFace([10],1, False))
|
||||
self.assertTrue(self.doTestFace([10],1, False, 'Tensor'))
|
||||
def test_FaceIP_2D_isotropic(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, False))
|
||||
self.assertTrue(self.doTestFace([10, 4],1, False, 'Tensor'))
|
||||
def test_FaceIP_3D_isotropic(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, False))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, False, 'Tensor'))
|
||||
def test_FaceIP_2D_anisotropic(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, False))
|
||||
self.assertTrue(self.doTestFace([10, 4],2, False, 'Tensor'))
|
||||
def test_FaceIP_3D_anisotropic(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, False))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, False, 'Tensor'))
|
||||
def test_FaceIP_2D_tensor(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],3, False))
|
||||
self.assertTrue(self.doTestFace([10, 4],3, False, 'Tensor'))
|
||||
def test_FaceIP_3D_tensor(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],6, False))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],6, False, 'Tensor'))
|
||||
|
||||
def test_FaceIP_1D_float_fast(self):
|
||||
self.assertTrue(self.doTestFace([10],0, True))
|
||||
self.assertTrue(self.doTestFace([10],0, True, 'Tensor'))
|
||||
def test_FaceIP_2D_float_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, True))
|
||||
self.assertTrue(self.doTestFace([10, 4],0, True, 'Tensor'))
|
||||
def test_FaceIP_3D_float_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, True))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, True, 'Tensor'))
|
||||
def test_FaceIP_1D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestFace([10],1, True))
|
||||
self.assertTrue(self.doTestFace([10],1, True, 'Tensor'))
|
||||
def test_FaceIP_2D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, True))
|
||||
self.assertTrue(self.doTestFace([10, 4],1, True, 'Tensor'))
|
||||
def test_FaceIP_3D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, True))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, True, 'Tensor'))
|
||||
def test_FaceIP_2D_anisotropic_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, True))
|
||||
self.assertTrue(self.doTestFace([10, 4],2, True, 'Tensor'))
|
||||
def test_FaceIP_3D_anisotropic_fast(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, True))
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, True, 'Tensor'))
|
||||
|
||||
def test_EdgeIP_1D_float(self):
|
||||
self.assertTrue(self.doTestEdge([10],0, False))
|
||||
self.assertTrue(self.doTestEdge([10],0, False, 'Tensor'))
|
||||
def test_EdgeIP_2D_float(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, False, 'Tensor'))
|
||||
def test_EdgeIP_3D_float(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, False, 'Tensor'))
|
||||
def test_EdgeIP_1D_isotropic(self):
|
||||
self.assertTrue(self.doTestEdge([10],1, False))
|
||||
self.assertTrue(self.doTestEdge([10],1, False, 'Tensor'))
|
||||
def test_EdgeIP_2D_isotropic(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, False, 'Tensor'))
|
||||
def test_EdgeIP_3D_isotropic(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, False, 'Tensor'))
|
||||
def test_EdgeIP_2D_anisotropic(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, False, 'Tensor'))
|
||||
def test_EdgeIP_3D_anisotropic(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, False, 'Tensor'))
|
||||
def test_EdgeIP_2D_tensor(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],3, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4],3, False, 'Tensor'))
|
||||
def test_EdgeIP_3D_tensor(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],6, False))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],6, False, 'Tensor'))
|
||||
|
||||
def test_EdgeIP_1D_float_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10],0, True))
|
||||
self.assertTrue(self.doTestEdge([10],0, True, 'Tensor'))
|
||||
def test_EdgeIP_2D_float_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, True, 'Tensor'))
|
||||
def test_EdgeIP_3D_float_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, True, 'Tensor'))
|
||||
def test_EdgeIP_1D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10],1, True))
|
||||
self.assertTrue(self.doTestEdge([10],1, True, 'Tensor'))
|
||||
def test_EdgeIP_2D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, True, 'Tensor'))
|
||||
def test_EdgeIP_3D_isotropic_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, True, 'Tensor'))
|
||||
def test_EdgeIP_2D_anisotropic_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, True, 'Tensor'))
|
||||
def test_EdgeIP_3D_anisotropic_fast(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, True))
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, True, 'Tensor'))
|
||||
|
||||
|
||||
|
||||
def test_FaceIP_2D_float_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, False, 'LRM'))
|
||||
def test_FaceIP_3D_float_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, False, 'LRM'))
|
||||
def test_FaceIP_2D_isotropic_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, False, 'LRM'))
|
||||
def test_FaceIP_3D_isotropic_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, False, 'LRM'))
|
||||
def test_FaceIP_2D_anisotropic_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, False, 'LRM'))
|
||||
def test_FaceIP_3D_anisotropic_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, False, 'LRM'))
|
||||
def test_FaceIP_2D_tensor_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],3, False, 'LRM'))
|
||||
def test_FaceIP_3D_tensor_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],6, False, 'LRM'))
|
||||
|
||||
def test_FaceIP_2D_float_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, True, 'LRM'))
|
||||
def test_FaceIP_3D_float_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, True, 'LRM'))
|
||||
def test_FaceIP_2D_isotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, True, 'LRM'))
|
||||
def test_FaceIP_3D_isotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, True, 'LRM'))
|
||||
def test_FaceIP_2D_anisotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, True, 'LRM'))
|
||||
def test_FaceIP_3D_anisotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, True, 'LRM'))
|
||||
|
||||
def test_EdgeIP_2D_float_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, False, 'LRM'))
|
||||
def test_EdgeIP_3D_float_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, False, 'LRM'))
|
||||
def test_EdgeIP_2D_isotropic_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, False, 'LRM'))
|
||||
def test_EdgeIP_3D_isotropic_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, False, 'LRM'))
|
||||
def test_EdgeIP_2D_anisotropic_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, False, 'LRM'))
|
||||
def test_EdgeIP_3D_anisotropic_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, False, 'LRM'))
|
||||
def test_EdgeIP_2D_tensor_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],3, False, 'LRM'))
|
||||
def test_EdgeIP_3D_tensor_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],6, False, 'LRM'))
|
||||
|
||||
def test_EdgeIP_2D_float_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, True, 'LRM'))
|
||||
def test_EdgeIP_3D_float_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, True, 'LRM'))
|
||||
def test_EdgeIP_2D_isotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, True, 'LRM'))
|
||||
def test_EdgeIP_3D_isotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, True, 'LRM'))
|
||||
def test_EdgeIP_2D_anisotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, True, 'LRM'))
|
||||
def test_EdgeIP_3D_anisotropic_fast_LRM(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, True, 'LRM'))
|
||||
|
||||
|
||||
|
||||
|
||||
def test_FaceIP_2D_float_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, False, 'Tree'))
|
||||
def test_FaceIP_3D_float_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, False, 'Tree'))
|
||||
def test_FaceIP_2D_isotropic_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, False, 'Tree'))
|
||||
def test_FaceIP_3D_isotropic_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, False, 'Tree'))
|
||||
def test_FaceIP_2D_anisotropic_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, False, 'Tree'))
|
||||
def test_FaceIP_3D_anisotropic_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, False, 'Tree'))
|
||||
def test_FaceIP_2D_tensor_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],3, False, 'Tree'))
|
||||
def test_FaceIP_3D_tensor_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],6, False, 'Tree'))
|
||||
|
||||
def test_FaceIP_2D_float_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],0, True, 'Tree'))
|
||||
def test_FaceIP_3D_float_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],0, True, 'Tree'))
|
||||
def test_FaceIP_2D_isotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],1, True, 'Tree'))
|
||||
def test_FaceIP_3D_isotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],1, True, 'Tree'))
|
||||
def test_FaceIP_2D_anisotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4],2, True, 'Tree'))
|
||||
def test_FaceIP_3D_anisotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestFace([10, 4, 5],3, True, 'Tree'))
|
||||
|
||||
def test_EdgeIP_2D_float_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, False, 'Tree'))
|
||||
def test_EdgeIP_3D_float_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, False, 'Tree'))
|
||||
def test_EdgeIP_2D_isotropic_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, False, 'Tree'))
|
||||
def test_EdgeIP_3D_isotropic_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, False, 'Tree'))
|
||||
def test_EdgeIP_2D_anisotropic_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, False, 'Tree'))
|
||||
def test_EdgeIP_3D_anisotropic_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, False, 'Tree'))
|
||||
def test_EdgeIP_2D_tensor_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],3, False, 'Tree'))
|
||||
def test_EdgeIP_3D_tensor_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],6, False, 'Tree'))
|
||||
|
||||
def test_EdgeIP_2D_float_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],0, True, 'Tree'))
|
||||
def test_EdgeIP_3D_float_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],0, True, 'Tree'))
|
||||
def test_EdgeIP_2D_isotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],1, True, 'Tree'))
|
||||
def test_EdgeIP_3D_isotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],1, True, 'Tree'))
|
||||
def test_EdgeIP_2D_anisotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4],2, True, 'Tree'))
|
||||
def test_EdgeIP_3D_anisotropic_fast_Tree(self):
|
||||
self.assertTrue(self.doTestEdge([10, 4, 5],3, True, 'Tree'))
|
||||
|
||||
if __name__ == '__main__':
|
||||
unittest.main()
|
||||
|
||||
@@ -349,10 +349,3 @@ class SimPEGLinearOperator(LinearOperator):
|
||||
@property
|
||||
def T(self):
|
||||
return self.__class__((self.shape[1],self.shape[0]),self.rmatvec,rmatvec=self.matvec,matmat=self.matmat)
|
||||
|
||||
|
||||
class DerivOperator(object):
|
||||
def __init__(self, f):
|
||||
self.f = f
|
||||
def __mul__(self, v):
|
||||
return self.f(v)
|
||||
|
||||
Reference in New Issue
Block a user