10.6μm几何相位超表面透镜的设计与仿真验证全流程

做红外光学系统的人对10.6μm都不会陌生这是CO2激光器最典型的工作波长工业切割、焊接、医疗整形里到处都能碰到。传统方案里聚焦光路基本靠ZnSe或者Ge材料的球面透镜组件厚、重量大热效应也让人头疼。超表面透镜提供了另一种思路在平面上用亚波长结构单元去调制相位厚度做到微米量级。而在众多超表面方案里几何相位也叫Pancharatnam-Berry相位PB相位又是实现起来最直观的一种——结构单元形状不变只靠旋转角度就能给出0到2π的完整相位覆盖。这篇文章我就拿10.6μm这个波长做例子把几何相位超表面透镜从原理、参数计算、代码生成到仿真验证的完整链路梳理一遍。文中的代码我实际跑过逻辑可以直接复用设计参数怎么定、单位怎么统一、量化怎么做、版图怎么导出这些容易出问题的细节我会单独拿出来说。适合正在做中远红外超表面、衍射透镜或者需要把“会聚焦的透镜”落到实际掩膜版图上的工程师参考。1. 设计思路与原理拆解1.1 为什么10.6μm几何相位超表面值得做先说说10.6μm这个波长在超表面设计里的特殊性。这个波长对应的光子能量只有约0.12eV所以很多在可见光波段不能用的材料在这里反而很合适。金属在红外波段依然有等离子体色散但欧姆损耗会拖累透射效率所以可见光超表面里常用的金属纳米天线到了10.6μm就不太划算。全介质超表面是更好的选择常用的柱体材料有锗Ge、高阻硅、硒化锌ZnSe等其中Ge的折射率约4.0透明度好是10.6μm附近非常经典的介质材料。再从加工角度看。可见光超表面的特征尺寸通常是100300nm对光刻和刻蚀要求很高10.6μm波长对应的单元周期可以放到几微米柱体尺寸也在微米量级普通的紫外光刻和干法刻蚀就能做深宽比压力小很多。这是中红外超表面落地相对容易的重要原因也是我觉得这个波长特别适合练手的原因。几何相位方案还有一个结构上的优势它不需要每个像素换一种柱体尺寸。共振相位超表面通常要扫描柱子的长、宽、高很难做到多个相位点都是高效率而几何相位超表面的所有柱体形状完全一样只是旋转角不同。这意味着单元结构优化只需要做一次找到最优的半波片参数后面成千上万个像素全是同样的柱体只是角度变化。对设计代码来说这一步简化非常多对工艺一致性也更友好。不过几何相位有一个硬前提入射光必须是圆偏振而且我们关心的是出射光中与入射旋性相反的那个偏振分量。换句话说一个LCP左旋圆偏振入射真正被聚焦的是转换后的RCP右旋圆偏振分量。系统里通常要加一枚四分之一波片来产生圆偏振这个在后面仿真和测试时都要记得。1.2 透镜相位分布从球面波到像素矩阵超表面透镜设计的第一步是确定每个位置需要补偿的相位。常规折射透镜靠曲面轮廓让光程连续变化超表面透镜则直接在平面上写相位分布。为了在焦距f处会聚理想相位分布是φ(x, y) - (2π / λ) * (√(x² y² f²) - f)这个公式的含义很直接透镜中心和边缘到焦点的光程差不同边缘光路更长所以需要提前补偿的相位也更多符号取负。如果要做发散透镜符号取正即可。设计代码里通常把相位包裹到02π因为超表面单元能提供的相位本身就覆盖一个周期。实际计算时要把连续口径离散成网格。每个网格大小等于一个超表面单元的周期p每个单元都有一个(x, y)坐标和对应的目标相位。网格量直接决定了透镜的像素数。比如口径5mm、周期5.3μm大概是943×943个单元接近90万个结构这个规模在版图里可以处理但已经不小了。这里还会遇到一个细节相位是按连续函数算出来后要量化的。制造时不可能每个像素都对应无限个相位值一般把2π分成8、16或者32层。量化层数越高聚焦效率损失越小但版图数据量会变大。对于几何相位超表面量化还可以直接作用到旋转角度上因为相位φ对应的是2θ的旋转所以角度量化等价于相位量化。1.3 单元结构设计的核心半波条件与转化效率几何相位的物理本质要稍微展开一点。一束圆偏振光经过一个各向异性的结构后出射的反转圆偏振分量会携带一个额外的相位大小为2θ其中θ是结构快轴或慢轴相对参考方向的旋转角。这个相位和波长没有显式关系所以几何相位本身是宽带适用的但单元结构的工作带宽会受到共振特性的限制。要让“旋转角度”这个机制生效单元结构必须近似一个半波片。也就是说在柱体的两个正交主轴上透射光复振幅的相位差要接近π同时两个主轴方向的透过率要尽量高、尽量接近。只有满足这个条件入射圆偏振才能高效地转换成反向圆偏振并且转换分量的相位确实随2θ线性变化。在10.6μm波段可以用椭圆柱、矩形柱或者十字形柱体。设计流程通常是先在FDTD或RCWA软件里搭一个单周期单元设置周期性边界入射平面波扫描柱体的长轴、短轴和高度。每次扫描都记录x偏振和y偏振下的复透射系数分别记为t_x和t_y。对于线偏振基LCP入射转换为RCP出射的效率可以写成η_cross |(t_x - t_y) / 2|²这个公式很实用。单元优化目标就是让|t_x - t_y|尽量大同时透射总能量尽量高。如果你用的工具能直接输出圆偏振转换效率那就更方便。实际典型参数范围上周期可以选约5.3μmλ/2柱体高度在25μm之间长轴在24μm短轴在12μm具体取值用扫描曲线定。2. 设计代码实现参数计算与版图生成2.1 参数怎么定一个5mm口径透镜的算例空讲公式很难有手感我直接给一个可以照抄的算例。假设设计一个用于CO2激光聚焦的几何相位超表面透镜波长λ10.6μm焦距f10mm口径D5mm单元周期p5.3μm。先算数值孔径NA D / (2f) 5 / (2×10) 0.25这个NA不算高适合第一版设计验证。理想情况下艾里斑半径约为1.22λf/D 25.9μm也就是说焦点的中心亮斑大致在50μm直径这个量级。焦深约为λ/NA² ≈ 170μm对后续装配公差还算友好。单元周期p为什么取5.3μm两个原因。第一p必须小于工作波长否则会产生高阶衍射周期性结构就会像光栅一样把能量分到多个级次上去。第二p也不宜太小否则柱体之间的近场耦合会变强以及刻蚀深宽比会变大。取λ/2左右是一个很稳的折中。还有一个约束值得注意对于聚焦透镜单元间隔还需要满足p ≤ λ/(2NA)这里p5.3μm远小于21.2μm完全满足。相位量化层数我习惯先取16层。量化层数带来的效率损失可以用近似公式η≈[sin(π/L)/(π/L)]²估算16层时大约是98.7%8层大约是94.9%。所以16层是一个“效率足够、数据量不大”的常见选择。对应到几何相位旋转角2π相位分为16层旋转角步长就是π/16 11.25°。2.2 Python主程序相位计算与坐标输出下面这段Python代码是我自己整理的版本输入参数后直接生成每个单元的坐标和旋转角输出CSV。代码逻辑很直接先生成网格再算理想相位再做量化和角度换算最后导出。import numpy as np import pandas as pd # ---------- 设计参数 ---------- LAMBDA 10.6e-6 # 波长, m FOCAL 10e-3 # 焦距, m DIAMETER 5e-3 # 透镜直径, m PERIOD 5.3e-6 # 单元周期, m N_LEVELS 16 # 相位量化层数 # ---------- 网格生成 ---------- N int(np.floor(DIAMETER / PERIOD)) if N % 2 0: N 1 # 确保奇数尺寸, 中心像素对准坐标原点 xv (np.arange(N) - (N - 1) / 2) * PERIOD XX, YY np.meshgrid(xv, xv) R np.sqrt(XX**2 YY**2) # ---------- 相位分布 ---------- phase -(2 * np.pi / LAMBDA) * (np.sqrt(R**2 FOCAL**2) - FOCAL) phase np.mod(phase, 2 * np.pi) # 包裹到 [0, 2π) mask R DIAMETER / 2 # 圆形口径 # ---------- 量化与旋转角 ---------- delta_phase 2 * np.pi / N_LEVELS phase_q np.round(phase / delta_phase) * delta_phase theta phase_q * 0.5 # 几何相位: φ 2θ theta_deg np.degrees(theta) # ---------- 导出有效单元 ---------- xx XX.ravel() yy YY.ravel() rr R.ravel() tt theta_deg.ravel() data np.column_stack([xx, yy, rr, tt]) data data[rr DIAMETER / 2] # 只保留圆形口径内单元 df pd.DataFrame(data, columns[x_m, y_m, r_m, theta_deg]) df.to_csv(metalens_10p6um.csv, indexFalse) print(f有效单元数: {len(df)}) print(f阵列尺寸: {N} x {N}, 物理口径: {N * PERIOD * 1e3:.3f} mm)有几个点需要特别提醒。第一np.mod之后相位范围是0到2π量化时我用的是四舍五入而不是向下取整这样每个相位区间内的最大误差只有半个量化步长。第二几何相位旋转角是φ/2不要忘了除以2这是新手最容易写错的地方。第三坐标用国际单位米导出CSV时也保留米后面转版图时再统一乘10^6换成μm省得单位混乱导致结构尺寸差千倍。如果你要用MATLAB实现核心思路完全一样只需注意MATLAB里meshgrid的维度顺序和Python有些差别计算时保持坐标一致即可。2.3 从坐标到GDSII版图CSV只是中间产物最终加工需要GDSII版图。以Python里的gdspy库为例写版图时有一个非常重要的经验不要在每个像素位置单独画一个旋转后的多边形文件会爆炸。90万个单元的GDS直接画多边形生成的文件可能几个GB甚至更大流片时传递文件都很痛苦。更聪明的做法是按旋转角度分cell再用CellReference做阵列引用。因为量化后的旋转角只有16种先建16个包含单个柱体的cell主cell里每个单元用一个引用并只带旋转参数。GDS里的reference记录的量很小文件体积会从几个GB降到几十MB这个差距在实际项目里非常关键。import gdspy import numpy as np # 柱体尺寸(μm)和周期(μm), 以仿真优化结果为准 w_um 1.8 l_um 3.4 period_um 5.3 lib gdspy.GdsLibrary() top lib.new_cell(TOP) # 按角度建立独立cell, 矩形本体不旋转, 旋转放到引用上 unique_angles np.unique(np.round(theta_deg, 2)) angle_cells {} for a in unique_angles: cell lib.new_cell(fUNIT_{a:.1f}) cell.add(gdspy.Rectangle((-w_um / 2, -l_um / 2), (w_um / 2, l_um / 2))) angle_cells[a] cell # 主阵列引用 for i in range(N): for j in range(N): if mask[i, j]: a round(theta_deg[i, j], 2) origin (YY[i, j] * 1e6, XX[i, j] * 1e6) top.add(gdspy.CellReference(angle_cells[a], originorigin, rotationa)) lib.write_gds(metalens_10p6um.gds)这里要注意gdspy的Rectangle默认中心在原点而CellReference的旋转是绕参考点旋转所以整个单元更像是在原点旋转后再平移到目标位置效果完全正确。还有一点旋转角度的参考方向要和你单元仿真时的慢轴方向保持一致否则一上来就是符号相反。3. 仿真验证与加工复盘3.1 不跑全波也能验证角谱传播法检查焦点90万个单元的透镜直接扔进FDTD做全波仿真是不现实的内存和时间都撑不住。实际项目里可以先做一个中间验证用角谱传播法把相位板产生的场传播到焦平面看焦点位置、焦斑尺寸和能量集中度。这个验证不需要解决Maxwell方程组只需要把每个像素近似为理想半波片出射场相位等于2θ即可。Python里角谱传播的实现并不复杂。核心是传入初始电场E0做二维傅里叶变换到频域乘上传播算子再反变换回空间域。传播算子的形式是exp(i·kz·z)其中kz由平面波色散关系决定。def angular_spectrum_propagate(E0, dx, wavelength, z): ny, nx E0.shape fx np.fft.fftfreq(nx, ddx) fy np.fft.fftfreq(ny, ddx) FX, FY np.meshgrid(fx, fy) k 2 * np.pi / wavelength kz2 k**2 - (2 * np.pi * FX)**2 - (2 * np.pi * FY)**2 kz np.sqrt(np.maximum(kz2, 0)) H np.exp(1j * kz * z) E0f np.fft.fft2(np.fft.ifftshift(E0)) E_out np.fft.fftshift(np.fft.ifft2(E0f * H)) return E_out # 初始场: 口径内相位为2θ, 即phase_q src np.exp(1j * phase_q) src[~mask] 0 E_focal angular_spectrum_propagate(src, PERIOD, LAMBDA, FOCAL) I_focal np.abs(E_focal)**2跑完之后在焦平面看强度分布理论上的艾里斑直径应该在50μm左右。如果焦点明显偏大先查相位分布是不是反了如果焦点位置不对查焦距符号和坐标原点。角谱法的好处是几十毫秒就出结果非常适合在写版图之前做逻辑自检。3.2 单元仿真和全场仿真的配合角谱验证通过之后还是要回到单元仿真去确定柱子的具体尺寸。FDTD单元仿真的设置有几个关键点模型里包含一个周期单元边界用周期边界波长范围覆盖10.6μm及周边入射光源用圆偏振平面波。如果你用的工具只支持线偏振源就分别跑x偏振和y偏振两次再在数据处理时合成圆偏振照样能得到t_x和t_y。仿真扫描变量包括柱体长轴、短轴和高度。扫描完成后画出转化效率η_cross随尺寸变化的二维图找全局最优。我的经验是不要只盯效率峰值还要看峰值附近有没有足够的工艺容忍度。如果效率曲线变化特别陡稍微偏一点尺寸效率就掉很多实际加工良率会很难看宁可选一个效率略低但是平缓的区域。如果条件允许可以做一个小口径子透镜的全波验证比如直径500μm左右只有不到一万个单元。这样FDTD能跑得动结果可以用来校正角谱法里“每个单元都是理想半波片”的假设。把子透镜的仿真聚焦效率做出来再按比例评估大口径透镜性能是当前比较稳妥的工程流程。3.3 制造公差与材料选择10.6μm波段超表面的工艺容差相对可见光要宽裕得多。100nm的线宽误差在可见光波段可能是致命的但在几微米的柱体尺寸面前只是几个百分点的变化。真正要关注的是刻蚀深度误差和侧壁倾角。柱体高度改变会影响两个方向偏振的相位差可能导致半波条件失配转化效率下降。材料选择上Ge柱体是常见起点折射率高、吸收低和标准半导体加工兼容。基片可以用ZnSe本身对10.6μm透明。如果用硅柱要注意高阻硅在10.6μm处仍有残余吸收低功率验证问题不大但如果用于高功率CO2激光聚焦吸收引起的温升是不可忽视的。仿真时不要把材料当成无损耗尽量用实际测量或数据库中的复折射率。还有热效应。超表面透镜的优势是轻薄坏处是热量不容易通过结构传导散掉。高功率应用场合要考虑基片冷却方案和材料损伤阈值。我见过不少设计在小功率下性能很好一上功率就焦距漂移甚至烧毁这个问题从设计阶段就要提前想清楚。4. 常见问题与调试心得4.1 焦点处始终有一个零级亮斑这个现象在用几何相位超表面时非常典型。几何相位透镜的出射场里真正聚焦的是反向圆偏振分量而同向圆偏振分量会形成未被调制的背景和中心零级亮斑。如果单元的半波条件不完美同向分量变大零级亮斑就会更加明显。排查路径就两条。第一确认入射光是纯正的圆偏振四分之一波片的快轴对准是否准确椭圆度对零级强度影响很大。第二检查单元结构是不是真正工作在最优半波点。一个很直接的判断方式是看单元仿真里η_cross的数值如果远低于理论最大值的75%以上那问题大概率出在单元设计上而不是透镜排布上。4.2 相位符号反了透镜变成发散镜这是几何相位超表面最容易遇到的逻辑错误之一。由于几何相位出射相位是±2θ符号取决于入射圆偏振手性和柱体慢轴定义。同一个版图LCP入射和RCP入射得到的聚焦/发散行为是相反的。我现在的习惯是设计前先做一个最简验证准备一个只有几十个柱子的线形相位梯度光栅周期远小于透镜口径用FDTD或者角谱法看它是否把光束偏向预期的一侧。如果偏转方向反了把所有旋转角度取负即可。这个验证只需要几分钟能避免整个版图做完之后才发现符号反了的悲剧。4.3 GDS文件太大打开卡死前面已经提到一定要用CellReference按角度构建阵列不要逐个画多边形。如果你接手的是别人已经写好的代码发现它直接给几百万个Ellipse或者Polygon第一时间改成reference结构。另一种常见做法是生成多个小GDS分块分别做DRC后再合并也能缓解版图软件的压力。还有一个小技巧在导出GDS之前先用代码检查相邻柱体是否重叠。柱体旋转后矩形对角线方向会占据更大空间相邻单元容易“打架”。检查方法是遍历每个柱体的四个角点做间距布尔判断一旦发现重叠优先调小柱体宽度而不是增大周期因为周期改了相位采样和衍射特性都会变。4.4 相位量化层数对焦斑质量的影响量化层数越低焦斑旁瓣越高聚焦效率越低。8层和16层之间的差距在仿真里一眼就能看出来8层时焦点附近会有明显的周期性格外强度分布16层时基本消失32层和16层之间的差距很小。所以我推荐第一次设计直接用16层不要一上来就在角度分辨率上省钱。量化还会影响焦点纯度。理想球面波在焦平面的能量集中在艾里斑里量化引入的周期性子相位会造成额外衍射级。这个效应本质上是周期性的相位误差和单元周期、量化层数都有关系。越是追求干净焦斑的应用量化层数越要放宽。如果版图数据量允许32层也没问题反正有CellReference方案撑着。我在实际项目里还有一个体会设计代码本身并不难难的是把每个环节的假设都验证清楚。几何相位超表面看起来只有一个旋转角度但因为符号、定义、单位这些小的不一致翻车的概率一点也不低。上面这套流程走下来10.6μm单波长的聚焦透镜基本能一次做对。最后分享一个我自己养成的习惯正式生成完整版图之前先在口径中心和边缘各取一块20×20的区域输出成小GDS或者图片检查一遍。检查相邻柱体间距、旋转方向连续性、圆形口径边缘是否有效封闭。10.6μm的柱体虽然已经是微米级但旋转后相邻柱子在斜对角方向上很容易碰到一起这个隐患在数据里看不见在版图上一查就会原形毕露。设计软件和代码可以帮我们把复杂逻辑变成数据但真正决定成品率的永远是这些不起眼的细节。