Edge inner products are derivatives w/ testing.

SimPEG/Mesh/InnerProducts.py

@@ -77,7 +77,6 @@ class InnerProducts(object):
             :rtype: scipy.csr_matrix
             :return: M, the inner product matrix (nF, nF)
         """
-
         fast = None
 
         if hasattr(self, '_fastFaceInnerProductDeriv'):
@@ -97,8 +96,6 @@ class InnerProducts(object):
             P[0].T * sp.identity(n) * P[0]
 
 
-
-
     def getEdgeInnerProduct(self, materialProperty=None, returnP=False, invertProperty=False):
         """
             :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
@@ -107,6 +104,14 @@ class InnerProducts(object):
             :rtype: scipy.csr_matrix
             :return: M, the inner product matrix (nE, nE)
         """
+        fast = None
+
+        if returnP is False and hasattr(self, '_fastEdgeInnerProduct'):
+            fast = self._fastEdgeInnerProduct(materialProperty=materialProperty, invertProperty=invertProperty)
+
+        if fast is not None:
+            return fast
+
         if invertProperty:
             materialProperty = invPropertyTensor(self, materialProperty)
 
@@ -146,6 +151,25 @@ class InnerProducts(object):
         else:
             return A
 
+
+    def getEdgeInnerProductDeriv(self, materialProperty=None, P=None):
+        """
+            :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
+            :rtype: scipy.csr_matrix
+            :return: M, the inner product matrix (nF, nF)
+        """
+
+        fast = None
+
+        if hasattr(self, '_fastEdgeInnerProductDeriv'):
+            fast = self._fastEdgeInnerProductDeriv(materialProperty=materialProperty)
+
+        if fast is not None:
+            return fast
+
+        raise NotImplementedError('Derivatives for the material property specified are not yet implemented.')
+
+
 # ------------------------ Geometries ------------------------------
 #
 #

SimPEG/Mesh/TensorMesh.py

@@ -503,35 +503,82 @@ class TensorMesh(BaseRectangularMesh, TensorView, DiffOperators, InnerProducts):
             :rtype: scipy.csr_matrix
             :return: M, the inner product matrix (nF, nF)
         """
+        return self._fastInnerProduct('F', materialProperty=materialProperty, invertProperty=invertProperty)
+
+
+    def _fastEdgeInnerProduct(self, materialProperty=None, invertProperty=False):
+        """
+            Fast version of getEdgeInnerProduct.
+            This does not handle the case of a full tensor materialProperty.
+
+            :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
+            :param bool returnP: returns the projection matrices
+            :param bool invertProperty: inverts the material property
+            :rtype: scipy.csr_matrix
+            :return: M, the inner product matrix (nE, nE)
+        """
+        return self._fastInnerProduct('E', materialProperty=materialProperty, invertProperty=invertProperty)
+
+
+    def _fastInnerProduct(self, AvType, materialProperty=None, invertProperty=False):
+        """
+            Fast version of getFaceInnerProduct.
+            This does not handle the case of a full tensor materialProperty.
+
+            :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
+            :param str AvType: 'E' or 'F'
+            :param bool returnP: returns the projection matrices
+            :param bool invertProperty: inverts the material property
+            :rtype: scipy.csr_matrix
+            :return: M, the inner product matrix (nF, nF)
+        """
         if materialProperty is None:
             materialProperty = np.ones(self.nC)
 
         if materialProperty.size == self.nC:
             if invertProperty:
                 materialProperty = 1./materialProperty
-            Av = self.aveF2CC
+            Av = getattr(self, 'ave'+AvType+'2CC')
             V = Utils.sdiag(self.vol)
             return self.dim * Utils.sdiag(Av.T * V * materialProperty)
         if materialProperty.size == self.nC*self.dim:
             if invertProperty:
                 materialProperty = 1./materialProperty
-            Av = self.aveF2CCV
+            Av = getattr(self, 'ave'+AvType+'2CCV')
             V = sp.kron(sp.identity(self.dim), Utils.sdiag(self.vol))
             return Utils.sdiag(Av.T * V * Utils.mkvc(materialProperty))
 
 
-
     def _fastFaceInnerProductDeriv(self, materialProperty=None):
         """
             :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
             :rtype: scipy.csr_matrix
             :return: M, the inner product matrix (nF, nF)
         """
+        return self._fastInnerProductDeriv('F', materialProperty=materialProperty)
+
+
+    def _fastEdgeInnerProductDeriv(self, materialProperty=None):
+        """
+            :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
+            :rtype: scipy.csr_matrix
+            :return: M, the inner product matrix (nE, nE)
+        """
+        return self._fastInnerProductDeriv('E', materialProperty=materialProperty)
+
+
+    def _fastInnerProductDeriv(self, AvType, materialProperty=None):
+        """
+            :param str AvType: 'E' or 'F'
+            :param numpy.array materialProperty: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
+            :rtype: scipy.csr_matrix
+            :return: M, the inner product matrix (nF, nF)
+        """
         if materialProperty is None or materialProperty.size == self.nC:
-            Av = self.aveF2CC
+            Av = getattr(self, 'ave'+AvType+'2CC')
             return self.dim * Av.T * Utils.sdiag(self.vol)
         if materialProperty.size == self.nC*self.dim: # anisotropic
-            Av = self.aveF2CCV
+            Av = getattr(self, 'ave'+AvType+'2CCV')
             V = sp.kron(sp.identity(self.dim), Utils.sdiag(self.vol))
             return Av.T * V
 

SimPEG/Tests/test_massMatricesDerivs.py

@@ -11,7 +11,6 @@ class TestInnerProductsDerivs(unittest.TestCase):
     def test_FaceIP_derivs_isotropic(self):
         for d in range(3):
             mesh = Mesh.TensorMesh([10,5,4][d:])
-            M,Ps = mesh.getFaceInnerProduct(returnP=True)
             v = np.random.rand(mesh.nF)
             def fun(sig):
                 M = mesh.getFaceInnerProduct(sig)
@@ -21,10 +20,22 @@ class TestInnerProductsDerivs(unittest.TestCase):
             passed = checkDerivative(fun, sig, plotIt=False)
             self.assertTrue(passed)
 
+
+    def test_EdgeIP_derivs_isotropic(self):
+        for h in [[10,5],[10,5,4]]:
+            mesh = Mesh.TensorMesh(h)
+            v = np.random.rand(mesh.nE)
+            def fun(sig):
+                M = mesh.getEdgeInnerProduct(sig)
+                Md = mesh.getEdgeInnerProductDeriv(sig)
+                return M*v, Utils.sdiag(v)*Md
+            sig = np.random.rand(mesh.nC)
+            passed = checkDerivative(fun, sig, plotIt=False)
+            self.assertTrue(passed)
+
     def test_FaceIP_derivs_anisotropic(self):
         for d in range(3):
             mesh = Mesh.TensorMesh([10,5,4][d:])
-            M,Ps = mesh.getFaceInnerProduct(returnP=True)
             v = np.random.rand(mesh.nF)
             def fun(sig):
                 M = mesh.getFaceInnerProduct(sig)