984阅读
0回复

[RFEDA原创]FDTD迷你代码-1D 连接边界

UID ：90389

注册：2012-03-16
登录：2014-09-11
发帖：17
等级：仿真新人

0楼发表于: 2014-06-11 08:52:38

clear;
clc;
c=3e8;
m0=1e-7*4*pi;
e0=1e-9/(36*pi);
eta0=120*pi;
nz=100;
nc=10;
ex=zeros(1,nz+1);
hy=zeros(1,nz);
lambda=1;%wavelength
f=c/lambda;
omega=2*pi*f;
nd=20;
d=lambda/nd;
s=1;
dt=s*d/c;
cb=dt/(d*e0);
cq=dt/(d*m0);
mur1=(s-1)/(s+1);
t0=(nz/2-nc+2)*d/c;
z=(-nz/2:nz/2)*d;
fh=@(t)(heaviside(t).*sin(t).*((heaviside(t)-heaviside(t-6*pi)).*(1-cos(pi*t/(6*pi)))/2+heaviside(t-6*pi)));
exi=@(t)(fh(omega.*t));
hyi=@(t)(exi(t)/eta0);
h=plot(ex);
set(gca,'ylim',2*[-1,1]);
for n=0:200
    ex0=ex([2,end-1]);
    ex(2:end-1)=ex(2:end-1)-cb*diff(hy);
    t=(n+0.5)*dt-t0-(z(1+nc)-0.5*d)/c;
    ex(1+nc)=ex(1+nc)+cb*hyi(t);%connection
    t=(n+0.5)*dt-t0-(z(1+nz-nc)+0.5*d)/c;%connection
    ex(end-nc)=ex(end-nc)-cb*hyi(t);
    ex([1,end])=ex0([1,end])+mur1*(ex([2,end-1])-ex([1,end]));%ABC
    hy=hy-cq*diff(ex);
    t=(n+1)*dt-t0-(z(1+nc))/c;
    hy(nc)=hy(nc)+cq*exi(t);
    t=(n+1)*dt-t0-(z(1+nz-nc))/c;
     ..

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

发帖回复