简介这份资源是面向激光物理、光学工程方向的学习者与科研人员整理的MATLAB仿真代码包聚焦高斯光束传输与谐振腔设计这一激光系统建模中的核心问题。包内共95个文件以67个m脚本为主体配合17个fig图形结果、5个mat数据文件以及docx说明、pdf文档、md与txt说明等压缩包约27.37MB覆盖从基础计算到结果可视化的完整流程。内容围绕高斯光束传播模型、热效应分析、非线性频率变换、光束调节与匹配算法、稳定驻波谐振腔与行波谐振腔仿真、ABCD矩阵计算等模块展开可帮助读者理解束宽变化、焦点位置、热透镜效应、倍频混频效率及模式匹配等关键环节。已有122人学习适合需要借助MATLAB快速搭建光路模型、验证谐振腔稳定性或优化激光输出参数的读者参考使用。1. 从一束光说起MATLAB 里把激光高斯光束传输与谐振腔仿真跑通很多做激光系统的人都有个习惯先在纸上画光路再在软件里把光路“跑”一遍。问题是纸上算高斯光束的束腰、发散角、ABCD 矩阵手推几页公式还行一旦涉及热效应、非线性频率变换、驻波腔和行波腔的联合仿真手算基本就翻车了。这套基于 MATLAB 的激光高斯光束传输与谐振腔仿真代码包解决的正是这个痛点用脚本把高斯光束传播模型、ABCD 矩阵计算、谐振腔稳定性判断、热效应分析和非线性频率变换串成一条可复现的链路。它适合做固体激光器、光纤激光器、倍频与和频系统、腔体设计的工程师和研究生也适合想用 MATLAB 把光学公式变成可调参数工具的人。下面我按自己搭这套东西的顺序把关键模型、参数和踩过的坑讲清楚。2. 高斯光束传播模型从束腰到 ABCD 矩阵的完整链路2.1 为什么先用 q 参数而不是直接追光线几何光学追光线在自由空间里是直线到了透镜就折一下看起来简单。但激光束不是一条线它是有束腰、有曲率半径、有发散角的复参数光束。高斯光束的复参数 q 把束腰位置、瑞利长度、波前曲率全打包在一起传播和经过光学元件时只需要做一次代数变换。常见做法是自由空间传播距离 Lq₂ q₁ L薄透镜焦距 f1/q₂ 1/q₁ − 1/f球面镜反射等效于焦距 f R/2 的透镜这套规则的好处是任何由自由空间、透镜、镜片组成的轴对称光路都能写成一个 2×2 的 ABCD 矩阵然后对 q 做分式线性变换。MATLAB 里实现起来就是矩阵乘法和复数运算不需要解波动方程。2.2 用 MATLAB 写一个最小可跑的高斯光束传播脚本下面这段代码是我一般用来验证光路是否搭对的最小示例。它定义了一个束腰在 z0 处的高斯光束然后让它经过一段自由空间、一个薄透镜再传播一段距离最后输出束腰半径和波前曲率。% 高斯光束 q 参数传播最小示例 lambda 1064e-9; % 波长单位 m w0 0.5e-3; % 束腰半径单位 m zR pi * w0^2 / lambda; % 瑞利长度 q0 1i * zR; % 束腰处 q 参数 % 光学系统自由空间 L1 - 薄透镜 f - 自由空间 L2 L1 0.2; % 第一段自由空间单位 m f 0.1; % 透镜焦距单位 m L2 0.3; % 第二段自由空间单位 m % 自由空间 ABCD 矩阵 M_free (L) [1, L; 0, 1]; % 薄透镜 ABCD 矩阵 M_lens (f) [1, 0; -1/f, 1]; M M_free(L2) * M_lens(f) * M_free(L1); q_out (M(1,1)*q0 M(1,2)) / (M(2,1)*q0 M(2,2)); % 从 q 参数反推束腰半径和曲率半径 inv_q 1 / q_out; R 1 / real(inv_q); % 波前曲率半径 w sqrt(lambda / (pi * abs(imag(inv_q)))); % 当前光束半径 fprintf(输出 q %.4f %.4fi\n, real(q_out), imag(q_out)); fprintf(光束半径 w %.4f mm\n, w*1e3); fprintf(曲率半径 R %.4f m\n, R);逻辑说明q0 是束腰处的复参数实部为 0虚部为瑞利长度。M_free 和 M_lens 分别生成自由空间和薄透镜的 ABCD 矩阵。矩阵乘法顺序是从右到左对应光先经过 L1再经过透镜最后经过 L2。q_out 用分式线性变换得到。最后从 1/q 的实部和虚部分别提取曲率半径和光束半径。参数说明lambda 换成 532 nm 就是倍频光w0 决定瑞利长度和发散角L1、f、L2 根据实际光路填写。如果 q_out 的虚部变成负数说明光束已经过了焦点这是正常现象不是错误。2.3 ABCD 矩阵计算里最容易搞混的符号和顺序ABCD 矩阵本身不难难的是符号约定和乘法顺序。我见过不少人把自由空间矩阵写成 [1, 0; L, 1]结果光束越传越收敛明显不对。自由空间是 [1, L; 0, 1]透镜是 [1, 0; -1/f, 1]球面镜反射等效焦距 f R/2其中 R 是曲率半径凹面镜 R 取正凸面镜 R 取负。乘法顺序是后经过的元件放在左边先经过的放在右边。如果光路里有反射要注意传播方向反转后 q 参数的符号处理常见做法是把反射后的光路重新定义坐标方向而不是直接套用透射矩阵。提示每次改完光路先用一个已知的简单例子验证比如束腰在透镜前焦面经过透镜后束腰应该在透镜后焦面且大小满足 w0 lambdaf/(piw0)。如果这个都对不上后面的谐振腔和热效应分析就不用看了。3. 稳定驻波谐振腔设计从 g 参数到自洽场3.1 驻波腔的稳定性判据和 g 参数驻波腔由两面镜子组成光在两面镜子之间来回反射。稳定腔的条件是 0 g1*g2 1其中 g1 1 - L/R1g2 1 - L/R2L 是腔长R1、R2 是两面镜子的曲率半径。这个判据的本质是光在腔内往返一次后光束参数必须回到自身否则光会跑出腔外。MATLAB 里判断稳定性只需要几行代码但实际设计时还要考虑镜面热透镜效应因为高功率激光的镜面会吸收部分光产生热变形等效于改变镜子的曲率半径。3.2 用往返矩阵求稳定腔的基模参数稳定腔的基模光束参数可以通过往返矩阵的自洽条件求出。往返矩阵 M_rt 是光在腔内走一圈的 ABCD 矩阵自洽条件是 q (Aq B)/(Cq D)。解这个方程可以得到 q 的解析表达式然后反推束腰位置和大小。% 驻波腔基模参数计算 L 0.5; % 腔长单位 m R1 1.0; % 镜 1 曲率半径单位 m R2 1.0; % 镜 2 曲率半径单位 m lambda 1064e-9; % 波长 % 单程矩阵从镜 1 到镜 2 M_half [1, L; 0, 1] * [1, 0; -2/R1, 1]; % 往返矩阵镜 1 - 镜 2 - 镜 1 M_rt [1, 0; -2/R2, 1] * M_half; A M_rt(1,1); B M_rt(1,2); C M_rt(2,1); D M_rt(2,2); % 稳定性判据 g1 1 - L/R1; g2 1 - L/R2; if g1*g2 0 g1*g2 1 fprintf(腔是稳定的g1*g2 %.4f\n, g1*g2); else fprintf(腔不稳定g1*g2 %.4f\n, g1*g2); end % 自洽解 q (A-D)/(2C) ± sqrt(((A-D)/(2C))^2 - B/C) disc ((A-D)/(2*C))^2 - B/C; q_sol (A-D)/(2*C) sqrt(disc); inv_q 1/q_sol; w sqrt(lambda / (pi * abs(imag(inv_q)))); fprintf(基模光束半径 w %.4f mm\n, w*1e3);逻辑说明M_half 是从镜 1 到镜 2 的单程矩阵包含自由空间传播和镜 1 的反射。M_rt 是往返矩阵在单程基础上再加镜 2 的反射和反向传播。自洽解用二次方程求根取物理上合理的那个根。最后从 q 的虚部提取光束半径。参数说明L 是腔长R1、R2 是镜面曲率半径lambda 是波长。如果 g1*g2 接近 0 或 1腔处于临界稳定实际中容易受扰动失稳一般会留一定余量。镜面热透镜效应可以通过把 R1、R2 替换为等效热透镜焦距来近似。3.3 热效应分析镜面变形和热透镜怎么进模型高功率激光腔里镜面吸收的光会让镜面局部升温产生热膨胀和折射率变化等效于一个热透镜。常见做法是把热透镜当作一个薄透镜放在镜面附近焦距由吸收功率、镜面材料和冷却条件决定。MATLAB 里可以先用一个经验公式估算热透镜焦距然后把它代入 ABCD 矩阵重新计算腔的稳定性和基模参数。% 热透镜焦距估算与腔参数修正 P_abs 10; % 镜面吸收功率单位 W kappa 1.2; % 材料热导率单位 W/(m*K) dn_dT 2.5e-6; % 折射率温度系数单位 1/K alpha 1e-6; % 热膨胀系数单位 1/K w_pump 1e-3; % 泵浦光半径单位 m % 经验公式热透镜焦距 f_th pi*kappa*w_pump^2 / (P_abs*(dn_dT alpha)) f_th pi * kappa * w_pump^2 / (P_abs * (dn_dT alpha)); fprintf(热透镜焦距 f_th %.4f m\n, f_th); % 把热透镜插入腔矩阵 M_th [1, 0; -1/f_th, 1]; M_rt_th M_th * M_rt; % 重新计算稳定性和基模参数步骤同上逻辑说明热透镜焦距的经验公式把吸收功率、材料热导率、折射率温度系数和热膨胀系数放在一起给出一个量级估算。实际中这个公式需要根据具体材料和冷却条件修正但作为初步设计足够。把热透镜矩阵乘到往返矩阵里就能看到热效应对腔稳定性的影响。参数说明P_abs 是镜面实际吸收的功率不是入射功率。kappa、dn_dT、alpha 是材料参数不同玻璃和晶体差别很大。w_pump 是泵浦光在镜面上的半径越小热透镜越强。如果算出来的 f_th 和腔长可比说明热效应不能忽略。注意热透镜焦距对吸收功率非常敏感P_abs 翻倍f_th 减半。实际中镜面吸收功率很难精确知道一般会留 20% 到 50% 的余量或者用可调热透镜模型做灵敏度分析。4. 行波谐振腔仿真环形腔里的方向性和模式匹配4.1 行波腔和驻波腔的本质区别驻波腔里光来回反射形成驻波腔内有波节和波腹。行波腔通常是环形结构光沿一个方向循环传播没有驻波效应。行波腔的好处是避免了空间烧孔效应适合需要单频运转的激光器。但行波腔的稳定性分析比驻波腔复杂因为光在腔内走一圈要经过多个镜片和可能的光学元件ABCD 矩阵的维度更高。MATLAB 里处理行波腔的思路和驻波腔类似只是往返矩阵变成绕腔一圈的单程矩阵。4.2 环形腔的 ABCD 矩阵和模式匹配环形腔一般由多个镜片组成光绕一圈回到起点。模式匹配的目标是让光在腔内循环时光束参数在每一点都保持合理尤其是增益介质和非线性晶体所在的位置光束半径要匹配泵浦光或满足相位匹配条件。% 环形腔模式匹配示例四镜环形腔 L1 0.3; L2 0.2; L3 0.3; L4 0.2; % 四段自由空间 R1 0.5; R2 0.5; R3 0.5; R4 0.5; % 四面镜曲率半径 lambda 1064e-9; % 绕腔一圈的 ABCD 矩阵 M [1, 0; -2/R4, 1] * [1, L4; 0, 1] * ... [1, 0; -2/R3, 1] * [1, L3; 0, 1] * ... [1, 0; -2/R2, 1] * [1, L2; 0, 1] * ... [1, 0; -2/R1, 1] * [1, L1; 0, 1]; % 自洽解 A M(1,1); B M(1,2); C M(2,1); D M(2,2); disc ((A-D)/(2*C))^2 - B/C; q_sol (A-D)/(2*C) sqrt(disc); inv_q 1/q_sol; w sqrt(lambda / (pi * abs(imag(inv_q)))); fprintf(环形腔基模光束半径 w %.4f mm\n, w*1e3); % 检查增益介质位置的光束半径 % 假设增益介质在 L1 中间传播 L1/2 后计算 M_gain [1, L1/2; 0, 1]; q_gain (M_gain(1,1)*q_sol M_gain(1,2)) / (M_gain(2,1)*q_sol M_gain(2,2)); w_gain sqrt(lambda / (pi * abs(imag(1/q_gain)))); fprintf(增益介质处光束半径 w_gain %.4f mm\n, w_gain*1e3);逻辑说明环形腔的 ABCD 矩阵按光路顺序从右到左乘每个镜片用薄透镜矩阵等效自由空间用传播矩阵。自洽解和驻波腔一样用二次方程求根。增益介质处的光束半径通过从参考面传播到增益介质位置再计算。参数说明L1 到 L4 是四段自由空间长度R1 到 R4 是四面镜的曲率半径。如果增益介质处光束半径太小功率密度过高可能损伤晶体太大则泵浦光重叠效率低。一般会调整镜片曲率半径和距离让增益介质处光束半径在 0.2 到 0.5 mm 之间。4.3 行波腔里的非线性频率变换位置选择非线性频率变换比如倍频和和频对光束半径和相位匹配非常敏感。倍频效率与基频光功率密度的平方成正比所以通常会把非线性晶体放在光束半径最小的位置也就是束腰附近。但束腰太小会导致晶体损伤太大则效率低。MATLAB 里可以扫描晶体位置计算不同位置的光束半径和对应的倍频效率找到最优位置。% 扫描非线性晶体位置找最优倍频效率 positions linspace(0, L1, 50); % 在 L1 段内扫描 eta zeros(size(positions)); for k 1:length(positions) M_pos [1, positions(k); 0, 1]; q_pos (M_pos(1,1)*q_sol M_pos(1,2)) / (M_pos(2,1)*q_sol M_pos(2,2)); w_pos sqrt(lambda / (pi * abs(imag(1/q_pos)))); % 简化倍频效率模型与 1/w_pos^2 成正比但受损伤阈值限制 eta(k) 1 / (w_pos^2 (0.3e-3)^2); % 0.3 mm 为损伤阈值对应的最小半径 end [~, idx] max(eta); fprintf(最优晶体位置 %.4f m光束半径 %.4f mm\n, ... positions(idx), w_pos*1e3);逻辑说明扫描晶体在 L1 段内的位置每次计算该位置的光束半径然后用一个简化模型估算倍频效率。模型里加了一个最小半径限制模拟损伤阈值。实际中倍频效率还和相位匹配角、晶体长度、走离效应有关但位置扫描是第一步。参数说明positions 是扫描范围一般覆盖整个自由空间段。损伤阈值对应的最小半径根据晶体材料不同常见 KTP 和 LBO 的损伤阈值在 GW/cm² 量级对应光束半径在 0.1 到 0.5 mm 之间。如果最优位置在镜片附近要考虑镜片像差和热效应的影响。5. 避坑与排查高斯光束和谐振腔仿真里最容易翻车的 5 个地方5.1 现象q 参数传播后光束半径变成 NaN 或无穷大原因ABCD 矩阵里出现了除零或者 q 参数在某个位置正好等于 −D/C导致分式线性变换的分母为零。常见于透镜焦距设为零或者腔长正好等于镜片曲率半径的两倍。解决检查所有光学元件的参数确保焦距和曲率半径不为零。在代码里加一个判断如果分母绝对值小于 1e-12就报错并输出当前矩阵和 q 值。另外腔长等于 2R 时腔处于临界稳定实际中应该避开这个点。5.2 现象驻波腔稳定性判据算出来是稳定的但实际光路不出光原因稳定性判据只考虑了几何光学没有考虑热透镜、镜面吸收和增益介质的折射率变化。高功率泵浦下增益介质本身就是一个强热透镜等效于改变腔的 g 参数。解决把增益介质的热透镜矩阵加入往返矩阵重新计算 g1*g2。如果热透镜焦距和腔长可比腔可能从稳定变成不稳定。常见做法是先用低功率泵浦测出冷腔参数再逐步增加功率观察输出功率和光束质量的变化。5.3 现象行波腔模式匹配算出来的光束半径和实验测的差很多原因环形腔的镜片不是理想薄透镜有像差和反射相位延迟。另外光在腔内循环时镜片的实际曲率半径可能因为安装应力或热变形而改变。解决先用一个已知曲率半径的标准镜片验证 ABCD 矩阵确认符号和乘法顺序正确。然后测量实际镜片的曲率半径或者用干涉仪测波前。如果差得不多可以在模型里加一个等效焦距修正项通过实验数据反推。5.4 现象非线性频率变换效率扫描结果和实验趋势相反原因简化模型只考虑了光束半径忽略了相位匹配和走离效应。如果晶体位置改变了相位匹配角效率可能反而下降。另外光束半径太小会导致晶体损伤实验上不敢用最优位置。解决在扫描模型里加入相位匹配角和晶体长度的修正或者直接用实验测得的效率曲线反推最优位置。如果损伤阈值限制严格就把最小半径设大一点牺牲一点理论效率换安全。5.5 现象MATLAB 脚本跑得慢扫描几百个点要等很久原因每次循环都重新计算矩阵乘法和自洽解没有向量化。MATLAB 的循环效率不高尤其是涉及复数运算时。解决把扫描参数写成向量用矩阵运算一次性计算所有位置。比如把 positions 写成行向量用 arrayfun 或者直接矩阵乘法批量处理。另外预分配数组避免在循环里动态增长。如果还是慢可以把自洽解写成函数用 parfor 并行计算。提示仿真结果和实验对不上时先检查单位。波长用米束腰用米焦距用米不要混用毫米和米。我见过有人把 w0 写成 0.5 毫米但 lambda 用 1064 纳米结果瑞利长度算出来差三个数量级。6. 进阶技巧用 MATLAB 把仿真和实验数据串起来6.1 用实测光束半径反推 q 参数和热透镜焦距仿真做得再好最终还是要和实验对上。我一般会先用 CCD 相机测几个位置的光束半径然后用最小二乘拟合高斯光束的传播公式反推束腰位置、束腰半径和热透镜焦距。MATLAB 的 lsqcurvefit 或者 fminsearch 都能做这件事。% 用实测光束半径拟合高斯光束参数 z_data [0.1, 0.2, 0.3, 0.4, 0.5]; % 测量位置单位 m w_data [0.4, 0.35, 0.32, 0.35, 0.42] * 1e-3; % 实测半径单位 m lambda 1064e-9; % 拟合模型w(z) w0 * sqrt(1 ((z - z0)/zR)^2) model (p, z) p(1) * sqrt(1 ((z - p(2))/(pi*p(1)^2/lambda)).^2); p0 [0.3e-3, 0.3]; % 初始猜测w0 和 z0 p_fit lsqcurvefit(model, p0, z_data, w_data); fprintf(拟合束腰半径 w0 %.4f mm\n, p_fit(1)*1e3); fprintf(拟合束腰位置 z0 %.4f m\n, p_fit(2));逻辑说明model 是高斯光束传播公式p(1) 是束腰半径p(2) 是束腰位置。lsqcurvefit 用实测数据拟合这两个参数。拟合完可以反推瑞利长度和发散角。参数说明z_data 和 w_data 是实测数据测量点至少 5 个覆盖束腰前后。如果拟合残差大可能是热透镜效应导致光束不是理想高斯需要在模型里加修正项。6.2 把热透镜焦距作为拟合参数如果实测数据偏离理想高斯可以把热透镜焦距也作为拟合参数用 ABCD 矩阵模型拟合。这样能同时得到冷腔参数和热透镜焦距。% 用 ABCD 矩阵模型拟合热透镜焦距 function w_model abcd_model(p, z_data, lambda, L, R1, R2) f_th p(1); % 热透镜焦距 w0 p(2); % 冷腔束腰 z0 p(3); % 束腰位置 q0 1i * pi * w0^2 / lambda; q0 q0 z0; % 把束腰移到 z0 w_model zeros(size(z_data)); for k 1:length(z_data) M [1, z_data(k); 0, 1] * [1, 0; -1/f_th, 1]; q (M(1,1)*q0 M(1,2)) / (M(2,1)*q0 M(2,2)); w_model(k) sqrt(lambda / (pi * abs(imag(1/q)))); end end逻辑说明这个模型把热透镜放在腔的起点然后传播到测量位置。拟合参数包括热透镜焦距、冷腔束腰和束腰位置。实际使用时可以把腔的 ABCD 矩阵完整写进去而不是只用一个薄透镜近似。参数说明f_th 初始猜测可以根据经验公式给w0 和 z0 根据冷腔测量给。如果拟合不收敛先固定 f_th拟合 w0 和 z0再放开 f_th。6.3 一个我常用的习惯每次改光路先跑一遍 ABCD 矩阵检查我自己的习惯是不管光路多复杂先在 MATLAB 里把 ABCD 矩阵写出来跑一遍 q 参数传播看看光束半径在关键位置是否合理。如果某个位置光束半径突然变得很大或很小说明矩阵符号或顺序有问题。这个检查只需要几秒钟但能省掉后面几个小时的调试。仿真和实验对不上时先怀疑单位再怀疑符号最后怀疑热效应。希望帮到你。本文还有配套的精品资源点击获取