diff --git a/SimPEG/PropMaps.py b/SimPEG/PropMaps.py index d7bd2ba2..995216f7 100644 --- a/SimPEG/PropMaps.py +++ b/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): @@ -119,6 +132,9 @@ class PropModel(object): self._nP = len(set(inds)) return self._nP + def __contains__(self, val): + return val in self.propMap + _PROPMAPCLASSREGISTRY = {} @@ -167,6 +183,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() @@ -238,3 +255,6 @@ class PropMap(object): def __call__(self, vec): return self.PropModel(self, vec) + def __contains__(self, val): + activeMaps = [name for name in self._properties if getattr(self, '%sMap'%name) is not None] + return val in activeMaps diff --git a/SimPEG/Tests/test_PropMaps.py b/SimPEG/Tests/test_PropMaps.py index c2726d71..b4c0eb8d 100644 --- a/SimPEG/Tests/test_PropMaps.py +++ b/SimPEG/Tests/test_PropMaps.py @@ -36,7 +36,16 @@ class TestPropMaps(unittest.TestCase): assert PM.muMap is None assert PM.muIndex is None + assert 'sigma' in PM + assert 'mu' not in PM + assert 'mui' not in PM + m = PM(np.r_[1.,2,3]) + + assert 'sigma' in m + assert 'mu' not in m + assert 'mui' not in m + assert m.mu == mu_0 assert m.muModel is None assert m.muMap is None @@ -67,6 +76,14 @@ class TestPropMaps(unittest.TestCase): assert pm.nP == 6 + assert 'sigma' in PM + assert 'mu' in PM + assert 'mui' not in PM + + assert 'sigma' in pm + assert 'mu' in pm + assert 'mui' not in pm + assert np.all(pm.sigmaModel == [1.,2,3]) assert np.all(pm.sigma == np.exp([1.,2,3])) assert np.all(pm.muModel == [4.,5,6]) @@ -83,11 +100,33 @@ class TestPropMaps(unittest.TestCase): assert pm.nP == 3 + assert 'sigma' in PM + assert 'mu' in PM + assert 'mui' not in PM + + assert 'sigma' in pm + assert 'mu' in pm + assert 'mui' not in pm + assert np.all(pm.sigmaModel == [1,2,3]) assert np.all(pm.sigma == np.exp([1,2,3])) 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) @@ -103,6 +142,17 @@ class TestPropMaps(unittest.TestCase): assert pm.sigmaDeriv is not None assert pm.rhoDeriv is not None + assert 'sigma' in PM + assert 'rho' not in PM + assert 'mu' not in PM + assert 'mui' not in PM + + + assert 'sigma' in pm + assert 'rho' not in pm + assert 'mu' not in pm + assert 'mui' not in pm + assert pm.mu == mu_0 assert pm.mui == 1.0/mu_0 assert pm.muMap is None @@ -121,6 +171,17 @@ class TestPropMaps(unittest.TestCase): assert pm.sigmaDeriv is not None assert pm.rhoDeriv is not None + assert 'sigma' not in PM + assert 'rho' in PM + assert 'mu' not in PM + assert 'mui' not in PM + + + assert 'sigma' not in pm + assert 'rho' in pm + assert 'mu' not in pm + assert 'mui' not in pm + self.assertRaises(AssertionError, MyReciprocalPropMap, [('rho', iMap), ('sigma', iMap)]) self.assertRaises(AssertionError, MyReciprocalPropMap, [('sigma', iMap), ('rho', iMap)])