随着微电子封装技术进步，电子器件不断向小型化和高度集成化方向发展。由于塑料焊球阵列(Plastic Ball Grid Array, PBGA)封装的集成度高、电性能好并且易于装配到印制电路板上，PBGA封装技术已被广泛地用到电子器件的封装中。PBGA封装产品在制造和使用中会受到温度、湿气、振动等环境应力的影响，易造成可靠性问题并影响产品正常使用^[1-2]。有关调查数据表明，温度是造成PBGA封装产品失效的主要原因。封装结构中不同材料之间存在热膨胀系数差异，PBGA封装在回流焊温变过程中会产生翘曲变形。结构的翘曲会影响封装结构的共面度，引发芯片断裂、界面分层和焊点装联缺陷等质量和可靠性问题^[3-4]。因而深入了解与研究回流焊过程中PBGA的翘曲变形机制与规律，对提高其装联质量及热可靠性具有相当重要的意义。

本文基于PBGA封装回流焊组装的质量控制需求，采用有限元法对模拟回流焊温变过程中PBGA封装的翘曲变形行为及应力应变进行仿真分析。并采用阴影云纹法对翘曲变形的模拟分析结果进行测试验证，分析结果可以为PBGA装联质量可靠性提供支撑。

1 翘曲变形分析理论

电子封装产品由于材料参数差异，在温度变化过程中会发生翘曲变形，对产品的装联质量及可靠性有着重要影响。对于装配体的变形分析，最早是由Timoshenko^[5]提出的，主要针对双层金属的变形问题；而针对SIP(System In a Package，系统级封装)、BGA (Ball Grid Array，焊球阵列封装)、晶圆片级封装等先进封装的翘曲问题，国内外也有大量的相关研究^[6-8]。

PBGA封装通常是由芯片、粘接层、塑封层、基板等多层结构组成，如图1所示，其结构可以被当成一个多层板。图2所示为简化多层板的结构模型，X轴是中面位置，其所在的层为中性层；Z₀~Z_n是以中面为基准，分别到各层板边界的距离^[9-10]。当温度发生变化ΔT时，对于第k层板，其应力−应变的关系可以表示为

$ \left\{ {{\sigma }} \right\} = {\left[ {{Q}} \right]_k}\left[ {\left\{ {{\varepsilon }} \right\} - {{\left\{ {\rm{\alpha }} \right\}}_k}\Delta T} \right]{\rm{}}\left({k = 1,2, \cdots , n} \right) $

(1)

式(1)中，σ是应力，[Q]_k是刚度矩阵，ε是应变，α是热膨胀系数，ΔT是温度变化值。其中，[Q]_k是与弹性模量(E)和泊松比(μ)有关的刚度矩阵，表示为

$ {\left[ {{Q}} \right]_k} = \frac{{{E_1}}}{{1 - {{\rm{\mu }}^2}}}\left[ {\begin{array}{*{20}{c}} 1&{\rm{\mu }}&0\\ {\rm{\mu }}&1&0\\ 0&0&{\frac{{1 - {\rm{\mu }}}}{2}} \end{array}} \right] $

(2)

根据多层板经典理论的基本假设^[11]，{ε}可表示为

$ \left\{ {{\varepsilon }} \right\} = \left\{ {{{{\varepsilon }}^0}} \right\} + Z\left\{ {{K}} \right\} $

(3)

$ \left\{ {{{{\varepsilon }}^0}} \right\} = \left\{ {\begin{array}{*{20}{c}} {{{\varepsilon }}_{\rm{x}}^0}\\ {{{\varepsilon }}_{\rm{y}}^0}\\ {{{\gamma }}_{{\rm{xy}}}^0} \end{array}} \right\} = \left\{ {\begin{array}{*{20}{c}} {\dfrac{{\partial {{{u}}^0}}}{{\partial {{x}}}}}\\ {\dfrac{{\partial {{{v}}^0}}}{{\partial {{y}}}}}\\ {\dfrac{{\partial {{{u}}^0}}}{{\partial {{x}}}} + \dfrac{{\partial {{{v}}^0}}}{{\partial {{y}}}}} \end{array}} \right\} $

(4)

