基膜谐振腔Matlab仿真

本文主要是介绍基膜谐振腔Matlab仿真，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！

设计目的

①探求基膜下激光强度分布特性

②探求基膜谐振腔稳定性条件

设计任务和要求

熟悉了解基膜谐振腔的设计流程及其稳定条件参数设置，观察仿真的现象，直观感受基膜谐振腔的设计过程即工作机制。

设计原理概述

激光二极管端面抽运全固态激光器中，激光晶体内的热量通过热传导方式到达温度相对恒定的周边。通过循环水冷的方式或半导体致冷方式对激光晶体周边进行冷却。具体方法是，首先依据激光晶体的尺寸设计紫铜夹块，将晶体周边涂抹导热硅脂后用铟包裹置于紫铜夹块中，循环冷却水或半导体模块对紫铜夹块进行冷却。

抽运光经过一定的耦合方式热射到晶体的端面，形成端面抽运方式。并且抽运光斑沿着晶体端面集合中心入射。分析端面抽运激光晶体内温场分布状态，需建立沿晶体端面几何中心抽运、边缘恒温长方形晶体热模型，耦合到激光晶体端面的抽运光具有理想的高斯分布，并且抽运光版沿晶体中心纵向抽运。

$\text{[math]}$

晶体通过光端面与空气接触，经两端面与空气通过热交流出的热量远小于传导出的热量，可设晶体两端满足绝热条件，得出长方形晶体具有的边界条件：

$\text{[math]}$ , $\text{[math]}$

温度场遵守Poisson方程：

$\text{[math]}$

[注]： $\text{[math]}$ 为激光晶体内部的热功率密度； $\text{[math]}$ 为各向异性热传导系数。

晶体吸收的光束能量转化为热辐射表述：

$\text{[math]}$

求解温度场：

$\text{[math]}$

[注]：式中β为激光晶体对抽运光的吸收系数，η为荧光量子效应和内损耗决定的热转换系数

（ $\text{[math]}$ ）, $\text{[math]}$ 为激光二极管抽运波长808nm， $\text{[math]}$ 为谐振腔激光振荡波长1064nm。下面选用晶体基本参数：a-cut，0.3at%,吸收系数为7.711/cm；热传导系数 $\text{[math]}$ ,

$\text{[math]}$ , $\text{[math]}$ ；晶体尺寸设为3*3*5mm3,泵浦光斑为0.4mm，在20W泵浦强度下

求温度场分布。

热透镜焦距估算公式：

$\text{[math]}$

[注]：其中 $\text{[math]}$ 0.24为热负荷比，v=0.33为泊松比， $\text{[math]}$ -2.2*10^{-12}为光弹系数， $\text{[math]}$ 为吸

收抽运功率, $\text{[math]}$ 为热导率， $\text{[math]}$ 为平均抽运光斑半径， $\text{[math]}$ 2.165为环境温度下的折射率，

$\text{[math]}$ 5.092*10^{-5}为折射率温度系数， $\text{[math]}$ 4.43*10^-6/K为热膨胀系数。

设计实现

求解温度场：

clc;clear;
P=20;
a=3;
b=3;
c=5;
w0=0.4;
F=@(x,y)exp(-2*((x-a/2).^2+(y-b/2).^2)/w0.^2);
I0=P/integral2(F,0,a,0,b);
[T,x,z]=temper(I0);
toc;
function [T,xx,zz]=temper(I0)
global n m a b c w0 P;
beta=3.4;%1/mm吸收系数
eta=0.32;
Ky=5.23e-3;%1/mm热导率
Kx=5.1e-3;% W/mm/K
Kz=Kx;
x=linspace(0,a,101);
y=b/2;
zb=2;%参杂2mm；zb=0%未参杂晶体长度
z=linspace(0,c,101);
[xx,zz]=meshgrid(x,z);
T=zeros(size(xx));
N=10;
M=10;
L=20;
for m=1:2:Mfor l=0:Lfor n=1:NAnml=8*I0*beta*eta*(beta*c*cos(l*pi*zb/c)-beta*c*...exp(-beta*(c-zb))*cos(l*pi)-l*pi*sin(l*pi*zb/c)).../(a*b*pi^2)/(beta^2*c^2+l^2*pi^2)/(Kx*n^2/a^2+...Ky*m^2/b^2+Kz*l^2/c^2);Anml=Anml*integral2(@tempt,0,a,0,b);T=T+Anml*sin(n*pi*xx/a).*sin(m*pi*y/b).*cos(l*pi*zz/c);endend
end
figure;
surf(zz,xx,T);
shading interp;
ylabel('x/mm');
xlabel('z/mm');
zlabel('T/a.u.');
title(['Temperature Distribution,\omega_P=',num2str(w0),',P=',...num2str(P),'W']);
end
function Anml=tempt(x,y)
global n m a b w0;
Anml=exp(-2*((x-a/2).^2+(y-b/2).^2)/w0.^2).*sin(n*pi*x/a).*sin(m*pi*y/b);
end

仿真结果：

未参杂长度为2mm

计算Nd:YVO4的热透镜焦距：

%计算Nd:YVO4的热透镜焦距
Px=linspace(5,20,500);
Kc=5.23;%热传导W/(K*m)
wp=0.6e-3;
n0=2.165;
eta=1*(1-808/1064);
dndt=5.092e-6;
absorp=0.34;%吸收系数
L=1e-2;%晶体长度
alpha=4.43e-6;%/K热膨胀
niu=0.33;
Cr=-2.2e-12;%光弹系数
wrate=0.85;%震荡泵浦光斑尺寸比
fth=2*pi*Kc*wp^2./eta./Px./(dndt+(n0-1)*alpha);
plot(Px,fth*1e3);
%另一种方法
hold on
fth=2*pi*Kc*wp^2./eta./Px./(dndt+(n0-1)*(1+niu)*alpha+2*n0^3*alpha*Cr);
plot(Px,fth*1e3,'r');

