留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

半导体多物理效应并行计算程序JEMS-CDS-Device设计与实现

李光荣 赵振国 王卫杰 游春光 周海京

引用本文:
Citation:

半导体多物理效应并行计算程序JEMS-CDS-Device设计与实现

    作者简介: 李光荣(1988—),男,助理研究员,从事高性能电磁计算研究,器件与电路仿真软件研发:li_guangrong@iapcm.ac.cn.
  • 基金项目: 科学挑战专题资助项目(TZ2018002);国家自然基金项目(11675022);中物院复杂电磁环境重点实验基金项目(FZ2019-001);中国工程物理研究院创新发展基金项目(2019034);国防基础科研计划项目(C1520110002)
  • 中图分类号: TN302

Design and implementation of semiconductor multi-physical parallel computing program JEMS-CDS-Device

  • CLC number: TN302

  • 摘要: 针对复杂电磁环境下器件多物理效应机理研究需求,研发了半导体多物理效应并行计算程序JEMS-CDS-Device。介绍了JEMS-CDS-Device的架构设计与实现技术。程序基于非结构网格并行框架JAUMIN实现,采用有限体积法(FVM)离散,使用牛顿法全耦合求解“电-漂移扩散输运”问题。程序采用“内核+算法库”形式架构,支持2维和3维非结构网格、千万自由度问题并行求解,支持物理方程、离散算法、材料物理模型等的扩展开发。
  • 图 1  插值网格与积分网格(控制体)

    Fig. 1  Interpolation grid and integral grid (control body)

    图 2  前向自动微分

    Fig. 2  Forward automatic differentiation

    图 3  积分网格(对偶网格)修正

    Fig. 3  Integral mesh (dual mesh) correction

    图 4  异质界面处的间断[15]与算法区域处理

    Fig. 4  Discontinuity heterogeneous interface and algorithm region processing (right figure, there are two logical nodes M1 and M2 in node A)

    图 5  网格层次结构与并行构件[17]

    Fig. 5  Grid hierarchy and parallel components

    图 6  JEMS-CDS-Device主要层次

    Fig. 6  JEMS-CDS-Device

    图 7  软件模块间关系

    Fig. 7  Relationships between software modules

    图 8  三维PN净掺杂密度分布与二极管阳极V-I曲线

    Fig. 8  PN diode ’s density distribution of net doping and Anode V-I curve

    图 9  二维NMOS管平衡态电势分布结果对比(电势单位:V,参考零点为无穷远)

    Fig. 9  2D NMOS equilibrium potential distribution comparison(potential unit: V, the potential reference zero is at infinity)

    图 10  二维PN管的净掺杂分布与载流子分布(正向偏压0.5 V)

    Fig. 10  Net doping and carrier distribution of 2D PN device (positive deviation 0.5 V)

    图 11  输入的波形与某时刻PIN的温度分布

    Fig. 11  Input waveform and temperature distribution of PIN (at 7.78 ns)

    表 1  DDM1非线性迭代中线性解法器典型收敛情况(时间单位:秒)

    Table 1  Typical convergence of linear solver in DDM1 nonlinear iteration (s)

    linear solverpreconditionmesh refinement
    0123
    degree of freedom (DOF)
    11 16543 925174 245694 085
    time/stime/stime/stime/s
    LU-0.209 6 (1)1.065 2 (1)7.018 4 (1)50.442 3 (1)
    BiCGSTABJacobi0.131 8 (106)0.868 3 (226)7.119 6 (455)68.205 1 (1 011)
    BiCGSTABASM0.126 5 (106)0.919 0 (226)7.666 1 (455)72.558 4 (1 011)
    BiCGSTABILU0.127 6 (106)0.840 7 (226)7.020 5 (455)69.624 0 (1 011)
    GMRESBJacobi0.263 5 (448)2.358 7 (940)26.995 4 (2 548)427.814 (9 450)
    GMRESASM0.291 4 (448)2.466 6 (940)27.758 8 (2 548)462.929 (9 450)
    GMRESILU0.283 8 (448)2.311 4 (940)27.126 1 (2 548)425.375 450)
    下载: 导出CSV

    表 2  弱扩展并行测试(Basic Newton,ASM+BiCGSTAB)

    Table 2  Weak extension parallel test (Basic Newton,ASM+BiCGSTAB)

    coresunknownsaverage iterationstime/sefficiency/%
    4 1.74×105 154 83.32 100
    16 6.94×105 353 151.70 54.9
    64 2.77×106 757 285.78 29.1
    256 1.11×107 1727 531.19 15.7
    下载: 导出CSV

    表 3  强扩展并行测试(Basic Newton,ASM+BiCGSTAB,2.771×106未知量)

    Table 3  Strongly extended parallel test (Basic Newton,ASM+BiCGSTAB,2.771×106 Unknowns)

    corestotal time/slinear solver time/sspeedupefficiency/%
    32 470.31 228.44 1.00 100.0
    64 284.78 113.67 1.65 82.6
    128 157.90 59.20 2.98 74.5
    256 99.56 32.31 4.72 59.1
    下载: 导出CSV
  • [1] 赵振国. PIN限幅器高功率微波效应机理仿真与实验研究[D]. 绵阳: 中国工程物理研究院, 2013: 1-53.

    Zhao Zhenguo. Simulation and experimental study of high power microwave effect of PIN limiter[D]. Mianyang: China Academy of Engineering Physics, 2013: 1-53
    [2] 孟凡宝, 杨周炳, 马弘舸, 等. 高功率微波超宽带电磁脉冲技术[M]. 北京: 国防工业出版社, 2011: 178-196.

    Meng Fanbao, Yang Zhoubing, Ma Hongge, et al. High-power microwave ultra-wideband electromagnetic pulse technology[M]. Beijing: National Defense Industry Press, 2011: 178-196
    [3] 马振洋, 柴常春, 任兴荣, 等. 双极晶体管微波损伤效应与机理[J]. 物理学报, 2012, 61:075011. (Ma Zhenyang, Chai Changchun, Ren Xinrong, et al. Microwave damage effect and mechanism of bipolar transistor[J]. Acta Physica Sinica, 2012, 61: 075011
    [4] 游海龙, 蓝建春, 范菊平, 等. 高功率微波作用下热载流子引起n型金属−氧化物−半导体场效应管特性退化研究[J]. 物理学报, 2012, 61:1085011. (You Hailong, Lan Jianchun, Fan Juping, et al. Research on characteristics degradation of n-metal-oxide-semiconductor field-effect transistor induced by hot carrier effect due to high power microwave[J]. Acta Physica Sinica, 2012, 61: 1085011
    [5] 陈曦, 杜正伟, 龚克. 脉冲宽度对PIN限幅器微波脉冲热效应的影响[J]. 强激光与粒子束, 2010, 22(7):1602-1606. (Chen Xi, Du Zhengwei, Gong Ke. Effect of pulse width on thermal effect of microwave pulse on PIN limiter[J]. High Power Laser and Particle Beams, 2010, 22(7): 1602-1606 doi: 10.3788/HPLPB20102207.1602
    [6] 周怀安. 电磁脉冲对电子设备中半导体有源器件作用的研究[D]. 北京: 清华大学电子工程系, 2005: 1-32.

    Zhou Huaian. Research on the effect of electromagnetic pulses on semiconductor active devices in electronic equipment[D]. Beijing: Department of Electronic Engineering, Tsinghua University, 2005: 1-32
    [7] 陈曦, 杜正伟, 龚克. 外电路在电磁脉冲对双极型晶体管作用过程中的影响[J]. 强激光与粒子束, 2007, 19(7):1197-1202. (Chen Xi, Du Zhengwei, Gong Ke. Influence of circuit during injection of EMP into bipolar junction transistor[J]. High Power Laser and Particle Beams, 2007, 19(7): 1197-1202
    [8] 韦源, 谢红刚, 贡顶, 等. 金属氧化物半导体场效应管长期辐射效应的数值模拟[J]. 强激光与粒子束, 2013, 25(4):1031-1034. (Wei Yuan, Xie Honggang, Gong Ding. Numerical simulation of long-term radiation effects for MOSFETs[J]. High Power Laser and Particle Beams, 2013, 25(4): 1031-1034 doi: 10.3788/HPLPB20132504.1031
    [9] Kramer K M, Hitchon W N G. Semiconductor devices: A simulation approach[M]. New Jersey: Prentice Hall PTR, 1997:20-54.
    [10] 贡顶, 张相华. 半导体器件的数值模拟: GSS软件用户手册[R]. 2008: 14-150.

    Gong Ding, Zhang Xianghua. Numerical simulation of semiconductor devices: GSS software user manual[R]. 2008: 14-150
    [11] 施敏, 伍国珏, 耿莉, 等. 半导体器件物理[M]. 西安: 西安交通大学出版社, 2008:5-193.

    Shi Min, Wu Guojue, Geng Li, et al. Physics of semiconductor devices[M]. Xi’an: Xi’an Jiaotong University Press, 2008:5-193
    [12] 李荣华. 偏微分方程数值解法[M]. 北京: 高等教育出版社, 2005: 45-213.

    Li Ronghua. Numerical solution of partial differential equations[M]. Beijing: Higher Education Press, 2005:45-213
    [13] Balay S, Abhyankar S, Adams M, et al. PETSC users manual revision 3.8[R]. No. ANL-95/11 Rev 3.8, 2017: 73-129.
    [14] Speelpenning B. Compiling fast partial derivatives of functions given by algorithms[R]. UIUCDCS-R-80-1002, 1980: 8-82.
    [15] Synopsys. Sentaurus device user guide[M]. Version H, 2013: 973-986
    [16] Lin P T. Improving multigrid performance for unstructured mesh drift-diffusion simulations on 147,000 cores[J]. International Journal for Numerical Methods in Engineering, 2012, 91(9): 971-989. doi: 10.1002/nme.4315
    [17] 领域编程框架团队. 并行自适用非结构网格应用支撑软件框架JAUMIN用户指南[M]. 北京: 北京应用物理与计算数学研究所, 2015.

    Domain Programming Famework Team. Parallel self-adaptive unstructured grid application support software framework JAUMIN user guide[M]. Beijing: Institute of Applied Physics and Computational Mathematics, 2015
    [18] Liu, Qingkai, Mo Zeyao, Zhang Aiqing, et al. JAUMIN: a programming framework for large-scale numerical simulation on unstructured meshes[J]. CCF Transactions on High Performance Computing, 2019, 1(1): 35-48. doi: 10.1007/s42514-019-00001-z
    [19] http://www.caep-scns.ac.cn/JAUMIN.php[EB/OL].
    [20] 肖丽, 曹小林, 王华维, 等. 激光聚变数值模拟中的大规模数据可视分析[J]. 计算机辅助设计与图形学学报, 2014, 26(5):675-686. (Xiao Li, Cao Xiaolin, Wang Huawei, et al. Large-scale data visual analysis for numerical simulation of laser fusion[J]. Chinese Journal of Computer-Aided Design and Computer Graphics, 2014, 26(5): 675-686
    [21] http://www.caep-scns.ac.cn/TeraVAP.php[EB/OL].
    [22] 郑澎, 方维, 徐权, 等. 面向JAUMIN的并行AFT四面体网格生成[J]. 计算机科学与探索, 2018, 12(4):567-574. (Zheng Peng, Fang Wei, Xu Quan, et al. Parallel AFT tetrahedral mesh generation for JAUMIN[J]. Journal of Frontiers of Computer Science and Technology, 2018, 12(4): 567-574 doi: 10.3778/j.issn.1673-9418.1611083
    [23] http://www.caep-scns.ac.cn/SuperMesh.php[EB/OL].
    [24] Lin P T, Shadid J N, Sala M, et al. Performance of a parallel algebraic multilevel preconditioner for stabilized finite element semiconductor device modeling[J]. Journal of Computational Physics, 2009, 228(17): 6250-6267. doi: 10.1016/j.jcp.2009.05.024
  • [1] 李勇, 贡顶, 宣春, 夏洪富, 谢海燕, 王建国. 器件效应数值模拟中漂移扩散模型的误差分析[J]. 强激光与粒子束, 2014, 26(06): 063204. doi: 10.11884/HPLPB201426.063204
    [2] 张家雷, 谭福利, 刘仓理, 王伟平. 激光熔蚀金属平板的三维有限体积法数值分析[J]. 强激光与粒子束, 2015, 27(08): 081002. doi: 10.11884/HPLPB201527.081002
    [3] 何定全, 李相强, 刘庆想. 高功率径向线阵列天线充气区域气体泄漏特性[J]. 强激光与粒子束, 2014, 26(06): 063035. doi: 10.11884/HPLPB201426.063035
    [4] 陈海林, 陈彬, 李正东, 易韵, 陆峰. 不同电磁脉冲作用下地面有限长电缆外导体感应电流的数值计算[J]. 强激光与粒子束, 2004, 16(10).
    [5] 李中杰, 李永东, 王洪广, 林舒, 刘纯亮. 半导体断路开关电路-流体耦合数值模拟[J]. 强激光与粒子束, 2010, 22(06).
    [6] 吴文斌, 于颖锐, 向宏志, 甯忠豪, 李庆. 基于大规模并行计算的三维多群中子扩散方程有限差分方法[J]. 强激光与粒子束, 2017, 29(08): 086001. doi: 10.11884/HPLPB201729.160328
    [7] 贡顶, 韩峰, 王建国. 基于流体动力学模型的2维砷化镓金属半导体场效应管数值模拟[J]. 强激光与粒子束, 2006, 18(07).
    [8] 孟雪松, 鲍献丰, 刘德赟, 周海京. 嵌入式薄片模型在时域有限差分算法中的应用[J]. 强激光与粒子束, 2017, 29(12): 123203. doi: 10.11884/HPLPB201729.170196
    [9] 林振, 梁昌洪. 负媒质模型的时域有限差分法分析[J]. 强激光与粒子束, 2006, 18(04).
    [10] 何锋, 苏建仓, 李永东, 刘纯亮, 孙鉴. 半导体断路开关数值模拟[J]. 强激光与粒子束, 2005, 17(12).
    [11] 李进玺, 程引会, 周辉, 吴伟. 用传输线和时域有限差分法计算电缆X射线响应[J]. 强激光与粒子束, 2007, 19(12).
    [12] 王强, 关永超, 王刚华, 谢龙, 蒋吉昊. C型固体电枢3维有限元电磁发射计算[J]. 强激光与粒子束, 2010, 22(04).
    [13] 李智慧, 唐靖宇, 张伦. 有限积分理论(FIT)及其在腔体计算中的应用[J]. 强激光与粒子束, 2002, 14(01).
    [14] 彭斓, 杨中海, 胡权, 黄桃, 李斌. 通电螺线管2维磁场有限元计算[J]. 强激光与粒子束, 2011, 23(08).
    [15] 朱湘琴, 王建国, 王玥, 王光强, 陈再高. 2.5维模拟中时域有限差分近-远场变换方法[J]. 强激光与粒子束, 2011, 23(08).
    [16] 袁志军, 楼祺洪, 周军, 董景星, 魏运荣, 王之江. 脉冲激光晶化非晶硅薄膜的有限差分模拟[J]. 强激光与粒子束, 2008, 20(07).
    [17] 林舒, 李永东, 王洪广, 刘纯亮. 半导体断路开关截断过程模拟的缩比模型[J]. 强激光与粒子束, 2013, 25(09): 2341-2345. doi: 10.3788/HPLPB20132509.2341
    [18] 张霖, 钟建, 陈玉成. 喷雾工艺制备的高性能有机光敏薄膜晶体管[J]. 强激光与粒子束, 2013, 25(03): 787-790. doi: 10.3788/HPLPB20132503.0787
    [19] 王文兵, 周辉, 刘逸飞, 马良, 程引会. 二维单步交替方向隐式时域有限差分法吸收边界性能分析[J]. 强激光与粒子束, 2018, 30(10): 103204. doi: 10.11884/HPLPB201830.180052
    [20] 张青洪, 廖成, 李瀚宇, 周海京, 刘强, 盛楠. 基于JASMIN框架的抛物方程有限差分解法并行计算及其应用[J]. 强激光与粒子束, 2015, 27(08): 083204. doi: 10.11884/HPLPB201527.083204
  • 加载中
图(11) / 表(3)
计量
  • 文章访问数:  15
  • HTML全文浏览量:  14
  • PDF下载量:  0
出版历程
  • 收稿日期:  2019-06-25
  • 录用日期:  2019-12-05
  • 网络出版日期:  2020-02-10

半导体多物理效应并行计算程序JEMS-CDS-Device设计与实现

    作者简介: 李光荣(1988—),男,助理研究员,从事高性能电磁计算研究,器件与电路仿真软件研发:li_guangrong@iapcm.ac.cn
  • 1. 中国工程物理研究院 高性能数值模拟软件中心,北京 100088
  • 2. 北京应用物理与计算数学研究所,北京 100094
  • 3. 中国工程物理研究院 复杂电磁环境科学与技术重点实验室,四川 绵阳 621900

摘要: 针对复杂电磁环境下器件多物理效应机理研究需求,研发了半导体多物理效应并行计算程序JEMS-CDS-Device。介绍了JEMS-CDS-Device的架构设计与实现技术。程序基于非结构网格并行框架JAUMIN实现,采用有限体积法(FVM)离散,使用牛顿法全耦合求解“电-漂移扩散输运”问题。程序采用“内核+算法库”形式架构,支持2维和3维非结构网格、千万自由度问题并行求解,支持物理方程、离散算法、材料物理模型等的扩展开发。

English Abstract

  • 半导体器件广泛存在于各种现代设备中,其组成的电子系统面临着复杂电磁环境如高强电磁辐射、电离辐射等的直接威胁。基于实验的器件效应分析昂贵、费时,且不易分析物理机理。半导体器件仿真软件(TCAD)是进行复杂环境下器件特性分析的一种重要手段,它是建立在半导体物理基础上的数值模拟仿真工具。国内学者采用TCAD对高功率微波作用下的器件特性退化与损伤、电离辐射对器件性能影响等问题进行了很多相关研究[1-8]。半导体器件模拟基于求解电磁-载流子输运耦合方程组,且须解决多介质、复杂边界、非线性多物理耦合等问题。JEMS-CDS-Device基于非结构网格并行框架JAUMIN开发,采用有限体积方法,支持2维和3维器件仿真和相关组件的扩展开发。本文介绍了JEMS-CDS-Device软件的主要架构与关键实现技术及其验证与应用。

    • 半导体仿真一般需耦合求解三个基本方程:电势泊松方程、电子连续性方程和空穴连续性方程,即

      $\begin{split} & \nabla \cdot \varepsilon \nabla \psi = - q\left( {p - n + N_D^ + - N_A^ - } \right) - {\rho _{\rm{s}}} \\ & \dfrac{{\partial n}}{{\partial t}} = \dfrac{1}{q}\nabla \cdot {{ J}_{\rm{n}}} - (U - G) \\ & \dfrac{{\partial p}}{{\partial t}} = - \dfrac{1}{q}\nabla \cdot {{ J}_{\rm{p}}} - (U - G) , \end{split}$

      (1)

      式中:$\varepsilon $是介电常数;$\psi $是电势;$q$是单位电荷量;$n$$p$分别是电子和空穴浓度;$N_A^ - $$N_D^ + $是杂质离子浓度;${\rho _s}$是固定电荷量;${ J_n}$${J_p}$分别是电子与空穴电流密度;UG分别是载流子的生成率与复合率。

      漂移−扩散输运模型(DDM)是波尔茨兹曼方程的近似特例,适用于描述较大尺度半导体器件中的载流子输运现象[9-11]

      $\begin{array}{l} {{{\vec J}_n} = q{\mu _{\rm{n}}}n{{ {{E}}}_n} + q{D_n{\rm{}}}\nabla n} \\ {{{\vec J}_p} = q{\mu _{\rm{p}}}p{{ {{E}}}_p} - q{D_{\rm{p}}}\nabla p} \end{array}$

      (2)

      式中:${\mu _{\rm{n}}}$${\mu _{\rm{p}}}$分别是电子与空穴的迁移速率;${D_{\rm{n}}}{\rm{ = }}\dfrac{{{k_{\rm{B}}}T}}{q}{\mu _{\rm{n}}},\;{D_{\rm{p}}}{\rm{ = }}\dfrac{{{k_{\rm{B}}}T}}{q}{\mu _{\rm{p}}}$kb是玻耳兹曼常数,T为温度;${ {{E}}_n}$${{{E}}_p}$是电子和空穴受到的有效电场强度,考虑到能带变化等效应[10],可表示为

      $\begin{split} & {{{{{E}}}_n} = \dfrac{1}{q}\nabla \left( { - q\psi - \chi - \Delta {E_{\rm{c}}}} \right) - \dfrac{{{k_{\rm{B}}}T}}{q}\nabla \left[ {\ln \left( {{N_{\rm{c}}}} \right) - \ln \left( {{T^{3/2}}} \right)} \right]} \\ & {{{{{E}}}_p} = \dfrac{1}{q}\nabla \left( {{E_{\rm{c}}} - {E_{\rm{g}}} + \Delta {E_{\rm{v}}}} \right) + \dfrac{{{k_{\rm{B}}}T}}{q}\nabla \left[ {\ln \left( {{N_{\rm{v}}}} \right) - \ln \left( {{T^{3/2}}} \right)} \right]} , \end{split}$

      (3)

      式中:$\chi $是电势;χ是电子亲和能;${E_{\rm{c}}}$是导带底能级;${E_{\rm{g}}}$是能隙;$\Delta {E_{\rm{c}}}$$\Delta {E_{\rm{v}}}$是重掺杂或压力引起的能带偏移;${N_{\rm{c}}}$${N_{\rm{v}}}$是态密度。

      如考虑温度等效应,还需与热传导方程等耦合求解[10]

      $\rho {c_{\rm{P}}}\frac{{\partial T}}{{\partial t}} = \nabla \cdot \kappa \nabla T + {{J}} \cdot {{E}} + \left( {{E_g} + 3{k_{\rm{B}}}T} \right) \cdot \left( {U - G} \right) ,$

      (4)

      式中:${c_{\rm{P}}}$$\rho $分别为材料的热容量和密度;$\kappa $为热导率。公式中右侧第二项为欧姆热,其中,$ {{J}}$为电流密度矢量;${{E}}$为电场;第三项为载流子复合与生成导致的放热与吸热,其中,U为载流子复合率,G为载流子电离率。考虑到热效应,载流子扩散电流${J_{{\rm{n}},{\rm{diff}}}}$${J_{{\rm{p}},{\rm{diff}}}}$应分别表达为[10]

      $\begin{split} & {{J_{{\rm{n}},{\rm{diff}}}} = {\mu _{\rm{n}}}{k_{\rm{b}}}\left( {T\nabla n + n\nabla T} \right)} \\ & {{J_{{\rm{p}},{\rm{diff}}}} = - {\mu _{\rm{p}}}{k_{\rm{b}}}\left( {T\nabla p + p\nabla T} \right) } \end{split} $

      (5)

      控制方程组改写为[10]

      $\left\{ \begin{array}{l} \dfrac{{\partial n}}{{\partial t}} = \nabla \cdot \left( {{\mu _{\rm{n}}}n{{{{E}}}_{\rm{n}}} + {\mu _{\rm{n}}}\dfrac{{{k_{\rm{b}}}T}}{q}\nabla n + {\mu _{\rm{n}}}\dfrac{{{k_{\rm{B}}}\nabla T}}{q}n} \right) - \left( {U - G} \right) \\ \dfrac{{\partial p}}{{\partial t}} = - \nabla \cdot \left( {{\mu _{\rm{p}}}p{{{{E}}}_{\rm{p}}} - {\mu _{\rm{p}}}\dfrac{{{k_{\rm{b}}}T}}{q}\nabla p - {\mu _{\rm{p}}}\dfrac{{{k_{\rm{B}}}\nabla T}}{q}p} \right) - \left( {U - G} \right) \\ \nabla \cdot \varepsilon \nabla \psi = - q(p - n + N_D^ + - N_A^ - ) - {\rho _{\rm{s}}} \\ \rho {c_{\rm{p}}}\dfrac{{\partial T}}{{\partial t}} = \nabla \cdot \kappa \nabla T + {{J}} \cdot {{E}} + \left( {{E_{\rm{g}}} + 3{k_{\rm{B}}}T} \right) \cdot \left( {U - G} \right) \end{array} \right..$

      (6)

      迁移率与载流子有效质量和被散射几率有关[11]。迁移率模型分为弱场模型、强场模型和表面模型,一般是载流子浓度、场强和温度等相关的非线性函数。常见的有Analytic模型、Philips模型和Lucent模型等[10]。载流子复合模型一般由SRH复合、Auger复合和直接复合三种机制组合而成[9-11]。载流子生成率模型一般考虑碰撞电离效应和能带间隧道跃迁两种机制[9-11]

    • 器件的三个控制方程皆可写为如下守恒形式

      $\frac{{\partial {Y_i}}}{{\partial t}} - \nabla \cdot {J_{{Y_i}}} + {R_{{Y_i}}} = 0 $

      (7)

      为保证守恒性,一般采用有限体积法离散此方程。程序采用顶点型有限体积法(CV-FVM),在原网格上进行场离散与通量插值,在对偶网格(外心型,如图1)上进行积分。方程两边对控制体Ci积分,并使用高斯公式将散度体积分转换为通量的表面积分,即[12]

      图  1  插值网格与积分网格(控制体)

      Figure 1.  Interpolation grid and integral grid (control body)

      $\int_{{C_i}} {\frac{{\partial {Y_i}}}{{\partial t}}{\rm{d}}V} - \int_{\partial {C_i}} {{J_{{Y_i}}} \cdot \hat n{\rm{d}}S{\rm{ + }}\int_{{C_i}} {{R_{{Y_i}}}{\rm{d}}V} = 0} $

      (8)

      载流子输运问题通常是一种对流占优问题,一般需要采用迎风离散格式保持计算稳定。半导体模拟通常使用Scharfetter-Gummel离散格式(SG)离散电流,SG格式是一种自适应权重迎风格式。不考虑温度变化时,基于漂移-扩散模型(DDM)的电流通量密度SG格式可以表示为[9-10]

      $\begin{split} & {\left. {{J_{\rm{n}}}} \right|_{{\rm{mid}}}} = \dfrac{{q{V_{\rm{T}}}}}{l}{\left. {{\mu _n}} \right|_{\rm{mid}}}\left[ {{n_j}B\left( {\dfrac{{{\psi _j} - {\psi _i}}}{{{V_{\rm{T}}}}}} \right) - {n_i}B\left( {\dfrac{{{\psi _i} - {\psi _j}}}{{{V_{\rm{T}}}}}} \right)} \right] \\ & {\left. {{J_{\rm{p}}}} \right|_{{\rm{mid}}}} = \dfrac{{q{V_T}}}{l}{\left. {{\mu _{\rm{p}}}} \right|_{{\rm{mid}}}}\left[ {{p_i}B\left( {\dfrac{{{\psi _j} - {\psi _i}}}{{{V_{\rm{T}}}}}} \right) - {p_j}B\left( {\dfrac{{{\psi _i} - {\psi _j}}}{{{V_{\rm{T}}}}}} \right)} \right] \\ & B(x) = \dfrac{x}{{{{\rm{e}}^x} - 1}} \end{split} $

      (9)

      式中:${\left. {{J_{\rm{n}}}} \right|_{{\rm{mid}}}}$${\left. {{J_{\rm{p}}}} \right|_{{\rm{mid}}}}$分别表示棱边中点的电子与空穴电流密度;${V_{\rm{T}}} = \dfrac{{{k_{\rm{B}}}T}}{q}$为热电压常数;$l$为棱边长度;np为节点上的电子、空穴,下标ij表示一条棱边上的两个节点编号。为考虑温度变化的SG电流密度离散公式较为复杂,可参考文献[10]。

      半导体方程离散后,采用全耦合牛顿法直接求解非线性方程组,对于$F\left( x \right) = 0$[13],有

      $ \begin{split} & {J\left( {{x_k}} \right)\Delta {x_k} = - F\left( {{x_k}} \right)}\\ & {{x_{k + 1}} = {x_k} + \Delta {x_k}} \end{split} $

      (10)

      式中:k表示迭代步。对于瞬态求解,半导体方程稳定时间步长一般远小于模拟总时间,故时间离散一般采用隐式欧拉法。

    • 牛顿法需建立Jacobia矩阵,其计算难度以平方比例增长。自动微分技术起源于Bill Gear的论文[14],其原理基于复合函数链式求导法则,可以获得偏导数的数值精确解,示例如图2所示。

      图  2  前向自动微分

      Figure 2.  Forward automatic differentiation

    • 在牛顿法迭代过程中,载流子密度可能出现负值(载流子密度物理上不小于0),负密度易引起牛顿迭代发散。需引入密度修正[10]

      ${\rho _{\rm{c}}} = \left\{ \begin{matrix} {\begin{array}{*{20}{l}} { {\rho _{\rm{c}}} },&{{\rho _c} > 0} \\ {\varepsilon \left( {{\rho _c}} \right)},&{{\rho _c} \leqslant 0} \end{array}} \end{matrix} \right.$

      (11)

      式中:ρc为电子或空穴载流子密度;$\varepsilon \left( {{\rho _{\rm{c}}}} \right)$为远小于ρc量级的一个正数。

    • CV-FVM方法对网格有一定要求,当原网格单元为含有钝角的非Delaunay网格或在边界上存在含有钝角的Delaunay网格时,其对偶网格会产生“负体积”,并引起积分错误,如图3所示。程序中需对积分网格按单元进行截断修正[15]

      图  3  积分网格(对偶网格)修正

      Figure 3.  Integral mesh (dual mesh) correction

    • 器件仿真中存在大量的间断型边界,间断型边界一般出现在两种不同材料的界面处,如金属−半导体界面、异质半导体界面、辐射问题中的SiO2-Si界面。物理量在边界处不连续,出现了“跳跃”,如图4(a)。程序中在界面物理节点处引入了额外的逻辑节点[15],用于区分不同算法区域,如图4(b)

      图  4  异质界面处的间断[15]与算法区域处理

      Figure 4.  Discontinuity heterogeneous interface and algorithm region processing (right figure, there are two logical nodes M1 and M2 in node A)

    • 半导体DD模型全耦合生成的Jacobia矩阵为稀疏带状,且一般为非对称、强多尺度。此类型问题求解一般需采用直接法(LU)或使用“ILU类预条件+Krylov迭代”[16]表1为某二维PN器件生成的线性方程串行求解时线性解法器收敛情况对比,每次加密后未知量增加到原来的4倍。计算规模较小的情况下,直接法速度较快。表1中时间后的括号中为线性解法器迭代次数(直接法LU迭代次数为1)。

      表 1  DDM1非线性迭代中线性解法器典型收敛情况(时间单位:秒)

      Table 1.  Typical convergence of linear solver in DDM1 nonlinear iteration (s)

      linear solverpreconditionmesh refinement
      0123
      degree of freedom (DOF)
      11 16543 925174 245694 085
      time/stime/stime/stime/s
      LU-0.209 6 (1)1.065 2 (1)7.018 4 (1)50.442 3 (1)
      BiCGSTABJacobi0.131 8 (106)0.868 3 (226)7.119 6 (455)68.205 1 (1 011)
      BiCGSTABASM0.126 5 (106)0.919 0 (226)7.666 1 (455)72.558 4 (1 011)
      BiCGSTABILU0.127 6 (106)0.840 7 (226)7.020 5 (455)69.624 0 (1 011)
      GMRESBJacobi0.263 5 (448)2.358 7 (940)26.995 4 (2 548)427.814 (9 450)
      GMRESASM0.291 4 (448)2.466 6 (940)27.758 8 (2 548)462.929 (9 450)
      GMRESILU0.283 8 (448)2.311 4 (940)27.126 1 (2 548)425.375 450)
    • JEMS-CDS-Device并行基于JAUMIN并行计算框架实现,该框架由北京应用物理与计算数学研究所研制。JAUMIN采用网格层-网格片-单元的分布式网格管理结构,使用并行构件搭建计算流程,构件遍历网格片并调用网格片数值子程序完成计算,用户需实现网格片策略类中特定函数来实现构件功能定制[17-19]。JEMS-CDS-Device计算结果的后处理绘图由TeraVAP[20-21]软件完成,软件输入网格生成使用到了SuperMesh[22-23]程序。

      图  5  网格层次结构与并行构件[17]

      Figure 5.  Grid hierarchy and parallel components

      JEMS-Device基于离散中台DDMS构建,主要层次关系如图6所示。DDMS主要由内核与特征功能库两部分组成。内核包含数据模型、算法框架及方程构件。方程构件集成了各组件,并对外提供了统一的接口。算法框架定义了时间积分、非线性FVM求解、流程控制等的算法框架(模板)及接口。数据模型提供了器件结构的抽象描述和相关数据的存储。物理模型库中提供了常用材料的能带、迁移率等的计算模型,并通过统一的接口被算法调用。方程模板中包含了特定输运模型方程组的数值离散与组装算法,通过加载不同的方程模板,DDMS可以生成用于解决特定输运模型的方程构件。各模块之间的关系如图7

      图  6  JEMS-CDS-Device主要层次

      Figure 6.  JEMS-CDS-Device

      图  7  软件模块间关系

      Figure 7.  Relationships between software modules

    • 基本算例1,三维二极管稳态问题,此算例用于验证仿真程序基本功能。算例器件为在掺杂浓度1015 cm−3的均匀N型Si衬底上注入1019 cm−3受主杂质形成。图8(a)为器件净掺杂密度分布图,器件尺寸为6 μm×3 μm×1 μm,图中器件上面中心区域为阳极,下面为阴极。本算例中采用了Analytic迁移率模型,仿真了器件在0~1 V范围内正偏时的状态,得到了器件阳极的电压-电流关系,并与Genius TCAD软件了对比,结果符合良好,如图8(b)所示。

      图  8  三维PN净掺杂密度分布与二极管阳极V-I曲线

      Figure 8.  PN diode ’s density distribution of net doping and Anode V-I curve

      基本算例2,2维NMOS管,此算例用于验证程序对异质间断型边界的处理能力。场效应晶体管(MOSFET)是现代集成电路中的重要器件,一般由金属、氧化物和半导体组成。算例器件为一个μm级工艺的NMOS管,由Si,SiO2和金属三种类型材料组成,仿真时Si-SiO2界面设定存在密度为1010 cm−2的界面电荷,程序中此界面作为间断型边界处理。本算例仿真得到了其平衡态电势分布,并与Genius TAD的仿真进行对比,电势分布结果相符,如图9所示。

      图  9  二维NMOS管平衡态电势分布结果对比(电势单位:V,参考零点为无穷远)

      Figure 9.  2D NMOS equilibrium potential distribution comparison(potential unit: V, the potential reference zero is at infinity)

    • 二维PN算例,此算例器件为一3 μm×3 μm的Si材质二极管,净掺杂分布如图10(a)中所示,测试中对此PN管施加了0.5 V的正向偏压,电子浓度分布如图10(b)所示。线性解法器采用BiCGSTAB,预条件使用ASM。ASM是一种加性区域分解预条件[13],每个分解区域上单独采用ILU[13, 16]生成子区域预条件。

      图  10  二维PN管的净掺杂分布与载流子分布(正向偏压0.5 V)

      Figure 10.  Net doping and carrier distribution of 2D PN device (positive deviation 0.5 V)

      表2表3展示了JEMS-CDS-Device的弱扩展测试(保持每核心计算规模近似不变)与强扩展测试(保持总计算规模不变)情况。由表3可知,程序在256核并行强扩展效率为59.1%(32核基准),效率较好。由表2中可知,程序在256核并行时相对4核并行的弱扩展效率为15.7%,效率下降的主要原因在于,随着每次未知量的增多(未知量4倍,并行核数4倍),牛顿迭代每步线性求解需要的平均迭代数也随之增加(约2倍,与未知量数量的平方根成正比[24])。

      表 2  弱扩展并行测试(Basic Newton,ASM+BiCGSTAB)

      Table 2.  Weak extension parallel test (Basic Newton,ASM+BiCGSTAB)

      coresunknownsaverage iterationstime/sefficiency/%
      4 1.74×105 154 83.32 100
      16 6.94×105 353 151.70 54.9
      64 2.77×106 757 285.78 29.1
      256 1.11×107 1727 531.19 15.7

      表 3  强扩展并行测试(Basic Newton,ASM+BiCGSTAB,2.771×106未知量)

      Table 3.  Strongly extended parallel test (Basic Newton,ASM+BiCGSTAB,2.771×106 Unknowns)

      corestotal time/slinear solver time/sspeedupefficiency/%
      32 470.31 228.44 1.00 100.0
      64 284.78 113.67 1.65 82.6
      128 157.90 59.20 2.98 74.5
      256 99.56 32.31 4.72 59.1
    • 电热耦合仿真通过联立求解基本DDM方程与晶格热传导方程,得到器件在外部激励源作用下的瞬态温度分布。本算例中PIN二极管模型为Si材质的PNN管,主体为直径约20 μm的柱体,净掺杂峰值绝对值约为9.17×1019 cm−3,采用其过中心截面进行二维仿真。器件两电极均设为欧姆接触型,阳极边界热传导系数取某金属3.17 W·cm−1·K−1,器件后部因截断了部分硅,故阴极导热系数取硅的导热系数1.56 W·cm−1·K−1。载流子迁移率模型采用Philips模型并使用了Caughey-Thomas模型进行了强场修正,复合率使用了SRH, Auger和直接复合机制模型计算,碰撞电离率采用了适用于高温的Valdinoci模型。器件端口激励的波形如图11(a)所示,图11(b)为此激励下7.78 ns时的温度分布情况,局部放大图中可见,PIN管PN结边缘处温度较高(最高707.5 K),边缘阳极接触部分亦有明显升温,这与实际经验相吻合。

      图  11  输入的波形与某时刻PIN的温度分布

      Figure 11.  Input waveform and temperature distribution of PIN (at 7.78 ns)

    • 本文介绍了JEMS-CDS-Device程序的架构与实现。JEMS-CDS-Device基于JAUMIN框架与DDMS中台构建,围绕多介质、多物理、非线性全耦合数值模拟这一需求,主要实现了3方面功能:(1)支持多介质,间断型边界,可用于半导体、导体、绝缘体组成的复杂器件模型,支持典型间断型边界问题的高效计算;(2)材料物理模型可扩展,支持物理模型库扩展开发;(3)基于Newton法与自动微分技术,实现非线性多物理方程全耦合求解。算例验证与并行性能测试的情况表明,JEMS-CDS-Device程序对非线性多物理全耦合问题具备千万未知量、数百CPU核的并行计算能力,可以应用于强电磁脉冲等复杂电磁环境下半导体器件多物理效应机理的研究。

参考文献 (24)

目录

    /

    返回文章
    返回