From b8bd011662e31d6e773f8b49306354add908816c Mon Sep 17 00:00:00 2001 From: GudniRos Date: Thu, 4 Feb 2016 02:28:22 -0800 Subject: [PATCH] Added a 3D forward example. Missing plots. --- SimPEG/Examples/MT_3D_Foward.py | 64 +++++++++++++++++++++ SimPEG/MT/Examples/simple3DfowardProblem.py | 41 ------------- 2 files changed, 64 insertions(+), 41 deletions(-) create mode 100644 SimPEG/Examples/MT_3D_Foward.py delete mode 100644 SimPEG/MT/Examples/simple3DfowardProblem.py diff --git a/SimPEG/Examples/MT_3D_Foward.py b/SimPEG/Examples/MT_3D_Foward.py new file mode 100644 index 00000000..fea5c95d --- /dev/null +++ b/SimPEG/Examples/MT_3D_Foward.py @@ -0,0 +1,64 @@ +# Test script to use SimPEG.MT platform to forward model synthetic data. + +# Import +import SimPEG as simpeg +from SimPEG import MT +import numpy as np +try: + from pymatsolver import MumpsSolver as Solver +except: + from SimPEG import Solver + +def run(plotIt=True): + """ + MT: 3D: Forward + ======================= + + Forward model 3D MT data. + + """ + + # Make a mesh + M = simpeg.Mesh.TensorMesh([[(100,5,-1.5),(100.,10),(100,5,1.5)],[(100,5,-1.5),(100.,10),(100,5,1.5)],[(100,5,1.6),(100.,10),(100,3,2)]], x0=['C','C',-3529.5360]) + # Setup the model + conds = [1e-2,1] + sig = simpeg.Utils.ModelBuilder.defineBlock(M.gridCC,[-1000,-1000,-400],[1000,1000,-200],conds) + sig[M.gridCC[:,2]>0] = 1e-8 + sig[M.gridCC[:,2]<-600] = 1e-1 + sigBG = np.zeros(M.nC) + conds[0] + sigBG[M.gridCC[:,2]>0] = 1e-8 + + ## Setup the the survey object + # Receiver locations + rx_x, rx_y = np.meshgrid(np.arange(-500,501,50),np.arange(-500,501,50)) + rx_loc = np.hstack((simpeg.Utils.mkvc(rx_x,2),simpeg.Utils.mkvc(rx_y,2),np.zeros((np.prod(rx_x.shape),1)))) + # Make a receiver list + rxList = [] + for loc in rx_loc: + # NOTE: loc has to be a (1,3) np.ndarray otherwise errors accure + for rxType in ['zxxr','zxxi','zxyr','zxyi','zyxr','zyxi','zyyr','zyyi','tzxr','tzxi','tzyr','tzyi']: + rxList.append(MT.Rx(simpeg.mkvc(loc,2).T,rxType)) + # Source list + srcList =[] + for freq in np.logspace(3,-3,7): + srcList.append(MT.SrcMT.polxy_1Dprimary(rxList,freq)) + # Survey MT + survey = MT.Survey(srcList) + + ## Setup the problem object + problem = MT.Problem3D.eForm_ps(M,sigmaPrimary = sigBG) + problem.pair(survey) + problem.Solver = Solver + + # Calculate the data + fields = problem.fields(sig) + dataVec = survey.projectFields(fields) + + # Make the data + mtData = MT.Data(survey,dataVec) + # Add plots + if plotIt: + pass + +if __name__ == '__main__': + run() diff --git a/SimPEG/MT/Examples/simple3DfowardProblem.py b/SimPEG/MT/Examples/simple3DfowardProblem.py deleted file mode 100644 index feabae28..00000000 --- a/SimPEG/MT/Examples/simple3DfowardProblem.py +++ /dev/null @@ -1,41 +0,0 @@ -# Test script to use SimPEG.MT platform to forward model synthetic data. - -# Import -import SimPEG as simpeg -from SimPEG import MT -import numpy as np - -# Make a mesh -M = simpeg.Mesh.TensorMesh([[(100,5,-1.5),(100.,10),(100,5,1.5)],[(100,5,-1.5),(100.,10),(100,5,1.5)],[(100,5,1.6),(100.,10),(100,3,2)]], x0=['C','C',-3529.5360]) -# Setup the model -conds = [1e-2,1] -sig = simpeg.Utils.ModelBuilder.defineBlock(M.gridCC,[-1000,-1000,-400],[1000,1000,-200],conds) -sig[M.gridCC[:,2]>0] = 1e-8 -sig[M.gridCC[:,2]<-600] = 1e-1 -sigBG = np.zeros(M.nC) + conds[0] -sigBG[M.gridCC[:,2]>0] = 1e-8 - -## Setup the the survey object -# Receiver locations -rx_x, rx_y = np.meshgrid(np.arange(-500,501,50),np.arange(-500,501,50)) -rx_loc = np.hstack((simpeg.Utils.mkvc(rx_x,2),simpeg.Utils.mkvc(rx_y,2),np.zeros((np.prod(rx_x.shape),1)))) -# Make a receiver list -rxList = [] -for loc in rx_loc: - # NOTE: loc has to be a (1,3) np.ndarray otherwise errors accure - for rxType in ['zxxr','zxxi','zxyr','zxyi','zyxr','zyxi','zyyr','zyyi']: - rxList.append(MT.RxMT(simpeg.mkvc(loc,2).T,rxType)) -# Source list -srcList =[] -for freq in np.logspace(3,-3,7): - srcList.append(MT.SurveyMT.srcMT(freq,rxList)) -# Survey MT -survey = MT.Survey(srcList) - -## Setup the problem object -problem = MT.ProblemMT.MTProblem(M) -problem.pair(survey) - -fields = problem.fields(sig,sigBG) -mtData = survey.projectFields(fields) -