基于MATLAB的2×2光纤定向耦合器设计

《基于MATLAB的2×2光纤定向耦合器设计》由会员分享，可在线阅读，更多相关《基于MATLAB的2×2光纤定向耦合器设计（22页珍藏版）》请在装配图网上搜索。

1、基于MATLAB的22光纤定向耦合器设计1设计原理1.1 单模光纤的传导场如图1，光纤的横截面有三层介质，分别是是芯层、包层和涂层，芯层折射率稍大于包层折射率，导波光由于全反射背包层约束在芯层中沿光纤延伸方向传播。假设光的传播方向为光纤中心轴方向。 图1 阶跃光纤横截面结构图为简化讨论，只考虑基模的耦合。已知光纤中传导场表达式为 （1-1）其中，为光纤中导波光场的场振幅，为光纤中导波光场的场分布，为基模场的传播常数，为角频率。某时刻在光纤中的传导场的空间分布就与，和为相关。1.1.1 单模光纤的场分布当给定波导（即光纤）的边界条件时，求解波导场方程可得本征解及相应的本征值，即模式。模式是波导结

2、构的固有电磁共振属性的表征。单模光纤中只能存在基模，其场分布是确定的，可由亥姆霍兹方程求得。在弱导光纤中的电磁波，其横向场分量Et、Ht远大于纵向场分量Ez、Hz ，而且横向场分量是线偏振的。于是我们把电场的横向分量取为y轴方向，即Et=Ey。亥姆霍兹方程为 （1-2）其中，为真空中的光波矢量。利用分离变量法，将方程1-2在圆柱坐标系中求解，并结合电磁场的边界条件，可以解出电场的横向分量Ey： （1-3）其中，r是点到光纤中心轴的距离，m是整数，Jm和Km分别是m阶贝塞尔函数和m阶变态汉克尔函数；a是光纤芯层半径，一般单模光纤的a=25m。q是极角。横向磁场只包含Hx分量，由于横向电场与横向磁

3、场的比等于波阻抗Z0=122，故可由Ey计算出Hx： （1-4）再由麦克斯韦方程组可求出场的纵向分量Ez、Hz： （1-5）（1-6）1.1.2 单模光纤的传播常数由于在无耦合的情况下，光纤中导波光场的场振幅。故当某时刻t的场分布只同传播常数有关系。由式1-3可以得出光纤的本征方程本征方程，即r =0连续。 （1-7）同理，由横向磁场的连续性同样可以得出本征方程: （1-8）由于特征方程是超越方程，只能靠计算机数值求解U(或W)。其中U与W定义为 （1-9） （1-10）U与W是场的横向传播常数，U反映了导模在芯区中的驻波场的横向振荡频率，W值反映了导模在包层中的消逝场的衰减速度。W越大衰减越

4、快。式1-6、1-7与相联系，因此本征方程实际是关于的一个超越方程。对于不同的U值，可求得相应的值。由于贝塞尔函数及其导数具有周期振荡性质，所以本征值方程可以有多个不同的解mn(m=0,1,2,3., n=1,2,3.),每一个mn都对应于一个导模。单模光纤的本征值为01 。光纤的结构参数可由其归一化频率V表征： （1-11）其中，其中相对折射率差，n1、n2（n1n2）分别是芯层和包层的折射率，在阶跃折射率光纤中均为常数。实际的光纤的D很小，一般小于1%，因此也称之为弱导光纤。归一化频率V描述了光纤的结构参数a，D，n1以及工作波长（包含在k0中），是一个重要的综合性参数。V值越大，允许存在

5、的导模数就越多。除了基模之外，其他导模都可能在某一个V值以下截止，这时导模转化为辐射模。当某一模式截止时，它就不能沿光纤传输了。当W2 0时，对于导波模，场在纤芯外是指数衰减的；当W2 0时，W为虚数，场在纤芯外不再衰减了，也就不被约束在芯层中传播了，这种波叫做辐射波。而W2 =0正好是临界状态，以此作为导波模被截止的条件，此时的W记做Wc，相应的U和V值记做Uc和Vc，并称Vc为归一化截止频率。截止条件(W=0)下，特征方程可以简化为 （1-12）对于m=0，n=1时，截止时由1-9式得Vc=Uc=0，说明这种模式永不截止，是光纤中的最低阶模，也称基模。第二个归一化截止频率较低的模式是当m=

