1. TMz极化
假设TMz极化均匀平面波垂直入射半径为a的无限长均匀介质圆柱,相对介电常数为εr,相对磁导率为μr,波的传播方向为+x,入射电场和入射磁场用柱面波展开,分别表示为 (1)
(2)
2014浙江高考语文作文散射场朝外传播,因此,散射电场和磁场用柱第二类Hankel函数展开,分别表示如下
(3)
(4)
透射场则由柱面基本波函数的线性组合表示,由于透射场在介质内部均为有限大,因此
(5)
(6)
根据介质表面的边界条件,切向电场和切向磁场连续,可以得到
(7)
(8)
求解方程组,从而得到展开项的系数为
(9)
(1 0)
令,将系数带入展开式,得到散射电场和磁场的表达式为
(11)
(12)
对于远区散射场,kρ → ∞,,则相应的电场和磁场为
(13)
(14)
c**********************************************************************
c Compute TMz Scattering from homogenerous lossless dielectric Circular
c Cylinder by Mie Series
c
c a INPUT, real(8)
c On entry, 'a' specifies the radius of the circular cylinder
c epsR INPUT, real(8)
c On entry, 'epsR' specifies the relative permittivity of
c the homogenerous dielectric circular cylinder
c MuR INPUT, real(8)
c On entry, 'muR' specifies the relative permeability of
c the homogenerous dielectric circular cylinder
c f INPUT, real(8)
c On entry, 'f' specifies the incident frequency
c r INPUT, real(8)
c On entry, 'r' specifies the distance between the observation
c point and the origin of coordinates
c ph INPUT, real(8)
c On entry, 'ph' specifies the observation angle
姜山实验中学
c Ez OUTPUT, complex(8)
c On exit, 'Ez' specifies the z component of the electric
c scattering field
c Hpho OUTPUT, complex(8)
c On exit, 'Hpho' specifies the pho component of the magnetic
c scattering field
c Hphi OUTPUT, complex(8)
c On exit, 'Hphi' specifies the phi component of the magnetic
c scattering field
华硕n81
c
c Programmed by Panda Brewmaster
c**********************************************************************
subroutine dSca_TM_DIE_Cir_Cyl_Mie(a, epsR, muR, f, r, ph, Ez, Hpho, Hphi)
c**********************************************************************
implicit none
c - Input Parameters
real(8) a, epsR, muR, f, r, ph
complex(8) Ez, Hpho, Hphi
c - Constant Numbers
real(8), parameter :: pi = 3.141592653589793
real(8), parameter :: eps0 = 8.854187817d-12
real(8), parameter :: mu0 = pi * 4.d-7
complex(8), parameter :: cj = dcmplx(0.d0, 1.d0)
c - Temporary Variables
integer k, nmax
real(8) eta0, wavek0, ka0, kr
杜邦分析法 real(8) eta1, wavek1, ka1
complex(8) coe1, coe2
懦艮 real(8), allocatable, dimension (:) :: Jnka0, Ynka0, DJnka0, Jnkr, Ynkr
real(8), allocatable, dimension (:) :: Jnka1, Ynka1, DJnka1
complex(8), allocatable, dimension (:) :: Hnka0, DHnka0, Hnkr, DHnkr
complex(8), allocatable, dimension (:) :: Hnka1, DHnka1
complex(8), allocatable, dimension (:) :: an
eta0 = dsqrt(mu0 / eps0)
wavek0 = 2.d0 * pi * f * dsqrt(mu0 * eps0)
ka0 = wavek0 * a
kr = wavek0 * r
eta1 = dsqrt(muR * mu0 / epsR / eps0)
wavek1 = 2.d0 * pi * f * dsqrt(muR * mu0 * epsR * eps0)
ka1 = wavek1 * a
nmax = ka1 + 10.d0 * ka1 ** (1.d0 / 3.d0) + 1
if(nmax < 5) nmax = 5
allocate(Jnka0(- 1 : nmax + 1), Ynka0(- 1 : nmax + 1), Hnka0(- 1 : nmax + 1), DJnka0(0 : nmax), DHnka0(0 : nmax))
call dBES(nmax + 1, ka0, Jnka0(0 : nmax + 1), Ynka0(0 : nmax + 1))
Jnka0(- 1) = - Jnka0(1); Ynka0(- 1) = - Ynka0(1)
Hnka0(-1 : nmax + 1) = dcmplx(Jnka0(-1 : nmax + 1), - Ynka0(-1 : nmax + 1))
do k = 0, nmax
DJnka0(k) = (Jnka0(k - 1) - Jnka0(k + 1)) / 2.d0
DHnka0(k) = (Hnka0(k - 1) - Hnka0(k + 1)) / 2.d0
enddo
allocate(Jnka1(- 1 : nmax + 1), Ynka1(- 1 : nmax + 1), Hnka1(- 1 : nmax + 1), DJnka1(0 : nmax), DHnka1(0 : nmax))
call dBES(nmax + 1, ka1, Jnka1(0 : nmax + 1), Ynka1(0 : nmax + 1))
iso17799 Jnka1(- 1) = - Jnka1(1); Ynka1(- 1) = - Ynka1(1)
Hnka1(-1 : nmax + 1) = dcmplx(Jnka1(-1 : nmax + 1), - Ynka1(-1 : nmax + 1))
do k = 0, nmax
DJnka1(k) = (Jnka1(k - 1) - Jnka1(k + 1)) / 2.0
DHnka1(k) = (Hnka1(k - 1) - Hnka1(k + 1)) / 2.0
enddo
allocate(an(0 : nmax))
do k = 0, nmax
coe1 = (eta1 * DJnka0(k) * Jnka1(k) - eta0 * Jnka0(k) * DJnka1(k))
coe2 = (eta0 * DJnka1(k) * Hnka0(k) - eta1 * Jnka1(k) * DHnka0(k))
an(k) = coe1 / coe2
enddo
if(r >= 0) then ! Finite Distance
allocate(Jnkr(- 1 : nmax + 1), Ynkr(- 1 : nmax + 1), Hnkr(- 1 : nmax + 1), DHnkr(0 : nmax))
call dBES(nmax + 2, kr, Jnkr(0 : nmax + 1), Ynkr(0 : nmax + 1))
Jnkr(- 1) = - Jnkr(1)
Ynkr(- 1) = - Ynkr(1)
Hnkr(- 1 : nmax + 1) = dcmplx(Jnkr(- 1 : nmax + 1), - Ynkr(- 1 : nmax + 1))
Ez = an(0) * Hnkr(0)
do k = 1, nmax
Ez = Ez + cj ** (- k) * an(k) * Hnkr(k) * cdexp(dcmplx(0.d0, k * ph))
Ez = Ez + cj ** k * an(k) * Hnkr(k) * cdexp(dcmplx(0.d0, - k * ph))
enddo
Hpho = 0.d0
do k = 1, nmax