We should add a sub module in downunder that provides an easy way to run basic geopysical models with
problem specific input, simple tests (susceptability >=0?) and description in the user's guide.
attached is the magnetic.py as an example. The script can be improved (use 3D, magnetization proportional to total field not background,
possibility to choose the faces to fix, e.g. self.fixTop(True), etc; getTotalField())
This would be for starters:
a) grav
b) mag
c) TE MT
d) TM MT
c) sonic wave in frequency domain with PML
from esys.escript import
from esys.escript.linearPDEs import LinearSinglePDE, SolverOptions
class MagneticModel2D(object):
"""
a simple wrapper for a magnetic forward model
"""
def init(self, domain):
assert domain.getDim() == 2
self.pde=LinearSinglePDE(domain)
optionsG=self.pde.getSolverOptions()
if hasFeature('trilinos'):
optionsG.setPackage(SolverOptions.TRILINOS)
optionsG.setPreconditioner(SolverOptions.AMG)
optionsG.setTrilinosParameter("problem:type", "Poisson-2D")
x=domain.getX()
self.pde.setValue(A=kronecker(domain), q=whereZero(x[0]-inf(x[0]))whereZero(x[1]-inf(x[1]) ))
self.setSusceptibility()
self.setBackgroundMagneticField()
def setSusceptibility(self, k=0):
"""
sets the susceptibility
"""
self.k=k
self.reset=True
def setBackgroundMagneticField(self, Bh=[ 0., 45000.0] ):
"""
sets background magnet field in nT
"""
self.Bh=Bh
self.reset=True
def getAnomalyPotential(self):
"""
get the potential of the the magnetic anomaly
"""
if self.reset:
self.pde.setValue(X=self.k*self.Bh)
return self.pde.getSolution()
def getMagneticFieldAnomaly(self):
"""
get the total Magnetic field
"""
return -grad(self.getAnomalyPotential(), ReducedFunction(self.pde.getDomain()))
LutzGross
commented
2 years ago
We should add a sub module in downunder that provides an easy way to run basic geopysical models with problem specific input, simple tests (susceptability >=0?) and description in the user's guide. attached is the magnetic.py as an example. The script can be improved (use 3D, magnetization proportional to total field not background, possibility to choose the faces to fix, e.g. self.fixTop(True), etc; getTotalField())
This would be for starters:
a) grav b) mag c) TE MT d) TM MT c) sonic wave in frequency domain with PML
This would be very useful for teaching!
------------------------------------------- magnetic.py
from esys.escript import from esys.escript.linearPDEs import LinearSinglePDE, SolverOptions class MagneticModel2D(object): """ a simple wrapper for a magnetic forward model """ def init(self, domain): assert domain.getDim() == 2 self.pde=LinearSinglePDE(domain) optionsG=self.pde.getSolverOptions() if hasFeature('trilinos'): optionsG.setPackage(SolverOptions.TRILINOS) optionsG.setPreconditioner(SolverOptions.AMG) optionsG.setTrilinosParameter("problem:type", "Poisson-2D") x=domain.getX() self.pde.setValue(A=kronecker(domain), q=whereZero(x[0]-inf(x[0]))whereZero(x[1]-inf(x[1]) )) self.setSusceptibility() self.setBackgroundMagneticField() def setSusceptibility(self, k=0): """ sets the susceptibility """ self.k=k self.reset=True def setBackgroundMagneticField(self, Bh=[ 0., 45000.0] ): """ sets background magnet field in nT """ self.Bh=Bh self.reset=True def getAnomalyPotential(self): """ get the potential of the the magnetic anomaly """ if self.reset: self.pde.setValue(X=self.k*self.Bh) return self.pde.getSolution() def getMagneticFieldAnomaly(self): """ get the total Magnetic field """ return -grad(self.getAnomalyPotential(), ReducedFunction(self.pde.getDomain()))