$ \left\{ {{K}} \right\} = \left\{ {\begin{array}{*{20}{c}} {{{{K}}_{\rm{x}}}}\\ {{{{K}}_{\rm{y}}}}\\ {{{{K}}_{{\rm{xy}}}}} \end{array}} \right\} = \left\{ {\begin{array}{*{20}{c}} { - \dfrac{{{\partial ^2}{{{w}}^0}}}{{\partial {{{x}}^2}}}}\\ { - \dfrac{{{\partial ^2}{{{w}}^0}}}{{\partial {{{y}}^2}}}}\\ { - 2\dfrac{{{\partial ^2}{{{w}}^0}}}{{\partial {{x}}\partial {{y}}}}} \end{array}} \right\} $

(5)

式(3)中，{ε⁰}是中面应变向量、{K}是中面曲率变化向量；u⁰、v⁰、w⁰为中面沿x、y、z方向的位移。将式(3)代入式(1)得

$ \left\{ {{\sigma }} \right\} = {\left[ {{Q}} \right]_k}\left[ {\left\{ {{{{\varepsilon }}^0}} \right\} + Z\left\{ {{K}} \right\} - {{\left\{ {\rm{\alpha }} \right\}}_k}\Delta T} \right]\left({k = 1,2, \cdots ,n} \right) $

(6)

多层板截面上单位宽度的内力{N}和内力矩{M}如式(7)所示；将式(6)代入式(7)的两式，做积分得到{N}、{M}与应变的关系式(8)：

$ \left. {\begin{array}{*{20}{c}} {\left\{ {{N}} \right\} = \displaystyle\int \nolimits_{ - h/2}^{h/2} \left\{ {{\sigma }} \right\}{\rm{d}}{{\textit{z}}} = \displaystyle\sum \nolimits_{k = 1}^n \left\{ {{\sigma }} \right\}\left({{Z_k} - {Z_{k - 1}}} \right)}\\ {\left\{ {{M}} \right\} = \displaystyle\int \nolimits_{ - h/2}^{h/2} \left\{ {{\sigma }} \right\}{{\textit{z}}}{\rm{d}}{{\textit{z}}} = \frac{1}{2} \displaystyle\sum \nolimits_{k = 1}^n \left\{ {{\sigma }} \right\}\left({Z_k^2 - Z_{k - 1}^2} \right)} \end{array}} \right\} $

(7)

$ \left. {\begin{array}{*{20}{c}} {\left\{ {{N}} \right\} = \left[ {{A}} \right]\left\{ {{{{\varepsilon }}^0}} \right\} + \left[ {{B}} \right]\left\{ {{K}} \right\} - \left\{ {{{{N}}^t}} \right\}}\\ {\left\{ {{M}} \right\} = \left[ {{B}} \right]\left\{ {{{{\varepsilon }}^0}} \right\} + \left[ {{D}} \right]\left\{ {{K}} \right\} - \left\{ {{{{M}}^t}} \right\}} \end{array}} \right\} $

(8)

其中，

$ \left[ {{A}} \right] = \displaystyle\int \nolimits_{ - h/2}^{h/2} {\left[ {{Q}} \right]_k}{\rm{d}}{{\textit{z}}} = \mathop \sum \nolimits_{k = 1}^n {\left[ {{Q}} \right]_k}\left({{Z_k} - {Z_{k - 1}}} \right) $

(9)

$ \left[ {{B}} \right] = \displaystyle\int \nolimits_{ - h/2}^{h/2} {\left[ {{Q}} \right]_k}{{\textit{z}}}{\rm{d}}{{\textit{z}}} = \frac{1}{2}\mathop \sum \nolimits_{k = 1}^n {\left[ {{Q}} \right]_k}\left({Z_k^2 - Z_{k - 1}^2} \right) $

(10)

$ \left[ {{D}} \right] = \displaystyle\int \nolimits_{ - h/2}^{h/2} {\left[ {{Q}} \right]_k}{{{\textit{z}}}^2}{\rm{d}}{{\textit{z}}} = \frac{1}{3}\mathop \sum \nolimits_{k = 1}^n {\left[ {{Q}} \right]_k}\left({Z_k^3 - Z_{k - 1}^3} \right) $

(11)

$ \left[ {{{{N}}^t}} \right] = \mathop \sum \nolimits_{k = 1}^n {\left[ {{Q}} \right]_k}{\left\{ {\rm{\alpha }} \right\}_k}\Delta T\left({{Z_k} - {Z_{k - 1}}} \right) $

(12)

