当前半导体市场上最先进的32 nm技术已经接近集成电路制造工艺的物理极限,依靠工艺技术本身的提升来进一步提高电子器件集成度已经很难,于是三维封装技术便应运而生.三维封装技术的主流是芯片的堆叠,各芯片间通过垂直互连实现电通讯,由此大大减少了芯片间引线的长度和电信号的延迟,使得电子器件集成度提高的同时电性能也得以改进.目前堆叠芯片之间的垂直互联主要是通过硅通孔(TSV,Through-Silicon-Vias)来实现的.
TSV的直径通常为数微米到数十微米量级,内壁上沉积有数纳米厚的扩散阻挡层和1~2 μm的SiO2绝缘层,并以电镀方式将孔用铜填充作为电信号通道.
由于TSV是三维电子封装的核心器件,准确掌握TSV的失效机制和影响因素,对于TSV的优化设计和提高整个三维封装的可靠性极为关键.最近两年,TSV铜填充的电迁移及其引发的可靠性问题开始受到重视[1, 2].电迁移是指导金属材料在通以密度较高的电流时,金属原子沿着电子运动方向或者是空位沿着电子运动反方向进行迁移的质流扩散现象[3].电迁移使得阴极和阳极处分别发生负的和正的体积应变,当变形受到约束时便会产生机械应力(即本文所述的电致应力),而且阳极的压应力和阴极的拉应力之间的应力梯度对原子迁移有抑制作用,这些均已得到了实验证实[4, 5, 6, 7, 8].电迁移引发的TSV铜填充的变形被周围坚硬的材料束缚,必然产生电致应力[9],而已有的关于TSV可靠性的研究中往往漏掉了这一项.按照目前通用的唯象学评价方法,电致应力项的遗漏有可能导致对TSV可靠性做出错误评价.
本文基于力-扩散耦合方程,详细地给出了弹性材料属性下电致应力、应变模型的有限元公式的推导过程,建立了弹性材料属性下电致应力、应变的有限元模型.此耦合模型准确地阐明了电迁移过程中扩散场和位移场的相互作用及其影响,并能计算电致应力场、应变场以及空位浓度场的演化过程.进一步,利用ABAQUS中的用户自定义单元(user defining element)接口,实现了有限元模型的计算,并对该有限元模型进行了验证性计算.最后,利用本文的电迁移有限元模型对TSV中铜填充的电致应力进行了数值仿真计算. 1 电致应力的有限元分析模型 1.1 基本方程
电迁移是质量输运的扩散过程,可采用以下空位运输方程来描述[10]：
式(1)将研究对象视为空位的封闭系统.忽略电流密度不均匀引起的温度梯度对质流的驱动,仅考虑浓度梯度、电迁移和应力梯度对空位通量的贡献(即驱动力),空位通量q的表达式如下所示[9]：
空位的产生或者消失率代表着应力和空位的关系,表达式如下[10]：
是正则化的平衡空位浓度;τs是空位产生或者消失的时间.式(3)仅适用于空位浓度偏离平衡空位浓度较小时的情况.
无体力的应力平衡方程为
考虑镀铜材料为线弹性情况,总应变包括了纯弹性应变分量εeij和电迁移引起的非弹性应变分量εelij,即
结合式(1),方程(6)可变为
建立应力-应变本构模型如下：
几何方程,即总应变和位移的关系如下：
利用加权余量法,式(1)可写为
利用分部积分和散度定理,式(11)的Δ·q变为
将式(12)代入式(11)中,均采用正则化参量,消去Cv0得到
利用空位通量的封闭边界条件,即q·n=0,同时将式(2)和式(3)代入上式可得
用伽辽金法,变量场用权函数插值得到
由此得到了扩散方程的离散形式,由于δcN的任意性,上式大括号中的项为零,即此时扩散控制方程为
利用虚功原理,可以得到平衡方程(忽略体力)的等效积分弱形式,其矩阵形式为[11]
同样,变量场用权函数插值得到
由于δuM的任意性,即此时平衡方程为
若给定的节点位移和空位浓度值是假设的,即它们不是精确解时,针对位移和空位浓度分别存在以下残差：
最佳近似计算的依据是使残差最小.为此,首先由方程(23)对cN求导可得
分别记
即有
由方程(23)再对uM求导可得
分别记
即有
由方程(22)对uM求导可得
由方程(22)再对cN求导可得
最后,得到的Jacobian矩阵为
目前一般的商业有限元软件不能直接求解电迁移问题,必须进行二次开发.为了简化问题,这里只考虑电致应力、应变的平面问题,并忽略电迁移引发的材料非线性.运用大型商业软件ABAQUS的用户自定义单元子程序接口来实现电迁移模型的有限元计算.针对平面电致应力模型,选择8节点矩形单元,节点自由度除位移u和v之外,还增加了一个自由度-正则化的空位浓度c.基于Fortran语言,子程序中定义了残差向量,Jacobian矩阵以及依赖于解的状态变量等,这些结果均已在1.3节中给出.
为了验证有限元计算模型的正确性和有效性,本文选择与Sarychev和Zhinikov的解析模型进行对比来验证有限元方法的正确性.计算的问题如图 1所示,有一个无限大平板,厚度为L,在x=0处的平面固定,另一边自由,两边均为封闭边界条件,即边界上无空位扩散,沿厚度方向的电流密度为j=－j1e1,其中j1>0.对应此问题,建立的有限元模型及网格如图 2所示,可将问题简化为平面应变问题,并设置宽度W远大于厚度L,这里取W=200 μm,L=10 μm,使得宽厚比为20∶1,左端固定,右端自由,四边为封闭边界条件.计算中以Al材料为例,在理论计算及有限元计算过程中采用相同参数,以便对两者的结果进行对比.采用的材料属性参数如表 1所示,通过的电流密度值为j1=1.0 MA/cm2.
 |
