发帖回复

12 3

5197阅读
23回复

[已解决]关于平面波源的加入问题

离线wq_463

UID ：20925

注册：2008-11-06
登录：2021-04-22
发帖：227
等级：仿真三级

0楼发表于: 2009-03-23 09:25:06

我再做一个盒体孔缝的耦合模型,现在想在总场中加入平面波源,应该如何加入(我用MATLAB编写程序,不知道大家的手头谁有类似的平面波加入程序,可以供我参考)?是否 ..

未注册仅能浏览部分内容，查看全部内容及附件请先登录或注册

离线changjian

UID ：22086

注册：2008-11-27
登录：2009-09-01
发帖：59
等级：仿真一级

1楼发表于: 2009-03-23 09:51:05

electromagnetic simulation using the FDTD Method
这本书上有，你看看吧只不过它是用C编写的。看你模拟的有多复杂了，在连接边界引入就可以了
主要考虑一下极化方向

共1条评分

cc81372365

rf币 +3

积极参与论坛交流，欢迎继续参与本贴交流！

2010-12-28

离线wq_463

UID ：20925

注册：2008-11-06
登录：2021-04-22
发帖：227
等级：仿真三级

2楼发表于: 2009-03-23 09:54:36

引用第1楼changjian于2009-03-23 09:51发表的 : E'g?44vyw
electromagnetic simulation using the FDTD Method Xtk3~@
这本书上有，你看看吧只不过它是用C编写的。看你模拟的有多复杂了，在连接边界引入就可以了 h/s8".\
主要考虑一下极化方向

我看了,不是很明白~~呵呵

离线gwzhao

方恨少

UID ：17098

注册：2008-08-24
登录：2019-01-09
发帖：1374
等级：荣誉管理员

3楼发表于: 2009-03-23 16:39:20

?
我记得葛德彪老师的那本书里面有边界条件处理的源代码啊。
如果你不做处理，直接加入源的话，那么空间都是总场了。
所以要看散射场，你必须在连接边界的地方把你加入的总场给去掉。

逆流而上

离线herosword

UID ：3621

注册：2007-07-01
登录：2019-11-03
发帖：30
等级：仿真二级

4楼发表于: 2009-03-23 20:07:58

不明白可以问啊，呵呵
葛德彪那本书上的是二维程序，写的很详细

离线wq_463

UID ：20925

注册：2008-11-06
登录：2021-04-22
发帖：227
等级：仿真三级

5楼发表于: 2009-03-24 10:57:47

