[求助]2D_Y型光波导形成图像问题（MATLAB程序）

大家帮我看下这段程序  为什么会形成不规则的Y型  有什么地方需要修改吗 ;As~TGiT  
%tic
Thggas,  
%Initial(原始的) parameters and other things. Igo`\JY  
W=0.36; %Normalized（规范的） frequency %xgP*%Sv2  
%The following parameters are control parameters. {Ydhplg{  
WaveGuide=1; %If this program is for Wave Guide? If so, please specify 1. $JFjR@j  
IsMovie=0; %If you want to play movie, please use 1. c~T{;  
IsFigure=1; %If it will plot the figures? If so, please specify 1. D2ggFxqe
 
WantToSeeEp=1; %Do you want to see the distrubution of Ep? If so, please specify 1. v^C\ GDH  
%End of defining control parameters \]
x`f3F  
~?/7:S  
MLatx=11; %11How many Lattice（格子） cell in x direction. 7ZUN;mr  
MLaty=11; %11How many Lattice cell in y direction. ktU:Uq  
{\$S585  
NMlat=21; %21The grid number in each Lattice Cell. 534pX7
dg  
          %SHOULD BE ODD INTEGER（奇整数）!!!!  4!!|P  
if mod(NMlat,2)==0 !@z9n\Yj  
   NMlat=NMlat+1; 4phCn5  
end       %Force it to be a odd(奇数) integer! >{^_]phlb  
lU1SN/'zx  
NTx=MLatx*NMlat+1; %It is the number of the Grid along x axis. Z!v,;MW  
NTy=MLaty*NMlat+1; %It is the number of the Grid along y axis. / 4Q=%n  
2zj` H9  
if WaveGuide==1 X_PzK'#m  
   Nrow=5; %4位置问题The row number of columns between the PML boundary and the waveguide. !&U75FpN}:  
end _c]}m3/  
2f(`HSC'  
NPML=12; %12大小就是时间长短问题How many PML layers will be used in our computation. f}c;s  
+k=*AQt^8  
NTimeSteps=2500; %2500Total number of Time Steps N<JI^%HBgP  
UC&$8^  
Meach=10; %20Define the interval(间隔) for plot figures if IsFigure==1. \9+,ynJH8z  
           %This also works for saving intervals. JQ+Mg&&Q  
           (u?s@/e:`/  
R=0.2; %刚好全反射The radius of dielectric columns r`"_D%kc  
ea=11; %11.4The dielectric constant of these columns. +:Zwo+\kSN  
29Uqdo  
Zmax=0.6;    %0.6The maximum value for z axis when plotting figures. whI4@#  
Colormax=0.6; %0.6越大越暗The maximum value for colormap when plotting figures. {V[xBL <  
|]kiH^Ap  
%Some constants cC8$oCR?  
mu0=4*pi*1.0e-7; %Epsilon Zero, if using Gauss Unit, it equals to 1. ,.Gp_BI  
e0=8.85*1e-12;   %Mu Zero, if using Gauss Unit, it equals to 1. }-[l)<F:  
c=1/sqrt(mu0*e0); %The light speed. br\3}  
factor=mu0/e0; %The factor between conductivity(导电率) and permeability（浸透性）. RY5e%/bg~U  
               %Permeability=Conductivity*factor, in PML. h),;j`PrC  
               5zS%F: 3  
a=1;%e-6;   % 1 The lattice constant. oypLE=H  
W=W*(2*pi*c/a); %frequency ~criZI/  
               B} qRz  
Dx=a/NMlat; %Delta x. v8-szW).  
Dy=Dx; %Delta y. ],c0nz^%BR  
Dt=1/sqrt(1/(Dx*Dx)+1/(Dy*Dy))/c; %Time interval <(E)M@2  
{JWixbA  
%tic }Z <I%GT  
%In the following partm we define the dielectric constants: t_/qd9Jv  
Ep=ones(NTx-1,NTy-1)*e0; mLq0;uGL|  
Ep_cell=ones(NMlat,NMlat)*e0; 9^Xndo]y  
x=-(NMlat-1)/2*Dx:Dx:(NMlat-1)/2*Dx; b8a(.}8*  
[X,Y]=meshgrid(x); !MSa -  
X=X'; L'y0$  
Y=Y'; k2Q[v  
flag=find(sqrt(X.^2+Y.^2)<R); R5sEQ| E  
Ep_cell(flag)=e0*ea; W@}@5,}f>  
Ep=repmat(Ep_cell,MLatx,MLaty); ,IiKe_B  
1XS~b-St  
t(}Y/'  
   if WaveGuide==1 x.gzsd  
   eb=1; %The dielectric constant in the waveguide 8z93ETv7`  
   Ep(  1:(MLatx-Nrow)*NMlat,                        Nrow*NMlat+1:(Nrow+1)*NMlat)=e0*ea; 4H;g"nWqO  
   Ep( (MLatx-Nrow-1)*NMlat+1:(MLatx-Nrow+1)*NMlat,(Nrow+1)*NMlat+1:(Nrow+2)*NMlat)=e0*ea;%e0*eb XzAXcxC6G  
   Ep( (MLatx-Nrow)*NMlat+1:(MLatx-Nrow+2)*NMlat,(Nrow+2)*NMlat+1:(Nrow+3)*NMlat)=e0*ea;%e0*eb >uCO=T,|  
   Ep( (MLatx-Nrow+1)*NMlat+1:(MLatx-Nrow+3)*NMlat,(Nrow+3)*NMlat+1:(Nrow+4)*NMlat)=e0*ea;%e0*eb nu^@}|UG  
   Ep( (MLatx-Nrow+2)*NMlat+1:(MLatx-Nrow+4)*NMlat,(Nrow+4)*NMlat+1:(Nrow+5)*NMlat)=e0*ea;%e0*eb d^'_H>x  
   Ep( (MLatx-Nrow+3)*NMlat+1:(MLatx-Nrow+5)*NMlat,(Nrow+5)*NMlat+1:(Nrow+6)*NMlat)=e0*ea;%e0*eb y95 #t  
   Ep( (MLatx-Nrow+3)*NMlat+1:(MLatx-Nrow+5)*NMlat,(Nrow-5)*NMlat+1:(Nrow-4)*NMlat)=e0*ea;%e0*eb ejg!1*H@n  
   Ep( (MLatx-Nrow+2)*NMlat+1:(MLatx-Nrow+4)*NMlat,(Nrow-4)*NMlat+1:(Nrow-3)*NMlat)=e0*ea;%e0*eb :/->m6C`0  
   Ep( (MLatx-Nrow+1)*NMlat+1:(MLatx-Nrow+3)*NMlat,(Nrow-3)*NMlat+1:(Nrow-2)*NMlat)=e0*ea;%e0*eb
)+FnwW  
   Ep( (MLatx-Nrow)*NMlat+1:(MLatx-Nrow+2)*NMlat,(Nrow-2)*NMlat+1:(Nrow-1)*NMlat)=e0*ea;%e0*eb r%: :q^b3  
   Ep( (MLatx-Nrow-1)*NMlat+1:(MLatx-Nrow+1)*NMlat,(Nrow-1)*NMlat+1:Nrow*NMlat)=e0*ea; H8HH) ^  
   `y'%dY}$n  
   end :{w3l O  
%toc ,5 ,r.  
|TJ gH<I  
if WantToSeeEp==1 A7QT4h&6  
   x=0:Dx:(NMlat*MLatx-1)*Dx; ;$eY#ypx  
   x=x-(NMlat*MLatx-1)* .. ufR|V-BWx  
IlEU6Rs  

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

楼主能把具体的问题阐述下吗？谢谢！

学习了。。。