| 图 2 无限大有限厚板的有限元模型及网格(1/2模型)Fig. 2 Finite element model and mesh of the plate shown in fig.1(1/2 model) |
文献[10]中给出的稳态应力、应变以及浓度的解析解为
图 3和图 4分别给出了稳态时球应力和球应变沿厚度方向的分布情况,从图中可以看出,稳态时球应力和球应变沿厚度方向是呈线性规律变化的,有限元计算出的规律与理论解一致,且有限元计算和理论结果吻合较好.
 |
| 图 3 球应力沿厚度方向的分布Fig. 3 Spheric stress distribution along the thickness direction |
 |
| 图 4 体积应变沿厚度方向的分布Fig. 4 Volumetric strain distribution along the thickness direction |
图 5给出了稳态时正则化空位浓度沿厚度方向的分布情况,从理论解式中可以看出空位浓度沿x方向是呈指数规律分布的,有限元计算出的分布规律与理论解基本一致.可以看出,左端空位浓度高,为空位集中区域,而右端空位浓度低,即为原子堆积区域.上述对比结果表明,本文的有限元分析模型和方法是正确、可靠的.
 |
| 图 5 空位浓度沿厚度方向的分布Fig. 5 Vacancy concentration distribution along the thickness direction |
本节考虑TSV在通以较高电流密度时的电迁移及电致应力现象[11].所建立的简化几何模型及其尺寸如图 6所示,其中,d=5 μm,h=30 μm,w=20 μm,t=2 μm,s=3 μm,记电流流入拐点为A点,流出拐点为B点.Cu连接四周为坚硬的Si及SiO2,故位移边界条件简化为四边完全约束,Cu连接周围有Cu原子扩散阻挡层,故对Cu原子/空位而言应为封闭系统.电流方向如图 6中所示,计算时取电流密度值为j=1.0×106 A/cm2,总时间t=2×105 s.Cu的计算参数如表 2所示[12].该结构与芯片内的双嵌入式Cu内连接的结构类似,研究结果具有参考意义.
 |
| 图 6 TSV结构模型及有限元网格Fig. 6 Model of TSV structure and finite element mesh |
首先,计算了TSV结构中的电流密度场,然后保持网格不变,将节点上的电流密度值应用于后续的电致应力的计算.图 7和图 8分别给出了2×105 s后TSV球应力和球应变的分布云图.此时电迁移已达到稳态,从图中可以看出,上端引线主要受压,从上到下逐渐过渡到拉应力区,由于空位由上向下迁移,阴极端附近空位集中,为拉应力区,阳极端附近空位浓度低,为原子集中区域,为压应力区,两端的应力最大,最大应力值为52.6 MPa.从图 8还可以看到,球应变主要集中在拐角处,球应变的最大值为1.9×10-5.
 |
| 图 7 稳态时TSV结构球应力的分布云图Fig. 7 Spheric stress distribution of TSV structure at steady state |
 |
| 图 8 稳态时TSV结构球应变(无量纲)的分布云图Fig. 8 Spheric strain distribution of TSV structure at steady state |
图 9给出的是稳态时正则化空位浓度的分布云图.可见,电迁移过程达到稳态时,整体空位浓度分布呈现明显的梯度分布形式(最高正则化空位浓度值为1.014,最低值为0.986 3).阳极空位浓度低,原子浓度就高,容易形成凸丘;阴极空位浓度高,容易形成孔洞等缺陷,这与有关文献报道[13, 14, 15]中双嵌入式Cu内连接结构的电迁移失效模式相符合.
 |
| 图 9 不同电流密度时TSV中的最大球应力Fig. 9 Maximum spheric stress at different current densities in TSV |
对TSV结构中电致应力随电流密度的变化也做了研究.图 9和图 10给出了电流密度分别为1,2,3,4，5 MA/cm2，并且迁移达到稳态的情况下,TSV结构中最大球应力和Mises应力随电流密度的变化.最大球应力出现在结构的两端(即引线处),而最大Mises应力出现在结构拐点处(即图 7中的A和B点).
 |
