登 录
註 冊
论坛
微波仿真网
注册
登录论坛可查看更多信息
微波仿真论坛
>
时域有限差分法 FDTD
>
紧急求助 2D FDTD 程序 &nbs ..
发帖
回复
1590
阅读
7
回复
紧急求助 2D FDTD 程序 会的指点一下 万分感谢
离线
flanky
UID :19721
注册:
2008-10-21
登录:
2011-11-15
发帖:
62
等级:
仿真一级
0楼
发表于: 2009-05-09 12:51:56
— 本帖被 destroyer 从 【互助速问速答有求速应】 移动到本区(2009-11-28) —
ic=Nx/4; % 源的X位置
[dUAb
jc=Ny/4; % 源的Y位置
!q$>6P
e'c3.sQ|?
for i=1:Nx+1;
Aq7`A^1t$
for j=1:Ny+1;
gg8)oc+w
dz(i,j)=0; %z方向电位移
#3eI4KJ4+l
ez(i,j)=0; % z方向电场
^eefR5^_w
hx(i,j)=0; % x方向磁场
<P&X0S`O
hy(i,j)=0; % y方向磁场
+2}Ar<elP
ihx(i,j)=0;
L<XX?I\p
ihy(i,j)=0;
Ey%[t
iz(i,j)=0; % z方向求和中间参量
KK]R@{ r
end;
I9rQX9#B
end;
4:733Q3oK
$}5M`p\&C
for i=2:Nx;
xDSiTp=)O
for j=2:Ny;
Ov4 [gHy&
ga(i,j)=1;
|2=w":2#
end;
4>H0a
end;
e=IbEm{|
%PML参数的设置
+ 65<|0
for i=1:Nx;
|u$*'EsP
gi2(i)=1;
Dq9f Fe
gi3(i)=1;
_ouZd.
fi1(i)=0;
yd'cLZd<}
fi2(i)=1;
/L2n ~/
fi3(i)=1;
1=E}X5
end
~k?t
G&MO(r}B
for j=1:Ny;
!=yO72dgLY
gj2(j)=1;
xB?!nd
gj3(j)=1;
s?nj@:4
fj1(j)=0;
K|6}g7&X
fj2(j)=1;
[nX{sM%
fj3(j)=1;
=,i?8Fuz
end
i.[k"(
;- Vs|X
for i=1:npml+1; %设置PML层中的参数
Wf%)::G*uR
xnum=npml+1-i;
IM[=]j.?
xn=0.33*(xnum/npml)^3; %这里的0.333并不是严格的计算,而是经验值
I\rjw$V#
gi2(i)=1.0/(1+xn);
^YfAsBs&
gi2(Nx-1-i)=1/(1+xn);
vl|3WYA
gi3(i)=(1-xn)/(1+xn);
#q9BU:
gi3(Nx-1-i)=(1-xn)/(1+xn);
C:bA:O
66+]D4(k
xn=0.25*((xnum-0.5)/npml)^3; %这里的0.25和0.333也是一个道理
upJy,|5
fi1(i)=xn;
;i+(Q%LO
fi1(Nx-2-i)=xn;
-XkjO$=!=
fi2(i)=1.0/(1+xn);
#>BC|/P}
fi2(Nx-2-i)=1/(1+xn);
LUz`P6
fi3(i)=(1-xn)/(1+xn);
e(`r"RrQ
fi3(Nx-2-i)=(1-xn)/(1+xn);
I/v#!`L
end
$[J\sokpY
,.gJ8p(0x
for i=1:npml+1;
KpS=oFX{}
xnum=npml+1-i;
c\tw#;\9
xn=0.33*(xnum/npml)^3;
9jl\H6JY|
gj2(i)=1.0/(1+xn);
Aqg$q* Y
gj2(Ny-1-i)=1/(1+xn);
!4t%\N6Ib
gj3(i)=(1-xn)/(1+xn);
`$3P@SO"
gj3(Ny-1-i)=(1-xn)/(1+xn);
,pkzNe`F
09G]t1!,
xn=0.25*((xnum-0.5)/npml)^3;
R\$6_
fj1(i)=xn;
HJ!)&xT
fj1(Ny-2-i)=xn;
X\<a|/{V A
|Z Cv>8?n
fj2(i)=1.0/(1+xn);
x0TE+rf5
fj2(Ny-2-i)=1/(1+xn);
=umF C[.W
fj3(i)=(1-xn)/(1+xn);
E!I4I'
fj3(Ny-2-i)=(1-xn)/(1+xn);
<!|=_W6
end
L9whgXD
+yHzp
%2.迭代求解电场和磁场
_F8-4
O&">%aU1I
for t=1:T;
kll,^A
for i=2:Nx; % 为了使每个电场周围都有磁场进行数组下标处理
1'B=JyR~K
for j=2:Ny;
eMH\]A~v"
dz(i,j)=gi3(i)*gj3(j)*dz(i,j)+gi2(i)*gj2(j)*0.5*(hy(i,j)-hy(i-1,j)-hx(i,j)+hx(i,j-1)); %
15DlD`QV
end;
)*o) iN 7l
end; % 电场循环结束
:4COPUBpPV
0nlh0u8#
pulse=sin(2*pi*f*t*dt); % 正弦波源
DFGgyFay
%pulse=exp(-4*pi*(t*dt-t0)^2/(2/f)^2);高斯脉冲情况
-OfAl~ 4
dz(ic,jc)=dz(ic,jc)+pulse; %软源,就是为了防止反射做的一个处理
2"Y=*s
BTE&7/i21
for i=1:Nx; %为了使每个电场周围都有磁场进行数组下标处理
FOcDBCrOe
for j=1:Ny;
,-Lv3
ez(i,j)=ga(i,j)* dz(i,j); %反映煤质的情况都是放到这里的
i l%9j
end;
$,T3vX]<
end; % 电荷密度循环结束
Io_7
b1ZHfe:
for j=1:Ny;
<)@^TRS
ez(1,j)=0;
uQWd`7
ez(Nx,j)=0;
{Z[kvXf"mZ
end
ACgWT
ZPY84)A_}
for i=1:Nx;
,z0E2
ez(i,1)=0;
dBW#PRg
ez(i,Ny)=0;
"n<u(m8E
end
9x9E+DG#(
uQW d1>
for i=1:Nx; % 为使每个磁场周围都有电场进行数组下标处理
H329P*P
for j=1:Ny-1;
l02aXxT)]
curl_e=ez(i,j)-ez(i,j+1);
I]ol[ X0S
ihx(i,j)=ihx(i,j)+fi1(i)*curl_e;
xtp55"g
hx(i,j)=fj3(j)*hx(i,j)+fj2(j)*0.5*(curl_e+ihx(i,j));
lt}U,p,S
end;
,k/<Nv;
end; % 磁场HX循环结束
?Ye%k
]JI A\|b6
for i=1:Nx-1; % 为了使每个磁场周围都有电场进行数组下标处理
sM%.=~AN
for j=1:Ny;
z7lbb*Xe
curl_e=ez(i+1,j)-ez(i,j);
/p@0Q[E
&n ..
]}AyDy6C
k${F7I(Tb
未注册仅能浏览
部分内容
,查看
全部内容及附件
请先
登录
或
注册
共
条评分
离线
flanky
UID :19721
注册:
2008-10-21
登录:
2011-11-15
发帖:
62
等级:
仿真一级
1楼
发表于: 2009-05-05 20:36:59
谁对FDTD比较擅长的,帮小弟一把,谢谢,这是一个用FDTD仿真2维TM波金属方柱散射
% FDTD_2D_RCS
.xuLvNyQr
% 本程序实现FDTD仿真2维TM波金属方柱散射
\aQBzEX
% TM波沿Y方向传播,X方向极化
`9Zoq=/
% 采用UPML吸收,脉冲源激励
jm'(t=Ze
Ps(3X@
KD*,u{v;
I2(5]85&]s
% 20050330
>Qz#;HI
% 模拟区域的范围以及编号的处理:
d>}pz
% 数字代表区域方向上的编号;
/i$E |[
% ------------------------------------------------------
IY?[ 0S
% | | BACK PML 2 | |
as- Z)h[B
% ------------------------------------------------------
[s"3g\L';
% |L | __________________________________ /| R|
at${^,&
% |E | | | (ib,jb) | I|
1 e1$x@\\
% |F | | | | G|
0:R}
% |T | | | | H|
8E"Ik~
% |3 | | MAIN GRID | | T|
f@T/^|`mh
% |P | | | |4 |
E$; =*0w
% |M | | | | P|
aR ao\Wp|
% |L | | | | M|
! R?r)G5E
% | |/ |__________________________________| | L|
<ZxxlJS)6
% -----------------------------------------------------
}x`W+r
% | | FRONT PML 1 | |
[WnX'R R
% -----------------------------------------------------
")"VQ|$y
?2_Oa%M
\,i9 m9;y
nVko]y
q [Rqy !,
clear;clc;
s\,F6c
% load('rcs1');
`Lb^!6`)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Wc ]BQn
%%*******************************初始化***********************************%%
{$)zC*l
tic;
\dTX%<5D
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
'{p/F $
v=3*10^8; % 波速
|g M|>
f=0.3*10^(9); % 频率
u;p{&\(]
lamda=v/f; % 波长
!*ct3{m
k=2*pi/lamda; % 波数
{v+,U}
b[<