cey_main.m

用matlab程序模拟了多层平板波导的导波情况 包括TE模和TM模 该程序可进一步扩展

Num_x=69;
Num_y=81;
N=Num_x*Num_y;
lembda=1.55;
k_=2*pi/lembda;


C_ey=sparse(N,N);
Cey_temp=zeros(1,5);
j=0;
for m=1:1:Num_x
    for n=1:1:Num_y
        boundary_b=0;   %是否是边界，0 不是；1 是普通边界;2 是角点 
        j=(m-1)*Num_y+n;%矩阵元素序号
        
        if m==1&n==1 %内移一个点的序号（考虑边界向内移一个点仍相同，因此直接将外部虚拟点的计算转化为内部点）
            boundary_b=2;
            Cey_temp=Cey_c(2,2);%(1,1)->(2,2)
            C_ey(j,j)=Cey_temp(1)+Cey_temp(2)+Cey_temp(3);%透明边界
            %C_ey(j,j)=Cey_temp(3);                       %消逝边界
            
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j+Num_y)=Cey_temp(5);          
        end
        if m==1&n==Num_y
            boundary_b=2;
            Cey_temp=Cey_c(2,Num_y-1);%(1,Num_y)->(2,Num_y-1)
            C_ey(j,j)=Cey_temp(1)+Cey_temp(4)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);            
            
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j+Num_y)=Cey_temp(5);          
        end  
        if m==Num_x&n==1
            boundary_b=2;
            Cey_temp=Cey_c(Num_x-1,2);%(Num_x,1)->(Num_x-1,2)
            C_ey(j,j)=Cey_temp(2)+Cey_temp(5)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);
            
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j-Num_y)=Cey_temp(1);          
        end      
        if m==Num_x&n==Num_y
            boundary_b=2;
            Cey_temp=Cey_c(Num_x-1,Num_y-1);%(Num_x,Num_y)->(Num_x-1,Num_y-1)
            C_ey(j,j)=Cey_temp(5)+Cey_temp(4)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);
            
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j-Num_y)=Cey_temp(1);          
        end      
           
           
        if m==1&boundary_b~=2%左边界
            boundary_b=1;
            Cey_temp=Cey_c(2,n);%向右移了一格
            C_ey(j,j)=Cey_temp(1)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);            
            
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j+Num_y)=Cey_temp(5);
        end
        if m==Num_x&boundary_b~=2%右边界
            boundary_b=1;
            Cey_temp=Cey_c(Num_x-1,n);%向左移了一格
            C_ey(j,j)=Cey_temp(5)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);
            
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j-Num_y)=Cey_temp(1);
        end
        if n==1&boundary_b~=2%下边界
            boundary_b=1;
            Cey_temp=Cey_c(m,2);%向上移了一格
            C_ey(j,j)=Cey_temp(2)+Cey_temp(3); 
            %C_ey(j,j)=Cey_temp(3);
            
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j-Num_y)=Cey_temp(1);
            C_ey(j,j+Num_y)=Cey_temp(5);            
        end
        if n==Num_y&boundary_b~=2%上边界
            boundary_b=1;
            Cey_temp=Cey_c(m,Num_y-1);%向下移了一格
            C_ey(j,j)=Cey_temp(4)+Cey_temp(3);
            %C_ey(j,j)=Cey_temp(3);
            
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j-Num_y)=Cey_temp(1);
            C_ey(j,j+Num_y)=Cey_temp(5);            
        end
        
        if boundary_b==0 %非边界
            Cey_temp=Cey_c(m,n);
            C_ey(j,j)=Cey_temp(3);
            C_ey(j,j-1)=Cey_temp(2);
            C_ey(j,j+1)=Cey_temp(4);
            C_ey(j,j-Num_y)=Cey_temp(1);
            C_ey(j,j+Num_y)=Cey_temp(5);            
        end           
    end
end


A_x=C_ey;
options.tol=1e-10;			    %为eigs设定计算精度。
options.maxit=10000;
options.disp=0;					%使中间迭代结果不显示。
[v,d]=eigs(A_x,1,'lr',options); %计算最大实部特征值和特征向量。

ss=zeros(Num_x*Num_y,1);
for i=1:Num_x*Num_y
    ss(i)=v(i);
end;

b2=sum(sum(ss*k_c)')

B_y=sqrt(d);   %传播常数
Neff_y=B_y/k_ %有效折射率

Ey=T_to_Emn(abs(ss));  
xx=zeros(Num_y,Num_x);          %列出x坐标矩阵
for n=1:Num_y
    for m=1:Num_x
        xx(n,m)=x_position(m,n);
    end;
end;

yy=zeros(Num_y,Num_x);          %列出y坐标矩阵
for n=1:Num_y
    for m=1:Num_x
        yy(n,m)=y_position(m,n);
    end;
end;

%figure(2);
%Z=contour(xx,yy,Ey,100);
%title('TM模的基模电场分布');   %绘制基模的场分布，为等高线图
%xlabel('x');
%ylabel('y');