dynamic00.m

分析了Tm

function x=dynamic00
%%%%%%%%%%%bbes 1 2(ap0.25) 3
tau2=334.7*10^(-6);tau4=14.2*10^(-6);tau42=0.72*tau4;tau6=567*10^(-6);tau8=1.4*10^(-6);tau9=0.022*10^(-6);tau41=0.28*tau4;
Taup=0.8964;Taus=0.7520;
h=6.63*10^(-34);c=3*10^8;A=1.39*10^(-11);Lp=790*10^(-9);Ls=2*10^(-6);
k4212=5*10^(-25);k2124=0.084*k4212;         % Tm交叉弛豫系数
k2156=3.3*10^(-22);k1265=k2156/12.54;       % Tm-Ho能量传递系数
k2168=1.0*10^(-24);k1286=0.1*k2168;         % 上转换1
k4169=1*10^(-24);k1496=0.1*k4169;           % 上转换2
bbap=0.25*10^(-25);bbep=5*10^(-25);bbas=2*10^(-25);bbes=3.0*10^(-25);  % 受激吸收与辐射截面

Ntm=6.18*10^(25);
%Nho=1.05*10^(25);
Nho=6.18*10^(24);
ap=0.00035;as=0.00035;L=1;
n=1.55;
P0=2*h*(c/n)^2/(Ls^3);

R4=0.9;R3=1;

t=0:1*10^(-8):3.5*10^(-4);
m=length(t);
dt=1*10^(-8);
%z=0.001;
z=2;
a_mirror=-log(R3*R4)/(2*L);

N1_0=Ntm;N2_0=0;N4_0=0;
N5_0=Nho;N6_0=0;N8_0=0;N9_0=0;
P_in=1;P_sf0=0.001;P_sb0=0.001;
P_s0=P_sf0+P_sb0;
%P_p0=P_p0*exp((-Taup*bbap*N1_0+ap)*z)
%P_p0=P_p0*(1-exp((-Taup*bbap*N1_0-ap))*z)/((Taup*bbap*N1_0+ap)*z)

for my=1:m
    P_p0(my)=P_in*(1-exp((-Taup*(bbap*N1_0-bbep*N4_0)-ap)*z))/((Taup*(bbap*N1_0)+ap)*z);
    %P_p0(my)=P_in*exp((-Taup*(bbap*N1_0-bbep*N4_0)+ap)*z);
    N1(my)=N1_0+(N4_0*tau41+N2_0/tau2-Taup*Lp/(h*c*A)*P_p0(my)*bbap*N1_0-k4212*N4_0*N1_0+k2124*N2_0^2+k2156*N2_0*N5_0-k1265*N6_0*N1_0+k2168*N2_0*N6_0-k1286*N8_0*N1_0+k4169*N4_0*N6_0-k1496*N9_0*N1_0)*dt;
    N1_0=N1(my);
    N4(my)=N4_0+(Taup*Lp/(h*c*A)*P_p0(my)*bbap*N1_0-N4_0/tau4-k4212*N4_0*N1_0+k2124*N2_0^2-k4169*N4_0*N6_0+k1496*N9_0*N1_0)*dt;
    N4_0=N4(my);
    N2(my)=N2_0+(N4_0/tau42-N2_0/tau2+2*k4212*N4_0*N1_0-2*k2124*N2_0^2-k2156*N2_0*N5_0+k1265*N6_0*N1_0-k2168*N2_0*N6_0+k1286*N8_0*N1_0)*dt;
    N2_0=N2(my);
    N1(my)=Ntm-N2_0-N4_0;
    N1_0=N1(my);
    N6(my)=N6_0+(-N6_0/tau6+k2156*N2_0*N5_0-k1265*N6_0*N1_0-k2168*N2_0*N6_0+k1286*N8_0*N1_0-k4169*N4_0*N6_0+k1496*N9_0*N1_0-Taus*Ls/(h*c*A)*(P_s0)*(bbes*N6_0-bbas*N5_0))*dt;
    N6_0=N6(my);
    N5(my)=N5_0+(N6_0/tau6+N8_0/tau8+N9_0/tau9-k2156*N2_0*N5_0+k1265*N6_0*N1_0+Taus*Ls/(h*c*A)*(P_s0)*(bbes*N6_0-bbas*N5_0))*dt;
    N8(my)=N8_0+(-N8_0/tau8+k2168*N2_0*N6_0-k1286*N8_0*N1_0)*dt;
    N8_0=N8(my);
    N9(my)=N9_0+(-N9_0/tau9+k4169*N4_0*N6_0-k1496*N9_0*N1_0)*dt;
    N9_0=N9(my);
    N6(my)=Nho-N5_0-N8_0-N9_0;
    N5_0=N5(my);
    if (Taus*(bbes*N6_0-bbas*N5_0)-as-a_mirror)<0
       
        P_s(my)=0;
        
    else
        g(my)=Taus*(bbes*N6(my)-bbas*N5(my));
        P_s(my)=P_s0+(Taus*(bbes*N6(my)-bbas*N5(my))-as-a_mirror)*P_s0*c*dt/n;
        %P_sf(my)=P_sf0+(Taus*(bbes*N6(my)-bbas*N5(my))-as)*P_sf0*c*dt/n;
        %P_sb(my)=P_sb0+(-Taus*(bbes*N6(my)-bbas*N5(my))+as)*P_sb0*c*dt/n;
        %P_s(my)=P_sf(my)+P_sb(my);
        P_s0=R4*P_s(my);
        %P_s0=2*R4*P_s(my)/(1-R4);
    end
    P_p(my)=P_in*exp((-Taup*(bbap*N1(my))-ap)*z);
     S(my)=2.5*P_s(my);
     P(my)=P_p(my);
end
t=t/10^-6;
%P_out=(1-R4)*P_s/(1+R4);
%P_out=(1-R4)*P_s;
P_out=(1-R4)*S;
%plot(t,P_p0);
figure(1);plot(t,P,t,P_out)   
figure(2);plot(t,N5,t,N6)
figure(3);plot(t,N1,t,N2,t,N4)
%plot(t,P_p0)
%legend('N_1','N_2','N_4')
%legend('N_5','N_6','N_8','N_9')