引用第3楼gwzhao于2009-03-23 16:39发表的 : V0a)9\x(\
? vGI?X#w3
我记得葛德彪老师的那本书里面有边界条件处理的源代码啊。 D?@e,e
如果你不做处理，直接加入源的话，那么空间都是总场了。 N Uq'96{Y
所以要看散射场，你必须在连接边界的地方把你加入的总场给去掉。

我要看电磁波通过缝隙进入盒体内部以后的影响,只要加个平面波源,让电磁波进入盒体就可以.至于盒体外面什么情况,我根本不用管的,因为我关心的是盒体内部场强的变化,所以我觉的散射场是否不用设置呢?直接加个平面波源就可以

离线wq_463

UID ：20925

注册：2008-11-06
登录：2021-04-22
发帖：227
等级：仿真三级

6楼发表于: 2009-03-24 11:02:11

引用第4楼herosword于2009-03-23 20:07发表的 : +oe%bk|A
不明白可以问啊，呵呵 .~yz1^ c
葛德彪那本书上的是二维程序，写的很详细

呵呵,好的,不懂我就贴出来,你有空给我回复一下,谢谢了

离线gwzhao

方恨少

UID ：17098

注册：2008-08-24
登录：2019-01-09
发帖：1374
等级：荣誉管理员

7楼发表于: 2009-03-24 12:17:58

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%                圆形金属               %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear;
%%%%%%%%%%%%%%%%%%%%%%%%定义变量
ob=ones(60,60);%%%%%%%%%存放望个截止的编号
obeps=ones(60,60);
Media=zeros(4,15);%个介质的电磁参数
WaveLength=1;
WL=40;%入射波长和没波长的望个数
Ismin=0;
Ismax=0;
Jsmin=0;
Jsmax=0;%目标区域界限
MediaNo=0;%目标区域所有截止数目
Record=0;
i=0;
j=0;
Iflag=1;
IflagI=1;

dt=WaveLength/(WL*2*3e8);
G1=0;
G2=0;
GH1=0;
GH2=0;
%%%%%%%%%%%%%%%%%%%变量的初始化
%WaveLength=1;
%WL=40;
%Ismin=1;
%Ismax=120;
%Jsmin=1;
%Jsmax=120;

%%%%%%%%%%%%%%%%%%%%%建立散射提
%for i=Ismin:Ismax
%    for j=Jsmin:Jsmax
%        ob(i,j)=0;
%    end
%end
%for i=Ismin:Ismax
%    for j=Jsmin:Jsmax
%        if abs(i-60)<WL & abs(j-60)<WL
%            ob(i,j)=1;
%        end
%    end
%end

%%%%%%%%%%%%%%%%%%%%%%定义电磁场变量
Ez=zeros(65,65);%FDTD采样点Ez
Dz=zeros(65,65);%FDTD采样点Dz
Ez_P=Ez;
Ez_PP=Ez;
Hx=zeros(65,65);%FDTD采样点Hx
Hy=zeros(65,65);%FDTD采样点Hy
SEz1=Ez;
SEz2=Ez;
SEz3=Ez;

Ein=zeros(1600);%入射波E采样点
Hin=zeros(1600);%入射波H采样点

IncidentStart=0;
IncidentEnd=0;
Isource=0;%一纬入射波数组的端点和原点位置
Isbias=0;
EBin=zeros(4);%入射波吸收边界条件临时变量

ECB=zeros(65,4);%连接边界处E入射波的分量
HCB=zeros(65,4);%连接边界处H入射波的分量

EAB=zeros(65,4);%吸收边界条件临时变量
EAC=zeros(4,2,2);%吸收边界条件角点临时变量

FE1=zeros(11,11);%FDTD地带系数
FE2=zeros(11,11);
FH1=zeros(11,11);
FH2=zeros(11,11);
FE=0;
FH=0;
TEMflag=1;%1代表TM，2代表TE

MU0=0;%自由空间电磁参数：磁到绿
EPS0=0;%接点常数
SMU0=0;%磁到率的平方跟
SEPS0=0;%接点常数的平方根

OMIGA=0;%数值角频率

Phi=0;%入射角
TimeStop=0;%总的时间补偿

MaxX=0;
MinX=0;
MaxY=0;
MinY=0;
Step=0;

XNum=0;
YNum=0;

Imax=0;%FDTD区域总边界（吸收边界）
Imin=0;
Jmax=0;
Jmin=0;

Itmax=0;%连接边界
Itmin=0;
Jtmax=0;
Jtmin=0;

Iomax=0;%输出边界，用于近远场的变换
Iomin=0;
Jomax=0;
Jomin=0;


Ismin=20;%目标区边界
Ismax=40;
Jsmin=20;
Jsmax=40;

T1=0;%村春临时变量
T2=0;
T3=0;
T4=0;
T5=0;
T6=0;
CT1=0;
CT2=0;
CT3=0;
CT4=0;
I=0;
J=0;
N=0;
II=0;
III=0;
IFlag=0;
NN=0;

MediaNumber=2;%目标所包含介质总数目

%%%%%%%%%%%%%%%%%%初始化变量

TEMFlag=1;

Isource=700;

Imin=1;
Imax=60;
Jmin=1;
Jmax=60;
Itmin=10;%连接边界
Itmax=50;
Jtmin=10;
Jtmax=50;
Phi=45;
TimeStop=200;
IncidentStart=1;
IncidentEnd=1420;
Isbias=2*max(Itmax,Jtmax);
Isbias=Isource;
%%%%%%%%%%%%%%常量初始化
MU0=pi*0.0000004;
EPS0=8.85e-12;
Z=sqrt(MU0/EPS0);
k=2*pi/WaveLength;

Media(1,1)=1;%背景为真空epsilon_r
Media(2,1)=1;%mu_r
Media(3,1)=0;%导电绿
Media(4,1)=0;%磁到绿

MediaNo=1;%只有等离子体一种介质
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%          等离子体介质先不定义，用真空代替，然后在具体换算的时候再重新处理            %%%%%%%%%%%%%%%%%%%%%
Media(1,2)=1;%epsilon_r
Media(2,2)=1;%mu_r
Media(3,2)=3.72e7;%导电绿
%Media(3,2)=0;%导电绿
Media(4,2)=0;%磁到绿

%ob(Ismin:Ismax,Jsmin:Jsmax)=2;
%%%%%%%%%%%%%%%%%%%%%建立散射提
%%%%%%%%%%%%%%%%%%%%%建立散射提
%for i=Ismin:Ismax
%    for j=Jsmin:Jsmax
%        ob(i,j)=0;
%    end
%end
for i=Imin:Imax
    for j=Jmin:Jmax
        %if abs(i-60)<WL & abs(j-60)<WL
        if ((i-(Imin+Imax)/2)^2+(j-(Jmin+Jmax)/2)^2)<100
            ob(i,j)=2;
        end
    end
end
%%%%%%%%5先初始化为真空

for i=1:MediaNumber
    Media(3,i)=-Media(3,i)*Z/k;
    Media(4,i)=-Media(4,i)/(Z*k);
end

Phi=pi*Phi/180;%换成角度
if TEMFlag==-1%TE
    EPS0=pi*0.0000004;
    MU0=8.85e-12;
    Z=sqrt(MU0/EPS0);
end

FE=0.5*TEMFlag;
FH=0.5*TEMFlag;
SMU0=sqrt(MU0);
SEPS0=sqrt(EPS0);
OMIGA=pi/WL;
cp=cos(Phi);
sp=sin(Phi);

if TEMFlag==1 %计算FDTD迭代公式系数
    for i=1:MediaNumber
        for j=1:MediaNumber
            FE1(i,j)=((Media(1,i)+Media(1,j))+pi*0.5*(Media(3,i)+Media(3,j))/WL)/((Media(1,i)+Media(1,j))-pi*0.5*(Media(3,i)+Media(3,j))/WL);
            FE2(i,j)=1./((Media(1,i)+Media(1,j))-pi*0.5*(Media(3,i)+Media(3,j))/WL);
            FH1(i,j)=((Media(2,i)+Media(2,j))+pi*0.5*(Media(4,i)+Media(4,j))/WL)/((Media(2,i)+Media(2,j))-pi*0.5*(Media(4,i)+Media(4,j))/WL);
            FH2(i,j)=1./((Media(2,i)+Media(2,j))-pi*0.5*(Media(4,i)+Media(4,j))/WL);
        end
    end
else
    for i=1:MediaNumber
        for j=1:MediaNumber
            FH1(i,j)=((Media(1,i)+Media(1,j))+pi*0.5*(Media(3,i)+Media(3,j))/WL)/((Media(1,i)+Media(1,j))-pi*0.5*(Media(3,i)+Media(3,j))/WL);
            FH2(i,j)=-1./((Media(1,i)+Media(1,j))-pi*0.5*(Media(3,i)+Media(3,j))/WL);
            FE1(i,j)=((Media(2,i)+Media(2,j))+pi*0.5*(Media(4,i)+Media(4,j))/WL)/((Media(2,i)+Media(2,j))-pi*0.5*(Media(4,i)+Media(4,j))/WL);
            FE2(i,j)=-1./((Media(2,i)+Media(2,j))-pi*0.5*(Media(4,i)+Media(4,j))/WL);
        end
    end
end
FH=FH2(1,1);
FE=FE2(1,1);
OMIGA=pi/WL;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%          主程序      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for N=1:TimeStop
    %%%%%%%%%%产生入射波
    for i=IncidentStart:IncidentEnd-1
        Hin(i)=Hin(i)-FH*(Ein(i+1)-Ein(i));
    end
    for i=IncidentStart+1:IncidentEnd-1
        Ein(i)=Ein(i)-FE*(Hin(i)-Hin(i-1));
    end
    %%%加入源
    Ein(Isource)=sin(OMIGA*N);
    if N<WL
        Ein(Isource)=Ein(Isource)*0.5*(1-cos(OMIGA*N));%开关函数，升余弦函数
    end
    Ein(IncidentStart)=EBin(4);
    EBin(4)=EBin(3);
    EBin(3)=Ein(IncidentStart+1);
    Ein(IncidentEnd)=EBin(2);
    EBin(2)=EBin(1);
    EBin(1)=Ein(IncidentEnd-1);
%    plot(Ein(600:800));
%    pause(0.1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%       连接边界处入射场各分量      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for I=Itmin:Itmax
    T1=I*cp+Jtmin*sp+Isbias;
    II=round(T1);
    III=II+1;
    T1=III-T1;
    ECB(I,1)=T1*(Ein(II)-Ein(III))+Ein(III);%%%%%%%%%底部连接边界处入射电场分量

    T1=I*cp+Jtmax*sp+Isbias;
    II=round(T1);
    III=II+1;
    T1=III-T1;
    ECB(I,3)=T1*(Ein(II)-Ein(III))+Ein(III);%%%%%%%%%%%%%顶部连接边界处入射电场分量
end
for J=Jtmin:Jtmax
    T1=J*sp+Itmin*cp+Isbias;
    II=round(T1);
    III=II+1;
    T1=III-T1;
    ECB(J,4)=T1*(Ein(II)-Ein(III))+Ein(III);%%%%%%%%%%%%%左边连接边界处入射电场分量

    T1=J*sp+Itmax*cp+Isbias;
    II=round(T1);
    III=II+1;
    T1=III-T1;
    ECB(J,2)=T1*(Ein(II)-Ein(III))+Ein(III);%%%%%%%%%%%%%右边连接边界处入射电场分量
end

for I=Itmin:Itmax
    T1=I*cp+(Jtmin-0.5)*sp+Isbias;
    II=round(T1-0.5);
    III=II+1;
    T1=III+0.5-T1;
    HCB(I,1)=T1*(Hin(II)-Hin(III))+Hin(III);
    HCB(I,1)=HCB(I,1)*sp;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%底部连接边界入射波磁场分量Hx

    T1=I*cp+(Jtmax+0.5)*sp+Isbias;
    II=round(T1-0.5);
    III=II+1;
    T1=III+0.5-T1;
    HCB(I,3)=T1*(Hin(II)-Hin(III))+Hin(III);
    HCB(I,3)=HCB(I,3)*sp;%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%顶部连接边界入射波磁场分量Hx
end

for J=Jtmin:Jtmax
    T1=J*sp+(Itmin-0.5)*cp+Isbias;
    II=round(T1-0.5);
    III=II+1;
    T1=III+0.5-T1;
    HCB(J,4)=T1*(Hin(II)-Hin(III))+Hin(III);
    HCB(J,4)=-HCB(J,4)*cp;%%%%%%%%%%%%%%%%%%%%左边连接边界入射波磁场分量

    T1=J*sp+(Itmax+0.5)*cp+Isbias;
    II=round(T1-0.5);
    III=II+1;
    T1=III+0.5-T1;
    HCB(J,2)=T1*(Hin(II)-Hin(III))+Hin(III);
    HCB(J,2)=-HCB(J,2)*cp;%%%%%%%%%%%%%%%%%%%%左边连接边界入射波磁场分量

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     Ez的FDTD迭代公式           %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Ez_PP=Ez;
for I=Imin+1:Imax-1
    for J=Jmin+1:Jmax-1
        T1=Hy(I,J)-Hy(I-1,J)-Hx(I,J)+Hx(I,J-1);
        Ez(I,J)=FE1(ob(I,J),ob(I,J))*Ez(I,J)+FE2(ob(I,J),ob(I,J))*T1;
        %Dz(I,J)=Dz(I,J)+(1/(2))*T1;
    end
end
Ez_P=Ez;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     连接边界处引入入射波Ez分量         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for I=Itmin+1:Itmax-1
    Ez(I,Jtmax)=Ez(I,Jtmax)-FE*HCB(I,3);
    Ez(I,Jtmin)=Ez(I,Jtmin)+FE*HCB(I,1);
end

for J=Jtmin+1:Jtmax-1
    Ez(Itmax,J)=Ez(Itmax,J)+FE*HCB(J,2);
    Ez(Itmin,J)=Ez(Itmin,J)-FE*HCB(J,4);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%         连接边界角点引入入射波         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
T1=HCB(Jtmax,2)-HCB(Itmax,3);
Ez(Itmax,Jtmax)=Ez(Itmax,Jtmax)+FE*T1;

T1=-HCB(Jtmin,4)+HCB(Itmin,1);
Ez(Itmin,Jtmin)=Ez(Itmin,Jtmin)+FE*T1;

T1=HCB(Jtmin,2)+HCB(Itmax,1);
Ez(Itmax,Jtmin)=Ez(Itmax,Jtmin)+FE*T1;

T1=-HCB(Itmin,3)-HCB(Jtmax,4);
Ez(Itmin,Jtmax)=Ez(Itmin,Jtmax)+FE*T1;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%         吸收边界条件Ez                 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for I=Imin+1:Imax-1
    T1=Hy(I,Jmin)-Hy(I-1,Jmin);
    T1=T1+Hy(I,Jmin+1)-Hy(I-1,Jmin+1);
    T1=T1*TEMFlag;
    T2=Ez(I,Jmin+1)-Ez(I,Jmin);
    Ez(I,Jmin)=EAB(I,1)-1/3*T2;
    Ez(I,Jmin)=Ez(I,Jmin)+1/6*T1;
    EAB(I,1)=Ez(I,Jmin+1);

    T1=Hy(I,Jmax)-Hy(I-1,Jmax);
    T1=T1+Hy(I,Jmax-1)-Hy(I-1,Jmax-1);
    T1=T1*TEMFlag;
    T2=Ez(I,Jmax-1)-Ez(I,Jmax);
    Ez(I,Jmax)=EAB(I,3)-1/3*T2;
    Ez(I,Jmax)=Ez(I,Jmax)+1/6*T1;
    EAB(I,3)=Ez(I,Jmax-1);
end

for J=Jmin+1:Jmax-1
    T1=Hx(Imin,J)-Hx(Imin,J-1);
    T1=T1+Hx(Imin+1,J)-Hx(Imin+1,J-1);
    T1=T1*TEMFlag;
    T2=Ez(Imin+1,J)-Ez(Imin,J);
    Ez(Imin,J)=EAB(J,4)-1/3*T2;
    Ez(Imin,J)=Ez(Imin,J)-1/6*T1;
    EAB(J,4)=Ez(Imin+1,J);

    T1=Hx(Imax,J)-Hx(Imax,J-1);
    T1=T1+Hx(Imax-1,J)-Hx(Imax-1,J-1);
    T1=T1*TEMFlag;
    T2=Ez(Imax-1,J)-Ez(Imax,J);
    Ez(Imax,J)=EAB(J,2)-1/3*T2;
    Ez(Imax,J)=Ez(Imax,J)-1/6*T1;
    EAB(J,2)=Ez(Imax-1,J);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%         吸收边界角点Ez                 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
IflagI=1-IflagI;%%交换存储标志
Iflag=IflagI+1;
Ez(Imin,Jmin)=0.2928932*EAC(1,1,Iflag)+0.7071068*EAC(1,2,Iflag);
EAC(1,1,Iflag)=Ez(Imin,Jmin);
EAC(1,2,Iflag)=Ez(Imin+1,Jmin+1);

Ez(Imax,Jmin)=.2928932*EAC(2,1,Iflag)+0.7071068*EAC(2,2,Iflag);
EAC(2,1,Iflag)=Ez(Imax,Jmin);
EAC(2,2,Iflag)=Ez(Imax-1,Jmin+1);

Ez(Imin,Jmax)=0.2928932*EAC(3,1,Iflag)+0.7071068*EAC(3,2,Iflag);
EAC(3,1,Iflag)=Ez(Imin,Jmax);
EAC(3,2,Iflag)=Ez(Imin+1,Jmax-1);

Ez(Imax,Jmax)=0.2928932*EAC(4,1,Iflag)+0.7071068*EAC(4,2,Iflag);
EAC(4,1,Iflag)=Ez(Imax,Jmax);
EAC(4,2,Iflag)=Ez(Imax-1,Jmax-1);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%      Hx Hy 的FDTD迭代公式              %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for I=Imin:Imax
    for J=Jmin:Jmax-1
        T1=Ez(I,J+1)-Ez(I,J);


        Hx(I,J)=FH1(ob(I,J),ob(I,J+1))*Hx(I,J)-FH2(ob(I,J),ob(I,J+1))*T1;
    end
end

for I=Imin:Imax-1
    for J=Jmin:Jmax
        T1=Ez(I+1,J)-Ez(I,J);
        Hy(I,J)=FH1(ob(I,J),ob(I+1,J))*Hy(I,J)+FH2(ob(I,J),ob(I+1,J))*T1;
    end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%         连接边界引入入射波Hx Hy分量               %%%%%%%%%%%%%%%%%%%%%%%%%
for I=Itmin:Itmax
    Hx(I,Jtmin-1)=Hx(I,Jtmin-1)+FH*ECB(I,1);
    Hx(I,Jtmax)=Hx(I,Jtmax)-FH*ECB(I,3);
end

for J=Jtmin:Jtmax
    Hy(Itmin-1,J)=Hy(Itmin-1,J)-FH*ECB(J,4);
    Hy(Itmax,J)=Hy(Itmax,J)+FH*ECB(J,2);
end

%%%%%%%%%%%%%%%
xx=1:60;
yy=1:60;

    [xxx,yyy]=meshgrid(xx,yy);
    subplot(2,1,1);
    pcolor(xxx,yyy,Ez(1:60,1:60));
    subplot(2,1,2);
    mesh(xxx,yyy,Ez(1:60,1:60));
    surf(xxx,yyy,Ez(1:60,1:60));
    %shading interp;
    %colormap(hot);
%    plot(Ez(round((Itmin+Itmax)/2),:));
    pause(0.1);

end

共2条评分

cc81372365	rf币 +10	感谢您提供的资料，期待您继续分享！	2010-12-28
cem-uestc	技术分 +2	感谢您的资料	2009-03-24

逆流而上

离线gwzhao

方恨少

UID ：17098

注册：2008-08-24
登录：2019-01-09
发帖：1374
等级：荣誉管理员

8楼发表于: 2009-03-24 12:19:24

这是二维的FDTD/Matlab代码，我把其中等离子体的东西删掉了，所以可能运行起来有点问题，但上面每一步都有注释，你可以参考一下。

逆流而上

离线wq_463

UID ：20925

注册：2008-11-06
登录：2021-04-22
发帖：227
等级：仿真三级

9楼发表于: 2009-03-24 13:28:44

fortran写的啊....

勉强着看吧,谢谢版主了

发帖回复