950阅读
0回复

[求助]请帮忙看下我做的fdtd程序

离线lylmy

UID ：58065

注册：2010-04-25
登录：2012-02-29
发帖：40
等级：仿真新人

0楼发表于: 2010-09-10 22:37:16

小弟刚学fdtd不久，照葛德彪那本书天线那张做了个二维柱坐标系下的TE波线电流源，激励源采用双指数函数，吸收边界采用一阶Mur
感觉激励源的加入还是有很大问题
另外步长改为1的话出来的波形明显震荡
请各位高手指点一下

clear all
close all
clc

%--------------------------------------------------
%FDTD法模拟柱坐标二维TE波的传播(Mur吸收边界条件)
%--------------------------------------------------

%初始化
c=3e8; %光速
delta=10;
% delta=1e-3; %空间步长
dt=delta/(2*c); %时间步长
u0=(4*pi)*1e-7; %磁导率
e0=(1/(36*pi))*1e-9; %介电常数
I0=15e3; %基电流
n0=200; %空间网格数
% n1=33; %源的位置
abi=4*dt/(pi*e0*delta^2); %激励源系数
CA=1;
CB=dt/(e0*delta); %Er与Ez中的系数
CP=1;
CQ=dt/(u0*delta); %Hphi中的系数
C1=(1/(2*delta)-1/(2*c*dt)-1/(8*delta*(n0)))/(1/(2*delta)+1/(2*c*dt)+1/(8*delta*(n0))); %一阶边界条件中的系数
C2=(-1/(2*delta)+1/(2*c*dt)-1/(8*delta*(n0)))/(1/(2*delta)+1/(2*c*dt)+1/(8*delta*(n0)));
C3=(1/(2*delta)+1/(2*c*dt)-1/(8*delta*(n0)))/(1/(2*delta)+1/(2*c*dt)+1/(8*delta*(n0)));
% C1=-(4*n0+1)/(12*n0+1); %一阶边界条件中的系数
% C2=-(-4*n0+1)/(12*n0+1);
% C3=-(-12*n0+1)/(12*n0+1);
% C1=-(4*(n0+1)+1)/(12*(n0+1)+1); %一阶边界条件中的系数
% C2=-(-4*(n0+1)+1)/(12*(n0+1)+1);
% C3=-(-4*(n0+1)+1)/(-12*(n0+1)+1);
C4=(1.3e8*dt-delta)/(1.3e8*dt+delta);
Er=zeros(n0+1,n0+1);
Ez=zeros(n0+1,n0+1);
Hphi=zeros(n0,n0);
Hphi_last=zeros(n0,n0); %大型矩阵维数的预先确定

for k=1:2000

   %迭代求电场

   %求Er
   for j=2:n0
   for i=1:n0
   Er(i,j)=CA*Er(i,j)-CB*(Hphi(i,j)-Hphi(i,j-1));
   end
   end %Er循环结束
   %求Ez
   for j=1:n0
   for i=2:n0
   Ez(i,j)=CA*Ez(i,j)+CB*(delta*i*Hphi(i,j)-delta*(i-1)*Hphi(i-1,j))/(delta*(i-0.5));
   end
   end %Ez循环结束
   %激励源
   t=k*dt;
   I=1.1*I0*(exp(-9.2e4*t)-exp(-0.5e7*t));
   for j=1:n0
   Ez(1,j)=CA*Ez(1,j)+4*CB*Hphi(1,j)/delta-abi*I;
   end

   %迭代求解磁场
   Hphi_last=Hphi;
   for j=1:n0
   for i=1:n0
   Hphi(i,j)=CP*Hphi(i,j)+CQ*(Ez(i+1,j)-Ez(i,j))-CQ*(Er(i,j+1)-Er(i,j));
   end
   end %磁场循环结束

   %Mur边界条件

   % r=h
   for j=1:n0
   Hphi(n0,j)=C1*Hphi(n0-1,j)+C2*Hphi_last(n0,j)+C3*Hphi_last(n0-1,j);
   end
   %z=0,h
   for i=1:n0
..

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

发帖回复