6、1，n=1时，称为二阶模，其归一化截止频率Vc=Uc=2.40483。当光纤中只有基模工作而其他模式均被截止时，就称为单模光纤。要保证单模传输，必须要求光纤的归一化频率小于二阶模LP11的截止频率Vc，即让二阶模截止： （1-13）此即单模传输条件。1.1.3LPmn模的功率密度以及功率平均坡印廷矢量Sav的大小就代表功率密度，可由以及场解1-31-6式，可以计算出光纤横截面上的功率密度(Sav)近似为： （1-14）平均坡印廷矢量对横截面积分便得到光纤中的功率P： （1-15）令P=1 (归一化功率)，可算出以上各式的常数A： （1-16）1.2 22光纤定向耦合器的场振幅在两根单模光纤耦合

7、的情况下，光纤中导波光场的场振幅，而是随着轴向距离z变化的。故某时刻的场传导方程就还需求解场振幅。22光纤定向耦合器的示意图如下图所示：图(a) 22光纤定向耦合器的示意图图(b) 22光纤耦合器横截面示意图图(c) 耦合器纵向剖面示意图(即rz面)如图所示，两平行直光纤的纤芯相互接近时，在其中传输的基模场分布就会互相渗透和交叠。这样1号光纤中传播的导模场将在2号光纤中产生极化作用，从而在2号光纤中激起传导模。设nc_01和 nc_02分别为1号光纤和2号光纤的芯层折射率，nc_clad为包层折射率。可将两光纤中的传导场写为： （1-17） （1-18）其中e1、e2是两根单模阶跃光纤中的横向

8、电场，其表达式由（1-3）式给出。a1(z)、a2(z)为复数振幅，可由耦合波方程解得： （1-19） （1-20）其中，Db=b1b2为失配相位常数。k12、k21为耦合系数，取决于波导结构以及波长，且满足： （1-21）耦合系数可由下式计算： （1-22） （1-23）其中，n1(r)、n2(r)分别代表1号光纤和2号光纤的折射率分布。P0（=1，已经归一化）为输入光功率。（1-19）式两边对z求导，并利用(1-20)式，可得到： （1-24）求解(1-24)式，解得 （1-25）再利用(1-19)式解得 （1-26）积分常数c1和c2由边界条件，给出则最终求解出 （1-27）1.3耦合光

9、纤中的耦合效率耦合光纤中的功率密度分布即为平均坡印廷矢量的大小。光耦合发生时，光横截面上功率密度还需要在(1-14)式上乘上复数振幅的绝对值平方（1号光纤）或（2号光纤）。这样在任意的耦合距离z处，直通臂（即1号光纤）中的相对功率 （1-28）耦合臂（2号光纤）中的相对功率为 （1-29）耦合效率为 （1-30）其中P0=1为归一化的输入光功率。2设计方案基于上述原理，现通过MATLAB进行设计仿真实现。MATLAB具有功能强大的数据可视化功能，可以给出数据的二维、三维的图形表现。本次课程设计可以通过对上述设计原理的MATLAB语言描述及求解，应用MATLAB中常用的二维、三维绘图函数以及和绘

10、图有关的命令，对两根单模光纤的耦合规律进行计算机仿真。在设计中，需要注意几个重要的问题，如如何求解单模光纤的本征值。从原理可以知道，传播常数跟光纤的本征值有着一一对应的关系，故而对传播常数的求解便是关键问题。同样对耦合系数k12、k21进行积分求解需要用到MATLAB中的二重积分，这就需要对函数变量进行取舍。显然对于柱形结构的光纤，选用极坐标更为合适。于是设计出如图2所示流程图。该流程图的关键部分为判断光纤是否为单模光纤，求解特征方程从而获得传播常数，对耦合系数k12、k21进行积分求解，求解光线中的场振幅和能流密度，以及最终的绘图输出。MATLAB语言设计流程框图如下求解1号光纤的传播常数1

11、求解1号光纤的常数A1定义光纤基本参数开始输入光纤1参数,纤芯半径、折射率判断1号光纤是否为单模光纤V12.40483输入光纤2参数,纤芯半径、折射率判断2号光纤是否为单模光纤V12.40483,disp(1号光纤不是单模光纤,请重新输入);%判断1号光纤是否满足单模条件r01_core=input(输入1号光纤的纤芯半径(单位为m)=)*1e-06; %定义1号光纤纤芯半径r01_clad=3*r01_core; %定义1号光纤的包层半径nc01=input(输入1号光纤的纤芯折射率=); %定义1号光纤的纤芯折射率delta1=(nc01-n_clad)/nc01; %相对折射率差V1=n

