光电成像器件CAD仿真实践:从扩散方程到量子效率计算
简介这份资源是一份光电成像器件计算机辅助设计实验报告面向光电成像、电子光学系统设计及微光夜视方向的学习者可用于课程设计或算法仿真项目参考。报告以电场三大定理为理论基础结合连续超张弛迭代算法与拉格朗日插值算法展示如何用C语言将物理模型转化为可计算圆柱形像管内任意点电位、绘制等势线并插值求对称轴电位的CAD程序完整覆盖从边界封闭处理、差分方程到迭代求解的思路。压缩包内为1个PDF文件大小约743KB虽仅1个文件但内容包含摘要、理论推导、程序实现说明及实验过程便于直接阅读、打印或对照梳理。该资源已有170人学习浏览适合需要理解像管电势数值计算、有限差分法及CAD程序设计流程的读者作为参考。1. 光电成像器件CAD实验先搞清楚这份报告在“仿真”什么第一次拿到“光电成像器件计算机辅助设计CAD实验”这个课题的人八成以为它是画版图、拉器件截面图、像机械制图那样把结构画漂亮就结束了。等真正跑起来才发现这个实验的核心动作不是“画”而是“算”——算光子进来之后有多少变成了可收集的电荷算这些电荷在电场和浓度梯度下怎么运动算最终从电极读出的电流信号对应多少入射光功率。换句话说CAD实验替你做了一次“虚拟流片”结构还没长出来你就先看到了它在400nm到1100nm波段到底能出多少量子效率。这个方向对两类人价值最大一类是刚进光电探测器课题组的在校生需要一套能落地的设计闭环而不是零散地调参数另一类是做成像芯片像素优化的工程师想在上线前把吸收层厚度、掺杂浓度、表面复合速率这些变量先扫一遍省掉一批试错晶圆。整篇笔记按“物理量怎么建模→最小可运行代码→二维扩展→参数扫描→踩坑清单→反推实验设计”的顺序展开每个环节都会给出能直接复算的命令和参数而不是停留在概念层。下面从建模仿真最核心的两个问题开始到底算哪些量以及用什么工具算。2. 光电成像器件CAD到底在算哪些量量子效率、光谱响应、暗电流与背后的半导体物理2.1 四个核心指标和它们的物理来源光电成像器件光电二极管、CCD的光敏单元、CMOS像素里的PD的CAD仿真最终都要落到四个可量化的指标上外量子效率EQE、光谱响应度、暗电流密度、响应时间。EQE描述的是“入射N个光子出来M个电子”的比值光谱响应度把它换算成A/W方便直接和测试仪器的读数对比暗电流决定动态范围的下限响应时间决定帧率上限。这四者的物理源头都在半导体结和输运方程里。EQE的损失主要来自三个地方表面的反射、短波在重掺杂表面层的复合、长波穿透吸收层后的漏吸收。暗电流则来自耗尽区产生复合电流和中性区扩散电流两者都和少子寿命、掺杂浓度强相关。仿真时如果你只在几何层面“画图”这些量一个也算不出来需要的是把入射光强分布、载流子产生率、连续性方程、边界条件耦合起来去解。2.2 为什么我选Python而非套件界面来做CAD常见做法是找某商用EDA套件或某实验室自研的光电器件仿真软件界面点点点就能看到电场分布和IV曲线。这类工具的优势是省去了数值方法的重复劳动但在这个实验里有两个问题一是参数非常多且很多藏在二级菜单里换一个结构就要到处找半天二是整个实验过程变成黑匣子——输入厚度、输出EQE中间怎么收敛的、边界条件默认是什么用户完全不知道一旦结果异常很难排查。我一般会用Python组装自己的CAD流程只依赖NumPy、SciPy和Matplotlib这三个通用科学计算包。原因很简单一维和准二维的器件仿真扩散方程加边界条件用有限差分法就能解代码量不超过200行参数扫描天然适合用循环和数组操作结果和实验曲线放在同一张图上对比也最简单。更重要的是每一个物理量都能追溯到代码里的某个矩阵或某个边界条件出了问题能一层层查下去这种“白盒”属性对理解器件物理的价值远大于套件界面。2.3 从“三层结构”到坐标系建模前的两个关键思考动手写代码前先做两件事。第一件把器件结构抽象成层叠。以最常见的垂直入射p-i-n光电二极管为例从光照面往下依次是重掺杂p型接触层、轻掺杂或本征i层耗尽区主体、重掺杂n型衬底层。每一层都要给出厚度、掺杂浓度、少子寿命、迁移率。第二件确定坐标系和光照方向。我习惯把入射面设为x0光沿x轴正向传播耗尽区边界记为x1和x2中性区在耗尽区两侧。这样吸收系数随波长的变化就只体现在产生率G(x)里方程形式统一。这两个决定直接决定了后面代码的复杂程度。如果你把坐标反着放或者把p层放在衬底侧边界条件的符号和区域划分会全部乱掉排查时非常痛苦。这也是为什么本篇所有的代码示例都固定用“上表面入射、p层朝上”的约定——先在一个约定下跑通再谈扩展。3. 跑通最小流程一维扩散方程求解量子效率的完整代码3.1 主脚本与结构从入射光到外量子效率最小可运行的CAD实验是做一条“EQE-波长”曲线。用有限差分法求解稳态少子扩散方程把p型中性区、耗尽区、n型中性区三段分别处理耗尽区内电场很强光生载流子全部扫出收集效率近似为1中性区靠少子扩散解稳态扩散方程得到向耗尽区边界输运的少子流。import numpy as np from scipy import constants as const # ---------- 材料与结构参数 ---------- q const.e # 元电荷 kbT const.k_B * 300 # 300K热电压对应能量 T 300.0 # 结构p层0.3um / i层1.0um / n层10um厚度单位cm xp 0.3e-4 # p型中性区实际还要去掉顶端复合层 xi 1.0e-4 # 耗尽区宽度 xn 10.0e-4 # n型中性区 Np 1e18 # p层掺杂浓度 cm^-3 Ni 1e13 # i层近似本征 Nn 1e16 # n层掺杂浓度 cm^-3 # 输运参数典型参考值实际可按材料库修正 tau_p 1e-6 # n区少子空穴寿命 s tau_n 1e-6 # p区少子电子寿命 s mu_n 800.0 # p区电子迁移率 cm^2/(V·s) mu_p 300.0 # n区空穴迁移率 cm^2/(V·s) Dn mu_n * kbT / q Dp mu_p * kbT / q Ln np.sqrt(Dn * tau_n) Lp np.sqrt(Dp * tau_p) S_front 1e3 # 前表面复合速率 cm/s # 吸收系数数据波长nm: alpha cm^-1典型参考值可用实测光谱替换 wl_list np.array([400, 450, 500, 550, 600, 650, 700, 750, 800, 850, 900, 950, 1000, 1050, 1100]) alpha_list np.array([6e4, 3e4, 1.2e4, 5e3, 4e3, 2.5e3, 2e3, 1.2e3, 8e2, 5e2, 3e2, 1.5e2, 8e1, 4e1, 2e1]) def R_front(wl_nm): 前表面反射率简化模型短波低反射长波略高。 return 0.08 if wl_nm 700 else 0.12 def generate_rate(alpha, phi0, R, xgrid): 光生载流子产生率G(x) phi0*(1-R)*alpha*exp(-alpha*x) return phi0 * (1 - R) * alpha * np.exp(-alpha * xgrid) # 各区间收集效率 def collection_p_side(alpha, xp_eff, x1): p型中性区少子电子收集概率解扩散方程得到边界电流占比 # 网格 N 500 x np.linspace(0, xp_eff, N) # 从表面到耗尽区边界 dx xp_eff / (N - 1) G generate_rate(alpha, 1.0, R_front(alpha), x) # 归一化光强 # 方程: Dn * d2n/dx2 - n/tau_n G 0 A np.zeros((N, N)) b np.zeros(N) for i in range(N): if i 0: # 表面边界Dn*dn/dx S*n差分近似 A[i,0] -(Dn/dx S_front) A[i,1] Dn/dx elif i N-1: # 耗尽区边界少子浓度近似0被强电场扫走 A[i,i] 1 else: A[i,i-1] Dn/dx**2 A[i,i] -(2*Dn/dx**2 1/tau_n) A[i,i1] Dn/dx**2 b[i] -G[i] n np.linalg.solve(A, b) # 边界扩散电流 Dn * dn/dx|_(xp_eff)归一化到入射光强 Jn Dn * (n[N-1] - n[N-2]) / dx # 注意这里n[N-1]0所以J近似负取绝对值 eta_p abs(Jn) return eta_p这段代码实现了p型中性区的收集概率计算。边界条件的物理含义是关键表面处i0用了复合边界即扩散流等于表面复合流耗尽区边界iN-1用了理想收集边界即少子浓度被强制归零模拟强电场把载流子立即扫走。中间网格点就是标准的稳态扩散方程离散形式右侧源项b由光生载流子产生率构成。参数说明N500表示p区内部网格数收敛性一般够用tau_n和mu_n分别控制少子寿命和迁移率两者共同决定扩散长度S_front对短波影响极大表面复合快会显著拉低蓝光响应alpha是波长相关的吸收系数它决定光子在哪一层被吸收——短波alpha大表层就基本吸收完长波alpha小光子要穿透很深才能被吸收所以p区厚度和扩散长度直接决定蓝紫光响应n区厚度决定近红外响应。3.2 让脚本跑起来三处最容易卡住的参数第一次跑这个脚本有三处几乎必卡。第一处单位混用。能带算出来的长度单位是cm但很多人查资料时习惯用微米xp 0.3直接带进去会让吸收深度全部错位。我的习惯是所有几何厚度统一转成cm波长相关的alpha保持cm^-1这样exp(-alpha*x)指数部分无量纲电流密度才是A/cm^2。第二处Kronecker符号缺失。扩散方程离散后每个内点的对角元是-(2*Dn/dx^2 1/tau_n)少子复合项只有在对角元里出现漏掉它矩阵就奇异。第三处边界条件的符号。表面复合边界处理成Dn*(n1-n0)/dx S*n0移项后写成代码里那种形式很多人会把S_front和Dn/dx的符号搞反导致算出来的收集概率大于1或出现负值。3.3 结果怎么读量子效率曲线的三个判断点跑完之后得到一串波长对应的EQE值三个判断点帮你快速确认曲线合理第一个是零偏时的峰值一般在500到700nm之间典型硅p-i-n结构峰值EQE在0.7到0.9之间第二个是短波端400nm的下坠幅度如果EQE不到长波峰值的六成基本可以断定表面复合S_front设得太高或者p层厚度超过了电子扩散长度第三个是长波端的滚降斜率900nm之后EQE应该平滑下降如果出现鼓包说明反映衬底反射的边界条件没处理干净这个后面避坑章节会展开。4. 步进到二维把串扰和边缘效应算出来并做参数扫描4.1 二维结构的离散化思路与边界条件一维结果给的是“无穷大平面”下的理想值像素阵列里真正影响成像质量的往往是相邻像素之间的串扰这是一维模型完全看不到的。二维建模通常在x-y平面做x是光的入射深度方向y是像素横向方向。光照在某个像素中心部分光生载流子会横向扩散到相邻像素的耗尽区形成串扰电流。离散化的思路跟一维一样只是把拉普拉斯算子展开成d2u/dx2 d2u/dy2矩阵从三对角变成块三对角。边界条件除了表面和耗尽区边界外多了像素边缘的对称边界相邻像素结构相同取对称面法向电流为零和底部衬底的欧姆接触边界。建议先用50×50网格跑通形状确认误差可接受再加密不要一上来就200×200矩阵求解时间会指数级增加调试期会非常痛苦。4.2 吸收层厚度扫描饱和厚度和最优工作点跑通单点之后CAD实验最有价值的部分就来了参数扫描。最常见的扫描对象是吸收层厚度i层厚度它直接影响长波响应和暗电流之间的平衡。下面这段代码对i层厚度做循环扫描画出不同厚度下的EQE曲线并自动找指定波长的最优值thickness_list np.linspace(0.5e-4, 3.0e-4, 6) # i层厚度从0.5um到3um target_wl 850 # 近红外成像关注的波长 eqe_at_target [] for xi_val in thickness_list: # 重新设置耗尽区宽度 xi xi_val # 这里调用3.1节中的单波长计算函数返回该波长下的EQE eqe_val calculate_eqe_at_wavelength(target_wl, xixi, structure_paramsparams) eqe_at_target.append(eqe_val) print(fxi{xi_val*1e4:.2f}um, EQE{target_wl}nm {eqe_val:.3f}) best_idx int(np.argmax(eqe_at_target)) best_xi thickness_list[best_idx] print(f建议吸收层厚度{best_xi*1e4:.2f}um)这段代码的意义在于把“凭经验猜厚度”变成“扫出来看趋势”。运行后你会发现一个规律总厚度从0.5μm增加到2μm时850nm处EQE显著上升再往后上升趋缓甚至回落。原因是i层加厚一方面让长波光子被吸收的概率更大另一方面耗尽区变宽后载流子渡越时间变长虽然稳态EQE不受渡越时间影响但若同时引入复合寿命较短的材料参数厚度增加反而让部分载流子在扫出前复合掉曲线出现饱和或回落的拐点这个拐点就是你要找的最优工作点。如果你扫描的是掺杂浓度而不是厚度要关注的趋势是耗尽区宽度随掺杂升高而变窄器件电容变大响应速度变快但EQE变低——高灵敏度、高速度是这个方向永远要权衡的一组矛盾。4.3 网格尺度选择什么时候用0.5μm什么时候用0.1μm二维网格步长直接决定收敛性和计算时间这个选择有规律可循。判断依据是物理量的最小特征长度在p-i-n结构里如果光吸收长度(1/alpha)和扩散长度都远大于网格步长那0.5μm步长就够长波端尤其如此。但在表层高掺杂区和结附近耗尽区边缘的电场变化剧烈少子浓度在很短距离内从平衡值降到接近零这时网格必须细到0.1μm甚至更小否则求解出的收集电流会出现锯齿状跳变。我通常在同一个结构里用非均匀网格表层前0.3μm用0.05μm步长结区和耗尽区用0.1μm中性大区用0.5μm。这样做的好处是计算量控制在单次几秒钟内同时边界层分辨率足够。判断收敛性的一个土办法是把网格步长减半看结果变化是否小于1%如果变化明显说明原网格太粗需要加密那个区域。这个方法虽然朴素但在工程调试中比任何理论判据都直接。5. 光电成像器件CAD实验的避坑清单五个反复出现的翻车现场5.1 现象曲线形状对但绝对值差了一截第一次跑完EQE曲线和文献对比会发现形状完全一致——峰值位置对、滚降斜率对——但绝对值整体偏低20%到30%。很多人的第一反应是调整吸收系数折腾半天也没用。原因几乎都出在材料输运参数上很多人从教科书上抄的少子寿命是1ms量级的高纯材料参数而实际外延层或离子注入层的寿命是0.1到10μs量级。寿命直接决定扩散长度进而决定中性区载流子能走多远。解决办法不是凭感觉调tau而是用实验测量的寿命值或者依据掺杂浓度查经验拟合曲线。5.2 现象加反向偏压后短路电流反而下降物理直觉告诉我们反向偏压加大耗尽区变宽收集应该更完全。但仿真结果显示电流略微下降看起来像见了鬼。排查后指向一个容易被忽略的点偏压增大导致耗尽区边界向两侧扩展原本在中性区离边界较近的载流子现在落入了耗尽区这部分载流子的产生复合电流不再被中性区扩散方程处理而你的模型里耗尽区复合可能没建模或建得不够准确。解决方式是在耗尽区加入SRH复合项让电场扫出和复合形成竞争关系否则偏压越高、耗尽的低电场区越宽复合损失越严重结果就出现反直觉的下降。5.3 现象量子效率在长波区出现“鼓包”正常EQE曲线随波长增加是单调下降或平滑滚降但有些模型在1000到1100nm处出现一个宽鼓包明显高于900nm处。这个现象的原因几乎总是同一个衬底反射没有处理好。光子穿过耗尽区进入n型衬底如果衬底足够厚且背表面有高反射结构或欧姆接触被误建模成反射镜部分长波光子会被反射回耗尽区再产生一次收集机会等效于增强了长波响应。真实器件里背反射确实存在但不能用默认全反射来建模。正确做法是给背表面设定一个反射率参数用实测或经验值而不是代码里无意间把边界设为诺依曼零梯度默认成全反射。5.4 现象网格缩小一倍结果却跳变网格从1μm加密到0.5μmEQE变化2%以内算正常但有人遇到结果是跳变5%以上的。排查后发现是网格和结构层厚不匹配i层厚度1.0μm如果网格步长取0.3μm层边界落在网格点的位置会对不齐等效结构厚度在加密前后改变了产生率的积分也跟着变。解决方法是把网格步长设计成能整除各层厚度或者用非均匀网格让边界正好落在网格线上。一个实用技巧是先用整数微米的厚度建模调试通了再改成连续参数。5.5 现象边界条件写错把欧姆接触当成了肖特基接触这是最隐蔽的坑。数值求解时重掺杂接触边界通常设为“欧姆接触”即少子浓度等于平衡值或接近零。但有时候为了数值稳定有人会把边界条件改成固定电流为零这实际上模拟的是绝缘边界或肖特基阻挡接触。两者的结果差异巨大欧姆接触让少子被边界快速复合等效于一个超级复合中心绝缘边界则把少子反射回来等效于表面复合速度为零。如果发现短波响应异常偏高先检查边界是不是被写成了dn/dx0正确的欧姆接触应该用nn0平衡浓度近零。6. 进阶玩法用光谱响应曲线反推吸收层厚度让CAD反过来指导实验设计模拟实验做完真正的价值在于把仿真结果倒着用在实验上——这步决定了CAD是交作业还是真的省晶圆批次。实操中我常用一个“反推法”先测出器件的相对光谱响应曲线拿到几个特征波段的曲线斜率然后回到仿真里做参数扫描寻找厚度值使模拟曲线与实测曲线在那些波段的斜率和交叉点重合。具体做法分三步一是在固定其它参数的条件下把厚度从0.5到3μm等间隔取20个值分别算出400到900nm的EQE曲线二是把每条模拟曲线的短波长积分400到600nm除以长波长积分700到900nm得到一个比值这个比值对厚度变化非常敏感三是把实测数据算出的同一个比值拿来比对插值出最接近的厚度。这个方法的优势是它对反射率、表面复合速度这些系统误差不敏感——系统误差会同样压缩短波和长波响应但不会显著改变两者的比值。所以即使在参数不够准确的情况下反推的厚度误差也能控制在20%以内对生长工艺的反馈已经足够用。类似的方法还能反推扩散长度用500nm附近的绝对EQE下降幅度去拟合虽然要标定反射率但趋势判断是可靠的。做这件事我沉淀出的习惯是每次仿真实验都保留一个“参数真值表”记录所有输入参数和对应的输出曲线同时留一列备注标注来源文献默认值、估算值、实测值。两个月后回头处理数据时不会面对一堆记不清含义的数字发呆。这个习惯让我在多个模拟项目里免于返工也算对“实验报告”这四个字最诚实的态度。还有一个我比较坚持的做法仿真完成不等于结束拿一条实测曲线来对比。哪怕只对比三个波长点也能逼着你认真检查边界条件和复合参数。没有实测数据时至少要和公开教材里的典型曲线做定性对比——形状差得远一定是模型或边界有问题宁可多查一步不要急着下结论。CAD的价值只有和实验数据碰撞才体现得出来希望这篇整理能帮你少走几次弯路。希望帮到你。本文还有配套的精品资源点击获取