3. 热结构耦合控制方程与弱/强耦合策略

在3D IC封装中，热与结构从来不是两个独立的问题。硅通孔（TSV）、微凸点、再布线层、底部填充胶和基板在温度变化下各自膨胀或收缩，而它们的约束关系又反过来改变热传导路径与界面接触状态。要可靠预测翘曲、应力与寿命，必须把温度场和位移场放在同一套方程体系里描述。本节从控制方程出发，说明耦合项从哪里来，再讨论顺序耦合与直接耦合的适用场景和计算代价。

一、热传导方程：温度场的控制方程

封装体内的瞬态温度场由傅里叶热传导方程描述。取温度 T(x,y,z,t)，各向同性材料的热传导方程为：

ρ·c_p·∂T/∂t = ∇·(k·∇T) + Q
其中：ρ 密度，c_p 比热容，k 导热系数，Q 体积热源（如TSV焦耳热）

稳态分析时左端为零，方程退化为 ∇·(k∇T) + Q = 0。边界条件通常有三类：第一类给定温度（如散热器底面恒温）；第二类给定热流密度（如芯片功耗面）；第三类对流换热 −k·∂T/∂n = h·(T − T∞)，h 为对流换热系数，T∞ 为环境温度。在3D IC中，TSV铜柱、微凸点、底部填充胶和硅基体的 k 值相差可达两个数量级（铜约400 W/(m·K)，硅约150 W/(m·K)，底部填充胶约0.3～1 W/(m·K)），因此材料界面的热流分配对温度梯度影响极大。

二、热弹性力学方程：位移场的控制方程

结构场由热弹性力学方程描述。以位移 u 为基本未知量，平衡方程为：

∇·σ + f = 0
σ = C : (ε − ε_th)
ε = ½(∇u + ∇uᵀ)
ε_th = α·(T − T_ref)·I

其中 σ 为应力张量，C 为弹性刚度张量，ε 为总应变，ε_th 为热应变，α 为热膨胀系数，T_ref 为无应力参考温度，I 为单位张量。这里的关键是热应变项 ε_th：温度变化直接产生应变，而结构约束（如基板对芯片的约束、TSV对周围硅的约束）阻止自由膨胀，从而产生热应力。对于各向异性材料（如某些基板或封装材料），α 和 C 需用张量形式表达。

三、耦合项来源：双向作用

热与结构的耦合不是单向的，而是双向的，耦合项来自两个方向：

• 热→结构：温度场通过热应变 ε_th = α(T − T_ref) 进入结构方程，这是最主要的耦合路径。温度升高使材料膨胀，受到约束时产生应力。
• 结构→热：变形改变几何构型，进而改变热边界条件与热传导路径。典型表现包括：界面接触热阻随接触压力变化；翘曲导致散热器与芯片之间间隙变化，改变对流或接触换热；裂纹或分层改变热流通道。此外，塑性耗散和摩擦生热也可作为热源反馈到温度场。

在多数封装热应力问题中，结构→热的影响相对较弱，但翘曲改变接触热阻这一项在散热路径中不可忽略，尤其当封装翘曲量达到数十微米时，界面间隙的热阻可能成倍增加。

四、顺序耦合与直接耦合的适用场景

顺序耦合（弱耦合）先求解热传导方程得到温度场，再把温度作为体载荷施加到结构模型上求解位移和应力。它假设结构变形对温度场的影响可以忽略，因此不需要迭代。这种策略的优点是两个物理场可以分别用各自成熟的求解器，网格可以不同，计算代价低，适合热→结构单向主导的问题，例如稳态功耗下芯片与基板的整体翘曲和焊点应力评估。在3D IC中，如果封装翘曲量相对散热路径尺寸很小（例如小于几微米），顺序耦合通常足够。

直接耦合（强耦合）把温度自由度和位移自由度放在同一方程组中同时求解，每一步迭代都更新耦合项。它能捕捉结构变形对热边界条件的反馈，例如接触热阻随接触压力变化、翘曲导致散热器与芯片之间间隙改变。代价是矩阵规模显著增大，非线性迭代次数增加，内存和计算时间可能比顺序耦合高一个数量级。适合翘曲显著、接触状态变化剧烈或存在热-力强反馈的场景，如大尺寸封装、散热器压接、以及需要精确预测界面热阻的场合。

五、教学示例：一个简化的顺序耦合计算

考虑一个简化模型：硅芯片尺寸 10 mm × 10 mm × 0.1 mm，底部有 100 μm 厚底部填充胶，再下方是 1 mm 厚基板。芯片功耗 2 W，环境温度 25 ℃，对流换热系数 h = 20 W/(m²·K) 作用于基板底面。假设材料参数为：硅 k=150 W/(m·K)，α=2.6×10⁻⁶ /K；底填胶 k=0.5 W/(m·K)，α=30×10⁻⁶ /K；基板 k=0.3 W/(m·K)，α=17×10⁻⁶ /K。参考温度 T_ref=25 ℃。

第一步热分析：稳态下，芯片产生的 2 W 热量主要通过底填胶和基板向下传导，再经对流散出。粗略估算热阻：底填胶层 R≈t/(kA)=100×10⁻⁶/(0.5×10⁻⁴)=2 K/W；基板 R≈1×10⁻³/(0.3×10⁻⁴)≈33 K/W；对流 R=1/(hA)=1/(20×10⁻⁴)=500 K/W。总热阻约 535 K/W，温升 ΔT≈2×535≈1070 K，显然对流是瓶颈，实际需要散热器或更高 h。若 h 提高到 5000 W/(m²·K)，对流热阻降为 2 K/W，总热阻约 37 K/W，温升约 74 K，芯片温度约 99 ℃。这个示例说明热边界条件对温度场的主导作用。

第二步结构分析：把温度场作为体载荷。由于底填胶 α 远大于硅和基板，且被上下层约束，底填胶层会产生较大剪切应变。用简化一维约束估算：若底填胶自由膨胀量 ΔL=α·ΔT·L=30×10⁻⁶×74×10×10⁻³≈2.2×10⁻⁵ m，而实际被硅和基板约束，应变约 2.2×10⁻³，对应应力约 E·ε（E 取 5 GPa）≈11 MPa。这个量级足以在界面处引发分层风险。实际仿真中需用三维有限元求解，并检查底填胶与硅、基板界面的剪应力和剥离应力。

六、选择策略与判断依据

实际工程中可按以下依据选择：若翘曲量远小于散热路径特征尺寸，且接触热阻变化可忽略，优先顺序耦合；若翘曲量与界面间隙同量级，或接触压力显著改变热阻，应使用直接耦合。也可以先用顺序耦合快速评估，再对关键区域做直接耦合细化。计算代价方面，顺序耦合通常可在几分钟到几十分钟完成，直接耦合可能需要数小时甚至更久，需根据项目阶段和精度要求权衡。

无论采用哪种策略，都要保证参考温度、材料参数和边界条件的一致性。热应变以 T_ref 为基准，若 T_ref 设置不当，会引入虚假应力。对于多层结构，各层 α 和 k 的差异是耦合效应的根源，也是3D IC热结构仿真中最需要关注的材料输入。