12、c01*k0*r01_core*sqrt(2*delta1); %归一化频率end; disp(1号光纤是单模光纤,输入正确) %1号光纤定义结束%由特征方程求1号光纤的归一化径向相位常数和归一化径向衰减常数fun1=(U1) U1.*besselj(m+1,U1)/besselj(m,U1)-sqrt(V12-U12).*besselk(m+1,sqrt(V12-U12)./besselk(m,sqrt(V12-U12); %特征方程U1=fzero(fun1,1);%求解归一化径向相位常数UW1=sqrt(V12-U12); %求解归一化径向衰减常数Wbeta01=sqrt(nc012*k

13、02-(U1/r01_core)2); %求解1号光纤的传播常数betaA1=U1/V1*besselk(m,W1)/sqrt(besselk(m+1,W1)*besselk(m+1,W1)*sqrt(4*z0/n_clad/pi/r01_core2); %求解1号光纤的A%求1号单模光纤中的横向场分量Ey1Ey01=(r,theta) A1*cos(m*theta).*besselj(m,U1*r/r01_core)./besselj(m,U1);Ey02=(r,theta) A1*cos(m*theta).*besselk(m,W1*r/r01_core)./besselk(m,W1);%

14、定义2号光纤的参数r02_core=input(输入2号光纤的纤芯半径(单位为m)=)*1e-06; %定义2号光纤纤芯半径r02_clad=3*r02_core; %定义2号光纤的包层半径nc02=input(输入2号光纤的纤芯折射率=); %定义2号光纤的纤芯折射率delta2=(nc02-n_clad)/nc02; %相对折射率差V2=nc02*k0*r02_core*sqrt(2*delta2); %归一化频率while V22.40483,disp(2号光纤不是单模光纤,请重新输入);%判断2号光纤是否满足单模条件r02_core=input(输入2号光纤的纤芯半径(单位为m)=)*

15、1e-06; %定义2号光纤纤芯半径r02_clad=3*r02_core; %定义2号光纤的包层半径nc02=input(输入2号光纤的纤芯折射率=); %定义2号光纤的纤芯折射率delta2=(nc02-n_clad)/nc02; %相对折射率差V2=nc02*k0*r02_core*sqrt(2*delta2); %归一化频率end;disp(2号光纤是单模光纤,输入正确)disp(开始计算，输出仿真图)%由特征方程求2号光纤的归一化径向相位常数和归一化径向衰减常数fun2=(U2) U2.*besselj(m+1,U2)/besselj(m,U2)-sqrt(V22-U22).*bes

16、selk(m+1,sqrt(V22-U22)./besselk(m,sqrt(V22-U22);U2=fzero(fun2,1); %求解归一化径向相位常数U W2=sqrt(V22-U22); %求解归一化径向衰减常数Wbeta02=sqrt(nc022*k02-(U2/r02_core)2); %求解2号光纤的传播常数beta%求解2号光纤的AA2=U2/V2*besselk(m,W2)/sqrt(besselk(m+1,W2)*besselk(m+1,W2)*sqrt(4*z0/n_clad/pi/r02_core2);%求2号单模光纤中的横向场分量Ey2Ey11=(r,theta) A

17、2*cos(m*theta).*besselj(m,U2*r/r02_core)./besselj(m,U2);Ey12=(r,theta) A2*cos(m*theta).*besselk(m,W2*r/r02_core)./besselk(m,W2);d=r01_clad+r02_clad;delta_beta=beta01-beta02;%定义失配相位常数F=(r,theta) conj(Ey01(r,theta).*Ey12(sqrt(d2+r.2-2*d*r.*cos(theta),theta)*(nc012-n_clad2).*r*omega*eption_0/4;K12=dblq

18、uad(F,0,r01_core,0,2*pi); 求解耦合系数K12K21=conj(K12); %求解耦合系数K21%求耦合光纤中的功率密度c2=(delta_beta+sqrt(delta_beta2+4*abs(K21)2)*a1_0+2*K12*a2_0)/(2*sqrt(delta_beta2+4*abs(K21)2);c1=a1_0-c2;%求解1号光纤中的耦合计算r_core1=linspace(-r01_core,r01_core,400); %rz平面上纤芯内某点的纵坐标z1=linspace(0,1,1000); %耦合距离Z1,R_core1=meshgrid(z1,r

19、_core1);theta0=0; %沿零度方位角取光纤的轴向剖面（rz平面）r_clad1=linspace(r01_core,(d+r02_clad),500); %rz包层内某点的纵坐标Z11,R_clad1=meshgrid(z1,r_clad1);a1=c1*exp(i*(delta_beta+sqrt(delta_beta2+4*abs(K21)2)/2*z1)+c2*exp(i*(delta_beta-sqrt(delta_beta2+4*abs(K21)2)/2*z1); %求解1号光纤中场振幅Sz_core1=A12/2/z0*(cos(m*theta0)2*nc01.*(b

20、esselj(m,U1*abs(r_core1)/r01_core)./besselj(m,U1).2*abs(a1).2; %求解纤芯功率密度Sz_clad1=A12/2/z0*(cos(m*theta0)2*n_clad.*(besselk(m,W1*r_clad1/r01_core)/besselk(m,W1).2*abs(a1).2; %求解包层功率密度%求解2号光纤中的耦合计算r_core2=linspace(d-r02_core,d+r02_core,500); %rz平面上纤芯内某点的纵坐标z2=linspace(0,1,1000); %耦合距离Z2,R_core2=meshgr

21、id(z2,r_core2);r_clad2=linspace(0,d-r02_core,500);Z22,R_clad2=meshgrid(z2,r_clad2);a2=-(c1*(delta_beta+sqrt(delta_beta2+4*abs(K21)2)/2*exp(i*(-delta_beta+sqrt(delta_beta2+4*abs(K21)2)/2*z2)+c2*(delta_beta-sqrt(delta_beta2+4*abs(K21)2)/2*exp(-i*(delta_beta+sqrt(delta_beta2+4*abs(K21)2)/2*z2)/K12; %求解

22、2号光纤中场振幅Sz_core2=A22/2/z0*(cos(m*theta0)2*nc02.*(besselj(m,U2*abs(r_core2-d)/r02_core)./besselj(m,U2).2*abs(a2).2; %求解2号光纤纤芯功率密度Sz_clad2=A22/2/z0*(cos(m*theta0)2*n_clad.*(besselk(m,W2*(d-r_clad2)/r02_core)/besselk(m,W2).2*abs(a2).2; %求解2号光纤包层功率密度%耦合功率密度绘图cmap=linspace(1,0,256);linspace(1,0,256);zero

23、s(1,256);colormap(cmap);subplot(1,1,1);surf(Z1,R_core1,Sz_core1),view(-10,60);hold on;surf(Z11,R_clad1,Sz_clad1),view(-10,60);subplot(1,1,1);surf(Z2,R_core2,Sz_core2),view(-10,60);hold on;surf(Z22,R_clad2,Sz_clad2),view(-10,60);shading flat;colorbar;axis tight;xlabel(耦合距离z（m）);ylabel(r,Fontsize,13,F

24、ontName,Times);title(1号和2号光纤的耦合功率密度分布);%绘制耦合光纤功率密度对比图figurecmap=linspace(1,0,256);linspace(1,0,256);zeros(1,256);colormap(cmap);subplot(2,1,2);mesh(Z1,R_core1,Sz_core1),view(0,0);hold on;mesh(Z11,R_clad1,Sz_clad1);shading flat;colorbar;axis tight;xlabel(耦合距离z(m);ylabel(r,Fontsize,13,Fontname,Times);

25、title(1号光纤的耦合功率密度分布);subplot(2,1,1);mesh(Z2,R_core2,Sz_core2);hold on;mesh(Z22,R_clad2,Sz_clad2),view(0,0);shading flat;colorbar;axis tight;xlabel(耦合距离z（m）);ylabel(r,Fontsize,13,FontName,Times);title(2号光纤的耦合功率密度分布);%绘制耦合效率曲线图figurePout_1=abs(a1).2;Pout_2=abs(a2).2;enta=Pout_2;subplot(3,1,1);plot(z1,Pout_1);xlabel(耦合距离z(m);ylabel(直通臂功率Pout_1);subplot(3,1,2);plot(z1,Pout_2);xlabel(耦合距离z(m);ylabel(耦合臂功率Pout_2);subplot(3,1,3);plot(z1,enta);xlabel(耦合距离z(m);ylabel(耦合效率);title(耦合效率耦合距离曲线图);