Compare commits

...
Author SHA1 Message Date
GudniRos 1a5e981dff Fixed bugs in meshutils, writing VTR files
Allow fileName to be None, which will output the VTKobject without save a file.
2015-12-17 22:52:41 -08:00
GudniRos 21c64cbe66 Added dpred to be saved in saveDict directive. 2015-12-15 19:36:31 -08:00
GudniRos eb24a70f31 Merge branch 'master' into pickleSupport 2015-10-26 17:30:16 -07:00
GudniRos 704776b8ba Commented out a reduce method. 2015-10-26 17:14:40 -07:00
GudniRos e4448c2f2e Progressing with pickling. Pickling of PropModels doesn't work which cause
many classes that use it not to pickle.
2015-08-14 11:57:02 -07:00
GudniRos 9d5db11b0e Added a new inversion derictive. 2015-08-13 10:19:45 -07:00
GudniRos e9957d7ec8 Fix difference in Regularization 2015-08-13 10:08:39 -07:00
GudniRos c74022a948 Fix differences in Regularization file. 2015-08-13 10:08:39 -07:00
Lindsey Heagy 7cd4ba7d61 updated inversion workflow image, discretization --> mesh and colors are same as in paper 2015-07-21 14:31:22 -07:00
12 changed files with 496 additions and 53 deletions
+28
View File
@@ -14,6 +14,34 @@ class BaseDataMisfit(object):
debug = False #: Print debugging information
counter = None #: Set this to a SimPEG.Utils.Counter() if you want to count things
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
def __init__(self, survey, **kwargs):
assert survey.ispaired, 'The survey must be paired to a problem.'
if isinstance(survey, Survey.BaseSurvey):
+98
View File
@@ -8,6 +8,34 @@ class InversionDirective(object):
def __init__(self, **kwargs):
Utils.setKwargs(self, **kwargs)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def inversion(self):
"""This is the inversion of the InversionDirective instance."""
@@ -207,6 +235,76 @@ class SaveOutputEveryIteration(_SaveEveryIteration):
f.close()
class SaveOutputDictEveryIteration(_SaveEveryIteration):
"""SaveOutputDictEveryIteration"""
def initialize(self):
print "SimPEG.SaveOutputDictEveryIteration will save your inversion progress as dictionary: '###-%s.npz'"%self.fileName
def endIter(self):
# Save the data.
ms = self.reg.Ws * ( self.reg.mapping * (self.invProb.curModel - self.reg.mref) )
phi_ms = 0.5*ms.dot(ms)
if self.reg.smoothModel == True:
mref = self.reg.mref
else:
mref = 0
mx = self.reg.Wx * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_mx = 0.5 * mx.dot(mx)
if self.prob.mesh.dim==2:
my = self.reg.Wy * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_my = 0.5 * my.dot(my)
else:
phi_my = 'NaN'
if self.prob.mesh.dim==3:
mz = self.reg.Wz * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_mz = 0.5 * mz.dot(mz)
else:
phi_mz = 'NaN'
# Save the file as a npz
np.savez('{:03d}-{:s}'.format(self.opt.iter,self.fileName), iter=self.opt.iter, beta=self.invProb.beta, phi_d=self.invProb.phi_d, phi_m=self.invProb.phi_m, phi_ms=phi_ms, phi_mx=phi_mx, phi_my=phi_my, phi_mz=phi_mz,f=self.opt.f, m=self.invProb.curModel)
class SaveOutputDictEveryIteration(_SaveEveryIteration):
"""SaveOutputDictEveryIteration
A directive that saves some relevant information from the inversion run to a numpy .npz dictionary file (see numpy.savez function for further info).
"""
def initialize(self):
print "SimPEG.SaveOutputDictEveryIteration will save your inversion progress as dictionary: '###-%s.npz'"%self.fileName
def endIter(self):
# Save the data.
ms = self.reg.Ws * ( self.reg.mapping * (self.invProb.curModel - self.reg.mref) )
phi_ms = 0.5*ms.dot(ms)
if self.reg.smoothModel == True:
mref = self.reg.mref
else:
mref = 0
mx = self.reg.Wx * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_mx = 0.5 * mx.dot(mx)
if self.prob.mesh.dim==2:
my = self.reg.Wy * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_my = 0.5 * my.dot(my)
else:
phi_my = 'NaN'
if self.prob.mesh.dim==3:
mz = self.reg.Wz * ( self.reg.mapping * (self.invProb.curModel - mref) )
phi_mz = 0.5 * mz.dot(mz)
else:
phi_mz = 'NaN'
# Save the file as a npz
np.savez('{:03d}-{:s}'.format(self.opt.iter,self.fileName), iter=self.opt.iter, beta=self.invProb.beta, phi_d=self.invProb.phi_d, phi_m=self.invProb.phi_m, phi_ms=phi_ms, phi_mx=phi_mx, phi_my=phi_my, phi_mz=phi_mz,f=self.opt.f, m=self.invProb.curModel,dpred=self.invProb.dpred)
# class UpdateReferenceModel(Parameter):
+28
View File
@@ -12,6 +12,34 @@ class Fields(object):
aliasFields = None #: Aliased fields, a dict with [alias, location, function], e.g. {"b":["e","F",lambda(F,e,ind)]}
dtype = float #: dtype is the type of the storage matrix. This can be a dictionary.
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
def __init__(self, mesh, survey, **kwargs):
self.survey = survey
self.mesh = mesh
+27 -3
View File
@@ -17,6 +17,30 @@ class IdentityMap(object):
Utils.setKwargs(self, **kwargs)
self.mesh = mesh
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
@property
def nP(self):
"""
@@ -289,7 +313,7 @@ class FullMap(IdentityMap):
"""
FullMap
Given a scalar, the FullMap maps the value to the
Given a scalar, the FullMap maps the value to the
full model space.
"""
@@ -314,8 +338,8 @@ class FullMap(IdentityMap):
:rtype: numpy.array
:return: derivative of transformed model
"""
return np.ones([self.mesh.nC,1])
return np.ones([self.mesh.nC,1])
class Vertical1DMap(IdentityMap):
"""Vertical1DMap
+1 -1
View File
@@ -33,7 +33,7 @@ class InnerProducts(object):
return self._getInnerProduct('E', prop=prop, invProp=invProp, invMat=invMat, doFast=doFast)
def _getInnerProduct(self, projType, prop=None, invProp=False, invMat=False, doFast=True):
"""
"""r
:param str projType: 'F' for faces 'E' for edges
:param numpy.array prop: material property (tensor properties are possible) at each cell center (nC, (1, 3, or 6))
:param bool invProp: inverts the material property
+28
View File
@@ -115,6 +115,34 @@ class Minimize(object):
Utils.setKwargs(self, **kwargs)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def callback(self):
return getattr(self, '_callback', None)
+30 -2
View File
@@ -22,6 +22,34 @@ class BaseProblem(object):
PropMap = None #: A SimPEG PropertyMap class.
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def mapping(self):
"A SimPEG.Map instance or a property map is PropMap is not None"
@@ -32,8 +60,8 @@ class BaseProblem(object):
val._assertMatchesPair(self.mapPair)
self._mapping = val
else:
self._mapping = self.PropMap(val)
self._mapping = self.PropMap(val)
def __init__(self, mesh, mapping=None, **kwargs):
Utils.setKwargs(self, **kwargs)
assert isinstance(mesh, Mesh.BaseMesh), "mesh must be a SimPEG.Mesh object."
+77 -1
View File
@@ -13,6 +13,34 @@ class Property(object):
self.doc = doc
Utils.setKwargs(self, **kwargs)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def propertyLink(self):
"Can be something like: ('sigma', Maps.ReciprocalMap)"
@@ -118,6 +146,36 @@ class PropModel(object):
self.vector = vector
assert len(self.vector) == self.nP
# Pickleing support methods
# def __reduce__(self):
# return (dict,{self.propMap,self.vector})
# def __getstate__(self):
# '''
# Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
# Used when doing:
# pickle.dump(pickleFile,object)
# '''
# self.__class__ = ProbModel
# odict = {}
# odict['vec'] = self.__dict__['vector']
# odict['pMap'] = self.__dict__['propMap']
# # Return the dict
# return odict
# def __setstate__(self,odict):
# '''
# Function that sets a pickle dictionary in to an object.
# Used when doing:
# object = pickle.load(pickleFile)
# '''
# # Update the dict
# # Re-hook the methods to the object
# self.propMap = odict['prMap']
# self.vector = odict['vec']
@property
def nP(self):
inds = []
@@ -207,6 +265,24 @@ class PropMap(object):
else:
raise Exception('mappings must be a dict, a mapping, or a list of tuples.')
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
pass
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
pass
def setup(self, maps, slices=None):
"""
@@ -239,7 +315,7 @@ class PropMap(object):
setattr(self, '%sMap'%name, mapping)
setattr(self, '%sIndex'%name, slices.get(name, slice(nP, nP + mapping.nP)))
nP += mapping.nP
self.nP = nP
self.nP = nP
@property
def defaultInvProp(self):
+29 -1
View File
@@ -27,6 +27,34 @@ class BaseRegularization(object):
self.mapping = mapping or Maps.IdentityMap(mesh)
self.mapping._assertMatchesPair(self.mapPair)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def parent(self):
"""This is the parent of the regularization."""
@@ -311,7 +339,7 @@ class Tikhonov(BaseRegularization):
if self.smoothModel == True:
mD1 = self.mapping.deriv(m)
mD2 = self.mapping.deriv(m - self.mref)
r1 = self.Wsmooth * ( self.mapping * (m))
r1 = self.Wsmooth * ( self.mapping * (m))
r2 = self.Ws * ( self.mapping * (m - self.mref) )
out1 = mD1.T * ( self.Wsmooth.T * r1 )
out2 = mD2.T * ( self.Ws.T * r2 )
+104 -1
View File
@@ -19,6 +19,35 @@ class BaseRx(object):
self._Ps = {}
Utils.setKwargs(self, **kwargs)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def rxType(self):
"""Receiver Type"""
@@ -129,6 +158,33 @@ class BaseSrc(object):
self.rxList = rxList
Utils.setKwargs(self, **kwargs)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def nD(self):
@@ -153,6 +209,25 @@ class Data(object):
if v is not None:
self.fromvec(v)
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
pass
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
pass
def _ensureCorrectKey(self, key):
if type(key) is tuple:
if len(key) is not 2:
@@ -210,11 +285,39 @@ class BaseSurvey(object):
mtrue = None #: True model, if data is synthetic
counter = None #: A SimPEG.Utils.Counter object
srcPair = BaseSrc #: Source Pair
def __init__(self, **kwargs):
Utils.setKwargs(self, **kwargs)
srcPair = BaseSrc #: Source Pair
# Pickleing support methods
def __getstate__(self):
'''
Method that makes the dictionary of the object pickleble, removes non-pickleble elements of the object.
Used when doing:
pickle.dump(pickleFile,object)
'''
odict = self.__dict__.copy()
# Remove fields that are not needed
del odict['hook']
del odict['setKwargs']
# Return the dict
return odict
def __setstate__(self,odict):
'''
Function that sets a pickle dictionary in to an object.
Used when doing:
object = pickle.load(pickleFile)
'''
# Update the dict
self.__dict__.update(odict)
# Re-hook the methods to the object
Utils.codeutils.hook(self,Utils.codeutils.hook)
Utils.codeutils.hook(self,Utils.codeutils.setKwargs)
@property
def srcList(self):
+46 -44
View File
@@ -149,7 +149,7 @@ def readUBCTensorModel(fileName, mesh):
Input:
:param fileName, path to the UBC GIF mesh file to read
:param mesh, TensorMesh object, mesh that coresponds to the model
:param mesh, TensorMesh object, mesh that coresponds to the model
Output:
:return numpy array, model with TensorMesh ordered
@@ -170,7 +170,7 @@ def writeUBCTensorMesh(fileName, mesh):
:param str fileName: File to write to
:param simpeg.Mesh.TensorMesh mesh: The mesh
"""
assert mesh.dim == 3
s = ''
@@ -205,7 +205,6 @@ def writeUBCTensorModel(fileName, mesh, model):
np.savetxt(fileName, modelMatTR.ravel())
def readVTRFile(fileName):
"""
Read VTK Rectilinear (vtr xml file) and return SimPEG Tensor mesh and model
@@ -216,7 +215,7 @@ def readVTRFile(fileName):
Output:
:return SimPEG TensorMesh object
:return SimPEG model dictionary
"""
# Import
from vtk import vtkXMLRectilinearGridReader as vtrFileReader
@@ -296,84 +295,87 @@ def writeVTRFile(fileName,mesh,model=None):
vtkObj.SetZCoordinates(numpy_to_vtk(vZ,deep=1))
# Assign the model('s) to the object
for item in model.iteritems():
# Convert numpy array
vtkDoubleArr = numpy_to_vtk(item[1],deep=1)
vtkDoubleArr.SetName(item[0])
vtkObj.GetCellData().AddArray(vtkDoubleArr)
# Set the active scalar
vtkObj.GetCellData().SetActiveScalars(model.keys()[0])
vtkObj.Update()
if model is not None:
for item in model.iteritems():
# Convert numpy array
vtkDoubleArr = numpy_to_vtk(item[1],deep=1)
vtkDoubleArr.SetName(item[0])
vtkObj.GetCellData().AddArray(vtkDoubleArr)
# Set the active scalar
vtkObj.GetCellData().SetActiveScalars(model.keys()[0])
# Check the extension of the fileName
ext = os.path.splitext(fileName)[1]
if ext is '':
fileName = fileName + '.vtr'
elif ext not in '.vtr':
raise IOError('{:s} is an incorrect extension, has to be .vtr')
# Write the file.
vtrWriteFilter = rectWriter()
vtrWriteFilter.SetInput(vtkObj)
vtrWriteFilter.SetFileName(fileName)
vtrWriteFilter.Update()
if fileName is not None:
ext = os.path.splitext(fileName)[1]
if ext is '':
fileName = fileName + '.vtr'
elif ext not in '.vtr':
raise IOError('{:s} is an incorrect extension, has to be .vtr')
# Write the file.
vtrWriteFilter = rectWriter()
vtrWriteFilter.SetInputData(vtkObj)
vtrWriteFilter.SetFileName(fileName)
vtrWriteFilter.Update()
else:
return vtkObj
def ExtractCoreMesh(xyzlim, mesh, meshType='tensor'):
"""
Extracts Core Mesh from Global mesh
xyzlim: 2D array [ndim x 2]
mesh: SimPEG mesh
This function ouputs:
This function ouputs:
- actind: corresponding boolean index from global to core
- meshcore: core SimPEG mesh
- meshcore: core SimPEG mesh
Warning: 1D and 2D has not been tested
"""
from SimPEG import Mesh
if mesh.dim ==1:
xyzlim = xyzlim.flatten()
xmin, xmax = xyzlim[0], xyzlim[1]
xind = np.logical_and(mesh.vectorCCx>xmin, mesh.vectorCCx<xmax)
xind = np.logical_and(mesh.vectorCCx>xmin, mesh.vectorCCx<xmax)
xc = mesh.vectorCCx[xind]
hx = mesh.hx[xind]
x0 = [xc[0]-hx[0]*0.5, yc[0]-hy[0]*0.5]
meshCore = Mesh.TensorMesh([hx, hy] ,x0=x0)
actind = (mesh.gridCC[:,0]>xmin) & (mesh.gridCC[:,0]<xmax)
elif mesh.dim ==2:
xmin, xmax = xyzlim[0,0], xyzlim[0,1]
ymin, ymax = xyzlim[1,0], xyzlim[1,1]
yind = np.logical_and(mesh.vectorCCy>ymin, mesh.vectorCCy<ymax)
zind = np.logical_and(mesh.vectorCCz>zmin, mesh.vectorCCz<zmax)
zind = np.logical_and(mesh.vectorCCz>zmin, mesh.vectorCCz<zmax)
xc = mesh.vectorCCx[xind]
yc = mesh.vectorCCy[yind]
hx = mesh.hx[xind]
hy = mesh.hy[yind]
x0 = [xc[0]-hx[0]*0.5, yc[0]-hy[0]*0.5]
meshCore = Mesh.TensorMesh([hx, hy] ,x0=x0)
actind = (mesh.gridCC[:,0]>xmin) & (mesh.gridCC[:,0]<xmax) \
& (mesh.gridCC[:,1]>ymin) & (mesh.gridCC[:,1]<ymax) \
elif mesh.dim==3:
xmin, xmax = xyzlim[0,0], xyzlim[0,1]
ymin, ymax = xyzlim[1,0], xyzlim[1,1]
zmin, zmax = xyzlim[2,0], xyzlim[2,1]
xind = np.logical_and(mesh.vectorCCx>xmin, mesh.vectorCCx<xmax)
yind = np.logical_and(mesh.vectorCCy>ymin, mesh.vectorCCy<ymax)
zind = np.logical_and(mesh.vectorCCz>zmin, mesh.vectorCCz<zmax)
zind = np.logical_and(mesh.vectorCCz>zmin, mesh.vectorCCz<zmax)
xc = mesh.vectorCCx[xind]
yc = mesh.vectorCCy[yind]
@@ -382,19 +384,19 @@ def ExtractCoreMesh(xyzlim, mesh, meshType='tensor'):
hx = mesh.hx[xind]
hy = mesh.hy[yind]
hz = mesh.hz[zind]
x0 = [xc[0]-hx[0]*0.5, yc[0]-hy[0]*0.5, zc[0]-hz[0]*0.5]
meshCore = Mesh.TensorMesh([hx, hy, hz] ,x0=x0)
actind = (mesh.gridCC[:,0]>xmin) & (mesh.gridCC[:,0]<xmax) \
& (mesh.gridCC[:,1]>ymin) & (mesh.gridCC[:,1]<ymax) \
& (mesh.gridCC[:,2]>zmin) & (mesh.gridCC[:,2]<zmax)
else:
raise(Exception("Not implemented!"))
return actind, meshCore
Binary file not shown.

Before

Width:  |  Height:  |  Size: 59 KiB

After

Width:  |  Height:  |  Size: 58 KiB