$ \left[ {{{{M}}^t}} \right] = \frac{1}{2}\mathop \sum \nolimits_{k = 1}^n {\left[ {{Q}} \right]_k}{\left\{ {\rm{\alpha }} \right\}_k}\Delta T\left({Z_k^2 - Z_{k - 1}^2} \right) $

(13)

式(7)中，h为层合板的总厚度；式(8)中，[A]是拉伸刚度矩阵，[B]是耦合刚度矩阵，[D]是弯曲刚度矩阵，[N^t]和[M^t]为温度变化ΔT时引起的热内力和热内力矩。将式(8)改写为

$ \left\{ {\begin{array}{*{20}{c}} {{N}}\\ {{M}} \end{array}} \right\} = \left[ {\begin{array}{*{20}{c}} {{A}}&{{B}}\\ {{B}}&{{D}} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {{{{\varepsilon }}^0}}\\ {{K}} \end{array}} \right\} - \left\{ {\begin{array}{*{20}{c}} {{{{N}}^t}}\\ {{{{M}}^t}} \end{array}} \right\} $

(14)

式(14)为温变条件下，层合板的内力与应变之间的关系式。利用矩阵运算，对式(14)进行求逆得出多层板的{ε⁰}、{K}与内力的关系式(15)：

$ \left\{ {\begin{array}{*{20}{c}} {{{{\varepsilon }}^0}}\\ {{M}} \end{array}} \right\} = \left[ {\begin{array}{*{20}{c}} {{{{A}}^ {\bullet } }}&{{{{B}}^ {\bullet } }}\\ {{{{C}}^ {\bullet } }}&{{{{D}}^ {\bullet } }} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {{N}}\\ {{K}} \end{array}} \right\} + \left[ {\begin{array}{*{20}{c}} {{{A}}^ {\bullet }}&0\\ {{{C}}^ {\bullet }}&{ - 1} \end{array}} \right]\left\{ {\begin{array}{*{20}{c}} {{{{N}}^{{t}}}}\\ {{{{M}}^{{t}}}} \end{array}} \right\} $

(15)

其中，