| 图 10 不同电流密度时TSV中的最大Mises应力Fig. 10 Maximum Mises stress at different current densities in TSV |
1) 本文给出了弹性材料属性下电致应力计算的有限元模型,此模型准确地阐明了电迁移过程中扩散场和位移场的相互作用及其影响,并能计算电致应力场、应变场以及空位浓度场的演化过程.
2) 利用ABAQUS中的用户自定义单元接口,实现了电致应力的有限元计算,并进行了验证性计算.
3) 针对TSV中铜填充的电迁移问题进行了计算,发现TSV结构中的电致球应力在电流密度为1~5 MA/cm2的情况下可达数十到数百兆帕,而Mises应力只有数十兆帕.阳极处过大的球应力可能是造成钝化层的破裂,而阴极处较大的球应力可能对界面的完整性造成不良影响.
上述计算结果为TSV可靠性的全面评估奠定了一定基础.文中忽略了铜材料本身的塑性以及电迁移引发的非线性问题,使得应力计算结果可能有所偏大,在后续的研究中会逐步解决.
致谢 在研究过程中同北京大学微电子所陈兢教授的研讨使作者受益匪浅,在此表示感谢！
| [1] | Tan Y C,Tan C M,Zhang X W,et al.Electromigration performance of through silicon via(TSV)-a modeling approach.[J].Microelectronics Reliability,2010,50(9-11):1336-1340 |
| Click to display the text | |
| [2] | Pak J,Pathak M,Lim S K,et al.Modeling of electromigration in through-silicon-via based 3D IC[C]//Electronic Components and Technology Conference (ECTC).Piscataway,NJ:Institute of Electrical and Electronics Engineers Inc,2011:1420-1427 |
| Click to display the text | |
| [3] | 杨卫.力电失效学[M].北京:清华大学出版社,2001 Yang Wei.Mechatronic reliability[M].Beijing:Tsinghua University Press,2001(in Chinese) |
| [4] | Blech I A.Electromigration in thin aluminum films on titanium nitride[J].Journal of Applied Physics,1976,47(4):1203-1208 |
| Click to display the text | |
| [5] | Blech I A,Herring C.Stress generation by electromigration[J].Applied Physics Letters,1976,29(3):131-133 |
| Click to display the text | |
| [6] | Wang P C,Cargill G S,Noyan I C,et al.Electromigration-induced stress in aluminum conductor lines measured by x-ray microdiffraction[J].Applied Physics Letters,1998,72(11):1296-1298 |
| Click to display the text | |
| [7] | Basaran C,Lin M.Electromigration induced strain field simulations for nanoelectronics lead-free solder joints [J].International Journal of Solids and Structures,2007,44(14/15):4909-4924 |
| Click to display the text | |
| [8] | Xu L H,Pang J H L.In-situ electromigration study on Sn-Ag-Cu solder joint by digital image speckle analysis[J].Journal of Electronic Materials,2006,35(11):1993-1999 |
| Click to display the text | |
| [9] | Ye H,Basaran C,Hopkins D C.Numerical simulation of stress evolution during electromigration in IC interconnect lines[J].IEEE Transactions on Component and Packaging Technologies,2003,26(3):673-681 |
| Click to display the text | |
| [10] | Sarychev M E,Zhinikov Y V,Borucki L,et al.General model for mechanical stress evolution during electromigration[J].Journal of Applied Physics,1999,86(6):3068-3075 |
| Click to display the text | |
| [11] | 王勖成.有限单元法[M].北京:清华大学出版社,2003 Wang Xucheng.Finite element method[M].Beijing:Tsinghua University Press,2003(in Chinese) |
| [12] | Basaran C,Lin M H,Li S D.Computational simulation of electromigration induced damage in copper interconnects[C]//Proceedings of the 2007 Summer Computer Simulation Conference.Vista,CA:the Society for Modeling and Simulation International,2007:261-268 |
| Click to display the text | |
| [13] | Oates A S,Lin M H.Electromigration failure distributions of Cu/low-k dual-damascene vias:impact of the critical current density and a newreliability extrapolation methodology[J].IEEE Trans Device Mater Reliab,2009,9(2):244-254 |
| Click to display the text | |
| [14] | Roy A,Tan C M.Probing into the asymmetric nature of electromigration performance of submicron interconnect via structure[J].Thin Solid Films,2007,515(7/8):3867-3874 |
| Click to display the text | |
| [15] | Chen Z H,Lv Z C,Wang X F,et al.Modeling of electromigration of the through silicon via interconnects[C]//11th International Conference on Electronic Packaging Technology & High Density Packaging.New York:Association for Computing Machinery,2010:1221-1225 |
| Click to display the text |