MATLAB激光谐振腔仿真代码

打印 上一主题 下一主题

主题 845|帖子 845|积分 2537



 引用激光光束传输与谐振腔仿真代码包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
    

图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

谐振腔的标志信息

这样,我们就可以快速计算大量腔型,并根据效果进行筛选!

如查看某一腔型的有效束宽随热焦距的变化情况

  1. next
  2. %% 腔参数
  3. lambda=1064*10^-9;
  4. %
  5. rho1=-1500*10^-3;
  6. d1=(10+0)*10^-3;
  7. % F1=500*10^-3;
  8. lenF1=30*10^-3;nF1=1.8;
  9. d3=15*10^-3;
  10. % F2=F1;
  11. lenF2=lenF1;nF2=nF1;
  12. d2=(40+0)*10^-3;
  13. rho2=-2500*10^-3;
  14. %
  15. limR=5e-3;
  16. LL=d1+lenF1+d3+lenF2+d2;
  17. %% 暴力求解,耗时短
  18. F1x=linspace(250,2500,250).*10^-3;
  19. www=zeros(size(F1x));
  20. for ii=1:length(F1x)
  21.     F1=F1x(ii);
  22.     F2=F1;
  23.     RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
  24.     [www(ii),~,~]=Rez4mThick00(lambda,RezPara);
  25. end
  26. % 剔除超过边界的数据
  27. www( (www>limR) )=0;
  28. %% 绘图
  29. plot(F1x,www.*10^3);
  30. grid on;hold on;
  31. xlabel('焦距/m');
  32. ylabel('有效束宽/mm')
  33. title('有效束宽随热焦距变化情况')
  34. %% 版本信息
  35. % 作者:                Quincy Howard
  36. % 联系方式:           quincy.hd@qq.com
  37. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  38. % 若使用请注明来源
  39. % 最后编辑于           2024 年 07 月 10 日
复制代码


如查看某一腔型的有效束宽随两腔镜曲率的变化情况

如查看某一腔型的稳定性随两腔镜曲率的变化情况

  1. next
  2. %%
  3. lambda=1064e-9;
  4. %
  5. % rho1=-1500*10^-3;
  6. d1=(10+0)*10^-3;
  7. F1=500*10^-3;lenF1=30*10^-3;nF1=1.8;
  8. d3=15*10^-3;
  9. F2=F1;lenF2=lenF1;nF2=nF1;
  10. d2=(40+0)*10^-3;
  11. % rho2=-2500*10^-3;
  12. %
  13. limR=5e-3;
  14. LL=d1+lenF1+d3+lenF2+d2;
  15. %% 暴力求解,耗时较短
  16. tic
  17. Rho=[5000,4500,4000,3500,3000,2500,2000,1500,1250,1000,800,600,500,450,400,350,300,250]*10^-3;
  18. Rhoall=[inf,Rho,flip(-Rho),-inf];
  19. www=zeros(length(Rhoall),length(Rhoall));
  20. FlagRez=zeros(length(Rhoall),length(Rhoall));
  21. for ii=1:length(Rhoall)
  22.     rho1=Rhoall(ii);
  23.     for jj=1:length(Rhoall)
  24.         rho2=Rhoall(jj);
  25.         RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
  26.         [www(ii,jj),~,FlagRez(ii,jj)]=Rez4mThick00(lambda,RezPara);
  27.     end
  28. end
  29. toc
  30. % 剔除超过边界的数据
  31. www( (www>limR) )=0;
  32. %% 绘图
  33. [XX,YY]=meshgrid(Rhoall,Rhoall);
  34. figure;
  35. mesh(XX.*10^3,YY.*10^3,www.*10^3);
  36. xlabel('rho1/mm');ylabel('rho2/mm');zlabel('www/mm');
  37. title('有效束宽随两腔镜曲率的变化情况');
  38. figure;
  39. mesh(XX.*10^3,YY.*10^3,FlagRez);
  40. xlabel('rho1/mm');ylabel('rho2/mm');zlabel('FlagRez');
  41. title('稳定性随两腔镜曲率的变化情况');
  42. %% 版本信息
  43. % 作者:                Quincy Howard
  44. % 联系方式:           quincy.hd@qq.com
  45. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  46. % 若使用请注明来源
  47. % 最后编辑于           2024 年 07 月 10 日