$ \left\{\begin{array}{l} {{{A}}^ {\bullet } } = {{{A}}^{ - 1}}\\ {{{B}}^ {\bullet } } = - {{{A}}^{ - 1}}{{B}}\\ {{{C}}^ {\bullet } } = {{B}}{{{A}}^{ - 1}} = -{\left[{{{{B}}^ {\bullet } }} \right]^{\rm{T}}}\\ {{{D}}^ {\bullet } } = - {{B}}{{{A}}^{ - 1}}{{B}} + {{D}} \end{array}\right. $

对式(15)进一步作逆转，得到

$ \left\{ {\begin{array}{*{20}{c}} {{{{\varepsilon }}^0}}\\ {{K}} \end{array}} \right\} = \left[ {\begin{array}{*{20}{c}} {{{{A}}^{\rm{*}}}}&{{{{B}}^{\rm{*}}}}\\ {{{{C}}^{\rm{*}}}}&{{{{D}}^{\rm{*}}}} \end{array}} \right]\left({\left\{ {\begin{array}{*{20}{c}} {{N}}\\ {{M}} \end{array}} \right\} + \left\{ {\begin{array}{*{20}{c}} {{{{N}}^{{t}}}}\\ {{{{M}}^{{t}}}} \end{array}} \right\}} \right) $

(16)

其中，

$ \left\{\begin{array}{l} {{{A}}^{{*}}} = {{{A}}^ {\bullet } } - {{{B}}^ {\bullet } }{{{D}}^{{ {\bullet - 1}}}}{{{C}}^ {\bullet } }\\ {{{B}}^{\rm{*}}} = {{{B}}^ {\bullet } }{{{D}}^{{ {\bullet - 1}}}}\\ {{{C}}^{\rm{*}}} = - {{{D}}^{{ {\bullet - 1}}}}{{{C}}^ {\bullet } } = {{{D}}^{{ {\bullet - 1}}}}{\left({{{{B}}^ {\bullet } }} \right)^{\rm{T}}} = {\left({{{{B}}^{\rm{*}}}} \right)^{\rm{T}}}\\ {{{D}}^{\rm{*}}} = {{{D}}^{{ {\bullet - 1}}}} \end{array}\right. $

由上述理论可以得出多层板的热弹性特性：当多层板的温度变化ΔT时，即使ΔT沿厚度h方向不变，也会产生曲率的变化而引起翘曲变形。

2 PBGA封装翘曲变形仿真

BGA封装是从引脚栅格阵列演变而来的一种高密度表面装配封装技术，而PBGA封装是其最常见一种封装形式。采用ANSYS软件对PBGA封装回流焊温变过程进行仿真模拟，可以得到其翘曲变形和热应力应变^[12]。本研究所用的PBGA样品如图3所示，内部结构如图4所示。PBGA样品由硅芯片、粘接层、塑封层和基板4个部分组成，基板为多层布线结构，粘接层材料为粘接胶。硅芯片通过粘接层和引线与基板实现连接。

2.1 PBGA封装仿真模型

针对PBGA样品的封装结构，本文对PBGA样品实体结构进行简化建模，实体模型包含硅芯片、芯片粘接层、塑封层和基板4个组成部分。PBGA样品的3D实体模型如图5所示。

PBGA样品所涉及的各种材料参数如表1所列，其中芯片、塑封层、粘接层采用线弹性材料，基板采用各向异性材料；且考虑材料的玻璃化转变温度(T_g)^[13]。

2.2 PBGA温变加载

PBGA封装翘曲变形属于热−结构耦合问题，模拟仿真分析用ANSYA软件中的热分析模块和结构分析模块。首先通过仿真模拟得到样品模型的热分析结果，再将热分析结果作为载荷加载到模型的结构物理场中，从而得出结构分析结果。采用手动网格划分对模型进行网格划分，其中粘结层高度部分相对细化；选用Solid70作为热分析单元类型，Solid185作为结构分析单元类型^[14]。初始温度和参考温度都设置为27 ℃，热对流系数为125 W/(m²·K)。选择热对流方式对模型外表面进行热加载，依据如图6所示的回流焊温变载荷时间历程，确定仿真模型的加热条件和温变过程，如表2所列。

PBGA封装翘曲变形仿真分析时依据芯片的对称性，将基板底面中心点附近9个节点进行全约束，相当于固定一个小平面，作为翘曲变形的参考平面，如图7所示；然后在此约束条件之下，分析整个封装体在模拟回流焊温变过程中的翘曲变形情况。初始时封装体基板下表面假设为绝对平面，对应实测中的参考面。

2.3 仿真结果与分析 2.3.1 翘曲变形仿真结果分析

模拟仿真的翘曲值为PBGA样品翘曲变形的等效位移，将基板底面中心点的坐标点作为翘曲值零点。t=600 s(回流焊峰值温度265 ℃)时，PBGA封装体的温度场、翘曲变形如图8所示。图8(a)结果表明：在265 ℃时PBGA封装的温度分布基本均匀，温度场的最大值与最小值仅相差0.4 ℃。图8(b)中：PBGA封装体呈现出四角向上翘起的翘曲变形，翘曲值从底面中心约束点向边角处逐渐增大；边角d处的翘曲值最大，为35.9 μm；封装体内芯片的最大翘曲值也在芯片的边角d处，其翘曲值为16.4 μm。PBGA封装回流焊过程中，边角d处翘曲值变化如图9所示：翘曲值随温度的变化发生相应的变化，在温度最高时达到最大值(35.9 μm)；随着温度的降低，恢复到初始值。

2.3.2 应力应变仿真结果分析

t=600 s(回流焊峰值温度265 ℃)时，PBGA封装体中心点X-Z剖面的等效应力和等效热应变结果如图10所示.

由图10的仿真结果得出：整个封装体中，等效热应力和等效热应变分布极不平衡，主要集中分布在粘接层及基板上半部分。等效热应力最大与最小值之比为87∶1；最大等效热应力在基板靠近粘接层的位置，其值为0.861×10⁸ Pa；最大等效热应变也在基板靠近粘接层的位置，其值为0.016。由于各材料之间的热膨胀系数差异较大，材料在温变过程中会相互制约，不能自由胀缩，从而各自产生不同的热应力应变，并造成应力应变分布极不平衡现象。

PBGA封装体在各个方向上的应力应变也存在着分布不平衡现象，以X方向的热应力为例，如图11所示：X方向热应力分布也极不平衡，主要最集中分布在粘接层及基板上半部分。粘接层上X方向热应力为负，大小在–0.848×10⁸~ –0.11×10 ⁸ Pa之间；由于粘接层的CTE最大，热膨胀时会受到其他材料的抑制作用。从纵向看基板靠近粘接层处，X方向热应力随着逐渐远离粘接层，从负应力转变为正应力。

2.3.3 不同温度点的仿真结果对比分析

回流焊温变过程中，随着温度的变化，PBGA封装的翘曲变形和热应力应变也会发生相应的变化。选取t=480 s(回流焊温度230 ℃)和t=760 s(回流焊温度217 ℃)时PBGA封装体的仿真结果，其中心点X-Z剖面的等效热应力和等效应变结果如图12和图13所示。

比较t=480 s(回流焊温度230 ℃)、t=600 s(回流焊峰值温度265 ℃)和t=760 s(回流焊温度217 ℃)3个时间点的仿真结果，可以看到：随着温度的变化，PBGA封装体的翘曲值和等效热应力应变也会发生相应的数值变化；但其翘曲变形最大位置和等效热应力应变集中分布的特点并没有变化。PBGA封装始终呈现出四角向上的翘曲变形，翘曲值从中心约束点到边角处逐渐增大，热应力应变集中区域始终在粘接层及基板接近粘接层的位置。

2.3.4 PBGA翘曲变形测试验证

为了验证PBGA封装回流焊翘曲变形的仿真结果，本研究采用阴影云纹(Shadow moire)法测量其翘曲值。阴影云纹法采用明暗相间光栅的相互重叠干涉生成moire条纹，并使用摄像头抓取图形并对图形进行分析得到样品的翘曲状况^[15-16]。阴影云纹法测试原理如图14所示，光源发出的白光以一定的角度照射在样品表面，光栅产生moire干涉条纹，再将干涉图样用计算机转换成3D形貌图，得出样品的变形图或翘曲量。

(1) 翘曲变形测试条件：PBGA封装翘曲变形试验选取3个相同型号的样品。测试过程模拟回流焊的温变过程并采用热风对流方式对样品进行加热，加热条件和回流焊温变过程如表2和图6所示。加热过程中的温变速率是变化的，尽量模拟实际的回流焊温变历程，每个时间点对应着相应的温度；整个温变过程是从室温(27 ℃)加热到回流焊的峰值温度(265 ℃)，最后再空冷到室温(27 ℃)。

(2) 测试结果：试验样品在不同温度−时间点的翘曲值(单位：μm)如表3所列：当回流焊温度达到峰值温度265 ℃时，3个样品的翘曲值也达到最大值，分别为36，32，47 μm。温变过程中，PBGA封装的翘曲值随着温度的升高逐渐增大，当温度降低时又会逐减小并恢复到初始状态。

3 仿真与试验结果对比分析

模拟回流焊温变过程中，PBGA封装体边角d处翘曲值的模拟值与实测值如图15所示。

由图15得出：PBGA封装体边角d处的翘曲值随回流焊温度的变化发生相应的变化，在温度最高时达到或接近最大值；模拟值与样品1的实测值在峰值温度265 ℃时最接近，分别为35.9，36 μm。回流焊温变过程中，样品翘曲值的模拟值与实测值的变化趋势具有一致性，表明模型的建立和仿真过程是合理的。根据样品的实测数据，初始时3个样品都有一定翘曲值，这与产品的生产制造因素有关；而模拟是一种理想的状态^[17]，无法考虑样品的初始翘曲状态。仿真结果与实测数据对比分析表明：回流焊温变过程中，PBGA封装体边角处的翘曲值往往最大，在回流焊峰值温度265 ℃时模拟翘曲值达到35.9 μm，样品的实测值最大达到47 μm。因此边角处的BGA焊点也最容易出现开路、虚焊等装联失效问题。

4 结论

本研究采用ANSYS有限元软件对PBGA封装在模拟回流焊温变过程中的翘曲变形行为及应力应变进行仿真分析，并用Shadow moire法进行了翘曲测试验证，通过分析可得出以下结论：

(1) 样品翘曲值的模拟值与实测值较为接近，整个温变过程中翘曲值的实测值与模拟值的变化趋势具有一致性，说明建立的模型和仿真过程是合理的。

(2) 回流焊过程中，PBGA封装体边角处的翘曲值最大，因此边角处的焊点最容易出现开路、虚焊等装联失效问题。

(3) PBGA回流焊过程中，热应力应变主要集中分布在粘接层及基板上半部分，最大值也都在基板靠近粘接层的位置，该位置是PBGA封装体中出现热可靠性失效问题的最大风险点。