二维TM波讨论平面波源(使用直接算方法)的加入
平面波公式
平面波公式
平面波公式是物理学中的一个概念，表达了质点在空间中随时间的变化。
平面波公式可以被用来描述振动，波浪等物理变化情况。
它为我们探索宇宙奥秘提供了重要的指导意义。
平面波公式的本质
平面波公式是由线性液动力学理论建立起来的，这里只考虑了液体在物体表面上的摩擦力，忽略了液体内部的变化现象，所以它可以简化为形如：
F=Acos（k（x-ut））
其中A是振幅，k是波数，u是速度，t是时间，x是空间变量。
该公式描述的是质点随时间的变化受空间变量影响的力学运动规律。
应用
以上的公式可以用来解释许多现象，比如波浪在大海里的传播；地震波在地底传播；声音在空气中传播等等。
平面波公式具有广泛的应用，可以应用于地面建筑和海洋结构物的动力设计，军事和航空航天技术中，将其用于测量和导航，以及电子技术中，将其用于电磁波的控制。
另外，平面波公式还可以用于液力学的研究。
液力学是研究液体力学运动的理论，其关键观点是，当液体运动时，液体中心会产生一个力，并且受空间变量影响，利用平面波公式，可以详细描述液体运动情况，为液力学学习者提供了重要的参考。
最后，平面波公式也可以用来模拟物理实验，比如微波加热实验
中，可以利用平面波公式来模拟微波在物质表面的传播，这从很大程度上提高了实验的可信度和准确性。
总结
总而言之，平面波公式是一个重要的概念，它是描述质点在不同空间位置随时间的变化的表达式，它为我们解开物理世界的奥秘提供了有价值的指导。
平面波公式的本质是线性液动力学，其应用非常广泛，可以用在液力学，航空航天，电子技术，建筑和海洋结构物等，还可以用来模拟实验。
二维波动方程
二维波动方程二维波动方程是物理学和数学中一个重要的概念。
它是一类非线性的偏微分方程，用来描述物理系统中的波动现象。
二维波动方程的应用非常广泛，可以用来研究电磁波、热传导、流体力学等领域中的波动现象。
其中最常用的二维波动方程之一是Korteweg-de Vries方程(KdV方程)，它用来描述一维流体中的纵波。
它可以用来描述潮汐、河流涨落、海浪等现象。
它的形式是：u_t + uu_x + u_{xxx} = 0。
其中u(x,t)表示流体的速度场，x和t分别表示空间坐标和时间。
另一个重要的二维波动方程是Sine-Gordon方程，它描述了非线性玄学领域中的波动现象，例如在线性电力线的电磁波传播，在非线性声学领域中的声波传播。
这个方程的形式是：u_{tt}-u_{xx}+sin(u)=0。
还有一种重要的二维波动方程是Nonlinear Schrödinger方程(NLS方程) ,它是研究光学波动、超流体波动等领域中的重要方程。
NLS方程的形式是：i∂u/∂t+1/2∂²u/∂x²+|u|²u=0。
在这些方程中，通常都采用数值方法和解析方法来解决它们。
数值方法包括有限差分法、有限元法、有限体积法等,而解析方法则有传统的小扰动分析和现代的非线性分析等。
总的来说，二维波动方程是物理学和数学领域中重要的概念，它们被广泛应用于研究物理系统中的波动现象。
使用数值方法和解析方法来解决这些方程是目前常用的方法。
综上所述，二维波动方程是一类重要的非线性偏微分方程，应用广泛，用来描述物理系统中的波动现象。
主要有Korteweg-de Vries方程，Sine-Gordon方程，Nonlinear Schrödinger方程等。
在研究中采用数值方法和解析方法来解决这些方程是常见的方法。
TM波FDTD讨论平面波源的加入
! TM波FDTD讨论平面波源的加入module data_moduleimplicit noneinteger,parameter::nx0=0,nx1=360,ny0=0,ny1=360,nz0=-100,nz1=1200integer,parameter::nxl1=nx0+80,nxl2=nx1-80,nyl1=ny0+80,nyl2=ny1-80 !连接边界real,parameter::f=2.0e8,c=3.0e8,delt=0.0177,deltt=delt/6.0e8,eps0=8.85e-12,miu0=1.2566e-6,pi= 3.14159real,parameter::w=2*pi*f,s=-0.477369real,parameter::p=-1.0/3.0,q=-miu0*c/6,r=-miu0*c/2,p1=1/(2*miu0*c),p2=1/(2*eps0*c)real,parameter::tal=2e-9,t0=0.8*tal,fai=pi/3.0real cez,chx,chyinteger,parameter::nt=2000,m0=200integer ncomplex Ez3(nx0:nx1,ny0:ny1)real Ez4(nx0:nx1,ny0:ny1),Ez2(nx0:nx1,ny0:ny1) !记录幅值提取时的实部和虚部real sita(nx0:nx1,ny0:ny1),Ez0(nx0:nx1,ny0:ny1) !记录幅值和相位real Ez(nx0:nx1,ny0:ny1),Hx(nx0:nx1,ny0:ny1),Hy(nx0:nx1,ny0:ny1),Ez1(nx0:nx1,ny0:ny1) real Ei0(nz0:nz1),Hi0(nz0:nz1),Ei1(nz0:nz1)real Ezi(nx0:nx1,ny0:ny1),Hxi(nx0:nx1,ny0:ny1),Hyi(nx0:nx1,ny0:ny1)end module data_module!/////////////////////////////////////////////////////////////////////////////////////////////////subroutine inc()use data_moduleimplicit noneinteger i,j,kreal t,dt=n*delttEi1=Ei0do k=nz0,nz1-1Hi0(k)=Hi0(k)-p1*(Ei1(k+1)-Ei1(k))end do!Ezido i=nxl1,nxl2do j=nyl1,nyl2d=real(i-nxl1)*cos(fai)+real(j-nyl1)*sin(fai)Ezi(i,j)=(d-int(d))*Ei0(m0+int(d)+1)+(1-(d-int(d)))*Ei0(m0+int(d))end doend dodo k=nz0+1,nz1-1Ei0(k)=Ei1(k)-p2*(Hi0(k)-Hi0(k-1)) ! 入射波的场量end doEi0(m0-2)=sin(w*t) !exp(-(4*pi*(t-t0)**2/Tal/tal))k=nz1Ei0(k)=Ei1(k-1)-1.0/3.0*(Ei0(k-1)-Ei1(k)) !入射波右吸收边界K=nz0Ei0(k)=Ei1(k+1)-1.0/3.0*(Ei0(k+1)-Ei1(k)) !入射波左吸收边界!Hxij=nyl1do i=nxl1,nxl2d=real(i-nxl1)*cos(fai)+(j-0.5-nyl1)*sin(fai)+0.5Hxi(i,j-1)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*sin(fai) end doj=nyl2do i=nxl1,nxl2d=real(i-nxl1)*cos(fai)+(j+0.5-nyl1)*sin(fai)+0.5Hxi(i,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*sin(fai) end do!Hyii=nxl1do j=nyl1,nyl2d=(i-0.5-nyl1)*cos(fai)+real(j-nyl1)*sin(fai)+0.5Hyi(i-1,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*(-cos(fai)) enddoi=nxl2do j=nyl1,nyl2d=(i+0.5-nyl1)*cos(fai)+real(j-nyl1)*sin(fai)+0.5Hyi(i,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*(-cos(fai)) enddowrite(13,*)n,Ei0(0),Ei0(20),Ei0(30),Ei0(100),sin(w*t)end subroutine inc!/////////////////////////////////////////////////////////////////////////////////////////////////program mainuse data_moduleimplicit noneinteger i,j,kreal topen(11,file='A.txt')open(12,file='fai.txt')open(13,file='Ez.txt',recl=300)Ez=0;Hx=0;Hy=0;Ez1=0Ei0=0;Hi0=0Ez4=0;Ez2=0;Ez3=cmplx(0.0,0.0)cez=1.0;chx=sin(fai);chy=-cos(fai)do n=0,ntcall inc()call FDTD()print *,'step',n,Ez(nx1/2,ny1/2)!write(13,*)n,Ei0(0),Ei0(20),Ei0(30),Ei0(100),sin(w*t)Ez1=Ezend doEz4=Ezdo n=nt+1,nt+42call inc()call FDTD()Ez1=Ezend doEz2=Ezcall AP()do j=ny0,ny1do i=nx0,nx1write(11,*)i,j,Ez0(i,j)write(12,*)i,j,sita(i,j)end doend doclose(11)close(12)close(13)end program main!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine FDTD()use data_moduleimplicit noneinteger i,jcall cal_H()call cal_E()end subroutine FDTD!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine cal_H()use data_moduleimplicit noneinteger i,jreal,external::Eido j=ny0,ny1-1Hx(i,j)=Hx(i,j)-p1*(Ez(i,j+1)-Ez(i,j))end doend dodo j=ny0,ny1do i=nx0,nx1-1Hy(i,j)=Hy(i,j)+p1*(Ez(i+1,j)-Ez(i,j))end doend do!连接边界i=nxl1-1do j=nyl1,nyl2Hy(i,j)=Hy(i,j)-p1*Ezi(i+1,j)end doi=nxl2do j=nyl1,nyl2Hy(i,j)=Hy(i,j)+p1*Ezi(i,j)end doj=nyl1-1do i=nxl1,nxl2Hx(i,j)=Hx(i,j)+p1*Ezi(i,j+1)end doj=nyl2do i=nxl1,nxl2Hx(i,j)=Hx(i,j)-p1*Ezi(i,j)end doend subroutine cal_H!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine cal_E()use data_moduleimplicit noneinteger i,jreal,external::Hido j=ny0+1,ny1-1do i=nx0+1,nx1-1Ez(i,j)=Ez(i,j)+p2*(Hy(i,j)-Hy(i-1,j)-Hx(i,j)+Hx(i,j-1))end doend do!Ez(nx1/2,ny1/2)=Ez(nx1/2,ny1/2)-p2/delt*sin(w*(n+0.5)*deltt) call mur()! 连接边界!四条边do j=nyl1+1,nyl2-1Ez(i,j)=Ez(i,j)-p2*Hyi(i-1,j)end doi=nxl2do j=nyl1+1,nyl2-1Ez(i,j)=Ez(i,j)+p2*Hyi(i,j)end doj=nyl1do i=nxl1+1,nxl2-1Ez(i,j)=Ez(i,j)+p2*Hxi(i,j-1)end doj=nyl2do i=nxl1+1,nxl2-1Ez(i,j)=Ez(i,j)-p2*Hxi(i,j)end do!四个角i=nxl1j=nyl1Ez(i,j)=Ez(i,j)+p2*(Hxi(i,j-1)-Hyi(i-1,j))i=nxl2Ez(i,j)=Ez(i,j)+p2*(Hxi(i,j-1)+Hyi(i,j))j=nyl2Ez(i,j)=Ez(i,j)+p2*(-Hxi(i,j)+Hyi(i,j))i=nxl1Ez(i,j)=Ez(i,j)+p2*(-Hxi(i,j)-Hyi(i-1,j))end subroutine cal_E!/////////////////////////////////////////////////////////////////////////////////////////////////subroutine mur()use data_moduleimplicit noneinteger i,j!四条边i=nx0do j=ny0+1,ny1-1Ez(i,j)=Ez1(i+1,j)+p*(Ez(i+1,j)-Ez1(i,j))+q*(Hx(i+1,j)+Hx(i,j)-Hx(i+1,j-1)-Hx(i,j-1)) end doi=nx1-1do j=ny0+1,ny1-1Ez(i+1,j)=Ez1(i,j)+p*(Ez(i,j)-Ez1(i+1,j))+q*(Hx(i+1,j)+Hx(i,j)-Hx(i+1,j-1)-Hx(i,j-1)) end doj=ny0do i=nx0+1,nx1-1Ez(i,j)=Ez1(i,j+1)+p*(Ez(i,j+1)-Ez1(i,j))-q*(Hy(i,j+1)+Hy(i,j)-Hy(i-1,j+1)-Hy(i-1,j)) end doj=ny1-1do i=nx0+1,nx1-1Ez(i,j+1)=Ez1(i,j)+p*(Ez(i,j)-Ez1(i,j+1))-q*(Hy(i,j+1)+Hy(i,j)-Hy(i-1,j+1)-Hy(i-1,j)) end do!角点i=nx0j=ny0Ez(i,j)=Ez1(i+1,j+1)+s*(Ez(i+1,j+1)-Ez1(i,j))i=nx1Ez(i,j)=Ez1(i-1,j+1)+s*(Ez(i-1,j+1)-Ez1(i,j))j=ny1Ez(i,j)=Ez1(i-1,j-1)+s*(Ez(i-1,j-1)-Ez1(i,j))i=nx0Ez(i,j)=Ez1(i+1,j-1)+s*(Ez(i+1,j-1)-Ez1(i,j))end subroutine mur!/////////////////////////////////////////////////////////////////////////////////////////////////////// subroutine AP()use data_moduleimplicit noneinteger i,jreal mreal tmp1,tmp2,tmp3complex tmptmp=cmplx(Ez2(nx1/2,ny1/2),Ez4(nx1/2,ny1/2))do j=ny0,ny1do i=nx0,nx1Ez3(i,j)=cmplx(Ez2(i,j),Ez4(i,j))Ez0(i,j)=sqrt(Ez2(i,j)**2+Ez4(i,j)**2)m=Imag(Ez3(i,j)/tmp)/real(Ez3(i,j)/tmp)if(real(Ez3(i,j)/tmp)==0 .and. Imag(Ez3(i,j)/tmp)>0) thensita(i,j)=pi/2elseif(real(Ez3(i,j)/tmp)==0 .and. Imag(Ez3(i,j)/tmp)<0) thensita(i,j)=-pi/2elsesita(I,j)=atan(m)end ifif (abs(sita(i,j))<=1e-6) thensita(i,j)=0end ifend doend doend subroutine AP!///////////////////////////////////////////////////////////////////////////////////////////////////////。
平面波的特性与波长的计算
平面波的特性与波长的计算波动现象是自然界中普遍存在的一种物理现象，而波长是用来描述波的特性的定量指标之一。
在物理学中，平面波是一种理想化的波模型，其特性与波长之间存在一定的关系。
本文将探讨平面波的特性以及如何计算波长。
1. 平面波的特性平面波是一种在无限远处传播、波前呈平面的波动。
与球面波和柱面波不同，平面波的波阵面是一个平面，不存在波源和波纹的聚焦或发散现象。
平面波特性的关键在于其波矢和传播方向。
平面波的波矢k定义为波矢的模长与传播方向的乘积，即k = |k| * n。
其中，|k|表示波矢的大小，n表示波的传播方向，即单位矢量。
平面波的表达式可以表示为Ae^(ik·r)，其中A为振幅，k为波矢，r为位置矢量。
平面波的传播速度与波矢的大小成正比，即v = ω/|k|，其中v表示传播速度，ω表示角频率。
2. 波长的计算波长是用来描述波动的空间周期性的物理量，指在一个完整波动周期内，波动传播的距离。
对于平面波来说，波长的计算可以通过波矢的大小求倒数得到。
波矢的模长|k|与波长λ之间存在以下关系：|k| = 2π/λ。
由此可得，波长λ = 2π/|k|。
这个公式说明了波矢的大小与波长之间的倒数关系，即波矢越大，波长越短，反之亦然。
以光波为例，光是一种电磁波，其波长范围很宽，从红外到紫外都有不同波长的光。
对于可见光来说，其波长范围约为380nm到780nm之间。
3. 示例计算现在我们以一个具体的例子来计算平面波的波长。
假设我们有一个平面波，其振幅A为2，波矢k为1.5。
根据前述公式，波长λ = 2π/|k| = 2π/1.5 ≈ 4.188。
根据计算结果可知，该平面波的波长约为4.188个长度单位，可以是米、厘米、纳米等等，具体的单位需根据具体的物理系统来确定。
4. 总结平面波是一种在无限远处传播且波前呈平面的波动，其特性由波矢和传播方向决定。
波长是描述波动空间周期性的物理量，对于平面波来说，波长可以通过波矢的大小求倒数得到。
PML吸收边界条件下TM波的时域有限差分法分析
会带来较强的数值反射，选取线性分布可以改善结果。通常，选取 n=2 或 3，可进一步改进 吸收效果。
4. TM 波实例分析
4. 1 问题提出
二 维 TM 波 在 自 由 空 间 传 播 ， 采 用 PML 吸 收 边 界 ， 激 励 源 使 用 脉 冲 源
为了让实际的计算空间与实际的无限的物理空间等效须对有限空间的周围边界做特殊处理使得向边界面行进的波在边界处保持向外行进特点无反射现象即好像被吸收了一样并且不会使内部空间的场产生畸变具有这种功能的边界条件称为吸收边界条件吸收边界条件处理成为了fdtd的一大重点和发展方向pml是近年来发展起来的吸收边界条件45二维情况下pml吸收层区分为图3所示的边和角顶两种区域
∆t ε (m)
− σ (m) 2
+ σ (m)
=
1− 1+
σ (m)∆t 2ε (m) σ (m)∆t
(5)
∆t 2
2ε (m)
∆t
CB(m)
=
ε(m)
1 +σ(m)
=
ε(m) 1+σ(m)∆t
(6)
∆t 2
2ε(m)
CP(m)
=
µ(m)
∆t µ(m)
− σm(m) 2
+ σm(m)
1− σm(m)∆t
外行波在四边的 PML 界面上是无反射的。
图 3 完全匹配层的设置
在四个角顶区域，介质参数为 PML( σx , σmx , σy , σmy ),而于之相邻的边上的参数为 PML(σx ,σmx ,0,0,)和 PML(0,0，σy ,σmy )。根据上面的讨论，满足匹配条件时，在角顶边和
二维超声速普朗特-迈耶系数波流场的数值解
课程设计题目：二维超声速普朗特-迈耶系数波流场的数值解学院：飞行器工程学院专业名称：飞行器设计与工程班级学号： 07034211学生姓名：李桂平指导教师：刘勇二O一O年十一月第一部分1．物理问题简介：普朗特——迈耶稀疏波的解析解图-1中，超声速流围绕着一个尖的扩张角膨胀，无数个无限弱的马赫波组成了稀疏波，在尖角处展开成扇形。
扇形稀疏波的波头与來流方向的夹角1，而2是其波尾与下游方向的夹角。
1和2称为马赫角，定义为：和Ma1和Ma2分别为上下游的马赫数。
通过稀疏波的流动是等熵流动。
当流体通过稀疏波后，马赫数增加，压力、温度和密度降低；图-1中标明了这些变化趋势。
在中心稀疏波前的流动是均匀的，马赫数为Ma1，而且流动平行与波前的壁面。
稀疏波后的流动是均匀的，马赫数为Ma2，并且平行于下游的壁面。
在稀疏波内，流动参数光滑变化，流线弯曲，如图-1所示。
稀疏波内的流动是二维的，唯一的例外是折角的定点，它是一个奇点，壁面流线的方向在此处有一个突然的变化，而且此处的流动参数也是不连续的。
这个奇点对流动的数值解会产生影响。
给定超声速来流条件和拐角处的偏转角，下游参数是唯一确定的。
对于完全气体，在稀疏波后的流动有精确的解析解，下面给出这个解。
流过中心稀疏波的流动，其解析解取决于简单的关系式……1 式中，f是普朗特——迈耶函数；是流动偏转角。
对于完全气体，普朗特——迈耶函数是Ma和γ的函数，定义为……2 解析解中如下依次得到。
对给定的Ma1，从式(2)计算函数f1。
然后，对给定的偏转角θ，从式(1)得到f2。
用这样得到f2的值，通过求解式(2)求出Ma2。
式(2)是关于Ma2的隐式关系式，需要用试凑法求解。
波后的压力温度和密度都可以由等熵流动关系式： (3) (4)和状态方程： (5)得到。
借助是式(1)～式(5)，中心稀疏波后的流动就完全确定了。
2.问题的提法考虑图-2所示的物理平面。
来流马赫数为2，来流的压力、密度、温度分别为：1.01x105N/m2、1.23kg/m3、286.1K。
平面波法计算能带的基本原理
平面波法计算能带的基本原理能带理论是描述固体中电子能量分布的理论，它对于理解固体的电子性质和材料的导电、绝缘以及半导体特性具有重要意义。
平面波法是计算能带的一种常用方法，本文将介绍平面波法计算能带的基本原理。
平面波法是一种基于量子力学的计算方法，它假设电子波函数能够用平面波表示。
平面波是一种特殊的波函数形式，它具有平直的波前和无限延伸的波长。
根据波动方程，平面波的波函数可以表示为Ae^(ik•r)的形式，其中A是振幅，k是波矢，r是位置矢量。
在能带计算中，我们通常采用周期性边界条件，即假设固体是无限重复的晶体结构。
根据布洛赫定理，电子的波函数可以表示为平面波和周期函数的乘积形式。
这样，我们可以通过计算晶体中的平面波的行为来研究电子的能量分布。
平面波法的基本思路是，首先选择一组适当的平面波作为基函数，然后将它们用作电子波函数的展开系数。
通常，我们选择的平面波基函数是满足布洛赫定理的，即满足晶格平移对称性的平面波。
这样，我们就可以通过求解薛定谔方程来确定电子的能量和波函数。
在平面波法中，薛定谔方程可以写为哈密顿算符作用在波函数上等于能量的形式。
由于平面波是哈密顿算符的本征函数，所以我们可以将波函数用平面波展开，并代入薛定谔方程中，从而得到一个简化的形式。
通过对展开系数的求解，我们可以得到电子的能量和波函数。
在实际计算中，我们通常采用周期性边界条件，将晶体划分为一个个无限重复的晶胞。
每个晶胞内的电子波函数可以用平面波展开，而不同晶胞之间的平面波之间存在相位差。
通过求解薛定谔方程，我们可以得到不同晶胞中的平面波的展开系数，从而确定电子的能量和波函数。
平面波法计算能带的基本原理就是通过选择适当的平面波基函数，并将其用作电子波函数的展开系数，从而求解薛定谔方程得到能量和波函数。
这种方法在实际应用中已经取得了很多成功，并对我们理解固体的电子性质提供了重要的帮助。
平面波法是一种常用的计算能带的方法，它基于平面波的特点和布洛赫定理，通过求解薛定谔方程来确定电子的能量和波函数。
二维核磁共振谱全解
O
2 C1
C
O CH2
CH3
H
5
6
5
46
1/2JCH F1 (HZ)
-100
0
100
δ C 150
100
50 F2
反式丙烯酸酯的13C,1H J分解2D谱。 F2轴为13C的化学位移。F1轴为1JCH偶合的多重峰。
21
10
7
9
5,6
16
2
8
O
3
5 4
由于DEPT等测定碳 原子级数的方法能代 替异核J 谱，且检测 速度快，操作方便， 因此异核J 谱较少应 用。
15
在同核1H,1H －2D J分解谱中，被测定的核为 1H核。
16
3
2
5
4
6
H
3
C CH3
4
O
C1
5
6
C
2
O CH2 CH3
H
F1
-10
J0 (HZ)
10
δH 6
4
2 F2
反式丙烯酸乙酯的1H－1H同核二维J分解谱
17
18
化学等价磁 不等价
化学等价磁 不等价
2位的a-H和e-H相互偶合， 裂分为双峰后进一步被3位 氢裂分而显示六重峰。同 样，3-H的a、e质子首先裂 分为双峰，再被邻位2-H分 裂为六重峰，和4-H质子进 一步裂分为十八重峰，但 由于有些峰相互重叠，在 2DJ分解谱中只显示9个点。 5位甲基没有受到偶合，因 此只在F1=0轴上显示单峰。
C6’,H6’ C6,H6
C1’,H1’ 36
F2
7,4
CH3CO
37
（二）1H检测的异核多量子相关谱（HMQC）
平面波公式
平面波公式平面波是声学中的一种波，它的特征是在某一平面上的传播，它的受力形式为前后两个振幅在某一平面上的波运动，此外还伴有一些横向的分部。
由此可以看出，平面波的定义中平面的概念是十分重要的，在物理学中，这一概念是指某一方向上能够看出常数的物理参量，就像抛物线形状的圆弧一样，当抛物线旋转着某个实际情况时，它形成的便就是一个平面。
在物理学中，任何一种波都有它的波公式，这是用来描述它特有的特性的。
其中，平面波公式是推导平面波波运动的基础，也是物理学中最为基本的公式。
平面波公式有时也和平面波方程有所区别，这是因为有时候平面波公式只描述了冲击的振幅及传播的距离，而平面波方程则描述了传播过程中的变化情况，它可以利用不同的变量来表达平面波的衰减情况。
平面波的运动也可以看作是一种波的传播过程，比如水面上的波，在水面上传播的过程其实就是平面波的运动。
平面波公式的形式如下：A = Acos(ωt - ks)其中A是波的振幅，A是初始振幅；ω是频率；t是时间；k是波动数；s是波行进的距离。
平面波的传播分为三种类型：静态传播、动态传播以及不完全传播。
静态传播是指输入的振幅以恒定的速度在一定距离内不变地传播；动态传播则指通过某种形式来改变传播的距离；不完全传播是指一开始在空间中曲率明显，经过一定时间后曲率逐渐减弱，变得平坦，最后又重新曲线上升。
这三种传播形式都可以看作是平面波传播的一种状态，由此可以推导出专有的平面波公式。
平面波公式是物理学中十分重要的一个公式，它可以用来描述静态的、动态的以及不完全传播的平面波运动情况。
它是物理学与声学领域中对平面波特性的理论推导和实验研究的基础。
因此，平面波公式在物理学领域中非常重要，学习平面波公式也是物理学中一项基础而关键的知识点。
总之，平面波公式是一种十分重要的物理学公式，平面波的特性及运动过程通过它的推导得到深入的分析，这也是平面波物理学研究的一个重要组成部分。
平面波的波动方程
探究平面波的波动方程
平面波是物理学中一种基本的波动形式，它在自然界中有着广泛的应用。
其中平面波的波动方程是研究平面波特性的重要基础。
本文将详细讲解平面波的波动方程，并深入探讨其物理意义。
平面波是指波的振动方向垂直于波传播方向的波动形式。
在数学上，平面波可以用以下公式表示：
A = A0sin(kx - ωt + φ)
其中A是波的振幅，k是波数，x是波的传播方向，t是时间，ω是角频率，φ是初相位。
根据波的定义，平面波的波动方程可以表示为：
∂2A/∂x2 = (1/v2)∂2A/∂t2
其中v是波的速度。
将平面波的表达式代入波动方程中，可得：
-k2A0sin(kx - ωt + φ) = (1/v2)(ω2A0sin(kx - ωt + φ)/∂t2)
整理后可得到平面波的波动方程：
∂2A/∂x2 + k2A = (ω2/v2)A
通过对波动方程的分析，我们可以得到以下结论：
1. 平面波的波动方程是二阶偏微分方程。
2. 平面波的波数k与角频率ω满足波速公式v = ω/k。
3. 平面波的波动方程可以描述出平面波的传播和振动状态。
通过以上的分析，我们可以进一步探讨平面波的物理意义。
首先，平面波的波动方程告诉我们，平面波的传播速度与波长和频率有关。
其次，平面波的波动方程将平面波的传播和振动状态联系在了一起，
揭示了平面波的本质特性。
总之，平面波的波动方程是研究平面波特性的重要基础。
通过对
波动方程的分析，我们可以深入探讨平面波的物理意义，并为平面波
的应用提供理论基础。
二维波动方程的解法
二维波动方程的解法二维波动方程是一种常见的偏微分方程形式，它描述了二维空间中的波动现象，如水波、电磁波等。
解决二维波动方程，对于理解复杂物理问题有很大帮助。
本文将介绍二维波动方程的解法，包括分离变量法、变换方法等。
二维波动方程是指一个偏微分方程，如下所示：$$\frac{\partial^2 u}{\partial x^2}+\frac{\partial^2 u}{\partialy^2}=\frac{1}{v^2}\frac{\partial^2 u}{\partial t^2}$$其中，$u$是波函数，$v$是波速，$t$是时间，$x$和$y$是二维空间的坐标。
我们假设边界条件为$u(x,0)=f(x)$和$u(x,y)|_{\partial D}=0$。
其中，$f(x)$是已知函数，$\partial D$是区域D的边界。
接下来，我们将介绍几种解决二维波动方程的方法。
分离变量法在这种方法中，我们假设解的形式为$u(x,y,t)=X(x)Y(y)T(t)$。
将$X(x)$，$Y(y)$和$T(t)$代入原方程，可以得到：$$\frac{X''(x)}{X(x)}+\frac{Y''(y)}{Y(y)}=\frac{1}{v^2}\frac{T''( t)}{T(t)}$$由于每个因子只与对应的变量相关，所以我们可以令每个因子等于一个常数，即$$\frac{X''(x)}{X(x)}=-k^2_x$$$$\frac{Y''(y)}{Y(y)}=-k^2_y$$$$\frac{T''(t)}{T(t)}=-v^2k^2$$这样，我们就得到了三个常微分方程，它们的解分别为：$$X(x)=A_x\sin(k_xx)+B_x\cos(k_xx)$$$$Y(y)=A_y\sin(k_yy)+B _y\cos(k_yy)$$$$T(t)=C\sin(vkx)$$其中，$A_x$，$B_x$，$A_y$，$B_y$和$C$是常数，$k_x$，$k_y$和$k$是相应常数。
分层背景2维FDTD中斜入射平面波的引入
图 １ 分 层背 景 ２维 Ｆ Ｄ 计 ， ＡＢ 为 内框 （ Ｂ 为连 接 边界 ， 总 Ｃ ） 即 场 （ ）散 射场 （ Ｆ 边 界 。在 图 中 ｉ ｌｉ ｊ Ｊ ６ 处 用 修 正 的 １维 Ｆ ＴＤ加 入 ２维 平 面入 射 波 。 ＴＦ 一 Ｓ） ＝ｂ，＝ｂ 和 ＝ｂ ， 一 ｄ Ｄ 沿 着 方 向的 １维 Ｍａ ｗｅ 方 程仅有 Ｅ ｘ ｌ ｌ 和 Ｈ 分量 ， 式 （ ） 而 １ 中含 有 Ｅ ， Ｈ 和 Ｈ， 。为 便 于加 入 入 射波 的 １
fig5犈狕distributionofdamagedcomponent图5带缺陷元件tm模式不同时刻犈狕近场分布fig6犎狕distributionofdamagedcomponent图6带缺陷元件te模式不同时刻犎狕近场分布fig7comparisonofnearzonescatterfieldbetweendamagedandpurifiedlayeredcomponent图7带缺陷元件和纯层状元件近区散射场比较表面覆盖薄膜的熔石英模型带气泡缺陷时tm模式下电场犈狕和te模式下磁场犎狕在250狋和350狋时刻30斜入射近场分布快照分别如图5和图6中a和b所示
（ ａ １）
（ ｂ １） （ ｃ １）
●
● 一
＿
Ｉ －—— —一 － －－ ——
Ｉ ＣＢ －— ・ － － －— — － －一
Ｔ Ｅ模式 频域 Ｍａ ｗ ｌ方 程为 ｘ ｅｌ
３ ／ ｙ— ｊ Ｈ ａ 旺Ｅ ３ ／ ｘ 一一 ｊ Ｈ＝ ａ 旺Ｅ ａ ／ ｘ一 ３ ／ ｙ一一 ｊ， Ｅ ａ Ｅ ａ ｃ Ｈ ｏ ｕ （ ａ ２） （ ｂ ２） （ ｃ ２）
１ 分 层 背 景 空 间 斜 入 射 平 面 波 引 入 理 论
关于平面波波动方程的推导
关于平面波波动方程的推导
一．平面波概念
1、平面波是指在某一领域内，波动运动的模式被限制在一个面上，即只具有一个
指向贯穿全域的波动方向的波动模式的称为平面波。
2、当一种介质中的介质振动，且振动波断面是平面时，就可以称为是平面波。
二．平面波波动方程
1、平面波波动方程即它的物理意义，指的是平面内找到各点上时间波速度的模式，当所有的点满足这个模式时，说明振动波断面是平面波。
2、它是一个非常重要的方程，用来描述特定领域内某种波运动模式的变化情况。
它也是一个包含了有关波速、张力系数等重要参数，受这些参数影响影响其形态及特性变化的方程。
三．平面波波动方程的推导
1、平面波波动方程是基于方程即事实来推导的，其主要思路是：先首先考虑可以
根据实验拟合出的线性模型来简化分析问题，再根据它的有限差值表达式，推导出平面波波动方程，就是Lamb导热方程。
2、根据常见维度理论中对Lamb导热方程的推导，由Lamb导热方程可推导出平
面波波动方程。
推导过程中，首先要分析定义波的空间偏微分方程中的压强U的
变化，其次，考虑在单位时间内产生的波的波动性质。
使用创新的振动模型，把压强U区分为正向和负向波，考虑受拉格朗日中值定理的影响，把正向和负向的压
强U的变化反应在波动方程中，通过不断的积分和相关变换，就可以得到平面波
波动方程。
3、最后，根据经验和实验数据，推导出量化参数，也就可以获得局部地区的平面
波波动方程了。
高效分裂平面波FDTD方法研究
Ｋ ｙｗ ｒ ｓ ｉｉ —ｉ ｅｅ ｃ ｉ —ｏ ｉ （ ＤＴ ；ｓｐ ｒｔ ｎｏ ｏａ－ ｅｄ ｓａｔｒ ｇｆ ；ｓ ｌｔ ｇ ｅ ｏｄ ：ｆ ｔ ｄｆ ｒｎ ｅｔ ｎｅ ｆ ｍｅｄ ｍａ Ｆ Ｄ） ｅ ａａｉ ｆ ｔｌ ｌ／ ｃｔｅｉ － ｄｄ ｐｉ ｉ ｎ ｏ ｔ 一 ｉ ｆ ｎ ｉ ｔｎ
第３ ３卷 第 ４期
２１ ０ ０年 ４月
合肥工业大学学报（ 自然科 学版）
Ｊ ＯＵＲＮＡＩ Ｉ ０Ｆ ＨＥＦＥ Ｉ ＩＵＮＩ ＶＥＲＳ ＴＹ Ｉ ＯＦ ＴＥＣＨＮＯＬ ， ＯＧＹ
Ｖｏ ． ３ Ｎｏ ４ １３ ．
Ａｐ ． ２ ０ ｒ ０１
高效 分 裂 平 面波 ｕ ｉｇ ｐａ ｅ ｗａ ｅ ｓ ｕ ｃ ｓ ｂ ｓ ｄ ｏ ｈ ｏ ｃｐ ｆｓ ｌ — ｅｄ ｎ ｍ ｅｙ ｔ ｅ ｓ ｌ ｔｎ ｔ ｏ ｏ ｎ ｒ ｄ ｃｎ ｌｎ ｖ ｏ ｒｅ ａ ｅ ｎ ｔ ｅ ｃ ｎ ｅ ｔｏ ｐ ｉ ｆ ｌ ， ａ ｌ ｈ ｐ ｉ ｉｇ ｖｉ ｔ ｐ ａ ｅ ｗａ ｅＦＤＴＤ （ Ｐ Ｆ ｌｎ － ｖ 一 Ｓ － ＤＴＤ）ｍｅｈ ｄ Ｔｈ ｓｍ ｅｈ ｄｆ ｎ ａ ｅ ｔｌｙｅｉ ｉａ ｅ ｈ ｒ ｏ ａ ｓ ｄｂ ｅ ｔｏ ． ｉ ｔ ｏ ｕ ｄ ｍ ｎ ａｌ ｌ ｎ ｔ ｓｔ ｅｅ ｒ ｒｃ ｕ ｅ ｙｔ ｍ ｈ
文献标识码 ： Ａ 文章 编 号 ：Ｏ ＆５ ６ （０ ０ ０ 一６４０ ｌ Ｏ 一０ ０ ２ １ ）４０ ０ —５ －
关键词 ： 时域有限差分 ； 总场／ 散射场分离 ； 分裂平面波 Ｆ Ｔ Ｄ Ｄ法
中 图分 类 号 ： Ｎ８ １ Ｔ ０
Ｒｅ ｅ ｒ ｈ ｏ ｎ ｅ ｆｃ ｅ ｔｓ ｌｔ ｉ ｇ ｐｌ ｎ － ｖ ｓ ａ ｃ ｎ ａ ｆ ｉ ｉ ｎ ｐ ｉｔｎ ａ ｅ ｗａ ｅ ＦＤＴＤ ｅ ｈ ｍ ｔ
二维普朗克公式推导
二维普朗克公式推导
二维普朗克公式是量子力学中一个重要的公式，它描述了光子在二维平面上的能量分布。
以下是其推导过程：
1. 首先，我们需要了解普朗克公式，它描述了单个光子在三维空间中的能量分布。
公式为：E=hν，其中E是光子的能量，h是普朗克常数，ν是光子的频率。
2. 在二维平面上，光子的能量分布可以通过将普朗克公式中的能量和频率分别除以面积来得到。
设二维平面的面积为A，则有：E/A=hν/A。
3. 进一步简化得到：E/A=hv/a，其中a为平面内的一个单位长度。
4. 由此可得，在二维平面上，光子的能量分布为：E=ahv。
这就是二维普朗克公式的形式。
通过以上步骤，我们可以推导出二维普朗克公式。
这个公式在量子力学和光学中有重要的应用，可以帮助我们理解光子在二维平面上的行为和性质。
二维波动方程初值问题
二维波动方程的初值问题可以描述为：
1.波的传播方向与坐标轴方向平行，所以方程中只含有二阶偏导数，不含有交叉偏导数。
2.波的传播具有周期性，即波前波后不断重复，因此方程具有周期性解。
3.波的传播具有散射性，即波在传播过程中会遇到障碍物或介质不均匀分布等情况，导致波的散射。
4.波的传播具有吸收性，即波在传播过程中会因为介质吸收能量而逐渐减弱。
对于初值问题，需要给出初始时刻波场的状态，即初始条件。
初始条件通常由波源的位置和强度决定，可以表示为：
1.初始位置：波源在初始时刻的位置。
2.初始速度：波源在初始时刻的速度。
3.初始位移：波源在初始时刻的位移。
4.初始能级：波源在初始时刻的能量级别。
根据不同的初始条件，可以求解出不同的二维波动方程的解。
求解方法通常采用分离变量法、傅里叶变换、拉普拉斯变换等数学方法。
工程电磁场 第7章 二维泊松方程的有限元法
式中 wi 为权函数， w1, w2 , , wk , 为权函数序列，
权函数之间要求线性无关。 权函数的不同选择导致不同的近似方法。
2、几种加权余量法
（1）配点法
在求解区域中选取 n 个点 P1, P2 , , Pn ， 让方程的余量在这 n 个点上为零。
即选权函数为
K (e) jm
K (e) mm
Kim
ui
bi
K mj
u
j
b
j
Kmm um bm
b(e) i
f (e)
e
N (e) i
即权函数为
wi
R ci
（ i 1,2, , n ）
（4）伽辽金法
选取权函数序列与基函数序列相同。
wi ui
ui Rd 0
m
ui L u jd ui fd
j 1
（ i 1,2, , n ）
在上述几种加权余量法中， 伽辽金法应用最广泛。 有限元法基于伽辽金法
场域离散
二维问题常使用三角形单元离散，便于处理复杂的场域形 状，容易实现。
单元：互不重叠，覆盖全部场域；每个单元内介质是 单一、均匀的。
节点：网格的交点，待求变量的设置点。 该步骤需要记录的信息： 节点编号、节点坐标 节点属性(激励源、是否边界等) 单元编号 单元节点编号 单元介质
三角形单元内的基函数 设三角形三个顶点处待求函数值 分别为u1, u2, u3。如果单元足够小， 可以采用线性近似，将单元内任 意p点的u(x,y)表示为
二维波动方程推导
0； y
H E ， t E H ， t
（1）
（2）
H 0， E 0，
由旋度公式可以得：
（3） （4）
E E E E E y i ( x z ) j y k z z x x H y H x H z H y H i ( )j k z z x x
同理可得以下亥姆霍兹方程：
2 Ex 1 2 Ex 2 Ex ( 2 ) t 2 z 2 x
2 Ey t 2
2 2 2H x 2H z 1 Ey Ey ( ) zt xt z 2 x 2
2 Ez 1 2 Ez 2 E z ( 2 ) t 2 x 2 z 2 H x 1 2 H x 2 H x ( ) t 2 z 2 x 2 2 H z 1 2 H z 2 H z ( ) t 2 x2 z 2
二维波动方程推导波动方程的推导三维波动方程的推导一维波动方程的推导波动方程推导二维波动方程波动方程波动方程表达式平面简谐波的波动方程一维波动方程
无源无损耗空间二维亥姆霍兹方程的表达形式
假设场分布在 x ， z 方向不均匀变化，在 y 方向均匀变化，即 无源无损耗空间，满足： v 0 ， 0 无外加电流的情况下，Maxwell 方程为：
（8）
（9）
E z t
（10）
由散度方程（3）（4）得： 、
E x E z x z H x H z x z
（11） （12）
（5）式两边分别对 t 求偏导：
2H y 2 Ez 2 Ex xt zt t 2 2H y z 2 2H y x 2
第四章-平面波
第四章 平面波本章从麦克斯韦方程及物质的本构关系出发，研究在均匀介质中平面波的传播及其主要特征。
首先讨论线性、均匀、各向同性介质中均匀平面波的传播，再推广到各向异性介质中的情况。
比平面波更复杂的电磁波也可用平面波展开，本章对此也作了讨论。
最后讨论平面波传播的传输线模型，为以后用传输线模型求解复杂的场问题打下基础。
4.1得出电场强度E 与磁场强度H 满足的波方程，4.2从波方程得到简单介质中的平面波解，4.3、4.4讨论平面波的极化特性以及平面波在有耗介质中的传播，4.5介绍色散与群速的基本概念，4.6与4.7分别研究电各向异性介质和磁各向异性介质中平面波的传播特征。
4.8讨论髙斯波束的平面波展开，4.9证明电磁波沿某一方向传播可与特定参数传输线上电压、电流波的传播等效，即电磁波传播的传输线模型。
4.1 波方程3.4已分析过，麦克斯韦方程组中两个旋度方程是独立的。
在两个旋度方程中电场强度E 与磁场强度H 耦合在一起。
从解方程角度看，先要将E 跟H “去耦”，即从两个旋度方程消去H （或E ），然后得到只关于E （或H ）的方程。
本节讨论无源、简单介质中麦克斯韦方程的解，所谓无源，就是指所研究的区域内不存在产生电磁场的源J 与ρv 。
对于简单介质，ε、μ是常量。
在这种特定情况下，将物质的本构关系（3.4.1）、（3.4.2）代入麦克斯韦方程（3.2.8）~（3.2.11），得到 ∇⨯E =–j ωμH （4.1.1） ∇⨯H = j ωεE （4.1.2） ∇⋅E = 0 （4.1.3） ∇⋅H = 0 （4.1.4） 式（4.1.1）、（4.1.2）两个方程中，只有E 和H 两个独立的场量，但E 和H 耦合在一起。
为了从这两个方程得到只关于E 或H 的方程，对式（4.1.1）取旋度，并将式（4.1.2）代入，得到 ()()()E E H E μεωωεωμωμ2=-=⨯∇-=⨯∇⨯∇j j j利用恒等关系()()E E E 2∇-⋅∇∇=⨯∇⨯∇，而根据式（4.1.3），0=⋅∇E ，所以上式成为022=+∇E E μεω（4.1.5）同样对式（4.1.2）取旋度，将式（4.1.1）代入，并利用式（4.1.4）及上面的矢量运算恒等关系，得到022=+∇H H μεω（4.1.6）式（4.1.5）、（4.1.6）可合并写成 ()022=⎩⎨⎧+∇HEk（4.1.7） 式中μεω22=k（4.1.8）在自由空间或真空中，μ = μ0，ε = ε0，k 记作k 000220εμω=k（4.1.9）式（4.1.5）、（4.1.6）或（4.1.7）叫做无源简单介质中的波方程，在这个方程中E 跟H 不再耦合在一起。
二维伊藤公式推导公式
二维伊藤公式推导公式二维伊藤公式是随机微积分中的基本公式之一，它描述了二维随机过程的演化规律。
具体而言，对于由$X_t$和$Y_t$组成的二维随机过程，其伊藤微分可以表示为$dX_t=a(X_t,Y_t)dt+b(X_t,Y_t)dW_t^1$和$dY_t=c(X_t,Y_t)dt+dW_t^2$的形式，其中$a$、$b$、$c$是确定性函数，$W_t^1$和$W_t^2$是标准布朗运动。
那么我们来推导二维伊藤公式的具体形式。
首先，我们考虑由$X_t$和$Y_t$组成的二维随机过程$f(X_t,Y_t)$的伊藤微分。
根据伊藤公式，其微分可以表示为：$$df=frac{partial f}{partial t}dt+frac{partial f}{partial X_t}dX_t+frac{partial f}{partialY_t}dY_t+frac{1}{2}frac{partial^2f}{partialX_t^2}(dX_t)^2+frac{1}{2}frac{partial^2f}{partialY_t^2}(dY_t)^2+frac{partial^2f}{partial X_tpartialY_t}dX_tdY_t$$将$dX_t$和$dY_t$的表达式带入上式中，得到：$$begin{aligned}df&=frac{partial f}{partial t}dt+frac{partial f}{partialX_t}(a(X_t,Y_t)dt+b(X_t,Y_t)dW_t^1)+frac{partial f}{partial Y_t}(c(X_t,Y_t)dt+dW_t^2)&+frac{1}{2}frac{partial^2f}{partialX_t^2}b^2(X_t,Y_t)dt+frac{1}{2}frac{partial^2f}{partialY_t^2}d^2(X_t,Y_t)dt+frac{partial^2f}{partial X_tpartialY_t}b(X_t,Y_t)dW_t^1dW_t^2end{aligned}$$根据布朗运动的性质，有$dW_t^1dW_t^2=rho dt$，其中$rho$是布朗运动的相关系数。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
! TM波FDTD讨论平面波源的加入module data_moduleimplicit noneinteger,parameter::nx0=0,nx1=360,ny0=0,ny1=360,nz0=-100,nz1=1200integer,parameter::nxl1=nx0+80,nxl2=nx1-80,nyl1=ny0+80,nyl2=ny1-80 !连接边界real,parameter::f=2.0e8,c=3.0e8,delt=0.0177,deltt=delt/6.0e8,eps0=8.85e-12,miu0=1.2566e-6,pi= 3.14159real,parameter::w=2*pi*f,s=-0.477369real,parameter::p=-1.0/3.0,q=-miu0*c/6,r=-miu0*c/2,p1=1/(2*miu0*c),p2=1/(2*eps0*c)real,parameter::tal=2e-9,t0=0.8*tal,fai=pi/3.0real cez,chx,chyinteger,parameter::nt=2000,m0=200integer ncomplex Ez3(nx0:nx1,ny0:ny1)real Ez4(nx0:nx1,ny0:ny1),Ez2(nx0:nx1,ny0:ny1) !记录幅值提取时的实部和虚部real sita(nx0:nx1,ny0:ny1),Ez0(nx0:nx1,ny0:ny1) !记录幅值和相位real Ez(nx0:nx1,ny0:ny1),Hx(nx0:nx1,ny0:ny1),Hy(nx0:nx1,ny0:ny1),Ez1(nx0:nx1,ny0:ny1) real Ei0(nz0:nz1),Hi0(nz0:nz1),Ei1(nz0:nz1)real Ezi(nx0:nx1,ny0:ny1),Hxi(nx0:nx1,ny0:ny1),Hyi(nx0:nx1,ny0:ny1)end module data_module!/////////////////////////////////////////////////////////////////////////////////////////////////subroutine inc()use data_moduleimplicit noneinteger i,j,kreal t,dt=n*delttEi1=Ei0do k=nz0,nz1-1Hi0(k)=Hi0(k)-p1*(Ei1(k+1)-Ei1(k))end do!Ezido i=nxl1,nxl2do j=nyl1,nyl2d=real(i-nxl1)*cos(fai)+real(j-nyl1)*sin(fai)Ezi(i,j)=(d-int(d))*Ei0(m0+int(d)+1)+(1-(d-int(d)))*Ei0(m0+int(d))end doend dodo k=nz0+1,nz1-1Ei0(k)=Ei1(k)-p2*(Hi0(k)-Hi0(k-1)) ! 入射波的场量end doEi0(m0-2)=sin(w*t) !exp(-(4*pi*(t-t0)**2/Tal/tal))k=nz1Ei0(k)=Ei1(k-1)-1.0/3.0*(Ei0(k-1)-Ei1(k)) !入射波右吸收边界K=nz0Ei0(k)=Ei1(k+1)-1.0/3.0*(Ei0(k+1)-Ei1(k)) !入射波左吸收边界!Hxij=nyl1do i=nxl1,nxl2d=real(i-nxl1)*cos(fai)+(j-0.5-nyl1)*sin(fai)+0.5Hxi(i,j-1)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*sin(fai) end doj=nyl2do i=nxl1,nxl2d=real(i-nxl1)*cos(fai)+(j+0.5-nyl1)*sin(fai)+0.5Hxi(i,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*sin(fai) end do!Hyii=nxl1do j=nyl1,nyl2d=(i-0.5-nyl1)*cos(fai)+real(j-nyl1)*sin(fai)+0.5Hyi(i-1,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*(-cos(fai)) enddoi=nxl2do j=nyl1,nyl2d=(i+0.5-nyl1)*cos(fai)+real(j-nyl1)*sin(fai)+0.5Hyi(i,j)=((d-int(d))*Hi0(m0-1+int(d)+1)+(1-(d-int(d)))*Hi0(m0-1+int(d)))*(-cos(fai)) enddowrite(13,*)n,Ei0(0),Ei0(20),Ei0(30),Ei0(100),sin(w*t)end subroutine inc!/////////////////////////////////////////////////////////////////////////////////////////////////program mainuse data_moduleimplicit noneinteger i,j,kreal topen(11,file='A.txt')open(12,file='fai.txt')open(13,file='Ez.txt',recl=300)Ez=0;Hx=0;Hy=0;Ez1=0Ei0=0;Hi0=0Ez4=0;Ez2=0;Ez3=cmplx(0.0,0.0)cez=1.0;chx=sin(fai);chy=-cos(fai)do n=0,ntcall inc()call FDTD()print *,'step',n,Ez(nx1/2,ny1/2)!write(13,*)n,Ei0(0),Ei0(20),Ei0(30),Ei0(100),sin(w*t)Ez1=Ezend doEz4=Ezdo n=nt+1,nt+42call inc()call FDTD()Ez1=Ezend doEz2=Ezcall AP()do j=ny0,ny1do i=nx0,nx1write(11,*)i,j,Ez0(i,j)write(12,*)i,j,sita(i,j)end doend doclose(11)close(12)close(13)end program main!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine FDTD()use data_moduleimplicit noneinteger i,jcall cal_H()call cal_E()end subroutine FDTD!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine cal_H()use data_moduleimplicit noneinteger i,jreal,external::Eido j=ny0,ny1-1Hx(i,j)=Hx(i,j)-p1*(Ez(i,j+1)-Ez(i,j))end doend dodo j=ny0,ny1do i=nx0,nx1-1Hy(i,j)=Hy(i,j)+p1*(Ez(i+1,j)-Ez(i,j))end doend do!连接边界i=nxl1-1do j=nyl1,nyl2Hy(i,j)=Hy(i,j)-p1*Ezi(i+1,j)end doi=nxl2do j=nyl1,nyl2Hy(i,j)=Hy(i,j)+p1*Ezi(i,j)end doj=nyl1-1do i=nxl1,nxl2Hx(i,j)=Hx(i,j)+p1*Ezi(i,j+1)end doj=nyl2do i=nxl1,nxl2Hx(i,j)=Hx(i,j)-p1*Ezi(i,j)end doend subroutine cal_H!///////////////////////////////////////////////////////////////////////////////////////////////// subroutine cal_E()use data_moduleimplicit noneinteger i,jreal,external::Hido j=ny0+1,ny1-1do i=nx0+1,nx1-1Ez(i,j)=Ez(i,j)+p2*(Hy(i,j)-Hy(i-1,j)-Hx(i,j)+Hx(i,j-1))end doend do!Ez(nx1/2,ny1/2)=Ez(nx1/2,ny1/2)-p2/delt*sin(w*(n+0.5)*deltt) call mur()! 连接边界!四条边do j=nyl1+1,nyl2-1Ez(i,j)=Ez(i,j)-p2*Hyi(i-1,j)end doi=nxl2do j=nyl1+1,nyl2-1Ez(i,j)=Ez(i,j)+p2*Hyi(i,j)end doj=nyl1do i=nxl1+1,nxl2-1Ez(i,j)=Ez(i,j)+p2*Hxi(i,j-1)end doj=nyl2do i=nxl1+1,nxl2-1Ez(i,j)=Ez(i,j)-p2*Hxi(i,j)end do!四个角i=nxl1j=nyl1Ez(i,j)=Ez(i,j)+p2*(Hxi(i,j-1)-Hyi(i-1,j))i=nxl2Ez(i,j)=Ez(i,j)+p2*(Hxi(i,j-1)+Hyi(i,j))j=nyl2Ez(i,j)=Ez(i,j)+p2*(-Hxi(i,j)+Hyi(i,j))i=nxl1Ez(i,j)=Ez(i,j)+p2*(-Hxi(i,j)-Hyi(i-1,j))end subroutine cal_E!/////////////////////////////////////////////////////////////////////////////////////////////////subroutine mur()use data_moduleimplicit noneinteger i,j!四条边i=nx0do j=ny0+1,ny1-1Ez(i,j)=Ez1(i+1,j)+p*(Ez(i+1,j)-Ez1(i,j))+q*(Hx(i+1,j)+Hx(i,j)-Hx(i+1,j-1)-Hx(i,j-1)) end doi=nx1-1do j=ny0+1,ny1-1Ez(i+1,j)=Ez1(i,j)+p*(Ez(i,j)-Ez1(i+1,j))+q*(Hx(i+1,j)+Hx(i,j)-Hx(i+1,j-1)-Hx(i,j-1)) end doj=ny0do i=nx0+1,nx1-1Ez(i,j)=Ez1(i,j+1)+p*(Ez(i,j+1)-Ez1(i,j))-q*(Hy(i,j+1)+Hy(i,j)-Hy(i-1,j+1)-Hy(i-1,j)) end doj=ny1-1do i=nx0+1,nx1-1Ez(i,j+1)=Ez1(i,j)+p*(Ez(i,j)-Ez1(i,j+1))-q*(Hy(i,j+1)+Hy(i,j)-Hy(i-1,j+1)-Hy(i-1,j)) end do!角点i=nx0j=ny0Ez(i,j)=Ez1(i+1,j+1)+s*(Ez(i+1,j+1)-Ez1(i,j))i=nx1Ez(i,j)=Ez1(i-1,j+1)+s*(Ez(i-1,j+1)-Ez1(i,j))j=ny1Ez(i,j)=Ez1(i-1,j-1)+s*(Ez(i-1,j-1)-Ez1(i,j))i=nx0Ez(i,j)=Ez1(i+1,j-1)+s*(Ez(i+1,j-1)-Ez1(i,j))end subroutine mur!/////////////////////////////////////////////////////////////////////////////////////////////////////// subroutine AP()use data_moduleimplicit noneinteger i,jreal mreal tmp1,tmp2,tmp3complex tmptmp=cmplx(Ez2(nx1/2,ny1/2),Ez4(nx1/2,ny1/2))do j=ny0,ny1do i=nx0,nx1Ez3(i,j)=cmplx(Ez2(i,j),Ez4(i,j))Ez0(i,j)=sqrt(Ez2(i,j)**2+Ez4(i,j)**2)m=Imag(Ez3(i,j)/tmp)/real(Ez3(i,j)/tmp)if(real(Ez3(i,j)/tmp)==0 .and. Imag(Ez3(i,j)/tmp)>0) thensita(i,j)=pi/2elseif(real(Ez3(i,j)/tmp)==0 .and. Imag(Ez3(i,j)/tmp)<0) thensita(i,j)=-pi/2elsesita(I,j)=atan(m)end ifif (abs(sita(i,j))<=1e-6) thensita(i,j)=0end ifend doend doend subroutine AP!///////////////////////////////////////////////////////////////////////////////////////////////////////。