diff --git a/SimPEG/Regularization.py b/SimPEG/Regularization.py index cdcfd59a..95bb40cc 100644 --- a/SimPEG/Regularization.py +++ b/SimPEG/Regularization.py @@ -285,22 +285,17 @@ class Tikhonov(BaseRegularization): class Simple(BaseRegularization): """ Only for tensor mesh - """ - + """ + smoothModel = True #: SMOOTH and SMOOTH_MOD_DIF options - alpha_s = Utils.dependentProperty('_alpha_s', 1.0, ['_W', '_Ws'], "Smallness weight") - alpha_x = Utils.dependentProperty('_alpha_x', 1.0, ['_W', '_Wx'], "Weight for the first derivative in the x direction") - alpha_y = Utils.dependentProperty('_alpha_y', 1.0, ['_W', '_Wy'], "Weight for the first derivative in the y direction") - alpha_z = Utils.dependentProperty('_alpha_z', 1.0, ['_W', '_Wz'], "Weight for the first derivative in the z direction") - alpha_xx = Utils.dependentProperty('_alpha_xx', 0.0, ['_W', '_Wxx'], "Weight for the second derivative in the x direction") - alpha_yy = Utils.dependentProperty('_alpha_yy', 0.0, ['_W', '_Wyy'], "Weight for the second derivative in the y direction") - alpha_zz = Utils.dependentProperty('_alpha_zz', 0.0, ['_W', '_Wzz'], "Weight for the second derivative in the z direction") + alpha_s = Utils.dependentProperty('_alpha_s', 1.0, ['_W', '_Ws'], "Smallness weight") + alpha_x = Utils.dependentProperty('_alpha_x', 1.0, ['_W', '_Wx'], "Weight for the first derivative in the x direction") + alpha_y = Utils.dependentProperty('_alpha_y', 1.0, ['_W', '_Wy'], "Weight for the first derivative in the y direction") + alpha_z = Utils.dependentProperty('_alpha_z', 1.0, ['_W', '_Wz'], "Weight for the first derivative in the z direction") def __init__(self, mesh, mapping=None, **kwargs): BaseRegularization.__init__(self, mesh, mapping=mapping, **kwargs) - - @property def Ws(self): """Regularization matrix Ws""" @@ -329,36 +324,15 @@ class Simple(BaseRegularization): self._Wz = Utils.sdiag((self.mesh.vol*self.alpha_z)**0.5)*self.mesh.unitCellGradz return self._Wz - @property - def Wxx(self): - """Regularization matrix Wxx""" - if getattr(self, '_Wxx', None) is None: - self._Wxx = Utils.sdiag((self.mesh.vol*self.alpha_xx)**0.5)*self.mesh.faceDivx*self.mesh.cellGradx - return self._Wxx - - @property - def Wyy(self): - """Regularization matrix Wyy""" - if getattr(self, '_Wyy', None) is None: - self._Wyy = Utils.sdiag((self.mesh.vol*self.alpha_yy)**0.5)*self.mesh.faceDivy*self.mesh.cellGrady - return self._Wyy - - @property - def Wzz(self): - """Regularization matrix Wzz""" - if getattr(self, '_Wzz', None) is None: - self._Wzz = Utils.sdiag((self.mesh.vol*self.alpha_zz)**0.5)*self.mesh.faceDivz*self.mesh.cellGradz - return self._Wzz - @property def Wsmooth(self): """Full smoothness regularization matrix W""" if getattr(self, '_Wsmooth', None) is None: - wlist = (self.Wx, self.Wxx) + wlist = (self.Wx,) if self.mesh.dim > 1: - wlist += (self.Wy, self.Wyy) + wlist += (self.Wy,) if self.mesh.dim > 2: - wlist += (self.Wz, self.Wzz) + wlist += (self.Wz,) self._Wsmooth = sp.vstack(wlist) return self._Wsmooth @@ -413,20 +387,18 @@ class Simple(BaseRegularization): return out class SparseRegularization(Simple): - - eps = 1e-1 - m = None + eps = 1e-1 + m = None gamma = 1. - p = 0. - qx = 2. - qy = 2. - qz = 2. - + p = 0. + qx = 2. + qy = 2. + qz = 2. + def __init__(self, mesh, mapping=None, **kwargs): Simple.__init__(self, mesh, mapping=mapping, **kwargs) - - + @property def Wsmooth(self): """Full smoothness regularization matrix W""" @@ -446,35 +418,35 @@ class SparseRegularization(Simple): wlist = (self.Ws, self.Wsmooth) self._W = sp.vstack(wlist) return self._W - + @property def Ws(self): """Regularization matrix Ws""" if getattr(self, 'm', None) is None: - self.Rs = Utils.speye(self.mesh.nC) - + self.Rs = Utils.speye(self.mesh.nC) + else: f_m = self.m self.rs = self.R(f_m , self.p, self.eps) #print "Min rs: " + str(np.max(self.rs)) + "Max rs: " + str(np.min(self.rs)) self.Rs = Utils.sdiag( self.rs ) - + self._Ws = Utils.sdiag((self.mesh.vol*self.alpha_s*self.gamma)**0.5)*self.Rs - + return self._Ws @property def Wx(self): """Regularization matrix Wx""" - + if getattr(self, 'm', None) is None: - self.Rx = Utils.speye(self.mesh.unitCellGradx.shape[0]) - + self.Rx = Utils.speye(self.mesh.unitCellGradx.shape[0]) + else: f_m = self.mesh.unitCellGradx * self.m self.rx = self.R( f_m , self.qx, self.eps) self.Rx = Utils.sdiag( self.rx ) - + if getattr(self, '_Wx', None) is None: self._Wx = Utils.sdiag((self.mesh.vol*self.alpha_x*self.gamma)**0.5)*self.Rx*self.mesh.unitCellGradx return self._Wx @@ -482,15 +454,15 @@ class SparseRegularization(Simple): @property def Wy(self): """Regularization matrix Wy""" - + if getattr(self, 'm', None) is None: - self.Ry = Utils.speye(self.mesh.unitCellGrady.shape[0]) - + self.Ry = Utils.speye(self.mesh.unitCellGrady.shape[0]) + else: f_m = self.mesh.unitCellGrady * self.m self.ry = self.R( f_m , self.qy, self.eps) self.Ry = Utils.sdiag( self.ry ) - + if getattr(self, '_Wy', None) is None: self._Wy = Utils.sdiag((self.mesh.vol*self.alpha_y*self.gamma)**0.5)*self.Ry*self.mesh.unitCellGrady return self._Wy @@ -500,18 +472,18 @@ class SparseRegularization(Simple): """Regularization matrix Wz""" if getattr(self, 'm', None) is None: - self.Rz = Utils.speye(self.mesh.unitCellGradz.shape[0]) - + self.Rz = Utils.speye(self.mesh.unitCellGradz.shape[0]) + else: f_m = self.mesh.unitCellGradz * self.m self.rz = self.R( f_m , self.qz, self.eps) self.Rz = Utils.sdiag( self.rz ) - + if getattr(self, '_Wz', None) is None: self._Wz = Utils.sdiag((self.mesh.vol*self.alpha_z*self.gamma)**0.5)*self.Rz*self.mesh.unitCellGradz - return self._Wz + return self._Wz + - def R(self, f_m , p, dec): eta = (self.eps**(1-p/2.))**0.5