3913阅读
5回复

[已解决]FDTD仿真PBG

离线febi

UID ：2537

注册：2007-05-11
登录：2016-05-05
发帖：123
等级：退休版主

0楼发表于: 2007-08-07 15:00:35

— 本帖被 wwwjjjsss 执行取消精华操作(2007-08-07) —

[ post]
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Finite-Difference-Time -Domain Method
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
v=sqrt(-1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Input Parameters

filename=input('output file = ? ','s');
q=1.3; %Incident wavelength(um)
k0=2*pi/q;
q0=1.3;
ratio1=0.35;
ratio2=0.2;
a=ratio1*q0;
RU=ratio2*a;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Calculation x-z Domain and Grid Size

dx=0.02;
dz=0.02;
p1=11;
p2=18;

Lx=p1*a;
Lz=p2*a;

x=-Lx/2:dx:Lx/2;
z=0:dz:Lz;
oo=size(x);
MU=oo(2);
oo=size(z);
NU=oo(2);

i=1:MU;
IU=2:MU;
j=1:NU;
JU=2:NU;

mz=3:NU;
cdt=dz/5;
P=(cdt-dz)/(cdt+dz);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Refractive Index Distribution PBG
% in air
n=sqrt(1)*ones(MU,NU);
r=zeros(MU,NU);
% in silicon
for n1=1:p1
for n2=1:p2
for n3=0:5:180
zr=-RU*sin(n3*pi/180):dz:RU*sin(n3*pi/180);
oz=size(zr);
xr=RU*cos(n3*pi/180);
x1=round((n1*a-a/2+xr)/dx);
z1=round((n2*a-a/2+zr)/dz);
if n1==4 & n2<=10
elseif n1>=4 & n2==10
else
n(x1,z1)=sqrt(11.4)*ones(1,oz(2));
end
end %n3
end %n2
end %n1

% Initial Condition

V1=n.^2+cdt*r/2*120*pi;
V2=n.^2-cdt*r/2*120*pi;
Hz=zeros(MU+1,NU);
Hx=zeros(MU,NU+1);
Ey=zeros(MU,NU);
x0=2*a;
E0=exp(-((x+x0)/a).^2)';

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Propagation Begin
Loop=round(Lz/cdt)*2;

for lp=1:Loop

ratio1=lp/Loop % 1st run, percent to complete

Hxi(i,1:2)=Hx(i,1:2);
Hxo(i,1:2)=Hx(i,NU:NU+1);
Hzi(1:2,j)=Hz(1:2,j);
Hzo(1:2,j)=Hz(MU:MU+1,j);

Ey(i,2)=Ey(i,2)+E0*exp(v*k0*lp*cdt);

Hx(i,JU)=Hx(i,JU)+cdt/dz*(Ey(i,JU)-Ey(i,JU-1));
Hz(IU,j)=Hz(IU,j)-cdt/dx*(Ey(IU,j)-Ey(IU-1,j));
Ey=cdt*((Hx(i,j+1)-Hx(i,j))/dz-(Hz(i+1,j)-Hz(i,j))/dx)./V1+Ey.*V2./V1;

Hx(i,1)=Hxi(i,2)+P*(Hx(i,2)-Hxi(i,1));
Hx(i,NU+1)=Hxo(i,1)+P*(Hx(i,NU)-Hxo(i,2));
Hz(1,j)=Hzi(2,j)+P*(Hz(2,j)-Hzi(1,j));
Hz(MU+1,j)=Hzo(1,j)+P*(Hz(MU,j)-Hzo(2,j));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
% Display
contour(z(mz),x,n(:,mz),1);
hold on;
pcolor(z(mz),x,real(Ey(:,mz)));shading flat;colormap(gray)
hold off
pause(1)
%%%%%%%%%

end&nb ..

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

共1条评分

wwwjjjsss

RFEDA币 +15

难得一见的好资料+5~15

2007-08-07

磁砖

离线森森

UID ：6275