Merge branch 'master' of https://github.com/simpeg/simpeg

README.rst

@@ -37,7 +37,7 @@ The vision is to create a package for finite volume simulation with applications
 
 
 Documentation:
-http://simpeg.3ptscience.com
+http://simpeg.rtfd.org
 
 
 Code:

SimPEG/Directives.py

@@ -88,7 +88,7 @@ class BetaEstimate_ByEig(InversionDirective):
     """BetaEstimate"""
 
     beta0 = None       #: The initial Beta (regularization parameter)
-    beta0_ratio = 0.1  #: estimateBeta0 is used with this ratio
+    beta0_ratio = 1e2  #: estimateBeta0 is used with this ratio
 
     def initialize(self):
         """
@@ -136,7 +136,7 @@ class BetaEstimate_ByEig(InversionDirective):
 class BetaSchedule(InversionDirective):
     """BetaSchedule"""
 
-    coolingFactor = 2.
+    coolingFactor = 8.
     coolingRate = 3
 
     def endIter(self):

SimPEG/Maps.py

@@ -12,7 +12,8 @@ class IdentityMap(object):
 
     mesh = None      #: A SimPEG Mesh
 
-    def __init__(self, mesh):
+    def __init__(self, mesh, **kwargs):
+        Utils.setKwargs(self, **kwargs)
         self.mesh = mesh
 
     @property
@@ -261,6 +262,62 @@ class Vertical1DMap(IdentityMap):
                     ), shape=(repNum, 1))
         return sp.kron(sp.identity(self.nP), repVec)
 
+
+class Map2Dto3D(IdentityMap):
+    """Map2Dto3D
+
+        Given a 2D vector, this will extend to the full
+        3D model space.
+    """
+
+    normal = 'Y' #: The normal
+
+    def __init__(self, mesh, **kwargs):
+        assert mesh.dim == 3, 'Only works for a 3D Mesh'
+        IdentityMap.__init__(self, mesh, **kwargs)
+        assert self.normal in ['X','Y','Z'], 'For now, only "Y" normal is supported'
+
+    @property
+    def nP(self):
+        """Number of model properties.
+
+           The number of cells in the
+           last dimension of the mesh."""
+        if self.normal == 'Z':
+            return self.mesh.nCx * self.mesh.nCy
+        elif self.normal == 'Y':
+            return self.mesh.nCx * self.mesh.nCz
+        elif self.normal == 'X':
+            return self.mesh.nCy * self.mesh.nCz
+
+    def _transform(self, m):
+        """
+            :param numpy.array m: model
+            :rtype: numpy.array
+            :return: transformed model
+        """
+        m = Utils.mkvc(m)
+        if self.normal == 'Z':
+            return Utils.mkvc(m.reshape(self.mesh.vnC[[0,1]], order='F')[:,:,np.newaxis].repeat(self.mesh.nCz,axis=2))
+        elif self.normal == 'Y':
+            return Utils.mkvc(m.reshape(self.mesh.vnC[[0,2]], order='F')[:,np.newaxis,:].repeat(self.mesh.nCy,axis=1))
+        elif self.normal == 'X':
+            return Utils.mkvc(m.reshape(self.mesh.vnC[[1,2]], order='F')[np.newaxis,:,:].repeat(self.mesh.nCx,axis=0))
+
+    def deriv(self, m):
+        """
+            :param numpy.array m: model
+            :rtype: scipy.csr_matrix
+            :return: derivative of transformed model
+        """
+        inds = self * np.arange(self.nP)
+        nC, nP = self.mesh.nC, self.nP
+        P = sp.csr_matrix(
+                    (np.ones(nC),
+                    (range(nC), inds)
+                ), shape=(nC, nP))
+        return P
+
 class Mesh2Mesh(IdentityMap):
     """
         Takes a model on one mesh are translates it to another mesh.

SimPEG/Tests/test_maps.py

@@ -13,9 +13,10 @@ class MapTests(unittest.TestCase):
         a = np.array([1, 1, 1])
         b = np.array([1, 2])
         self.mesh2 = Mesh.TensorMesh([a, b], x0=np.array([3, 5]))
+        self.mesh3 = Mesh.TensorMesh([a, b, [3,4]], x0=np.array([3, 5, 2]))
         self.mesh22 = Mesh.TensorMesh([b, a], x0=np.array([3, 5]))
 
