引用激光光束传输与谐振腔仿真代码包ABCDRez的内容
文档阅读:
激光光束传输与谐振腔仿真ABCDRez_abcdrez csdn-CSDN博客
文件下载:
Laser Beam Propagation and Resonator Simulation ABCDRez - File Exchange - MATLAB Central
https://download.csdn.net/download/qq_42712244/89533054
3.2 驻波谐振腔
可使用矩阵光学方法对多元件的稳定腔进行分析,采用“G参数等价腔”来描述腔内束宽。参见吕百达教授著《激光光学》第十章第一节。
3.2.1 两个厚透镜晶体加曲面腔镜的谐振腔
使用“G参数等价腔”相关内容,编写具体函数Rez4mThick。
![](https://i-blog.csdnimg.cn/direct/398e0740c4254ae292df89e28d0d9355.png)
图3.2.1.1 谐振腔示意图
图3.2.1.2 腔内束宽示意图
[www,wthetaL0,FlagRez]=Rez4mThick(lambda,RezPara,str)
变量
意义
lambda
波长
RezPara
谐振腔参数
str
若str为’plot’则绘图;否则,不绘制
www
有效的束宽(这里用来描述模式体积)
wthetaL0
靠近输出镜的光束的束腰、发散角、光腰位置信息
FlagRez
谐振腔的标志信息
其中:有效束宽www为振荡光斑在晶体四个端面的平均值。
略去绘图部分的代码,重新改写函数为
[www,wthetaL0,FlagRez]=Rez4mThick00(lambda,RezPara)
变量
意义
lambda
波长
RezPara
谐振腔参数
www
有效的束宽(这里用来描述模式体积)
wthetaL0
靠近输出镜的光束的束腰、发散角、光腰位置信息
FlagRez
谐振腔的标志信息
这样,我们就可以快速计算大量腔型,并根据结果进行筛选!
如查看某一腔型的有效束宽随热焦距的变化情况
next
%% 腔参数
lambda=1064*10^-9;
%
rho1=-1500*10^-3;
d1=(10+0)*10^-3;
% F1=500*10^-3;
lenF1=30*10^-3;nF1=1.8;
d3=15*10^-3;
% F2=F1;
lenF2=lenF1;nF2=nF1;
d2=(40+0)*10^-3;
rho2=-2500*10^-3;
%
limR=5e-3;
LL=d1+lenF1+d3+lenF2+d2;
%% 暴力求解,耗时短
F1x=linspace(250,2500,250).*10^-3;
www=zeros(size(F1x));
for ii=1:length(F1x)
F1=F1x(ii);
F2=F1;
RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
[www(ii),~,~]=Rez4mThick00(lambda,RezPara);
end
% 剔除超过边界的数据
www( (www>limR) )=0;
%% 绘图
plot(F1x,www.*10^3);
grid on;hold on;
xlabel('焦距/m');
ylabel('有效束宽/mm')
title('有效束宽随热焦距变化情况')
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
如查看某一腔型的有效束宽随两腔镜曲率的变化情况
如查看某一腔型的稳定性随两腔镜曲率的变化情况
next
%%
lambda=1064e-9;
%
% rho1=-1500*10^-3;
d1=(10+0)*10^-3;
F1=500*10^-3;lenF1=30*10^-3;nF1=1.8;
d3=15*10^-3;
F2=F1;lenF2=lenF1;nF2=nF1;
d2=(40+0)*10^-3;
% rho2=-2500*10^-3;
%
limR=5e-3;
LL=d1+lenF1+d3+lenF2+d2;
%% 暴力求解,耗时较短
tic
Rho=[5000,4500,4000,3500,3000,2500,2000,1500,1250,1000,800,600,500,450,400,350,300,250]*10^-3;
Rhoall=[inf,Rho,flip(-Rho),-inf];
www=zeros(length(Rhoall),length(Rhoall));
FlagRez=zeros(length(Rhoall),length(Rhoall));
for ii=1:length(Rhoall)
rho1=Rhoall(ii);
for jj=1:length(Rhoall)
rho2=Rhoall(jj);
RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
[www(ii,jj),~,FlagRez(ii,jj)]=Rez4mThick00(lambda,RezPara);
end
end
toc
% 剔除超过边界的数据
www( (www>limR) )=0;
%% 绘图
[XX,YY]=meshgrid(Rhoall,Rhoall);
figure;
mesh(XX.*10^3,YY.*10^3,www.*10^3);
xlabel('rho1/mm');ylabel('rho2/mm');zlabel('www/mm');
title('有效束宽随两腔镜曲率的变化情况');
figure;
mesh(XX.*10^3,YY.*10^3,FlagRez);
xlabel('rho1/mm');ylabel('rho2/mm');zlabel('FlagRez');
title('稳定性随两腔镜曲率的变化情况');
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
有效束宽随两腔镜曲率的变化情况:
稳定性随两腔镜曲率的变化情况:
当然,也可调用ABCDRez中其他谐振腔模拟函数查看相关变化情况,这里就不在赘述。
其中引用的函数如下:
function [www,wthetaL0,FlagRez]=Rez4mThick00(lambda,RezPara)
% RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
rho1= RezPara(01);d1= RezPara(02);F1= RezPara(03);
lenF1= RezPara(04);nF1= RezPara(05);d3= RezPara(06);
F2= RezPara(07);lenF2= RezPara(08);nF2= RezPara(09);
d2= RezPara(10);rho2= RezPara(11);
%% 不绘图,用于大量计算
syms z
%总腔长
LL=d1+lenF1+d3+lenF2+d2;
[~,MF1,~,~]=deF2Rho(F1,lenF1,nF1);
if (F1==F2 && lenF1==lenF2 && nF1==nF2)
MF2=MF1;
else
[~,MF2,~,~]=deF2Rho(F2,lenF2,nF2);
end
%
MM=[1,d2;0,1]*MF2*[1,d3;0,1]*MF1*[1,d1;0,1];
a=MM(1,1);b=MM(1,2);c=MM(2,1);d=MM(2,2);
%令
G1=a-b/rho1;
G2=d-b/rho2;
FlagRez(1)=G1*G2;
%% 是否满足0<G1*G2<1稳定条件
% figure(22);
if ( FlagRez(1)<0 || FlagRez(1)>1 )
www=0;
wthetaL0=zeros(1,3);
else
%成立稳定性条件
%0<G1*G2<1
%% 计算
%靠近两腔镜束腰大小
w001=abs(sqrt(+(lambda*b/pi)*(sqrt(G1*G2*(1-G1*G2))/abs(G1+a*a*G2-2*a*G1*G2))));
w002=abs(sqrt(+(lambda*b/pi)*(sqrt(G1*G2*(1-G1*G2))/abs(G2+d*d*G1-2*d*G1*G2))));
%束腰位置
L001=(b*G2*(a-G1))/(G1+a*a*G2-2*a*G1*G2);
L002=LL-( (b*G1*(d-G2))/(G2+d*d*G1-2*d*G1*G2) );
%基模发散角
theta001=sqrt((lambda/(pi*w001))^2);
theta002=sqrt((lambda/(pi*w002))^2);
%% 腔内束宽
%腔内的束宽变换
%由镜1全反镜向输出镜变换
[w003,theta003,L003]=fLRMm(w001,theta001,L001,lambda,MF1,lenF1,d1);
%
w001_z=fwz(w001,theta001,L001);
w003_z=fwz(w003,theta003,L003);
w002_z=fwz(w002,theta002,L002);
www=1/4*double(subs(w001_z,z,d1)+subs(w003_z,z,d1+lenF1)+...
subs(w003_z,z,d1+lenF1+d3)+subs(w002_z,z,d1+lenF1+d3+lenF2));
wthetaL0=double([w002,theta002,L002]);
end
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
function [Rho,MF,MFront,MBack]=deF2Rho(F,lenF,nF)
% F 厚透镜有效焦距
% lenF 厚透镜长度
% nF 厚透镜介质折射率
n1=1;
% 曲率表示
syms phox
MFront=[1,0;(nF-n1)/(nF*-phox),n1/nF];
Mmedian=[1,lenF;0,1];
MBack=[1,0;(n1-nF)/(n1*phox),nF/n1];
Mf1=MBack*Mmedian*MFront;
%求解代换
phoxx=vpasolve(Mf1(2,1)==-1/F);
[~,nn]=max(abs(phoxx));
Rho=phoxx(nn);
% % % % % % % %
% 多个解??!!
% % % % % % % %
MFront=double([1,0;(nF-n1)/(nF*-Rho),n1/nF]);
Mmedian=double([1,lenF;0,1]);
MBack=double([1,0;(n1-nF)/(n1*Rho),nF/n1]);
MF=MBack*Mmedian*MFront;
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
%% 高阶高斯光束由左向右传输
%% [w02m,theta02m,L02m]=fLRMm(w01m,theta01m,L01m,lambda,cM,lencM,dcM)
function [w02m,theta02m,L02m]=fLRMm(w01m,theta01m,L01m,lambda,cM,lencM,dcM)
% 高阶转基模
[w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda);
% 基模由左向右传输
[w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM);
% 基模转高阶
[w02m,theta02m,L02m,~]=fB2M(w02,theta02,L02,M2);
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
%% 高阶模转基模
%% [w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda)
function [w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda)
% w01m为高阶光束光腰 半径
% theta01m为高阶光束 半角发散角
% L01m为高阶光束光腰位置
% M2为光束质量
% w01为基模输出光腰 半径
% theta01为基模光束 半角发散角
% L01为基模光束光腰位置
%默认使用1064nm
% lambda=1064*10^-9;
%默认在空气中折射率nn=1
nn=1;
M2=(w01m*theta01m*pi*nn)/(lambda);
w01=(w01m)/sqrt(M2);
theta01=(theta01m)/sqrt(M2);
L01=L01m;
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
%% 基模由左向右传输
%% [w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM)
function [w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM)
% 默认使用基模计算,由左向右传输
% w01,theta01,L01 为入射光束参数
% lambda 为波长
% cM 为光学系统的传输矩阵
% lencM 为光学系统长度
% dcM 为 入射光束与光学系统的第一接触面 相对于 原点所处位置
% w02,theta02,L02 为出射光束参数
%% 物方光束
Z0=w01/theta01;
syms z;
w_1_z=w01*sqrt(1+((z-L01)/Z0)^2);
R_z1=Z0*((z-L01)/Z0+Z0/(z-L01));
R01=subs(R_z1,z,eps);
w_001=subs(w_1_z,z,0);
%% 传输矩阵相关定义及计算
syms X1 Y1 z
% syms R_1
% q参数定义
% syms s1 s2 a b c d
% 取n1=n2=1;
n1=1;n2=1;
%传输矩阵元素
s2=z;
a=cM(1,1);b=cM(1,2);c=cM(2,1);d=cM(2,2);
s1=dcM;
%其他特征参量
M=[1 s2;0,1]*[a b;c d]*[1,s1;0,1];
A=M(1,1);B=M(1,2);
% C=M(2,1);D=M(2,2);
Y2=(n1/n2)*Y1/...
(A^2+2*X1*A*B+(X1^2+Y1^2)*B^2);
%% 像方光束
w_2=sqrt(lambda/(pi*Y2));
w_2=subs(w_2,[X1 Y1],[1/R01 lambda/(pi*w_001^2)]);
%光腰参考点平移
w_2=subs(w_2,(z),(z-s1));
%% 输出参数
% disp('w_2*w_2=');pretty(w_2*w_2);
T1=vpa((taylor((w_2*w_2),z,'Order',1)));
T2=vpa((taylor((w_2*w_2),z,'Order',2)));
T3=vpa((taylor((w_2*w_2),z,'Order',3)));
C=T1;
B=(T2-T1)/z;
A=(T3-T2)/(z*z);
w02=vpa(sqrt(C-B*B/(4*A)));
theta02=vpa(sqrt(A));
L02=vpa(-B/(2*A));
L02=L02+lencM;
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
%% 基模转高阶模
%% [w01m,theta01m,L01m,M2]=fB2M(w01,theta01,L01,M2)
function [w01m,theta01m,L01m,M2]=fB2M(w01,theta01,L01,M2)
% w01m为高阶光束光腰 半径
% theta01m为高阶光束 半角发散角
% L01m为高阶光束光腰位置
% M2为光束质量
% w01为基模输出光腰 半径
% theta01为基模光束 半角发散角
% L01为基模光束光腰位置
% M2=(w01m*theta01m*pi*nn)/(lambda)
% 默认在空气中折射率nn=1
% nn=1;
% lambda=(w01m*theta01m*pi*nn)/(M2);
w01m=(w01)*sqrt(M2);
theta01m=(theta01)*sqrt(M2);
L01m=L01;
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
%% 光束随距离变化
%% w_z = fwz(w0,theta0,L0,z)
function w_z = fwz(w0,theta0,L0,z)
if nargin==3
syms z;
end
w_z=sqrt(w0^2+theta0^2*(z-L0)^2);
end
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
** next脚本:**
close all;
clear;
clc;
%% 版本信息
% 作者: Quincy Howard
% 联系方式: [email protected]
% 文件信息更新平台为 matlab云文件、GitHub、CSDN博客、知乎
% 若使用请注明来源
% 最后编辑于 2024 年 07 月 10 日
版权归原作者 QuincyHoward 所有, 如有侵权,请联系我们删除。