复制代码
有效束宽随两腔镜曲率的变化情况:

稳定性随两腔镜曲率的变化情况:


固然,也可调用ABCDRez中其他谐振腔模拟函数查看相关变化情况,这里就不在赘述。

其中引用的函数如下:
  1. function [www,wthetaL0,FlagRez]=Rez4mThick00(lambda,RezPara)
  2. % RezPara=[rho1,d1,F1,lenF1,nF1,d3,F2,lenF2,nF2,d2,rho2];
  3. rho1=    RezPara(01);d1=    RezPara(02);F1=     RezPara(03);
  4. lenF1=   RezPara(04);nF1=   RezPara(05);d3=     RezPara(06);
  5. F2=      RezPara(07);lenF2= RezPara(08);nF2=    RezPara(09);
  6. d2=      RezPara(10);rho2=  RezPara(11);
  7. %% 不绘图,用于大量计算
  8. syms z
  9. %总腔长
  10. LL=d1+lenF1+d3+lenF2+d2;
  11. [~,MF1,~,~]=deF2Rho(F1,lenF1,nF1);
  12. if (F1==F2 && lenF1==lenF2 && nF1==nF2)
  13.     MF2=MF1;
  14. else
  15.     [~,MF2,~,~]=deF2Rho(F2,lenF2,nF2);
  16. end
  17. %
  18. MM=[1,d2;0,1]*MF2*[1,d3;0,1]*MF1*[1,d1;0,1];
  19. a=MM(1,1);b=MM(1,2);c=MM(2,1);d=MM(2,2);
  20. %令
  21. G1=a-b/rho1;
  22. G2=d-b/rho2;
  23. FlagRez(1)=G1*G2;
  24. %% 是否满足0<G1*G2<1稳定条件
  25. % figure(22);
  26. if ( FlagRez(1)<0  ||  FlagRez(1)>1 )
  27.     www=0;
  28.     wthetaL0=zeros(1,3);
  29. else
  30.     %成立稳定性条件
  31.     %0<G1*G2<1
  32.     %% 计算
  33.     %靠近两腔镜束腰大小
  34.     w001=abs(sqrt(+(lambda*b/pi)*(sqrt(G1*G2*(1-G1*G2))/abs(G1+a*a*G2-2*a*G1*G2))));
  35.     w002=abs(sqrt(+(lambda*b/pi)*(sqrt(G1*G2*(1-G1*G2))/abs(G2+d*d*G1-2*d*G1*G2))));
  36.     %束腰位置
  37.     L001=(b*G2*(a-G1))/(G1+a*a*G2-2*a*G1*G2);
  38.     L002=LL-(  (b*G1*(d-G2))/(G2+d*d*G1-2*d*G1*G2)  );
  39.     %基模发散角
  40.     theta001=sqrt((lambda/(pi*w001))^2);
  41.     theta002=sqrt((lambda/(pi*w002))^2);
  42.     %% 腔内束宽
  43.     %腔内的束宽变换
  44.     %由镜1全反镜向输出镜变换
  45.     [w003,theta003,L003]=fLRMm(w001,theta001,L001,lambda,MF1,lenF1,d1);
  46.     %
  47.     w001_z=fwz(w001,theta001,L001);
  48.     w003_z=fwz(w003,theta003,L003);
  49.     w002_z=fwz(w002,theta002,L002);
  50.     www=1/4*double(subs(w001_z,z,d1)+subs(w003_z,z,d1+lenF1)+...
  51.         subs(w003_z,z,d1+lenF1+d3)+subs(w002_z,z,d1+lenF1+d3+lenF2));
  52.     wthetaL0=double([w002,theta002,L002]);
  53.    
  54. end
  55. end
  56. %% 版本信息
  57. % 作者:                Quincy Howard
  58. % 联系方式:           quincy.hd@qq.com
  59. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  60. % 若使用请注明来源
  61. % 最后编辑于           2024 年 07 月 10 日
