From 9c9944c1a2621bb110b62aedf54b3333b5c831d2 Mon Sep 17 00:00:00 2001 From: rowanc1 Date: Sat, 12 Jul 2014 15:24:11 -0500 Subject: [PATCH 1/5] change defaults in directives --- SimPEG/Directives.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/SimPEG/Directives.py b/SimPEG/Directives.py index 226fb860..6a492a50 100644 --- a/SimPEG/Directives.py +++ b/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): From cf5be08d96ff92dedb68786fc75ab9dba50f20aa Mon Sep 17 00:00:00 2001 From: rowanc1 Date: Sat, 12 Jul 2014 15:24:29 -0500 Subject: [PATCH 2/5] add a 2d to 3d map and test it! --- SimPEG/Maps.py | 59 ++++++++++++++++++++++++++++++++++++++- SimPEG/Tests/test_maps.py | 41 ++++++++++++++++++++++++++- 2 files changed, 98 insertions(+), 2 deletions(-) diff --git a/SimPEG/Maps.py b/SimPEG/Maps.py index b6af6f13..706e83f2 100644 --- a/SimPEG/Maps.py +++ b/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. diff --git a/SimPEG/Tests/test_maps.py b/SimPEG/Tests/test_maps.py index 986fa500..6e595732 100644 --- a/SimPEG/Tests/test_maps.py +++ b/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() From c570b294e47985923cb94e00e406446cdf2eaa9b Mon Sep 17 00:00:00 2001 From: rowanc1 Date: Sat, 12 Jul 2014 18:33:22 -0500 Subject: [PATCH 3/5] update documentation for the maps --- docs/api_Maps.rst | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/docs/api_Maps.rst b/docs/api_Maps.rst index 606f1282..cf33b7ad 100644 --- a/docs/api_Maps.rst +++ b/docs/api_Maps.rst @@ -154,6 +154,14 @@ Vertical 1D Map :undoc-members: +Map 2D Cross-Section to 3D Model +-------------------------------- + +.. autoclass:: SimPEG.Maps.Map2Dto3D + :members: + :undoc-members: + + Mesh to Mesh Map ---------------- @@ -196,7 +204,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. From 5a667090222121a7c9e1ef3075318b9be0018086 Mon Sep 17 00:00:00 2001 From: rowanc1 Date: Sat, 12 Jul 2014 18:38:53 -0500 Subject: [PATCH 4/5] updates to maps docs --- docs/api_Maps.rst | 12 ++---------- 1 file changed, 2 insertions(+), 10 deletions(-) diff --git a/docs/api_Maps.rst b/docs/api_Maps.rst index cf33b7ad..0ab32024 100644 --- a/docs/api_Maps.rst +++ b/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 ======= @@ -178,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)') @@ -212,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: From 619c2a2caef395bca2ca7e60b8fd8ce5ec6d6148 Mon Sep 17 00:00:00 2001 From: Rowan Cockett Date: Tue, 19 Aug 2014 19:58:35 -0700 Subject: [PATCH 5/5] Update README.rst --- README.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/README.rst b/README.rst index b99b6eaf..9878a7b5 100644 --- a/README.rst +++ b/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: