www.pudn.com > Mutilayer_FD.rar > FDmethod.m, change:2007-05-22,size:3049b


clear,close all 
na=3.1693; 
nb=3.3884; 
nc=3.1693; 
nd=1; 
ne=3.1693; 
nf=1; 
ng=1; 
n=[na nb nc nd ne nf ng] 
lengh=1.3; 
k0=2*pi/lengh; 
%%%%%%%%%%%%%%%%%%%%% 网格划分 %%%%%%%%%%%%%%%%%%%%%%% 
%%%%%%%% 非等间距划分 横向 31点、纵向81点 
deltx=0.2;               % 横向 31 点  
delty1=0.24;            % 纵向 81点 
delty2=0.04; 
delty3=0.015; 
delty4=0.01; 
delty5=0.01; 
delty6=0.14; 
delty7=0.01; 
delty8=0.095 
              
y1=1/delty1; 
y2=1/delty2; 
y3=1/delty3; 
y4=1/delty4; 
y5=1/delty5; 
y6=1/delty6; 
y7=1/delty7; 
y8=1/delty8; 
x=1/deltx; 
for i=1:31              % 给每个格点的折射率赋值 
    for j=1:81 
        if j<=16 
            N(j,i)=na; 
        elseif j>16&j<=36 
            N(j,i)=nb; 
        elseif j>36&j<=46 
            N(j,i)=nc; 
        elseif j>46&j<=66&i>1&i<=11 
            N(j,i)=nd; 
        elseif j>46&j<=66&i>11&i<=21 
            N(j,i)=ne; 
        elseif j>46&j<=66&i>21&i<=31 
            N(j,i)=nf; 
        else N(j,i)=ng 
        end 
    end 
end 
  
          
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 构造本征方程的矩阵 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 
CEX=zeros(31*81,31*81);        % 创建一个矩阵,就是系数矩阵 
M=81 
k=1 
for m=1:31 
   for n=1:81 
                         
% x=(n-1)*M+m 格点数与坐标的转换关系 即N(n,m)=N((n-1)*M+m) 
% 周围四个点的转换 
 
a2=(m-1)*M+(n-1);     % m,n-1 
a1=((m-1)-1)*M+n;     % m-1,n 
a3=(m-1)*M+n;         % m,n 
a5=((m+1)-1)*M+n;     % m+1,n 
a4=(m-1)*M+(n+1);     % m,n+1 
    if (n-1)<=0 
        a2=0; 
          end 
    if (n+1)>81 
        a4=0; 
          end 
    if (m-1)<=0 
        a1=0; 
          end    
    if (m+1)>31 
        a5=0; 
          end 
 
if a1~=0     
Tered=2*N(a1)/(N(a1)+N(a3)); 
else Tered=0; 
end 
if a5~=0 
Teadd=2*N(a5)/(N(a5)+N(a3)); 
else Teadd=0; 
end 
 
if n<=11 
    y=y1; 
elseif n>11&n<=21 
    y=y2; 
elseif n>21&n<=31 
    y=y3; 
elseif n>31&n<=41 
    y=y4; 
elseif n>41&n<=51 
    y=y5; 
elseif n>51&n<=61 
    y=y6; 
elseif n>61&n<=71 
    y=y7; 
else 
    y=y8; 
end 
Cx=[x^2*Tered y^2 -2*y^2-(4-Teadd-Tered)*x^2+k0^2*na^2 y^2 x^2*Teadd ]; 
  
%%%%%%%%%%%给矩阵赋值     
    if a1~=0 
        CEX(k,a1)=Cx(1); 
         end     
    if a2~=0 
         CEX(k,a2)=Cx(2); 
         end     
    if a3~=0 
             CEX(k,a3)=Cx(3); 
         end     
    if a4~=0 
                 CEX(k,a4)=Cx(4); 
         end     
     if a5~=0 
                     CEX(k,a5)=Cx(5); 
         end 
         
k=k+1; 
end 
end 
 
 
%%%%%%%%%%%%%%%%%%%%%%%%  求本征值和本征向量 %%%%%%%%%%%%%%%%%%%%%%%%%% 
[V,D]=eig(CEX);      % V 代表本征向量构成的矩阵,D 代表本征值构成的对角矩阵 
%  取出本征值,求出其中的最大数,作为所求的准TE模的本征值,也即传播常数 
B=diag(D);            
b=B(1); 
for i=1:2511 
    y=B(i)-b; 
    if y>0 
        b=B(i); 
        k=i 
    end 
end 
neff=sqrt(b/k0^2) 
E=V(:,k); 
k=1 
for i=1:31 
for j=1:81 
EX(i,j)=abs(E(k)); 
k=k+1 
end 
end 
% 构造本征向量矩阵,矩阵中的点的值对应与分割网格上的电场的值 
%%%%%%% 绘图 %%%%%%%%% 
 
 
 
 
 
%%%%%%%%%  %%%%%%%%% 
%Cex*Ex=b^2*Ex