1. 项目概述:为什么这个光学仿真值得做
马赫-曾德干涉仪是我在光学工程里打交道最多的结构之一。它原理不复杂:一束光被分束器分成两路,经过不同的光程后再合束,形成干涉条纹。可一旦涉及到实际应用——无论是测量折射率变化、检测微小位移,还是做光学相干断层扫描,你需要的不再只是“有条纹”,而是条纹背后的相位信息。这就要求把波动光学模型吃透,用仿真把电场传播、叠加、相位调制这些过程完整跑一遍,才能在实验前预判结构、光路参数和误差来源。
这个项目就是基于上面这个需求来的:用波动的视角,把马赫-曾德干涉仪从光束生成、分束、相位延迟到合束干涉的完整流程在计算机里复现,并做参数扫描和误差分析。我当时的直接目标是两个:第一,得到一套干净的干涉条纹图,能直观看到光程差和条纹周期的关系;第二,能通过仿真反推相位变化,验证可见度计算公式的适用条件,为后来做实验装置校准提供参考依据。
如果你正准备做光学仿真,或者手里有干涉测量实验却拿不准参数,这篇文章应该能帮你省掉不少试错时间。我假设你懂一点光的干涉概念,但不需要提前熟悉仿真软件,所有操作我会按步骤拆开讲,包括选型思路和踩坑记录。
1.1 需求拆解:不只是画两条干涉条纹
很多人听到“马赫-曾德干涉仪仿真”第一反应是:这不就是两束光叠加画个图吗?实际动手之后会发现,问题远没有这么简单。
首先是物理模型的取舍:用几何光学还是波动光学?几何光学可以告诉你光线怎么走,但解释不了干涉条纹的对比度、光强分布和衍射效应,所以做干涉仪仿真必须回到波动光学,直接对复振幅场做叠加。
其次是数值实现的问题:连续的光场在计算机里要用离散网格表示,网格大小决定了能不能分辨出条纹细节;采样太粗条纹直接糊掉,采样太细计算量爆炸。光束需要建在二维网格上,还是用一维剖面?束腰半径、波长、网格间距、光程差这些参数互相牵扯,任何一个不合适,仿真结果都会失真。
最后是分析层的需求:干涉条纹只是原始数据,我们需要从中提取可见度、相位变化、条纹间距等量化指标。这些指标和相位调制深度、分束器分束比、光源相干性之间的关系,都是这个项目要回答的问题。
所以这个项目本质上是三件事:建模、仿真、分析。建模决定物理上对不对,仿真决定数值上稳不稳,分析决定结论上有没有用。
1.2 谁适合学,能解决什么问题
这个项目适合三类人。
第一类是光学工程专业的学生,课程里学到干涉仪,但实验室设备贵、调试周期长,可以先在仿真里把干涉现象和参数关系摸清楚。
第二类是搞精密测量的工程师,需要利用马赫-曾德干涉仪做折射率、温度、振动测量,仿真可以帮助确定光程差范围、波长选择和探测器分辨率需求。
第三类是刚入门光学仿真但手里没有现成模板的人。因为这个项目不依赖商业软件,用Python加几个常用库就能跑通,很多数值处理技巧在别的光学仿真里同样适用。
我自己在这个项目里最大的收获是理解了一个道理:仿真不是“画个图就完事”,所有肉眼可见的现象背后都有一个数值原因,而这个数值原因又对应一个物理量。把这条因果链捋顺,实验里遇到诡异现象就不会慌。
2. 波动光学核心原理:先把干涉这件事说透
2.1 相干叠加与相位差的数学表达
马赫-曾德干涉仪的本质是两个相干光场的叠加。设分束后两路光的复振幅分别为E1和E2,合束处的总场强用E1加E2来算,光强I等于总场强的模平方。
平时我们看到的干涉公式长这样:I = I1 + I2 + 2·sqrt(I1·I2)·cos(Δφ)。其中Δφ是两路光的相位差。如果分束器是50:50,那么I1和I2理论上都是入射光强的四分之一?不,这里有个细节容易搞混。分束器是在振幅域分光的,功率按比例分配,但后续合束时两束光的振幅接近,干涉项强度才能最大。
相位差从哪来?主要来自几何光程差ΔL,那么Δφ = 2π·ΔL/λ。如果一路还经过电光调制器或者声光调制器,还要加上动态相位项。仿真里最容易控制的就是在某一臂上直接乘一个复指数因子exp(i·φ),这样扫相位非常方便。
从这个公式能看到几个关键结论:
- 当Δφ = 2mπ时,出现亮条纹,光强最大。
- 当Δφ = (2m+1)π时,出现暗条纹,光强最小。
- 两路光的偏振方向必须是平行的,投影分量才可能完全干涉,否则干涉对比度会下降。
- 光源的相干长度要远大于光程差,否则条纹可见度衰减,这一点在仿真里经常被忽略,尤其是在后续把宽谱光源纳入分析时。
关于相位差还有一个细节:如果两路光不经过任何器件,只是纯几何路径不同,ΔL就是两臂物理长度差。但如果在某一臂放了玻璃片、液体池或光纤,光程还要乘上对应介质的折射率。我就在这个点上吃过亏,仿真里光程差改了但忘了折射率项,结果条纹周期对不上。
2.2 从理想平面波到真实高斯光束
很多教材上的干涉分析默认光源是无限大平面波,但这个假设在仿真里非常危险。
平面波的波前是平的,叠加后条纹在空间里是等间距的直线,干净、容易分析。但真实激光器输出的是高斯光束,振幅在横截面上呈高斯分布,波前在近场和远场都有曲率变化。如果直接用平面波建模,你会得到一个过于理想化的干涉图,实验里看到的边缘光强衰减根本没有体现。
高斯光束的电场表达式是E(r, z) = E0·(w0/w(z))·exp(-r²/w(z)²)·exp(-i·k·r²/(2·R(z)))·exp(-i·k·z)……这里w(z)是光斑半径,R(z)是波前曲率半径,z是传播距离。
这个表达式看着吓人,但仿真时可以简化。如果干涉仪两臂的光束传播距离较短,且不考虑严格的发散效应,可以用一个二维高斯函数作为初始复振幅分布,再乘上传播相位exp(-i·k·z)来表示。这么做有一个好处:计算快,参数直观,而且能看到光束尺寸对条纹覆盖范围的影响。
不过要提醒一句:如果你仿真的目标是研究长距离传播或者需要精确的衍射行为,那就不能用这种简化模型,最好还是走角谱传播或者菲涅尔衍射积分。我的建议是先用简单模型把物理逻辑跑通,再按需升级模型复杂度。
3. 仿真工具选型与算法设计
3.1 为什么我选了Python+NumPy而不是商业软件
选工具的第一原则是看仿真规模和可扩展性。有的团队用COMSOL做光波传播,有的用Lumerical FDTD做严格电磁场仿真,这些软件精度高,但学习成本和学习门槛都比较高,而且对单次干涉模拟来说有点“杀鸡用牛刀”。
马赫-曾德干涉仪属于线性光学系统,不涉及强非线性效应,也不涉及亚波长结构散射,所以在“波动光学”框架下用标量复振幅来描述已经完全够用。基于这个判断,我选了Python加NumPy,理由有三条:
第一,代码透明,每一步物理操作都能对应到一行代码上。用商业软件经常遇到“黑盒”过程,仿真结果异常时很难定位是物理模型错了还是软件设置错了。而自己写代码,网格、分束、相位、合束每一步都能单独调试。
第二,参数扫描方便。用Python写个for循环就能批量跑不同相位差、不同分束比、不同波长,数据直接存成NumPy数组,后面画图、统计一步到位。
第三,可扩展性好。这套代码后续加偏振、加宽谱、加噪声都很方便,只需要在原来的复振幅模型上做文章,不需要重新搭平台。
我不是说商业软件不能碰。如果你要做微纳结构附近的严格场分布,FDTD类工具是必要的。但纯干涉分析阶段,标量衍射和Nelder-Mead,啊不,和数值叠加已经能给出非常好的预测。
3.2 波前构建和分束合束的建模思路
任何一个光场都可以表示成复振幅数组,这是整个仿真最底层的逻辑。比如二维网格上,每个点的值是一个复数,幅度代表电场强度,辐角代表相位。光在介质里传播一段距离,就是给复振幅乘一个相位因子。分束器操作本质上是在振幅域对光场做线性变换。
马赫-曾德干涉仪里的50:50分束器可以建模成:E_ref = E_in / sqrt(2),E_trans = E_in / sqrt(2)。这里的sqrt(2)来自能量守恒。分束器本质上是把入射场分成两个正交的输出端口,如果是理想无损50:50分束,每个端口的振幅是原来的1/√2,功率就是原来的1/2。
如果你要研究分束比不均衡的情况,就引入一个分束系数α,比如反射臂振幅系数为sqrt(α),透射臂为sqrt(1-α)。这个参数后面做误差分析时用处很大。
合束器同理,就是两个输入场E1和E2相加,再乘一个输出系数。由于合束器也是分束器结构,实际中两路光在合束器里会有一半往回反射,一半输出到探测器。仿真时只关心输出到探测器的那一路,所以总输出场可以记为E_out = (E1 + E2) / sqrt(2)。
这里有个非常关键的物理细节:合束时两个输入场必须来自同一个物理端口方向关系,否则会引入附加相位。仿真里最稳妥的做法是在两个分束器之间用坐标变换,保持E1和E2的传播方向一致,只让光程相位有差异。我在第一次建模时没注意轴向坐标,导致两束光的相位在不同位置对不上,条纹图出现了奇怪的“等倾干涉”变形,后来加上坐标对齐才解决。
4. 仿真实操:完整跑通马赫-曾德干涉仪
4.1 参数设置:波长、束腰、网格、光程差
仿真参数是决定结果质量的首要因素,我建议先从一个标准参考状态开始,逐个调参,不要一上来就扫描大范围。
我选用的参考参数是氦氖激光的波长632.8 nm,高斯光束束腰半径w0=500 μm,二维计算区域边长为5 mm,网格数1024×1024。这样每个网格对应的物理尺寸dx大约5 μm。干涉条纹的周期取决于光程差和波长,如果光程差ΔL=1 mm,那么在横向拍频下条纹周期大概是λ·L_eff/ΔL量级,这个公式不同模型下不太一样,但可以估计条纹频率是否低于采样频率。
网格划分有一条铁律:一个条纹周期内至少要采5个点,否则条纹形状和可见度会发生严重畸变。由于我们是在二维面上模拟相干叠加,如果两束光有横向夹角,条纹空间频率会很高。这时要保证dx小于条纹周期的一半,甚至更好是1/5周期。如果发现条纹过密,优先考虑减小光程差,或者增加网格数。
相位差的设置可以在参考臂上做,我给两臂设定一个基准光程差,然后动态相位phi从0变化到2π,步距取π/16。这样做的好处是能在屏幕上直观看到条纹“移动”的过程,并且能记录每个点的光强随时间(相位)变化曲线。
4.2 分束与合束的组合场计算
核心代码可以写成下面这样。这里我剔除了无关的可视化代码,只保留物理逻辑:
import numpy as np wavelength = 632.8e-9 k = 2 * np.pi / wavelength N = 1024 size = 5e-3 x = (np.arange(N) - N / 2) * (size / N) X, Y = np.meshgrid(x, x) w0 = 500e-6 E_in = np.exp(-(X**2 + Y**2) / w0**2) alpha = 0.5 E1 = np.sqrt(alpha) * E_in E2 = np.sqrt(1 - alpha) * E_in * np.exp(1j * 0.1) # 本底相位差 # 动态扫描相位 phi_values = np.linspace(0, 2 * np.pi, 32) I_record = [] for phi in phi_values: E2_mod = E2 * np.exp(1j * phi) E_out = (E1 + E2_mod) / np.sqrt(2) I_out = np.abs(E_out)**2 I_record.append(I_out) # 最后画图可以存I_out,分析用I_record这是一个基础版本,里面有两个近似需要注意。第一,E1和E2都基于同一个E_in,但实际分束前光束经过不同路径,光斑大小会有微小变化,仿真中为对比方便暂时忽略。第二,传播相位exp(-i·k·z)没有显式写出来,而是通过φ这个变量体现两路相位差。
如果你要仿真两束光以一定夹角合束,需要在其中一个场乘一个横向相位因子,比如exp(1j * kx * x),其中kx与夹角相关。这样出来的条纹就不再是整体明暗变化,而是空间上周期分布的条纹,更像是真正的干涉仪输出。
我之前在这个步骤里犯过一个低级错误:把sqrt(1-alpha)写成了1-alpha,结果能量不守恒,干涉条纹的最大值远低于预期。调试了好久才发现是分束系数的问题。分束器建模一定要时刻盯住能量守恒,功率总和应该等于入射功率,这个检查只需几行代码,建议每个版本都跑一遍。
4.3 相位扫描与干涉条纹提取
相位扫描的核心目的是看到光强随相位差的正弦变化,这组数据可以拟合出干涉仪的灵敏度。我在仿真时除了保存最后一张干涉图,还把中心点处的光强序列单独提取出来,存成数组。
中心点光强序列理论上满足I(phi) = I0·(1 + V·cos(phi - phi0)),其中V就是可见度。用这个公式拟合,可以精确提取可见度和初相位。我实际用scipy.optimize.curve_fit拟合,效果很好,标准误差大概在10的负三次方量级。
如果还要研究二维条纹图,比如检测条纹间距,可以在固定相位下提取一条横向光强剖面曲线,用find_peaks找到峰值位置,然后计算相邻峰间距。实测下来,只要网格够细,峰值定位精度能到亚像素水平。
提取相位这个环节更讲究。仿真里相位可以直接从复振幅数组算出来:phase_map = np.angle(E_out)。但生成的相位图拿过来直接用会有问题。因为相位是缠绕在[-π, π]区间的,如果真实相位变化超过2π,图上会出现相位跳变。这时你需要做相位解包裹,也就是把每相邻点之间的跳变加上或减去2π,让相位连续起来。
二维解包裹算法从原理上不复杂,但实现时边界条件很容易出错。我建议先用scikit-image的unwrap_phase接口,它对一般干涉图处理效果足够好。我测试时发现,如果在网格太稀疏的地方强行解包裹,容易产生“拉链”状伪影,这时优先提高网格分辨率,比算法调参更有效。
5. 参数扫描与误差分析:仿真不能只看条纹
5.1 光程差与波长对条纹周期的影响
参数扫描是这个项目里最有价值的部分,因为它能帮你建立直觉。
首先是光程差ΔL的影响。其他参数不变,只改变ΔL,你会发现条纹空间频率随ΔL近似线性增大。也就是说,光程差越大,输出端的干涉条纹越密。这个结论在光纤干涉仪设计里非常关键:如果要实现高灵敏度测量,我们希望光程差大,但同时探测器分辨率必须足够高,否则条纹直接混叠,信号就丢了。
其次是波长的影响。测量温度或折射率时,不同的激光波长会给出不同的相位灵敏度,因为相位变化Δφ = 2π·ΔL/λ,波长越短,相同光程差对应的相位变化越大。我在仿真里对比了632.8nm和1550nm两个波长,前者的灵敏度差不多是后者的2.45倍。这个倍数反过来印证了为什么很多高精度干涉测量系统选择可见光波段,而远距离光纤传感又偏偏用1550nm——那是因为光纤损耗和光源成本以及探测器技术更成熟,灵敏度反而排在后面了。
你可以直接扫描波长,然后记录固定光程差下中心点强度变化。随着波长连续变化,I(λ)会呈现振荡。这种曲线对“色散测量”很有用,我建议大家自己跑一下,感受相位变化的非线性。
5.2 分束比不均衡时可见度怎么退化
可见度V的定义是V = (I_max - I_min) / (I_max + I_min)。理想情况下,如果两路光振幅完全相等,V等于1,但实际分束器不可能完美,总会差一点。
我在仿真里把分束系数alpha从0.5扫到0.7,得到一组可见度变化数据。结果是:可见度随|alpha - 0.5|增大而单调下降,但下降速度不是线性的。在alpha=0.6时,可见度已经降到0.98左右;alpha=0.7时,可见度只有0.84。这个趋势满足公式V = 2·sqrt(alpha·(1-alpha))。注意这里如果alpha和1-alpha分别代表两臂的振幅系数,那么可见度等于2√(αβ),前提是两路光没有其他振幅损耗。
这个公式非常好用。实验里如果发现干涉条纹对比度不够,就可以先怀疑分束比是否偏离了设计值,再怀疑偏振对齐是否出了问题。仿真里验证公式,到实验现场就能省掉盲猜的过程。
关于振幅系数再补充一点:分束比误差不只会降低可见度,还会让干涉信号的极大值和极小值不对称吗?其实不会,只要两个振幅都是常数,光强仍然是cos函数的线性变换,极大值和极小值关于中线对称。只有噪声存在时,不对称才会被观察到。这点仿真和实验是一致的。
6. 常见问题与调试经验速查
6.1 问题对照表与排查思路
我在整个仿真过程中整理了一张问题速查表,分享给你。
| 异常现象 | 可能原因 | 解决办法 |
|---|---|---|
| 条纹过于密集,出现摩尔纹 | 光程差太大或网格数不够 | 增大网格数或减小光程差 |
| 条纹边缘光强突然衰减 | 高斯光束边缘超出计算区域 | 增大计算区域或减小束腰半径 |
| 可见度始终不到1 | 分束比不是50:50 | 检查alpha系数,改用0.5 |
| 相位图出现黑白跳变 | 相位包裹效应 | 用相位解包裹算法处理 |
| 干涉图有倾斜条纹 | 两路光没有共线合束 | 检查横向相位因子是否包含夹角 |
| 中心光强序列不是标准正弦 | 复振幅没有归一化或存在数值误差 | 检查能量守恒和I_min是否为0 |
这些现象单看都不算大事,但它们叠加起来会让人非常头疼。我的排查习惯是:每跑一个版本先在中心点看光强序列,如果序列是标准余弦,说明整体物理过程是稳定叠加,再去分析空间条纹。如果连中心点都乱,就先别管条纹,回到分束和合束的建模上去查。
6.2 独家避坑技巧:从网格到相位解包裹
避坑这件事,说多了都是泪,但几个关键点值得你记下来。
第一个技巧是网格分辨率检查。不要等到画图发现条纹异常才回去调网格,建议在仿真前做个快速估算:先大略计算条纹周期,然后用周期除以5得到建议的dx上限,再对比当前网格尺寸。这个步骤两分钟能完成,省下的调试时间以小时计。
第二个技巧是能量守恒自检。在分束和合束这两个环节后,分别计算总功率:sum(np.abs(E_out)**2)。理想情况是不论相位怎么变化,总功率应该恒定。这个自检能同时发现分束系数错误和数组坐标错位问题。
第三个技巧是相位解包裹的边界处理。如果只对中心区域做分析,可以先用掩膜屏蔽边缘低光强区域,再解包裹。低光强区域的相位本身噪声很大,强行解包裹会把噪声扩散到有效区域内。我在做相位提取时就吃过这个亏,加上掩膜之后结果稳定了很多。
第四个技巧是保存中间变量。不要只存最终干涉图,建议把分束后的E1、E2、合束后的E_out、光强I都存成npz文件。调试时可以直接加载任意一步的数据重放,不用重复跑仿真。这套仿真代码重跑一次虽然只要几秒,但如果你跑的是参数扫描几千组的情况,保存中间变量几乎就是必须的了。
最后一个建议是版本控制。别觉得一个光学仿真脚本没必要用git,我就经历过“昨天跑通今天全乱”的窘境。其实加个git commit花不了两分钟,后面改参数、加功能都安心得多。
这次仿真跑通之后,我自己最大的体会是:波动光学仿真的难点从来不是物理公式,而是数值实现过程中的各种隐性假设。每当你对仿真结果产生怀疑,先回到网格和相位这两个基本参数上检查一遍,大多数问题都能浮出水面。这套马赫-曾德干涉仪仿真的经验完全可以迁移到其他干涉系统,比如迈克尔逊干涉仪、法布里-珀罗腔,核心的复振幅叠加逻辑是完全通用的。工具只是手段,把物理模型翻译到底层代码里的每一步都看得清、算得准,才是这个项目真正留给我的收获。