-    def test_transforms(self):
+    def test_transforms2D(self):
         for M in dir(Maps):
             try:
                 maps = getattr(Maps, M)(self.mesh2)
@@ -24,6 +25,15 @@ class MapTests(unittest.TestCase):
                 continue
             self.assertTrue(maps.test())
 
+    def test_transforms3D(self):
+        for M in dir(Maps):
+            try:
+                maps = getattr(Maps, M)(self.mesh3)
+                assert isinstance(maps, Maps.IdentityMap)
+            except Exception, e:
+                continue
+            self.assertTrue(maps.test())
+
     def test_Mesh2MeshMap(self):
         maps = Maps.Mesh2Mesh([self.mesh22, self.mesh2])
         self.assertTrue(maps.test())
@@ -90,5 +100,34 @@ class MapTests(unittest.TestCase):
         self.assertRaises(ValueError, lambda: actMap * vertMap * expMap )
 
 
+    def test_map2Dto3D_x(self):
+        M2 = Mesh.TensorMesh([2,4])
+        M3 = Mesh.TensorMesh([3,2,4])
+        m = np.random.rand(M2.nC)
+        m2to3 = Maps.Map2Dto3D(M3, normal='X')
+        m = np.arange(m2to3.nP)
+        self.assertTrue(m2to3.test())
+        self.assertTrue(np.all(Utils.mkvc( (m2to3 * m).reshape(M3.vnC,order='F')[0,:,:] ) == m))
+
+
+    def test_map2Dto3D_y(self):
+        M2 = Mesh.TensorMesh([3,4])
+        M3 = Mesh.TensorMesh([3,2,4])
+        m = np.random.rand(M2.nC)
+        m2to3 = Maps.Map2Dto3D(M3, normal='Y')
+        m = np.arange(m2to3.nP)
+        self.assertTrue(m2to3.test())
+        self.assertTrue(np.all(Utils.mkvc( (m2to3 * m).reshape(M3.vnC,order='F')[:,0,:] ) == m))
+
+    def test_map2Dto3D_z(self):
+        M2 = Mesh.TensorMesh([3,2])
+        M3 = Mesh.TensorMesh([3,2,4])
+        m = np.random.rand(M2.nC)
+        m2to3 = Maps.Map2Dto3D(M3, normal='Z')
+        m = np.arange(m2to3.nP)
+        self.assertTrue(m2to3.test())
+        self.assertTrue(np.all(Utils.mkvc( (m2to3 * m).reshape(M3.vnC,order='F')[:,:,0] ) == m))
+
+
 if __name__ == '__main__':
     unittest.main()

docs/api_Maps.rst

@@ -119,7 +119,6 @@ When these are used in the inverse problem, this is extremely important!!
     expMap.test(m, plotIt=True)
 
 
-
 The API
 =======
 
@@ -154,6 +153,14 @@ Vertical 1D Map
     :undoc-members:
 
 
+Map 2D Cross-Section to 3D Model
+--------------------------------
+
+.. autoclass:: SimPEG.Maps.Map2Dto3D
+    :members:
+    :undoc-members:
+
+
 Mesh to Mesh Map
 ----------------
 
@@ -170,8 +177,8 @@ Mesh to Mesh Map
     v = Utils.mkvc(V)
     modh = Maps.Mesh2Mesh([M,M2])
     modH = Maps.Mesh2Mesh([M2,M])
-    H = modH.transform(v)
-    h = modh.transform(H)
+    H = modH * v
+    h = modh * H
     ax = plt.subplot(131)
     M.plotImage(v, ax=ax)
     ax.set_title('Fine Mesh (Original)')
@@ -196,7 +203,7 @@ Combo Map
 ---------
 
 The ComboMap holds the information for multiplying and combining
-maps. It also uses the chain rule create the derivative.
+maps. It also uses the chain rule to create the derivative.
 Remember, any time that you make your own combination of mappings
 be sure to test that the derivative is correct.
 
@@ -204,10 +211,3 @@ be sure to test that the derivative is correct.
     :members:
     :undoc-members:
 
-
-Non Linear Map
---------------
-
-.. autoclass:: SimPEG.Maps.NonLinearMap
-    :members:
-    :undoc-members: