当局部不再平衡：氢脆相场建模中热力学一致性的缺失与恢复

一、引言

氢被广泛视为未来脱碳能源体系的重要载体，然而溶解氢与微观组织特征的交互作用会导致材料性能的灾难性退化并引发突发断裂，即氢脆现象。氢脆涉及多种机制的耦合，目前已被广泛接受的机理包括氢增强脱聚（HEDE）、氢增强局部塑性（HELP）与氢增强应变诱导空位（HESIV）等。由于氢脆具有多物理场、多尺度特征，仅靠原子尺度或宏观尺度建模均难以完整预测，在原子与宏观尺度之间建立一个显式包含关键机制的介观模型成为现实途径。

相场方法是在介观尺度模拟复杂系统演化的参考性方法，已成功应用于凝固、相变、电沉积与裂纹扩展等问题。近年来，研究者通过脆性断裂力学与氢扩散耦合的方式，发展了一批面向氢辅助开裂的相场模型。然而这些模型存在两方面不足：其一，它们并未显式再现氢偏聚这一导致裂纹表面能下降、进而触发HEDE机制的核心物理过程；其二，它们均依赖Langmuir-McLean等温吸附关系由局部氢浓度计算裂纹表面覆盖度，而该等温线本质上是均匀体相与表面相之间的宏观平衡关系，不能作为非均匀微结构（尤其是存在弹性交互时）中的局部定律加以合理使用。

G. F. Bouobda-Moladje、 A. Ruffini、 Y. Le Bouar和 A. Finel在“Materials Science"期刊发表了文章“Thermodynamically consistent phase field model for hydrogen-assisted cracking"，该文基于Kim-Kim-Suzuki（KKS）多相场格式，提出了一种用于模拟多晶材料氢辅助开裂的热力学一致性相场模型。模型在统一变分框架内同时描述裂纹扩展、氢在裂纹表面与晶界上的偏聚，以及偏聚引起的界面能下降，从而在不引入Langmuir-McLean局部定律的前提下再现HEDE机制。模型首先在严格热力学框架下推导，随后解析计算了氢偏聚引起的裂纹有效表面能退化，最后通过多晶镍算例展示了模型捕捉穿晶断裂向氢辅助沿晶断裂转变的能力。

二、相场模型理论框架

2.1 状态变量与总自由能

模型考虑含预制裂纹的多晶金属材料，晶界假设保持静态。氢原子在固溶体中占据间隙位置，这些位置分布于体相、裂纹表面、完整晶界与断裂晶界四类微观结构特征中。系统状态由三组场变量描述：表征材料完整性的损伤场η（完好时为1，断裂时为0）、间隙位置的氢占据概率c，以及一组分别代表不同取向晶粒的序参数（ϕ1,…,ϕq）。与连续相场建模惯例一致，所有界面均以弥散区域表示。

系统总自由能由四部分构成：氢自由条件下的裂纹表面断裂能（含各向同性表面能γ0、内聚函数w(η)与正则化长度ξ）、线弹性框架下的弹性能（含刚度退化函数g(η)与氢致特征应变ε0(c)）、裂纹与晶界的耦合项（通过参数ζij控制晶间断裂能Γij）以及化学自由能密度fch(c,η)。其中化学自由能的构造是本文的核心。

2.2 KKS化学自由能与氢偏聚

化学自由能采用最初发展于凝固领域的KKS格式建立。该格式的突出优势是将界面性质与化学组元浓度场解耦：定义体相间隙占据率cb与裂纹表面相占据率cck两套浓度变量，化学自由能按插值函数h(η)在两相能量之间分配，总氢浓度c为两相浓度的加权之和。由于两套浓度对应同一种化学组元，二者之间通过化学势相等约束关联，实际只有一个自由度。

体相与裂纹相的化学自由能均在理想固溶体假设下写出，包含溶液能与构型熵两项。相别通过GB相关的氢溶液能E引入：体相溶液能Eb与裂纹相溶液能Eck分别在纯体相与纯裂纹相溶液能的基础上叠加序参数梯度平方的修正项，修正系数αij与βij分别控制完整晶界与断裂晶界处的氢偏聚能。由此，单晶中裂纹表面的氢偏聚能定义为∆Eck=Eck0−Eb0，数值为负时氢自发向裂纹表面富集。

2.3 控制方程与数值方法

系统演化假设为准静态过程，控制方程包括力学平衡、氢场平衡与损伤场平衡三组方程。应力由退化刚度与总应变（含氢致特征应变）给出。数值上，氢场与损伤场的平衡方程采用Fourier谱方法积分，力学平衡则采用一种快速鲁棒的离散FFT求解器获得。模型采用KKL断裂模型的内聚与退化函数，插值函数p(η)取五阶多项式以保证其一阶与二阶导数在η=0和η=1处消失，从而避免晶界区域的非物理软化。

三、模型解析性质

3.1 Langmuir-McLean等温线的自然涌现

针对含氢条件下的一维平面裂纹构型（弹性松弛），模型可解析导出裂纹表面覆盖度与体相浓度之间的经典Langmuir-McLean关系式。值得强调的是，该等温线在模型中并非作为局部定律被先验假设，而是在体相与裂纹相均为均匀相这一特定条件下作为平衡解自然涌现。这正说明：Langmuir-McLean关系只适用于宏观均匀相间的热力学平衡，将其作为非均匀微结构中的普适局部定律缺乏物理依据，而本文模型则从根本上回避了这一假设。

3.2 氢致表面能变化的解析表达

氢偏聚引起的表面能变化∆γch可以在张开的平面裂纹构型下解析计算：利用平衡态化学势的均匀性以及损伤剖面方程，最终得到∆γch=Chw·ξ·(fch_ck,eq−fch_b,eq)，其中Chw为仅依赖于插值函数的积分常数（本文取值0.5749）。在低氢浓度极限下，该式可进一步简化为关于偏聚能∆Eck的指数表达式。由于体相化学能密度与偏聚能均为负值，∆γch小于零，即氢偏聚必然降低裂纹表面能。该关系还确定了正则化长度ξ与给定表面能—氢覆盖度曲线之间的联系，使模型参数可由原子尺度数据标定。

四、数值验证与应用

4.1 单晶I型裂纹扩展验证

首先在二维单晶中验证模型：均匀外加应力下测量I型裂纹扩展的临界应力σc，材料为各向同性模型材料（E=80 GPa，ν=1/3），氢性能参数取自铝的相关数据（偏聚能∆Eck=−0.6 eV）。无HEDE机制的相场结果与Griffith理论阈值良好吻合，其中裂纹半长因网格离散而带有量级为网格间距的不确定度∆a。

图1  单晶中I型裂纹扩展临界应力σc随初始裂纹半长a的变化：无HEDE与含HEDE相场结果均遵循相应的Griffith阈值

开启HEDE机制后，临界应力明显降低（图1），其物理来源是氢致表面能变化∆γch参与裂纹扩展的能量平衡。含氢条件下相场结果仍遵循Griffith准则的形式，但有效表面能降低约4 J/m²，与解析公式预测的∆γch=−4.29 J/m²高度一致，验证了模型定量再现HEDE效应的能力。

4.2 多晶中穿晶—沿晶断裂转变

随后将模型应用于多晶材料的氢辅助裂纹扩展。模拟采用镍的弹性常数（E=209 GPa，ν=0.30）与原子尺度计算的氢偏聚能参数（∆Eck=−0.5 eV、∆EGB=−0.25 eV、∆Eck|GB=−0.9 eV），温度700 K，晶界/体相断裂能比初始设为ΓGB/G0c=0.5。氢致沿晶开裂通常被解释为氢引起的晶界韧性损失，原子模拟已证实氢偏聚可使镍晶界强度在慢断裂极限下降低约50%；而线弹性断裂力学进一步表明，裂纹在弱区域的偏转取决于断裂能比与偏转角。