仿真结果：

热透镜焦距随泵浦功率变化趋势(泵浦光斑半径0.6mm)

分析计算Z字型激光腔特性：

%计算"Z"字型激光腔的特性
clc;clear;
global lambda;
lambda=1.064e-3;
num=3000;
L1=20;%M1到热透镜的焦距
L23=20+250;%M3到等效薄透镜的距离
L4=115;%M4到M3的距离
R1=300e3;
R2=300e3;
R3=100;
R4=50;
z=linspace(0,L1+L23+L4,num);
wp=zeros(size(z));
st=wp;
%工作点基膜参数分析
thermallens=linspace(50,800,num);
for gk=1:nummatrix=lens(thermallens(gk))*space(L23)*mirror(R3)*space(L4)...*mirror(R4)*space(L4)*mirror(R3)*space(L23)*...lens(thermallens(gk))*space(L1)*mirror(R1)*space(L1);[spotsize,Rcurve,stable]=cavitystate(matrix,0,1,[],lambda);st(gk)=stable;
end
plot(thermallens,st,'k');
flens=100;%热透镜焦距
for gk=1:numif z(gk)<=L1matrix=space(L1-z(gk))*lens(flens)*space(L23)*mirror(R3)...*space(L4)*mirror(R4)*space(L4)*mirror(R3)*space(L23)*...lens(flens)*space(L1)*mirror(R1)*space(z(gk));[spotsize,Rcurve,stable]=cavitystate(matrix,0,1,[],lambda);wp(gk)=spotsize;elseif z(gk)<=L1+L23matrix=space(L23+L1-z(gk))*mirror(R3)*space(L4)*mirror(R4)...*space(L4)*mirror(R3)*space(L23)*lens(flens)*space(L1)*...mirror(R1)*space(L1)*lens(flens)*space(z(gk)-L1);[spotsize,Rcurve,stable]=cavitystate(matrix,0,1,[],lambda);wp(gk)=spotsize;elseif z(gk)<=L1+L23+L4matrix=space(L4+L23+L1-z(gk))*mirror(R4)*space(L4)*mirror(R3)...*space(L23)*lens(flens)*space(L1)*mirror(R1)*space(L1)*...lens(flens)*space(L23)*mirror(R3)*space(z(gk)-L1-L23);[spotsize,Rcurve,stable]=cavitystate(matrix,0,1,[],lambda);wp(gk)=spotsize;endwp(gk)=spotsize;
end
figure
plot(z,wp,'k',z,-wp,'k');
hold on
function matrix=space(L,n)
%用于产生一个长度为L的自由空间传输矩阵
if nargin<2n=1;
end
matrix=[1,L/n;0,1];
endfunction matrix=lens(f)
%用于产生焦距为f的透镜
matrix=[1,0;-1./f,1];
endfunction matrix=mirror(R,theta,s)
%用于产生一个曲率半径为R的反射镜,R>0为凹面镜
%theta为反射镜倾斜角，s标注子午面和弧矢面，默认为子午面
%s为1时为子午面，s为2时为弧矢面
if nargin<3s=1;
end
if nargin<2theta=0;
end
if s==1matrix=[1,0;-2./R./cos(theta),1];
elseif s==2matrix=[1,0;-2./R.*cos(theta),1];
elseerror('子午面和弧矢面光束用1和2区分，如mirror(R,theta,1).');
end
end
function varargout=cavitystate(matrix,s,MM2,q,lambda)
%计算ABCD传输矩阵的稳定性或者高斯参数的传输
%matrix为传输矩阵，s指明谐振腔(0)或外部光束(1)
%指明为谐振腔时，将返回腔的稳定性参数和基模高斯光束参数
%指明为外部光束时，只返回高斯光束的传输参数，这时需要输入传输光束的q参数
%MM2为光束的光束质量因子
A=matrix(1,1);
B=matrix(1,2);
C=matrix(2,1);
D=matrix(2,2);
if nargin<3MM2=1;
end
if nargin<4q=1/(-j*lambda/pi/0.5^2);%默认束腰宽度为0.5mm的光束
end
if nargin<2s=0;
end
if s==0stab=(A+D)^2/4;if stab>=1disp('warning：谐振腔不满足稳定性条件');wspot=0;Rcurv=0;elsewspot=sqrt(MM2*lambda*abs(B)/pi/sqrt(1-stab));Rcurv=2*B/(D-A);end
elseif s==1if isempty(q)q=input('请输入外部传输光束的q参数');MM2=input('请输入光束的M2因子：');endq=(A*q+B)/(C*q+D);wspot=sqrt(lambda*MM2/pi/abs(imag(1/q)));Rcurv=1/real(1/q);
elseerror('腔内光束和外部光束用0和1区分，如cavitystate(matrix,0).');
end
if nargout==1varargout{1}=wspot;
elseif nargout==2varargout{1}=wspot;varargout{2}=Rcurv;
elseif nargout==3&&s==0varargout{1}=wspot;varargout{2}=Rcurv;varargout{3}=stab;
elseif nargout==3&&s==1varargout{1}=wspot;varargout{2}=Rcurv;varargout{3}=q;
elseif nargout>3error('输出参数不正确!');
end
end

仿真结果：