复制代码

  1. function [Rho,MF,MFront,MBack]=deF2Rho(F,lenF,nF)
  2.     % F 厚透镜有效焦距
  3.     % lenF 厚透镜长度
  4.     % nF 厚透镜介质折射率
  5.     n1=1;
  6.     % 曲率表示
  7.     syms phox
  8.     MFront=[1,0;(nF-n1)/(nF*-phox),n1/nF];
  9.     Mmedian=[1,lenF;0,1];
  10.     MBack=[1,0;(n1-nF)/(n1*phox),nF/n1];
  11.     Mf1=MBack*Mmedian*MFront;
  12.     %求解代换
  13.     phoxx=vpasolve(Mf1(2,1)==-1/F);
  14.     [~,nn]=max(abs(phoxx));
  15.     Rho=phoxx(nn);
  16.     % % % % % % % %
  17.     % 多个解??!!
  18.     % % % % % % % %
  19.    
  20.     MFront=double([1,0;(nF-n1)/(nF*-Rho),n1/nF]);
  21.     Mmedian=double([1,lenF;0,1]);
  22.     MBack=double([1,0;(n1-nF)/(n1*Rho),nF/n1]);
  23.     MF=MBack*Mmedian*MFront;
  24. end
  25. %% 版本信息
  26. % 作者:                Quincy Howard
  27. % 联系方式:           quincy.hd@qq.com
  28. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  29. % 若使用请注明来源
  30. % 最后编辑于           2024 年 07 月 10 日
复制代码
  1. %% 高阶高斯光束由左向右传输
  2. %% [w02m,theta02m,L02m]=fLRMm(w01m,theta01m,L01m,lambda,cM,lencM,dcM)
  3. function [w02m,theta02m,L02m]=fLRMm(w01m,theta01m,L01m,lambda,cM,lencM,dcM)
  4. %         高阶转基模
  5.         [w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda);
  6. %         基模由左向右传输
  7.         [w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM);
  8. %         基模转高阶
  9.         [w02m,theta02m,L02m,~]=fB2M(w02,theta02,L02,M2);
  10. end
  11. %% 版本信息
  12. % 作者:                Quincy Howard
  13. % 联系方式:           quincy.hd@qq.com
  14. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  15. % 若使用请注明来源
  16. % 最后编辑于           2024 年 07 月 10 日
复制代码
  1. %% 高阶模转基模
  2. %% [w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda)
  3. function [w01,theta01,L01,M2]=fM2B(w01m,theta01m,L01m,lambda)
  4.     % w01m为高阶光束光腰 半径
  5.     % theta01m为高阶光束  半角发散角
  6.     % L01m为高阶光束光腰位置
  7.     % M2为光束质量
  8.     % w01为基模输出光腰 半径
  9.     % theta01为基模光束  半角发散角
  10.     % L01为基模光束光腰位置
  11.    
  12.     %默认使用1064nm
  13.     % lambda=1064*10^-9;
  14.     %默认在空气中折射率nn=1
  15.     nn=1;
  16.    
  17.     M2=(w01m*theta01m*pi*nn)/(lambda);
  18.     w01=(w01m)/sqrt(M2);
  19.     theta01=(theta01m)/sqrt(M2);
  20.     L01=L01m;
  21.    
  22. end
  23. %% 版本信息
  24. % 作者:                Quincy Howard
  25. % 联系方式:           quincy.hd@qq.com
  26. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  27. % 若使用请注明来源
  28. % 最后编辑于           2024 年 07 月 10 日
