wassname/simpeg
1
0
Fork 0
mirror of https://github.com/wassname/simpeg.git synced 2026-07-11 03:02:29 +08:00
Code Issues Packages Projects Releases Wiki Activity

made mu a vector in the casing soln, added E fields from an electric dipole

Browse Source
This commit is contained in:
Lindsey
2015-06-05 18:07:35 -07:00
parent 678d8611d8
commit 611e930925
2 changed files with 44 additions and 17 deletions
Download Patch File Download Diff File Expand all files Collapse all files
simpegEM/Analytics/FDEM.py
+27
View File
@@ -111,3 +111,30 @@ def AnalyticMagDipoleWholeSpace(XYZ, srcLoc, sig, f, m=1., orientation='X'):
return Bx, By, Bz
def ElectricDipoleWholeSpace(XYZ, srcLoc, sig, f, m=1., orientation='X'):
XYZ = Utils.asArray_N_x_Dim(XYZ, 3)
dx = XYZ[:,0]-srcLoc[0]
dy = XYZ[:,1]-srcLoc[1]
dz = XYZ[:,2]-srcLoc[2]
r = np.sqrt( dx**2. + dy**2. + dz**2.)
k = np.sqrt( -1j*2.*np.pi*f*mu_0*sig )
kr = k*r
front = moment / (4. * np.pi * sig * r**3) * exp(-1j*k*r)
mid = -k**2 * r**2 + 3*1j*k*r + 3
Ex = front*((dx**2 / r**2)*mid + (k**2 * r**2 -1j*k*r))
Ey = front*(dx*dy / r**2)*mid
Ez = front*(dx*dz / r**2)*mid
if orientation.upper() == 'X':
return Ex, Ey, Ez
elif orientation.upper() == 'Y':
return Ez, Ex, Ey
elif orientation.upper() == 'Z':
return Ey, Ez, Ex
simpegEM/Analytics/FDEMcasing.py
+17 -17
View File
@@ -10,20 +10,20 @@ def getKc(freq,sigma,a,b,mu=mu_0,eps=epsilon_0):
def _r2(xyz):
return np.sum(xyz**2,1)
def _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1.):
Kc1 = getKc(freq,sigma[1],a,b,mu,eps)
def _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
Kc1 = getKc(freq,sigma[1],a,b,mu[1],eps)
nobs = obsloc.shape[0]
dxyz = obsloc - np.c_[np.ones(nobs)]*np.r_[srcloc]
r2 = _r2(dxyz[:,:2])
sqrtr2z2 = np.sqrt(r2 + dxyz[:,2]**2)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
return Kc1 * moment / (4.*np.pi) *np.exp(-1j*k2*sqrtr2z2) / sqrtr2z2
def _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1.):
def _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
@@ -31,12 +31,12 @@ def _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=eps
r2 = _r2(dxyz[:,:2])
sqrtr2z2 = np.sqrt(r2 + dxyz[:,2]**2)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
return -HertzZ * np.sqrt(r2) / sqrtr2z2 * (1j*k2 + 1./ sqrtr2z2)
def _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1.):
def _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
nobs = obsloc.shape[0]
@@ -44,11 +44,11 @@ def _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=eps
r2z2 = _r2(dxyz)
sqrtr2z2 = np.sqrt(r2z2)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
return -HertzZ*dxyz[:,2] /sqrtr2z2 * (1j*k2 + 1./sqrtr2z2)
def _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1.):
def _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
dHertzZdr = _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
@@ -59,11 +59,11 @@ def _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=
r = np.sqrt(r2)
z = dxyz[:,2]
sqrtr2z2 = np.sqrt(r2 + z**2)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
return dHertzZdr*(-z/sqrtr2z2)*(1j*k2+1./sqrtr2z2) + HertzZ*(z*r/sqrtr2z2**3)*(1j*k2 + 2./sqrtr2z2)
def _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1.):
def _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1.):
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
dHertzZdz = _getCasingHertzMagDipoleDeriv_z(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
@@ -74,24 +74,24 @@ def _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=
r = np.sqrt(r2)
z = dxyz[:,2]
sqrtr2z2 = np.sqrt(r2 + z**2)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
return (dHertzZdz*z + HertzZ)/sqrtr2z2*(-1j*k2 - 1./sqrtr2z2) + HertzZ*z/sqrtr2z2**3*(1j*k2*z + 2.*z/sqrtr2z2)
def getCasingEphiMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1):
def getCasingEphiMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1):
return 1j * omega(freq) * mu * _getCasingHertzMagDipoleDeriv_r(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
def getCasingHrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1):
def getCasingHrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1):
return _getCasingHertzMagDipole2Deriv_z_r(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
def getCasingHzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1):
def getCasingHzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1):
d2HertzZdz2 = _getCasingHertzMagDipole2Deriv_z_z(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
k2 = k(freq,sigma[2],mu,eps)
k2 = k(freq,sigma[2],mu[2],eps)
HertzZ = _getCasingHertzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
return d2HertzZdz2 + k2**2 * HertzZ
def getCasingBrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1):
def getCasingBrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1):
return mu_0 * getCasingHrMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)
def getCasingBzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0,eps=epsilon_0,moment=1):
def getCasingBzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu=mu_0*np.ones(3),eps=epsilon_0,moment=1):
return mu_0 * getCasingHzMagDipole(srcloc,obsloc,freq,sigma,a,b,mu,eps,moment)