updates to adding projections to the models

SimPEG/PropMaps.py

@@ -1,4 +1,4 @@
-import Utils, Maps, numpy as np
+import Utils, Maps, numpy as np, scipy.sparse as sp
 
 class Property(object):
 
@@ -87,10 +87,23 @@ class Property(object):
             mapping = getattr(self, '%sMap'%prop.name)
             if mapping is None:
                 return None
-            index = getattr(self.propMap, '_%sIndex'%prop.name, slice(None))
+            index = getattr(self.propMap, '%sIndex'%prop.name)
             return self.vector[index]
         return property(fget=fget)
 
+    def _getModelProjProperty(self):
+        prop = self
+        def fget(self):
+            mapping = getattr(self, '%sMap'%prop.name)
+            if mapping is None:
+                return None
+            inds = getattr(self.propMap, '%sIndex'%prop.name)
+            if type(inds) is slice:
+                inds = range(*inds.indices(self.nP))
+            nI, nP = len(inds),self.nP
+            return sp.csr_matrix((np.ones(nI), (range(nI), inds) ), shape=(nI, nP))
+        return property(fget=fget)
+
     def _getModelMapProperty(self):
         prop = self
         def fget(self):
@@ -167,6 +180,7 @@ class _PropMapMetaClass(type):
 
             attrs[attr          ] = prop._getProperty()
             attrs[attr + 'Map'  ] = prop._getModelMapProperty()
+            attrs[attr + 'Proj' ] = prop._getModelProjProperty()
             attrs[attr + 'Model'] = prop._getModelProperty()
             attrs[attr + 'Deriv'] = prop._getModelDerivProperty()
 

SimPEG/Tests/test_PropMaps.py

@@ -100,6 +100,20 @@ class TestPropMaps(unittest.TestCase):
         assert np.all(pm.muModel == [1,2,3])
         assert np.all(pm.mu == [1,2,3])
 
+    def test_Projections(self):
+        m = Mesh.TensorMesh((3,))
+        iMap = Maps.IdentityMap(m)
+        PM = MyReciprocalPropMap([('sigma', iMap)])
+        v = np.r_[1,2.,3]
+        pm = PM(v)
+
+        assert pm.sigmaProj is not None
+        assert pm.rhoProj   is None
+        assert pm.muProj    is None
+        assert pm.muiProj   is None
+
+        assert np.all(pm.sigmaProj * v == pm.sigmaModel)
+
     def test_Links(self):
         m = Mesh.TensorMesh((3,))
         expMap = Maps.ExpMap(m)