复制代码
  1. %% 基模由左向右传输
  2. %% [w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM)
  3. function [w02,theta02,L02]=fLRM(w01,theta01,L01,lambda,cM,lencM,dcM)
  4. % 默认使用基模计算,由左向右传输
  5. % w01,theta01,L01 为入射光束参数
  6. % lambda 为波长
  7. % cM 为光学系统的传输矩阵
  8. % lencM 为光学系统长度
  9. % dcM 为 入射光束与光学系统的第一接触面 相对于 原点所处位置
  10. % w02,theta02,L02 为出射光束参数
  11. %% 物方光束
  12. Z0=w01/theta01;
  13. syms z;
  14. w_1_z=w01*sqrt(1+((z-L01)/Z0)^2);
  15. R_z1=Z0*((z-L01)/Z0+Z0/(z-L01));
  16. R01=subs(R_z1,z,eps);
  17. w_001=subs(w_1_z,z,0);
  18. %% 传输矩阵相关定义及计算
  19. syms X1 Y1 z
  20. % syms R_1
  21. % q参数定义
  22. % syms s1 s2 a b c d
  23. % 取n1=n2=1;
  24. n1=1;n2=1;
  25. %传输矩阵元素
  26. s2=z;
  27. a=cM(1,1);b=cM(1,2);c=cM(2,1);d=cM(2,2);
  28. s1=dcM;
  29. %其他特征参量
  30. M=[1 s2;0,1]*[a b;c d]*[1,s1;0,1];
  31. A=M(1,1);B=M(1,2);
  32. % C=M(2,1);D=M(2,2);
  33. Y2=(n1/n2)*Y1/...
  34.     (A^2+2*X1*A*B+(X1^2+Y1^2)*B^2);
  35. %% 像方光束
  36. w_2=sqrt(lambda/(pi*Y2));
  37. w_2=subs(w_2,[X1 Y1],[1/R01 lambda/(pi*w_001^2)]);
  38. %光腰参考点平移
  39. w_2=subs(w_2,(z),(z-s1));
  40. %% 输出参数
  41. % disp('w_2*w_2=');pretty(w_2*w_2);
  42. T1=vpa((taylor((w_2*w_2),z,'Order',1)));
  43. T2=vpa((taylor((w_2*w_2),z,'Order',2)));
  44. T3=vpa((taylor((w_2*w_2),z,'Order',3)));
  45. C=T1;
  46. B=(T2-T1)/z;
  47. A=(T3-T2)/(z*z);
  48. w02=vpa(sqrt(C-B*B/(4*A)));
  49. theta02=vpa(sqrt(A));
  50. L02=vpa(-B/(2*A));
  51. L02=L02+lencM;
  52. end
  53. %% 版本信息
  54. % 作者:                Quincy Howard
  55. % 联系方式:           quincy.hd@qq.com
  56. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  57. % 若使用请注明来源
  58. % 最后编辑于           2024 年 07 月 10 日
复制代码
  1. %% 基模转高阶模
  2. %% [w01m,theta01m,L01m,M2]=fB2M(w01,theta01,L01,M2)
  3. function [w01m,theta01m,L01m,M2]=fB2M(w01,theta01,L01,M2)
  4.     % w01m为高阶光束光腰 半径
  5.     % theta01m为高阶光束  半角发散角
  6.     % L01m为高阶光束光腰位置
  7.     % M2为光束质量
  8.     % w01为基模输出光腰 半径
  9.     % theta01为基模光束  半角发散角
  10.     % L01为基模光束光腰位置
  11.     %     M2=(w01m*theta01m*pi*nn)/(lambda)
  12.     % 默认在空气中折射率nn=1
  13.     %     nn=1;
  14.     %     lambda=(w01m*theta01m*pi*nn)/(M2);
  15.     w01m=(w01)*sqrt(M2);
  16.     theta01m=(theta01)*sqrt(M2);
  17.     L01m=L01;
  18.    
  19. end
  20. %% 版本信息
  21. % 作者:                Quincy Howard
  22. % 联系方式:           quincy.hd@qq.com
  23. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  24. % 若使用请注明来源
  25. % 最后编辑于           2024 年 07 月 10 日
复制代码
  1. %% 光束随距离变化
  2. %% w_z = fwz(w0,theta0,L0,z)
  3. function w_z = fwz(w0,theta0,L0,z)
  4.     if nargin==3
  5.         syms z;
  6.     end
  7.    
  8.     w_z=sqrt(w0^2+theta0^2*(z-L0)^2);
  9. end
  10. %% 版本信息
  11. % 作者:                Quincy Howard
  12. % 联系方式:           quincy.hd@qq.com
  13. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  14. % 若使用请注明来源
  15. % 最后编辑于           2024 年 07 月 10 日
复制代码
 next脚本:
  1. close all;
  2. clear;
  3. clc;
  4. %% 版本信息
  5. % 作者:                Quincy Howard
  6. % 联系方式:           quincy.hd@qq.com
  7. % 文件信息更新平台为   matlab云文件、GitHub、CSDN博客、知乎
  8. % 若使用请注明来源
  9. % 最后编辑于           2024 年 07 月 10 日
复制代码


免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!更多信息从访问主页:qidao123.com:ToB企服之家,中国第一个企服评测及商务社交产业平台。
回复

使用道具 举报

0 个回复

倒序浏览

快速回复

您需要登录后才可以回帖 登录 or 立即注册

本版积分规则

拉不拉稀肚拉稀

金牌会员
这个人很懒什么都没写!

标签云

快速回复 返回顶部 返回列表