图2  多晶系统中的裂纹扩展：(a)初始裂纹构型；(b)无HEDE机制的穿晶扩展；(c)含HEDE机制的沿晶扩展；(d)归一化氢浓度场

无氢参考工况下，裂纹仅以穿晶模式扩展（高韧性晶界）。引入氢后，裂纹路径由穿晶转变为沿晶模式（图2b、c）。相场测量给出的有效晶界/体相断裂能比由0.5骤降至0.07：由于断裂晶界GB型表面的氢偏聚能（−0.9 eV）高于晶内裂纹体相型表面（−0.5 eV），晶界断裂能的下降幅度远大于体相，导致断裂阻力比对氢含量急剧敏感，从而表现为沿晶开裂。氢浓度场（图2d）清晰呈现出晶界、晶内裂纹表面与断裂晶界表面三处偏聚能的差异所对应的浓度对比。

五、模型参数标定

5.1 插值函数与晶界断裂能标定

插值函数p(η)的选择对晶界行为至关重要。若取p(η)=1−4η³+3η⁴，晶界区域会出现非物理软化，导致裂纹沿晶界的非物理路径扩展；而取五阶形式p(η)=1−10η³+15η⁴−6η⁵时，因其一阶与二阶导数在两端点处消失，软化现象消失（图3）。

图3  不同插值函数p(η)下的裂纹扩展路径：(a)初始状态；(b)(c)低阶p(η)引起的晶界非物理软化；(d)(e)五阶p(η)消除软化

晶界断裂能通过参数ζ标定。对平面晶界的解析推导表明，ΓGB/G0c与ζ/(γ0d)呈线性关系，斜率由依赖于p(η)、损伤剖面与晶粒剖面的积分给出，连续场估计值为0.16。考虑空间离散后的一维数值模拟给出的斜率约为0.20（图4），本文即采用该数值估计，取ζ/(γ0d)=−2.5以实现ΓGB/G0c=0.5。

图4  晶界/体相断裂能比ΓGB/G0c随参数ζ的线性演化

5.2 氢偏聚能标定

完整晶界与断裂晶界处的氢偏聚能分别由参数α与β控制。解析推导表明，∆EGB与∆Eck|GB均与α、β呈线性关系。空间离散模拟给出的线性斜率比连续解低约20%（图5），本文采用数值估计，取α/(kBTd²)=−30、β/(kBTd²)=−45，分别对应∆EGB=−0.25 eV与∆Eck|GB=−0.9 eV。

图5  完整晶界与断裂晶界中氢偏聚能随参数α、β的线性变化

氢致表面能变化与表面氢浓度均依赖于偏聚能的取值。图6给出了表面能变化随有效偏聚能与氢覆盖度的演化关系，以及断裂晶界处的相应结果，二者共同验证了模型在参数—物理量映射上的自洽性，也为实际应用中的参数选择提供了依据。

图6  有效表面能随氢偏聚能（等价于氢覆盖度）的依赖关系：(a)晶内裂纹表面；(b)断裂晶界

六、结论与展望

该文在KKS格式下建立了耦合线弹性断裂力学与HEDE机制的新相场模型。模型在变分一致的热力学框架内再现了裂纹表面的氢偏聚，无需借助仅适用于宏观均匀相平衡的Langmuir-McLean等温线；该等温线被证明仅在静态无应力的一维极限下因体相与裂纹相均匀平衡而依然成立。氢含量对表面能的降低效应因此得以严格而普适地再现，使模型特别适合在微观组织尺度上描述HEDE机制。

模型在镍中的应用利用真实原子尺度输入数据，成功捕捉了HEDE机制诱导的穿晶向沿晶开裂模式转变。本文忽略了可能影响裂纹扩展的表面能各向异性，后续可通过在能量泛函中引入更高阶项予以考虑。未来的工作将致力于纳入其他氢脆机制（如HELP）并研究多机制对裂纹扩展路径的协同效应；就HELP而言，第一步可借助位错相场模型在系统中引入离散位错，进而发展位错、裂纹与氢场之间一致的热弹性与动力学耦合。