changdu.m

利用龙格库塔法解光波导激光器的重叠因子的传输方程

function [ output_args ] = Untitled1( input_args )
%UNTITLED1 Summary of this function goes here
%  Detailed explanation goes here
clear
clc
er12=6.5e-25;%%%截面
er13=2.58e-25;
er21=9e-25;
yb56=1e-24;
yb65=1e-24;
A21=1/7.9e-3;
A65=1/2e-3;
h=6.636e-34;
Cup=5e-23;
Ccr=2e-22;
NER=2.5e26;
NYB=2.5e27;
vs=3e8/1540e-9;%%%%%%%信号中心频率
vp=3e8/980e-9;%%%%%%%泵浦频率
A=2E-12;%%%%%%%%波导宽度
Ts=0.6;%%%%%%重叠因子
Tp=0.6;
dz=2e-4;%%%%%%%步长 0.1mm
pp(1)=80e-3;%%%%%%%泵浦功率
ps(1)=1e-3;%%%%%%%%信号功率
n1=0;
n2=0;
n5=0;
n6=0;
i=0;
ii=0;
i2=0;
l(1)=0;
p2=zeros(1,17);
for i2=0:16
    L=1+i2*0.1;
    i=0;
for z=0:dz:L*100e-4
i=i+1;
l(i+1)=l(i)+dz;
%%%%%%%%%%求能级粒子数%%%%%%%%%%%%%%%%
W1=-er12*ps(i)*Ts/h/vs/A-er13*pp(i)*Tp/h/vp/A;
W2=A21+er21*ps(i)*Ts/h/vs/A;
W22=Cup;
W16=-Ccr;
R5=-yb56*pp(i)*Tp/h/vp/A;
R6=yb65*pp(i)*Tp/h/vp/A+A65;
P=[W22*W16
W22*(R5-R6)+W16*(W1-W2-2*W22*NER)
(W1-W2-2*W22*NER)*(R5-R6)+W16*(W2*NER+W22*NER^2)+W16*R5*NYB
(W2*NER+W22*NER^2)*(R5-R6)];
x=roots(P);
for ii=1:3
if (isreal(x(ii)))
if (x(ii)>0)
if(x(ii)<=NER)
n1=x(ii);
end
end
end
end
n2=NER-n1;
n6=R5*NYB/(R5-R6+W16*n1);
n5=NYB-n6;
%%%%%R-K 法解泵浦光传输方程 dpp/dz=cp*pp%%%%%%%%%
cp=-(er12*n1+yb56*n5-er21*n6)*Tp;
k1=cp*pp(i);
k2=cp*(pp(i)+dz/2*k1);
k3=cp*(pp(i)+dz/2*k2);
k4=cp*(pp(i)+dz*k3);
pp(i+1)=pp(i)+dz/6*(k1+2*k2+2*k3+k4);
%%%%%%%%R-K 法解输入光传输方程 dps/dz=cs*ps%%%%%%%%%
cs=(er21*n2-er12*n1)*Ts;
k1=cs*ps(i);
k2=cs*(ps(i)+dz/2*k1);
k3=cs*(ps(i)+dz/2*k2);
k4=cs*(ps(i)+dz*k3);
ps(i+1)=ps(i)+dz/6*(k1+2*k2+2*k3+k4);
end
p1=pp(end);
p2(i2+1)=ps(end)*1e+3;
end
pp=zeros(1,17);
for i=0:16
    pp(i+1)=1+i*0.1;
end
 plot(pp(:),p2(:))
 xlabel('波导长度/cm')
 ylabel('输出功率/mw')
 grid on
 hold on