力-化耦合断裂相场模型研究进展

申强 ,  张世一 ,  马立峰

科学技术与工程 ›› 2026, Vol. 26 ›› Issue (3) :879-900

PDF (19082KB)
科学技术与工程 ›› 2026, Vol. 26 ›› Issue (3) :879-900 DOI: 10.12404/j.issn.1671-1815.2407797
综述:力学

力-化耦合断裂相场模型研究进展

作者信息 +

Review of Coupling Mechano-chemical Fracture Phase-field Model

Author information +
文章历史 +

摘要

化学场所诱导的失效是常见材料失效模式之一。由于其复杂的多物理场耦合失效机理,传统的有限元方法在分析这种复杂的力-化耦合失效时具有一定的局限性。近年来,相场断裂模型在分析裂纹萌生、扩展和分支方面具有独特的优势。相场断裂模型基于断裂力学的变分原理,该模型可以有效地处理多场耦合断裂问题。详细介绍了力-化耦合断裂相场模型的理论框架、控制方程及有限元实现,并简述了该模型在锂离子电池颗粒充放电过程中和金属材料氢脆中的应用。最后,阐述了该模型在工程应用中存在的主要问题,并探讨了未来研究方向。

Abstract

Chemical field-induced failure is one of the common failure modes of materials. Due to its complex multi-physics coupling failure mechanism, the traditional finite element method has certain limitations in analyzing this complex mechanical-chemical coupling failure. In recent years, the phase-field fracture model has unique advantages in analyzing crack initiation, propagation and branching. The phase-field fracture model is based on the variational principle of fracture mechanics, which can effectively deal with multi-field coupling fracture problems. The theoretical framework, control equations and finite element implementation of the mechanical-chemical coupling fracture phase field model was introduced in detail, and the application of the model in the charging and discharging process of lithium-ion battery particles and the hydrogen embrittlement of metal materials was briefly described. Finally, the main problems of the model in engineering application were expounded, and the future research direction was discussed.

Graphical abstract

关键词

断裂 / 损伤 / 相场模型 / 锂离子电池 / 氢脆

Key words

fracture / damage / phase-field model / lithium-ion batteries / hydrogen embrittlement

引用本文

引用格式 ▾
申强, 张世一, 马立峰. 力-化耦合断裂相场模型研究进展[J]. 科学技术与工程, 2026, 26(3): 879-900 DOI:10.12404/j.issn.1671-1815.2407797
SHEN Qiang, ZHANG Shi-yi, MA Li-feng. Review of Coupling Mechano-chemical Fracture Phase-field Model[J]. Science Technology and Engineering, 2026, 26(3): 879-900 DOI:10.12404/j.issn.1671-1815.2407797

登录浏览全文

4963

注册一个新账户 忘记密码

断裂是工程材料和结构中常见的破坏形式之一,其中由于化学浓度变化引起的断裂尤为常见,防止开裂引起的破坏是工程设计中的一个主要问题[1-2]。处于深海环境中的金属结构往往受到高浓度氢的影响,氢原子通过吸附、渗透等方式进入金属材料中。这种氢扩散会导致局部应力集中加剧和微裂纹的形成,会大幅度降低材料的延展性、断裂韧性,使得材料在远远低于断裂韧性的条件下发生破坏[3-5]。锂离子电池作为一种重要的能量存储设备,在电动汽车、便携电子设备领域发挥着关键的作用。然而,锂离子电池在长期循环使用过程中会面临各种问题,其中之一是电极材料中锂离子的脱出与嵌入过程会引起浓度的变化,会对电池颗粒的表面造成损伤,严重影响电池的性能和使用寿命[6-9]。然而,设计实验来研究这些问题不仅昂贵,而且很难考虑到所有可能影响裂纹扩展的因素。因此,考虑由化学场所诱导断裂的数值计算模型来预测材料的断裂过程是十分必要的。
数值计算模型能够通过数学建模和计算仿真方法,在计算机中模拟和预测材料的断裂行为,将物理过程和力学现象转化为数学表达式,并结合数值算法进行数值求解。数值计算模型的优势在于可以对多物理场作用下的复杂断裂问题进行全面、精确的模拟和分析。此外,数值计算模型可以有效节约成本和时间,探索和评估多种断裂模式和方案,为实际工程应用提供指导和优化选择。目前,比较常见的处理力-化耦合断裂问题的数值计算模型主要分为两大类:离散型裂纹断裂模型和弥散型裂纹断裂模型。
离散型裂纹断裂模型应用较为广泛的主要有扩展有限元法[10-14]和内聚力模型[15-16]等。虽然,扩展有限元法能在一定程度上模拟力-化耦合机制下的断裂问题,但模拟的结果很大程度上依赖于网格单元的大小和质量,选取不恰当的网格单元会导致模拟的结果与实际相差很大。对于一些复杂的裂纹形态和非线性载荷加载情况下具有一定的局限性。内聚力模型在模拟氢脆等问题时有一定的优势,但是内聚力模型通常基于一些简化假设。例如,假设金属晶格是均匀连续的,忽略了晶格缺陷、晶界等复杂的微观结构因素。这种简化可能导致对氢脆现象的描述不够准确。近年来,随着对力-化耦合断裂过程的不断研究,学者们开始意识到传统的离散型裂纹断裂模型在这方面存在不足。为弥补这一问题,越来越多的研究者开始转向发展弥散型裂纹断裂模型。这种模型考虑了裂纹周围的局部损伤和变形情况,能更有效地模拟裂纹扩展过程中的力学行为与断裂过程之间的相互作用,从而提高了对材料断裂行为的理解和预测能力。应用较为广泛的为近场动力学[17-18]和相场模型[19-21]。袁训锋等[19]基于相场模型,采用有限差分法研究了初始晶核条件、扰动强度和时间步长对多晶硅凝固界面形态演化的影响。王辉明等[20]基于提出的统一相场理论和内聚力模型,建立由钢纤维、骨料、砂浆及界面过渡区构成的钢纤维混凝土(steel fiber reinforced concrete, SFRC)细观有限元模型并结合试验予以验证。结果表明:基于相场损伤模型和内聚力模型建立的SFRC细观有限元模型具有较好的准确性和可靠性。李明耀等[21]基于相场法研究了岩石中裂纹的起裂、扩展和贯通问题,克服了传统断裂力学理论及数值计算方法显式追踪裂纹面的挑战,能很好地模拟裂纹起裂、扩展、分叉和多裂纹相互贯通。力-化耦合断裂相场模型在模拟氢脆现象[22-26]、电池充放电的疲劳损伤[27-29]等力-化耦合断裂领域也展现出了非常广泛的应用前景。
断裂相场模型是以位移场和裂纹相场为自变量的耦合断裂模型,该模型基于变分原理[30-31]可以较为容易地将化学势能添加到总势能函数中,从而进一步考虑关于浓度运输、离子扩散等化学行为对裂纹扩展的影响。在利用相场法求解力-化耦合断裂问题时,相场变量的时间演化均为耦合的偏微分方程,可根据热力学定律对每个相场变量进行推导求解。力-化耦合断裂相场模型的本质是关于位移场、化学势场和相场的耦合偏微分方程组求解能量最小化的过程。因此力-化耦合断裂相场模型被广泛应用求解化学势场诱导裂纹开裂的问题中。
在复杂的结构和加载条件下,传统的数值方法无法准确地捕捉材料在化学浓度变化过程中发生的失效行为。力-化耦合断裂相场模型是一种基于相场理论的损伤模型,该模型是一种通过将相场变量与物理场耦合的方法来描述材料断裂现象的理论框架,该模型的变分性质使得它可以很好地处理断裂与化学势场的耦合问题。现详细介绍力-化耦合断裂相场模型的理论基础,公式推导、化学扩散、力学变形和裂缝演化等耦合过程,并简要总结力-化耦合断裂相场模型在锂离子电池颗粒疲劳充放电与金属材料氢脆中的应用,最后对力-化耦合断裂相场模型的未来研究方向和工程应用中需要克服的主要问题进行讨论。

1 脆性断裂相场模型

1.1 断裂裂纹的相场近似

相场断裂模型是在Griffith断裂理论基础上基于Bourdin等[30]传统脆性断裂相场模型扩展而来的。该模型在处理断裂问题上不再把裂纹表面离散化,而是引入一个标量场,即相场[图1(a)],来近似地描述裂纹,将尖锐的裂纹表面弥散为正则化弥散裂纹,如图1(b)所示。

相场数值d(x)为1和0分别代表材料完全破坏和完好两种极限状态,而它们之间的数值代表一种损伤状态,并且裂纹的弥散程度由相场裂纹宽度l0来控制,数值越大表示弥散宽度越大,反之则越小。Amor等[31]将系统的总势能定义为由于变形引起的变形能与断裂能之和,将结构的位移场与裂纹面作为自变量。但在数值计算中,将离散的裂纹面作为未知量求解是十分困难的。因此,Bourdin等[30]和Miehe等[32]在此基础上对裂纹表面进行处理,通过引入一个与相场相关的裂纹面密度函数对系统的断裂能进行重构。

在一维情况下,不连续尖锐裂纹可以用一个强间断的标量场函数d(x)∈[0,1]表示为

d(x)= 1 ,   x = 0 0 , x ≠ 0

离散裂纹模型虽然能够描述裂纹的几何位置,但是当裂纹扩展、分支演化过程中,往往面临裂纹扩展路径的跟踪和网格重新划分等问题,这种强间断的标量场会在有限元数值计算中带来诸多不便。断裂相场模型很好地解决了这两个棘手的问题,通过引入一个连续的场函数将离散的裂纹面转换为分布在整个模型上的连续裂纹场,如图1(a)所示。Bourdin等[30]和Miehe等[32]将连续裂纹的分布情况用一个指数函数的形式表示为

d(x)=exp - x l 0

式(2)中:l0为裂纹宽度,表示裂纹弥散的程度。l0越大,裂纹的弥散程度越大;l0越小,裂纹的弥散程度就越小。

Bourdin等[30]和Miehe等[32]对连续裂纹场的分布情况表示为指数形式外,Pham等[33]在AT1模型中,连续裂纹场的分布情况表示为分段函数,Wu[34-35]的脆性断裂与准脆性断裂相场WN模型中,连续裂纹场的分布情况表示为正弦函数。

1.2 能量泛函与控制方程

基于Griffith断裂理论和Francfort-Marigo变分原理[31]可得,系统的总势能由4部分组成,表达式为

Ψtotal=Ψelas+Ψfrac+Wc-Wext

式(3)中:Ψelas 、 Ψfrac 分别为系统的弹性应变能、产生新裂纹面所耗散的断裂能;Wc和Wext分别表示由于浓度变化所引起的化学功和外力所做的功。

在考虑体力存在的情况下,系统外力的势能可以表示为

Wext=∫ΩbudV+∫∂ΩhudA

式(4)中:b为体力;u为位移场;h为面力;Ω为单元体;dV体微分;dA为面微分。系统的弹性应变能可表示为

Ψelas=∫Ωψ(ε)dV

式(5)中:ψ(ε)为弹性应变能密度。

系统的断裂能为临界能量释放率与裂纹表面的乘积,可表示为

$ \begin{aligned}\Psi_{\text {frac }} & =G_{\mathrm{c}} \Gamma_{1} \approx \int_{\Omega} G_{\mathrm{c}} \gamma(d, \nabla d) \mathrm{d} V \\& =\int_{\Omega} G_{\mathrm{c}}\left(\frac{1}{2 l_{0}} d^{2}+\frac{l_{0}}{2}|\nabla d|^{2}\right) \mathrm{d} V\end{aligned}$

式(6)中:Gc为材料的临界能量释放率;Γl为裂纹面;$ \gamma(d, \nabla d)=\frac{1}{2 l_{0}} d^{2}+\frac{l_{0}}{2}|\nabla d|^{2}$为裂纹表面密度函数。除上述裂纹表面密度函数的表达,Wu[34-35]也提出了一种新的用于描述裂纹分布的裂纹表面密度函数,即

$ \gamma(d, \nabla d)=\frac{1}{c_{0}}\left[\frac{1}{l_{0}} \alpha(d)+l_{0}|\nabla d|^{2}\right]$

式(7)中:c0为归一化参数;α(d)∈[0,1]为反映局部耗散能的几何函数。Wu[34-37]将几何函数定义为

α(d)=2d-d2

其他学者对于几何函数的研究如表1所示。

在线弹性材料中,系统的应变势能就是弹性应变能。由于相场模型是连续损伤模型,相场数值不为0的区域,会存在不同程度的损伤。损伤使得结构储存的弹性能会耗散掉一部分,因此需要引入一个退化函数w(d)来反映由于损伤的存在而导致弹性能的耗散。退化函数的选择与材料的损伤本构有关,尽管不同相场模型的退化函数有所差异,但通常都满足如下4个条件:①w(d)=0表示相场数值为0的区域,材料没有损伤,弹性应变能没有变化;②w(d)=1表示相场数值为1的区域,材料完全损伤,不具备承载能力;③w(1)=0表示材料完全破坏后相场数值不再继续演化;④ w ·(d)<0表示退化函数应该是一个严格递减的函数。

如表2所示为近些年来不同学者对能量退化函数的贡献。表2中第一项Bourdin等[30]提出二次多项式,主要应用于脆性材料中准静态的裂纹演化,至今还在沿用;第二项为Karma等[39]提出的三次多项式,将准静态断裂扩展到动态断裂;第三项为Borden等[40] 提出了在裂纹萌生前具有弹性域的附加参数s三次多项式,该函数提供了裂纹萌生之前的应力-应变响应,更接近于线性弹性行为,应用于韧性材料的裂纹演化;表2中后三项提出的能量退化函数主要应用于内聚断裂和准脆性断裂,也可适用于脆性断裂模型[34-37,41-44]。

将能量退化函数考虑到弹性应变能中,考虑损伤后的真实状态下的弹性应变能为

Ψelas(u,d)=∫Ωw(d)ψ(ε)dV=∫Ωψ[ε(u),d]dV

式(9)中:ε为应变张量。考虑各向同性固体的线弹性理论,弹性应变能密度可表示为

ψ(ε)= 1 2εT:C0:ε

式(10)中:C0为线弹性刚度矩阵;:表示双点积运算。在小变形理论中,应变张量ε可以用位移梯度表示为

$ \boldsymbol{\varepsilon}=\frac{1}{2}\left[\nabla \boldsymbol{u}+(\nabla \boldsymbol{u})^{\mathrm{T}}\right]$

式(11)中:u和$ \nabla$u分别为位移向量和位移梯度。将考虑损伤后的各部分势能代入式(3)中,系统的总势能函数为

$ \begin{aligned}\Psi_{\text {total }}(\boldsymbol{u}, d)= & \int_{\Omega} \psi[\boldsymbol{\varepsilon}(\boldsymbol{u}), d] \mathrm{d} V+ \\& \int_{\Omega} G_{\mathrm{c}} \gamma(d, \nabla d) \mathrm{d} V- \\& \int_{\Omega} \boldsymbol{b} \boldsymbol{u} \mathrm{d} V-\int_{\partial \Omega} \boldsymbol{h} \boldsymbol{u} \mathrm{d} A\end{aligned}$

在断裂相场模型中,系统的总势能变为以位移场u和相场d两个连续场为自变函数的泛函。当总势能达到某一临界值时,裂纹将开始扩展。根据Francfort-Marigo变分原理,位移场u和相场d使得系统的总势能最小,即

(u,d)=Arg{ i n f u , d ∈ B C sΨtotal(u,d)}

式(13)中:BCs为裂纹弥散边界。对总势能函数进行位移场u和相场d的变分,得

$ \begin{array}{l}\delta \Psi_{\text {total }}(\boldsymbol{u}, d)=\int_{\Omega} w(d) \boldsymbol{\sigma}: \delta \boldsymbol{\varepsilon} \mathrm{d} V+ \\\quad \int_{\partial \Omega} \frac{\partial \psi[\boldsymbol{u}(\boldsymbol{\varepsilon}), d]}{\partial d} \delta d \mathrm{~d} V+\int_{\Omega} G_{\mathrm{c}}\left[\frac{\partial \gamma(d, \nabla d)}{\partial d} \delta d+\right. \\\left.\quad \frac{\partial \gamma(d, \nabla d)}{\partial \nabla d} \delta \nabla d\right] \mathrm{d} V-\int_{\Omega} \boldsymbol{b} \delta \boldsymbol{u} \mathrm{d} V-\int_{\partial \Omega} \boldsymbol{h} \delta \boldsymbol{u} \mathrm{d} A\end{array}$

式(14):σ为应力张量。

对式(14)应用散度定理可得

$ \begin{array}{c}\delta \Psi_{\text {total }}(\boldsymbol{u}, d)=-\int_{\Omega}[w(d) \nabla \cdot \boldsymbol{\sigma}+\boldsymbol{b}] \delta \boldsymbol{u} \mathrm{d} V+ \\\int_{\partial \Omega}(\boldsymbol{\sigma} \boldsymbol{n}-\boldsymbol{h}) \delta \boldsymbol{u} \mathrm{d} A+\int_{B}\left[\frac{\partial \psi[\boldsymbol{u}(\boldsymbol{\varepsilon}), d]}{\partial d}+\right. \\\left.G_{\mathrm{c}} \delta d \gamma(d, \nabla d)\right] \delta d \mathrm{~d} V+ \\\int_{\partial B} G_{\mathrm{c}}\left[\frac{\partial \gamma(d, \nabla d)}{\partial \nabla d} \boldsymbol{n}_{B}\right] \mathrm{d} A\end{array}$

式(15)中:n为系统边界上的单位法向量;nΒ为弥散裂纹带B边界∂Β的单位法向量。耦合场(u,d)应满足d∈[0,1]及可逆边界条件 d ·≥0,相场边界条件$ \frac{\partial \gamma(d, \nabla d)}{\partial \nabla d} \boldsymbol{n}_{B}=\mathbf{0}$。δd和δu取任意值时式(15)都满足,可以得到断裂相场模型的控制方程

$ \left\{\begin{array}{ll}w(d) \nabla \cdot \boldsymbol{\sigma}+\boldsymbol{b}=0, & \text { in } \Omega \\G_{\mathrm{c}}\left(\frac{1}{l_{0}} d-l_{0} \Delta d\right)-\dot{w}(d) \psi[\boldsymbol{\varepsilon}(\boldsymbol{u}), d], & \text { in } \Omega \\w(d) \boldsymbol{\sigma} \boldsymbol{n}-\boldsymbol{h}=0, & \text { in } \partial \Omega\end{array}\right.$

裂纹相场演化过程是不可逆的,为了确保材料卸载或者压缩时出现相场数值减少,裂纹愈合等情况发生,应在断裂相场的控制方程中施加不可逆条件来防止裂纹的愈合。Miehe等[32]提出了历史应变能,表达式为

H= m a x τ ∈ [ 0 , t ]{ψ[ε(x,τ)]}

式(17)中:H表示在[0,t]时间段内;点x处的弹性应变能密度函数取最大值;τ为时间变量;t为当前时刻。相场的控制方程则被修改为

Gc 1 l 0 d - l 0 Δ d- w ·(d)H

对各向同性线弹性材料[45],弹性应变能密度函数可用主应变分量表示为

ψ(ε)= 1 2εT:C0:ε = 1 2λ0tr2(ε)+μ0( ε 1 2+ ε 2 2+ ε 3 2)

式(19)中:λ0和μ0为Lame'常数;tr(ε)为应变张量的迹;ε1、ε2、ε3为主应变。

为了防止各向同性线弹性材料在压缩状态下的失效,将弹性能量密度分解为拉伸部分ψ+(ε)和压缩部分ψ-(ε)来保证裂纹仅由拉伸部分ψ+(ε)驱动。这里介绍一下Miehe等[32]基于应变张量普分解的应变能分解方式,即

ψ±(ε)= 1 2λ0 < ε 1 + ε 2 + ε 3 > ± 2±( < ε 1 > ± 2+ < ε 2 > ± 2+ < ε 3 > ± 2)

式(20)中:<·>±为Macaulay符号。分解的应变能中,假设只有受拉伸部分应变能驱动相场演化,而受压缩部分并不会引起相场演化,这种分解方式能够有效避免裂纹面闭合时材料的互相侵入问题。

在各向异性模型中,只考虑拉应力引起的损伤,并且仅考虑弹性能的拉伸分量通过能量退化函数退化,弹性应变能密度ψ(ε,d)分解为

ψ(ε,d)=w(d)ψ+(ε)+ψ-(ε)

此模型最早由Amor等[31]提出。此模型中,拉伸和压缩分量可表示为

ψ+(ε)= 1 2K < t r ( ε ) > + 2+μεdev:εdev
ψ-(ε)= 1 2K < t r ( ε ) > - 2

式(23)中:K=λ+2μ/3,为体积模量,其中λ为材料的拉梅常数;μ为剪切模量;εdev为应变的偏张量;<x>±被定义为(x± x)/2。

Miehe能量分解和Amor能量分解的相场模型,虽然应变能定义的数学表达式相同,但应变能的分解选择方式是不同的。

除了上述能量分解的相场模型外,还有许多方式对拉伸和压缩的应变能进行分解。对于剪切断裂,Lancioni等[46]假设,仅弹性能量的偏分量退化。基于Amor等[31]开发的模型,Zhang等[47]将谱分解应用于偏应变,并开发了改进的体积-偏分裂模型。Wu等[48]将柯西应力张量分解为拉伸和压缩分量,并将开发的模型应用于脆性和准脆性断裂。在Ambati等[45]混合相场模型(各向同性-各向异性),仅对相场演化驱动力中的应变能进行分解,以保证只有与损伤有关的能量才会驱动相场演化。

2 力-化耦合断裂相场模型

前面重点介绍断裂相场理论的基础理论,主要探讨相场模型的发展过程以及相关的概述和演化方程。断裂相场模型是一种通过将相场变量与物理场耦合的方法来描述脆性断裂现象的理论框架。该模型的变分性质使得它可以很好地处理断裂与其他物理现象的耦合问题,尤其在一些涉及力-化耦合断裂问题中。

在一些处于深海腐蚀的环境中,金属材料常常面临应力腐蚀和氢脆断裂机制所导致的损伤问题。深海环境中的金属结构常常承受着高压力和高浓度离子等腐蚀介质的侵蚀,引发应力集中和晶粒腐蚀等问题,导致金属材料的应力腐蚀损伤,使材料强度、韧性下降,严重会导致金属结构的疲劳断裂和失效[3-5]。与此同时,氢脆机制也是深海金属材料损伤的另一重要原因,氢脆的断裂机制如图2所示。氢原子通过吸附、渗透、离子输运等方式进入金属的晶格结构中,当氢原子在晶格结构中移动、扩散、汇聚时,会加剧应力集中和晶格之间的相互作用,会大幅度降低材料的强度、断裂韧性,使材料在远远低于断裂韧性的情况下发生断裂。Anand等[49]对氢气含量对金属的延展性进行了研究,结果如图3所示。可以发现,碳钢在空气中的延展性要远远高于暴露于氢气中,与在空气中碳钢延展性相比,下降幅度为60%。

2.1 质量运输与反应做功

由于浓度变化所引起的化学功Wc,实际上是一个变形与扩散的耦合问题。在材料的晶格内,氢的传输受到浓度梯度和静水应力的共同驱动。根据质量守恒原理,要求单位时间内通过单位面积的扩散通量与该面积处的浓度梯度成正比,可表示为

∫Ω ∂ C ∂ tdV+∫∂ΩJndA=0

式(24)中:C为浓度;J为离子通量。

对质量守恒方程应用散度定理,可以得出此平衡方程的强形式为

$ \frac{\partial C}{\partial t}+\nabla \cdot \boldsymbol{J}=0$

利用式(25)进行关于浓度C的变分处理,并对处理后的方程应用散度定理,可以得到此平衡方程的弱形式为

$ \int_{\Omega}\left[\delta C\left(\frac{\mathrm{~d} C}{\mathrm{~d} t}\right)-\boldsymbol{J} \nabla \delta C\right] \mathrm{d} V+\int_{\partial \Omega_{q}} \delta C q \mathrm{~d} A=0$

式(26)中:q=Jn代表通过扩散离开系统表面∂Ωq的浓度通量。扩散由系统化学势的梯度$\nabla$μ驱动,氢在晶格中的化学势可表示为

μ=μ0+RTln θ L 1 - θ L-VHσH

式(27)中:μ为氢在晶格中的化学势;R为理想气体常数;T为绝对温度;θL为晶格位置占有率,与浓度C和晶格位置的数量N相关,可表示为θL=C/N,通常情况下采用低占有率(θL≪1);μ0为不受浓度变化下的标准化学势;σH为静水应力;VH为关于氢的偏摩尔体积。离子通量J遵循与化学势梯度 μ的线性Onsager关系,可表示为

$ \boldsymbol{J}=-\frac{D C}{R T} \nabla \mu$

式(28)中:D为扩散系数;C为浓度。将式(28)代入式(26),扩散方程进一步表示为

$ \begin{array}{c}\int_{\Omega}\left[\delta C\left(\frac{1}{D} \frac{\mathrm{~d} C}{\mathrm{~d} t}\right)+\nabla \delta C \nabla C-\nabla \delta C\left(\frac{V_{\mathrm{H}} C}{R T} \nabla \sigma_{\mathrm{H}}\right)\right] \mathrm{d} V+ \\\int_{\partial \Omega_{q}} \frac{1}{D} \delta C q \mathrm{~d} A=0\end{array}$

式(29)中: ∂Ωq为通量边界条件;q为边界上的氢通量。

氢通过降低金属原子之间的键能使金属材料变脆,从而降低了断裂的阻力[50-52]。Alvaro等[50]于第一性原理计算给出了氢覆盖θ的二次退化函数形式

φ(θ)=1-1.046 7θ+0.168 7θ2

但针对此问题,更多学者采用密度泛函理论[50,53-54]来研究具有不同氢覆盖率对断裂能的影响。如图4所示,Alvaro等[50]计算了镍金属在Σ3和Σ5晶界处存在氢原子时断裂能的变化情况,用归一化表面能表示与氢气覆盖率的关系,并与数据进行了线性拟合。Jiang等[54]还研究了铁和铝金属存在氢原子时断裂能的变化情况。

Martínez-pañeda等[55]认为材料的临界能量释放率同样受到氢覆盖率的影响,且满足关系式

G c ( θ ) G c 0=1-χθ

式(31)中:Gc(0)为没有氢气存在情况下的临界能量释放率;Gc(θ)为不同氢覆盖率下材料的临界能量释放率;χ为不同材料的损伤系数,可以通过拟合DFT数据来估计其他材料的损伤系数;θ为氢气覆盖率。

在均相相场问题中,Martínez-pañeda等[55]研究了不同氢覆盖率下铁的本构关系,如图5所示。

从图5可以看出,对于铁基类材料,氢覆盖率θ会显著降低材料强度。为了建立氢浓度扩散与断裂相场之间的关系,需要根据氢的浓度来计算氢的表面覆盖率。根据Serebrinsky[52]的研究,这种关系可由Langmuir-McLean等温线[52,56]表示为

θ= C C + e x p - Δ G b 0 R T

式(32)中:Δ${G}_{b}^{0}$为吉布斯自由能。

2.2 变形-扩散-断裂耦合问题

虑了浓度扩散的相场断裂问题是一个三场耦合问题。如图6所示,位移场u,外部边界分为∂Ωu和∂Ωh两部分,分别为位移u引导的Dirichlet型边界条件,拉力h所引导的Neumann型边界条件。断裂相场由位移场驱动,内部弥散裂纹表面Γl满足Dirichlet型边界条件。对于质量运输,外部边界分为∂ΩC和∂ΩJ两部分,分别为浓度C引导的Dirichlet型边界条件和通量J引导的Neumann型边界条件。

在考虑浓度C变化和体力b存在情况下,根据虚功原理可得

$ \begin{aligned}W & =\int_{\Omega}(\boldsymbol{\sigma}: \delta \boldsymbol{\varepsilon}+\omega \delta d+\boldsymbol{\zeta} \nabla \delta d-\dot{C} \delta \mu+\boldsymbol{J} \nabla \delta \mu) \mathrm{d} V \\& =\int_{\partial \Omega}(\boldsymbol{h} \delta \boldsymbol{u}+q \delta \mu+f \delta \mathrm{~d}) \mathrm{d} A+\int_{\Omega} \boldsymbol{b} \delta \boldsymbol{u} \mathrm{d} V\end{aligned}$

式(33)中: ζ和ω分别为相场微应力矢量和相场微应力标量,其分别与相场梯度$\nabla$d和相场d共轭;f为边界损伤通量。

将散度定理应用于虚功原理方程,可得

$ \begin{aligned}\delta W= & \int_{\partial \Omega}(\boldsymbol{\sigma} \boldsymbol{n} \delta \boldsymbol{u}+\boldsymbol{\zeta} \boldsymbol{n} \delta d+\boldsymbol{J} \boldsymbol{n} \delta \mu) \mathrm{d} A- \\& \int_{\Omega}[(\nabla \cdot \boldsymbol{\sigma}+\boldsymbol{b}) \delta \boldsymbol{u}+(\nabla \cdot \boldsymbol{\zeta}-\omega) \delta d+ \\& (\nabla \cdot \boldsymbol{J}+\dot{C}) \delta \mu] \mathrm{d} V\end{aligned}$

式(34) 中:要求对任意变量(δu,δd,δμ)都成立,可以得到两组平衡方程

$ \left\{\begin{array}{ll}\nabla \cdot \boldsymbol{\sigma}+\boldsymbol{b}=0, & \text { in } \Omega \\\nabla \cdot \boldsymbol{\zeta}-\omega=0, & \text { in } \Omega \\\frac{\mathrm{d} C}{\mathrm{~d} t}+\nabla \cdot \boldsymbol{J}=0, & \text { in } \Omega\end{array}\right.$
h = σ n   , i n   ∂ Ω h f = ζ n , i n   ∂ Ω f q = J n   , i n   ∂ Ω q

2.3 热力学定律与本构理论

热力学第一定律和第二定律是热力学中两个基本的定律,用来解释连续体内的能量和总熵的动态转化规律。在发生热量交换的封闭系统内,只考虑机械能和热能,单位时间内内能的变化率等于系统与外界交换的热量和对外界所做功的总和。这种连续体的内能变化率可以表示为

d d t∫ΩUdV= W · e x t-∮cndA+∫ΩρzdV

式(37)中:U为比内能;“=”后两项为连续体内单位时间所获得的总热量;c为单位时间单位面积上的热通量(热流矢量);ρ为质量密度;z为单位时间单位质量获得的热辐射。

根据热力学第二定律,总熵S随时间的变化率等于或大于通过连续体表面的熵通量与内部产生的熵总和。在数学上,这个熵的原理可以表示为一个积分形式——Clausius-Duhem不等式(耗散函数),即

d d t∫ΩSdV≥-∮ c n tdA∫Ω ρ z tdV

在考虑存在质量运输的热力学定律中,要求任一部分体积增加的自由能随时间的变化率要小于或等于外界所做功的功率加上由于扩散物质通过其边界流入的能量通量,表达式为

d d t∫ΩΨtotaldV≤∫∂Ω W · e x tdA+∫∂ΩμJndA

将式(39)代入等式(37),应用散度定理可以得到耗散不等式的局部形式为

$ \begin{array}{c}\int_{\Omega}\left[\Psi_{\text {total }}-(\boldsymbol{\sigma}: \dot{\boldsymbol{\varepsilon}}+\omega \dot{d}+\boldsymbol{\xi} \nabla \dot{d})+\right. \\\mu \dot{C}-\boldsymbol{J} \cdot \nabla \mu] \mathrm{d} V \leqslant 0\end{array}$

式(40):ξ为扩散系数。

自由能函数Ψtotal是关于应变ε、相场d、相场梯度$\nabla$d和化学势μ的函数。耗散不等式的局部形式可进一步写为

$ \begin{array}{c}\left(\boldsymbol{\sigma}-\frac{\partial \Psi_{\text {total }}}{\partial \boldsymbol{\varepsilon}}\right): \dot{\boldsymbol{\varepsilon}}+\left(\omega-\frac{\partial \Psi_{\text {total }}}{\partial d}\right) \dot{d}+\left(\boldsymbol{\xi}-\frac{\partial \Psi_{\text {total }}}{\partial \nabla d}\right) \nabla \dot{d}- \\{\left[\left(\mu-\frac{\partial \Psi_{\text {total }}}{\partial C}\right) \dot{C}-\boldsymbol{J} \nabla \mu\right] \geqslant 0}\end{array}$

要求变量(δu,δd,δC)在任意条件下都成立,可得出关于浓度变化下的本构方程

$ \left\{\begin{array}{l}\boldsymbol{\sigma}=\frac{\partial \Psi_{\text {total }}}{\partial \boldsymbol{\varepsilon}} \\\omega=\frac{\partial \Psi_{\text {total }}}{\partial d} \\\boldsymbol{\xi}=\frac{\partial \Psi_{\text {total }}}{\partial \nabla d} \\\mu=\frac{\partial \Psi_{\text {total }}}{\partial C}\end{array}\right.$

在许多关于氢脆断裂[50-52]的力-化耦合断裂相场模型中,均采用了如下的自由能函数[55],即

$ \begin{aligned}\Psi_{\text {total }} & (\varepsilon, d, \nabla d, C)= \\& \underbrace{w(d) \psi_{0}(\boldsymbol{\varepsilon})-K V_{\mathrm{H}}\left(C-C_{0}\right) \operatorname{tr}(\boldsymbol{\varepsilon})}_{\Psi_{\text {elas }}}+ \\& \underbrace{G_{\mathrm{c}}(\theta)\left(\frac{1}{2 l_{0}} d^{2}+\frac{l_{0}}{2}|\nabla d|^{2}\right)}_{\Psi_{\text {frac }}} \mu_{0} C+ \\& \underbrace{R T N\left[\theta_{\mathrm{L}} \ln \theta_{\mathrm{L}}+\left(1-\theta_{\mathrm{L}}\right) \ln \left(1-\theta_{\mathrm{L}}\right)\right]}_{\Psi_{C}}\end{aligned}$

式(43)中:ΨC为浓度相关项;C0为参考晶格氢浓度;R为气体常数;T为绝对温度;N为粒子总数。

根据自由能函数Ψtotal(ε,d, $\nabla$d,C)和本构方程,柯西应力张量可表示为

σ= ∂ Ψ t o t a l ∂ ε=w(d)(C0:ε)-KVH(C-C0)I

微应力标量ω与相场d共轭,微应力矢量ζ与相场梯度$\nabla$d共轭

ω= ∂ Ψ t o t a l ∂ d= w ·(d)ψ0(ε)+Gc(θ) 1 l 0d
$ \zeta=\frac{\partial \Psi_{\text {total }}}{\partial \nabla d}=G_{\mathrm{c}}(\theta) l_{0} \nabla d$

化学势μ表示为

$ \begin{array}{c}\mu=\frac{\partial \Psi_{\text {total }}}{\partial C}=u_{0}+R T \ln \frac{\theta_{\mathrm{L}}}{1-\theta_{\mathrm{L}}}-V_{\mathrm{H}} \sigma_{\mathrm{H}}+ \\G_{\mathrm{c}}(\theta) \dot{\theta}\left(\frac{1}{2 l_{0}} d^{2}+\frac{l_{0}}{2}|\nabla d|^{2}\right)\end{array}$

通量J可进一步表示为

$ \boldsymbol{J}=-\frac{D C}{R T} \nabla \mu=-\frac{D C}{\theta_{\mathrm{L}}\left(1-\theta_{\mathrm{L}}\right)} \nabla \theta_{\mathrm{L}}+\frac{D C}{R T} V_{\mathrm{H}} \nabla \sigma_{\mathrm{H}}$

因为晶格占有率θL≪1,因此通量J可进一步表示为

$ \boldsymbol{J}=-D \nabla C+\frac{D C}{R T} V_{\mathrm{H}} \nabla \sigma_{\mathrm{H}}$

2.4 有限元实现

力-化耦合断裂相场模型作为一个典型的多物理场问题。在进行有限元空间离散时,可将所有区域内的单元节点赋予位移u、相场d和浓度C的自由度。在空间上,对于每个单元节点、位移场u、相场d和浓度C可以离散为

u= ∑ i = 1 m N'iui, d= ∑ i = 1 m N'idi, C= ∑ i = 1 m N'iCi

式(50)中:m为节点数量;ui、d和Ci分别表示单元在第i个节点的位移场、相场和浓度场;N'i为插值函数矩阵分量,为位移场u、相场d和浓度C的形函数。

基于此离散方式,相应梯度向量(ε, $\nabla$d, $\nabla$C)可以被离散为

$ \boldsymbol{\varepsilon}=\sum_{i=1}^{m} \boldsymbol{B}_{i}^{u} \boldsymbol{u}_{i}, \nabla d=\sum_{i=1}^{m} \boldsymbol{B}_{i} d_{i}, \nabla C=\sum_{i=1}^{m} \boldsymbol{B}_{i} C_{i}$

式(51)中:Bi为形函数的空间导数的向量;Bui为标准的应变-位移矩阵。变形问题中平衡方程的弱形式一般可写为

∫Ω[w(d)σ:δε-bδu]dV-∫∂Ω(hδu)dA=0

扩散方程的弱形式可表示为

$ \begin{aligned}\int_{\Omega} & {\left[\delta C\left(\frac{1}{D} \frac{\mathrm{~d} C}{\mathrm{~d} t}\right)+\nabla \delta C \nabla C-\nabla \delta C\left(\frac{\bar{V}_{\mathrm{H}} C}{R T} \nabla \sigma_{\mathrm{H}}\right)\right] \mathrm{d} V=} \\& -\frac{1}{D} \int_{\partial \Omega_{q}} \delta C q \mathrm{~d} A\end{aligned}$

考虑到变分δu和δC在任意条件下都成立,对于此平衡方程和扩散方程的离散形式可以表示为关于位移场u和浓度C的残差向量,即

rui=∫Ωw(d)$({B}_{i}^{u}{)}^{T}$σdV-∫ΩNTibdV-${\int }_{\partial {\Omega }_{h}}{{N}^{T}}_{i}$hdA
$ \begin{aligned}\boldsymbol{r}_{i}^{C}= & \int_{\partial \Omega_{q}}\left[\boldsymbol{N}_{i}^{\mathrm{T}}\left(\frac{1}{D} \frac{\mathrm{~d} C}{\mathrm{~d} t}\right)+\boldsymbol{B}_{i}^{\mathrm{T}} \nabla C-\boldsymbol{B}_{i}^{\mathrm{T}}\left(\frac{\bar{V}_{\mathrm{H}} C}{R T} \nabla \sigma_{\mathrm{H}}\right)\right] \mathrm{d} V+ \\& \frac{1}{D} \int_{\partial \Omega q} \boldsymbol{N}_{i}^{\mathrm{T}} q \mathrm{~d} A\end{aligned}$

式中:rui为位移场的残差向量;Ni为用于描述位移场分布为形函数矩阵(有限元方法中);${r}_{i}^{C}$为浓度的残差向量。

相场d的残差向量可表示为

$ \begin{aligned}\boldsymbol{r}_{i}^{d}= & \int_{\Omega}\left[\dot{w}(d) \boldsymbol{N}_{i} H+G_{\mathrm{c}}(\theta)\left(\frac{1}{l_{0}} \boldsymbol{N}_{i} d+l_{0} \boldsymbol{B}_{i}^{\mathrm{T}} \nabla d\right)\right] \mathrm{d} V- \\& \int_{\partial \Omega_{f}} \boldsymbol{N}_{i} f \mathrm{~d} A\end{aligned}$

式(56)中:${r}_{i}^{d}$为相场的残差向量。

在力-化耦合断裂相场模型中,三场(位移场、相场和化学势场)之间互相耦合,一起求解会大幅度增加计算难度。为了解决这个难题,通常需要将位移场,相场和化学势场进行解耦处理,Miehe等[32]提出了一种交错迭代格式进行求解,即在计算位移场,相场或者化学势场时,仅计算其中一个,将另外两个场作为定值。值得注意的是,在Miehe等[32]提出的交错迭代格式算法中,任意改变所选取被固定场的顺序,其计算结果并不会受到影响,这也就大大的提高了收敛的稳定性。

由于残差方程式(54)~式(56)为非线性代数方程组,为了获得此数值解,采用Newton-Raphson迭代法来求解非线性代数方程组,迭代格式为

K u 0 0 0 K d 0 0 0 K C u d C+ 0 0 0 0 0 0 0 0 M u · d · C ·= r u r d r C

式(57)中:Ku、Kd和KC为整体切线刚度矩阵;M为浓度矩阵,可表示为

Mij=∫Ω N T i 1 DNjdV
${K}_{ij}^{u}$=$\frac{\partial {r}_{i}^{u}}{\partial {u}_{j}}$=∫Ωw(d)$({B}_{i}^{u}{)}^{T}$C0${B}_{j}^{u}$dV
$ \begin{aligned}\boldsymbol{K}_{i j}^{d}= & \frac{\partial \boldsymbol{r}_{i}^{d}}{\partial d_{j}}=\int_{\Omega}\left\{\left[2 H+\frac{\boldsymbol{G}_{\mathrm{c}}(\theta)}{l_{0}}\right] N_{i} N_{j}+\right. \\& \left.G_{\mathrm{c}}(\theta) l_{0} \boldsymbol{B}_{i}^{\mathrm{T}} \boldsymbol{B}_{j}\right\} \mathrm{d} V\end{aligned}$
$ \boldsymbol{K}_{i j}^{C}=\frac{\partial \boldsymbol{r}_{i}^{C}}{\partial C_{j}}=\int_{\Omega}\left(\boldsymbol{B}_{i}^{\mathrm{T}} \boldsymbol{B}_{j}-\boldsymbol{B}_{i}^{\mathrm{T}} \frac{\bar{V}_{H}}{R T} \nabla \sigma_{H} N_{j}\right) \mathrm{d} V$

以上就是采用Newton-Raphson迭代法来求解关于力-化耦合断裂相场模型的基本流程。

3 力-化耦合断裂相场模型在氢脆中的应用

金属结构常常暴露在含氢的环境中,环境中的氢原子通过吸附、电镀、淬火或腐蚀等方式进入金属的晶格结构中,通常会使金属的延展性、断裂韧性和抗疲劳性方面显著降低,会在远低于无氢环境下可承受的负载水平下失效。当氢存在时,韧性金属以脆性方式断裂,裂纹通常沿晶界成核扩展,严重影响金属结构的使用寿命。

Dinachandra等[57]基于自适应网格的氢脆相场模型模拟了单边裂纹方板在氢环境下受拉伸下的裂纹扩展,研究了不同氢浓度环境下对裂纹扩展行为的影响。模拟结果如图7所示,结果显示,随着氢浓度的增加,裂纹扩展的速度显著加快,这表明氢浓度是影响裂纹扩展行为的关键因素。研究还给出了在不同氢浓度下载荷随时间变化的曲线,如图8所示。随着氢浓度的增加,峰值载荷会出现明显的降低,并且达到峰值载荷所需的位移也会减少。与之前的模拟结果相比,自适应网格下的氢脆相场模型所得到的曲线与实际情况吻合度较好,进一步验证了该模型的准确性和可靠性。

Si等[58]基于氢脆相场模型模拟了双缺口不对称板和X形板在不同氢浓度环境中的裂纹扩展过程,研究发现,裂纹尖端区域氢浓度在静水应力的影响下总是达到最大值。

Wu等[59]也针对单边缺口板在氢环境下受拉伸裂纹扩展进行了研究,提出了两种氢辅助断裂的相场正则化内聚力模型(phase-field regularized cohesive zone model,PF-CZM),分别为模型Ⅰ和模型Ⅱ。模型Ⅰ主要指材料的断裂能与拉伸强度之比Gf/ft不变,而模型Ⅱ主要是指Irwin’s内部长度不变。研究的重点是探究氢辅助断裂模型中长度尺度参数l0的敏感性。研究中,使用不同的相场尺度参数b对模型Ⅰ和模型Ⅱ进行了模拟。在浓度C=0.5 mg/kg的条件下,绘制了模型Ⅰ和模型Ⅱ的载荷-位移曲线,如图9所示。可以发现,不同情况下的相场尺度参数l0对于载荷-位移曲线的吻合程度都非常高,这意味着相场尺度参数对该模型的敏感性较低。在研究中,还预测了裂纹扩展过程中的损伤演化、浓度分布和静水应力的分布情况,结果显示静水应力的分布情况和氢浓度的扩散始终在缺口尖端附近达到最大值。

Mandal等[4]针对此单边缺口板在氢环境下受拉伸裂纹扩展进行了更加充分的研究,考虑了氢浓度的两种分布情况;一种是稳态情况下的均匀浓度C = 0.5 mg/kg(无扩散),另一种是瞬态情况下的浓度C=0.5 mg/kg(有扩散)。如图10所示为均匀浓度C=0.5 mg/kg情况下的裂纹扩展过程,与没有氢的情况类似,裂纹在缺口尖端处成核并且以脆性方式水平地扩展到右边缘。如图11所示为这两类氢浓度分布下不同相场模型中载荷随位移变化的曲线,可以发现,无扩散情况下的峰值载荷要高于有扩散情况下的峰值载荷。如图12所示为在4个时刻下两种情况下的裂纹扩展情况,可以看到有扩散情况下的裂纹扩展速率要快于无扩散情况下。这是因为在扩散作用下,氢能够通过材料更快地扩散到裂纹尖端区域,在裂纹尖端附近的浓度升高,形成局部浓度梯度,导致材料的脆性特性显著增强。这使得裂纹能够以更快的速率扩展,从而导致裂纹扩展速率在有扩散情况下明显增加。如图13所示为在扩散情况下不同时刻下的裂纹扩展、浓度分布和应力分布情况,可以看出,浓度分布和应力分布情况总是在裂纹缺口附近达到最大,这与许多研究结论[4,5,60-61]是相同的。

Huang等[62]基于相场理论框架建立了氢脆数值计算模型来研究氢在材料中引起的韧性到脆性的转变。该模型引入了氢增强局部塑性效应(hydrogen-enhanced localized plasticity, HELP)和氢增强脱黏效应(hydrogen enhanced decohesion, HEDE)来降低材料的屈服应力和临界释放能,通过模拟分析了在氢浓度存在下双缺口拉伸试样中的裂纹扩展,讨论了有无HELP和HEDE机制对裂纹扩展的影响。

如图14所示,该模型很好地展示了裂纹萌生前试样的裂纹扩展和氢浓度的分布情况。在氢存在的情况下,裂纹尖端区域的正静水应力有助于氢浓度的扩散,随着裂纹的不断扩展,更多的氢在静水应力的影响下向裂纹尖端前方的区域扩散,导致氢浓度最高点在裂纹尖端区域。值得注意的是,裂纹尖端前方的区域会发生一定的塑性变形,由于氢的总量是固定的,裂纹尖端区域的氢浓度不会随着裂纹扩展而持续升高,在达到一定峰值后会出现下降的趋势。

HELP和HEDE机制降低了材料的屈服应力和临界能量释放率。该研究将氢浓度与HELP和HEDE机制耦合,分别讨论了无HELP和HEDE机制、仅有HELP机制、仅有HEDE机制和HELP和HEDE机制同时存在对裂纹扩展行为的影响。在两类机制单独作用时,无论是HELP效应还是HEDE效应,都会促进裂纹的扩展,加速材料的破坏。然而,这两种机制在单独作用时产生的最大塑性应变并没有太大的差异。但是,在两种机制同时作用下,材料表现出最小的塑性变形和峰值载荷。这表明HELP和HEDE机制之间存在协同效应,即它们共同作用时能够加剧材料的破坏。通过结合这两种机制,可以更好地模拟金属材料的氢脆和预测氢的存在对金属材料所造成的由韧性到脆性的转变。

4 力-化耦合断裂相场模型在锂离子电池中的应用

力-化耦合断裂相场模型可以用于模拟各种材料的断裂行为,包括锂离子电池颗粒的疲劳行为。关于裂纹生长的锂离子电池材料颗粒的研究往往只考虑了材料颗粒的单一(脱)锂化过程,而忽略了材料颗粒由于反复充放电而引起的疲劳效应。实际上,电极颗粒在使用过程中,会不断经历应力变化和损伤累积,这些变化会导致颗粒内部的应力集中,从而促进颗粒开裂。因此,在研究电池颗粒开裂行为时,应该考虑损伤累积等其他因素对颗粒开裂行为的影响,从而更准确地预测电池容量的衰减。有不少研究[6,63-67]将疲劳效应应用到力-化耦合相场模型中,通过引入一个与历史变量相关的疲劳退化函数在加载期间来降低临界能量释放率Gc,该方法源于Alessi等[68]和Carrara等[69]的研究。通过引入疲劳退化函数对断裂能重新改写为

$ \Psi_{\text {frac }}=\int_{\Omega} \int_{0}^{t} f[\bar{\alpha}(\tau)] G_{\mathrm{c}} \gamma(d, \nabla d) \mathrm{d} \tau \mathrm{~d} V$

式(62)中: α -为历史变量;f[ α -(τ)]为疲劳退化函数。

在针对锂离子电池的研究[70]中,化学扩散通常采用的通量边界条件为

J= V SCmax c 3   600

式(63)中:Cmax为电极材料中锂的最大理论浓度;V和S分别为活性材料颗粒的体积和表面积;c为表征锂电池充放电快慢的参数;c越大,充放电频率越快。

以锂离子电池颗粒为例,基于力-化耦合断裂相场模型研究了球形电池颗粒锂离子脱出过程,考虑了不含有缺陷与含有缺陷的锂离子电池颗粒的断裂行为,材料属性如表3所示。

在已有的研究中,通常使用的假设是将活性材料电池颗粒建模为球形、椭圆形或片状体,并假设裂纹按照Ⅰ型裂纹(即法向张开型)的方式扩展。以球形截面的锂离子电池颗粒为研究对象,为了方便计算,仅仅考虑了球形锂离子电池颗粒的截面的裂纹扩展。因此,此问题可以简化为一个二维平面问题。

电池颗粒的几何形状和边界条件如图15(a)所示。由于此球形锂离子电池颗粒为中心轴对称模型,为了方便计算,在本模拟中,采用1/4模型即可。其中通量边界条件与整个球形电池颗粒的中心轴对称,如图15(b)所示。在锂离子脱出的过程中,由于电池颗粒表面的浓度低于锂离子电池颗粒的中心区域,会形成浓度差,红色区域表示拉伸区域,蓝色区域表示压缩区域,电池颗粒表面红色区域会发生收缩而产生环向的拉应力,电池颗粒内部蓝色区域会产生环向压应力,受此循环往复加载,当应变产生的应变能驱动裂纹扩展时,电池颗粒会发生损伤破坏,如图15(c)所示。

4.1 有缺陷的锂离子电池颗粒损伤演化

将针对有缺陷的锂离子电池颗粒的损伤演化过程进行数值模拟分析。首先,考虑了半径r=4 μm的无缺陷锂离子电极颗粒,设置充放电频率Crate=2,设置加载时间为2 500 s,分析其疲劳损伤演化过程。

如图16所示为该模拟条件下的最大主应力随时间的变化曲线。可以看到,该情况下的最大主应力仅为0.5 GPa,显著低于锂离子电池颗粒的破坏强度1.3 GPa。随着脱锂过程的结束,最大主应力仍然不能推动裂纹的扩展,此时可以认为锂离子电池颗粒处于安全状态。

有缺陷的电池颗粒会导致应力集中加剧,会过早的造成电池颗粒的损伤。如图17所示,考虑0.1 μm和0.5 μm的两种预置裂纹的锂离子电池颗粒,半径r=3 μm。在较低的充放电频率Crate=3的情况下,对锂离子电池的疲劳损伤演化进行模拟分析,其中最下面的红色区域代表初始裂纹长度的大小。可以看出,虽然有初始裂纹的存在,但裂纹并没有向前扩展。如图18所示为不同初始裂纹长度下,最大主应力随时间的变化曲线。可以看出,两种缺陷下的最大主应力均已超过了材料的破坏强度1.3 GPa,但这并不能驱动裂纹的进一步扩展。这是由于裂纹的扩展准则是基于能量准则,由于电池颗粒的损伤要同时受到半径r、充放电频率Crate和预制裂纹长度的影响,该条件下所造成的应变能仍然不能满足裂纹扩展的条件,所以此时电池颗粒仍是处于一个安全的状态。然而,与无预制裂纹的电池颗粒相比,带有预制裂纹的电池颗粒达到最大主应力所需时间明显减少。总之,考虑到预制裂纹的存在,与没有预制裂纹的电池颗粒相比,裂纹的扩展时间会大幅缩短。这意味着有预制裂纹的电池颗粒会更容易发生疲劳损伤。

下面针对有缺陷的电池颗粒进行进一步的研究,保持颗粒半径r=3 μm,预制裂纹长度为0.5 μm,继续增大电池颗粒的充放频率Crate=15,模拟该条件下的锂离子电池颗粒的疲劳损伤演化。如图19所示为该条件下的裂纹扩展的全过程。

为了更清晰地观察到裂纹尖端区域的相场演化和裂纹扩展过程,对裂纹扩展路径进行了更加详细的描述,如图20所示。由于增大了锂离子电池颗粒的充放电频率Crate,会大幅度提高锂离子脱出的速率,因而在锂离子电池颗粒表面与内部中心区域产生较大的浓度差,锂离子电池颗粒的拉伸区域会产生更大的环向拉应力,在此情况下所造成的应变能能够推进裂纹的进一步扩展。但是,由于锂离子电池颗粒的中心区域还有环向压应力的存在,随着脱锂过程的进行,当应变能不足以抵抗裂纹扩展所需的能量时,裂纹最终会停止扩展。在此条件下,最终裂纹向前扩展了0.12 μm左右。

4.2 无缺陷的锂离子电池颗粒损伤演化

前面的研究主要针对有初始缺陷对锂离子电池颗粒损伤的影响,下面将重点放在无初始缺陷对锂离子电池颗粒的影响。考虑无缺陷颗粒半径r=3 μm的锂离子电池颗粒,将充放电频率提升Crate= 30,加载时间设置为2 500 s,此情况下的电池颗粒损伤演化如图21所示。经历2 500 s的持续加载,锂离子电池颗粒表面并没有出现微裂纹,裂纹相场数值在达到0.183左右时就停止了增长。可以看出,虽然显著提升了充放电频率Crate,但由于较小的颗粒半径r和没有预置裂纹的存在,导致了裂纹并没有向内发生扩展。这也说明了预置裂纹的存在对锂离子电池颗粒的损伤会造成较大的影响。

4.3 锂离子电池参数对损伤演化的影响

锂离子电池颗粒的疲劳损伤不仅与是否存在初始缺陷和初始缺陷的长短有关。还与充放电频率、材料的临界能量释放率和电池颗粒的大小有关。频繁的充放电循环会增加电池颗粒所受到的应力和应变,导致裂纹的扩展速度增加,从而加速疲劳损伤的发生。电池颗粒材料的临界能量释放率是裂纹扩展每单位面积过程中释放的应变能。较低的临界能量释放率意味着裂纹扩展所需的能量较少,裂纹更容易扩展,从而加速疲劳损伤的进程。颗粒的大小会影响电池内部的应力分布和裂纹扩展路径。下面将针对这些因素来考虑对锂离子电池充放电造成的疲劳损伤。

4.3.1 充放电频率对锂离子电池颗粒的影响

3种不同的锂离子充放电频率,分别为Crate、2Crate和3Crate。4种不同的临界能量释放率Gc分别为0.5、1.0、1.5、2.0 N/m。无缺陷的锂离子电池颗粒的初始半径设置为r=20 μm。如图22所示为在相同临界能量释放率下,不同充放电频率,最大主应力随时间的变化曲线。从曲线中可以发现,随着充放电频率Crate的不断增大,最大主应力会呈现下降的趋势,这会导致锂离子电池颗粒表面出现更大的收缩变形,使电池颗粒更早进入破坏阶段,更容易裂纹的产生和扩展。

为了进一步研究充放电频率对锂离子电池颗粒疲劳损伤的影响,以临界能量释放率0.5 N/m为研究对象,同一时刻下不同充放电频率影响下的相场分布情况、静水应力分布和浓度分布情况如图23所示。从图23(a)可以明显看出,随着充放电频率的增加,锂离子电池颗粒的外表面产生更大的收缩变形,导致裂纹扩展更为容易,电池颗粒的使用寿命会大幅度下降。同时,图23(b)展示了同一时刻下的静水应力分布情况,随着充放电频率的增大,电池颗粒内部压缩区域和拉伸区域的应力集中程度更加明显,裂纹尖端附近形成的拉应力区域,可以极大地促进浓度的扩散。这意味着静水应力的存在加速了锂离子的浓度扩散,导致浓度的最大值总是出现在裂纹尖端区域,如图23(c)所示。

如图24所示为在相同临界能量释放不同充放电频率,锂离子电池颗粒表面损伤程度随时间的变化曲线。锂离子电池颗粒表面的损伤程度随着充放电频率的增大而增大,充放电频率越大,疲劳损伤程度越高。在相同临界能量释放率下,较高的充放电频率会缩减达到电池颗粒损伤程度的时间。这是因为随着充放电频率的增大,锂离子扩散的速度会显著加快,电池颗粒内部与外部形成的浓度梯度增大,引起的表面拉伸应力越大。

4.3.2 临界能量释放率对锂离子电池颗粒的影响

锂离子电池颗粒的不均匀膨胀和收缩引起的变形不仅受到充放电频率和预置裂纹长度的影响,而且还与材料本身的临界能量释放率有关。本节将研究临界能量释放率对锂离子电池颗粒疲劳损伤的影响。

如图25所示为相同充放电频率下不同临界能量释放率,最大主应力随时间的变化曲线。可以观察到随着临界能量释放率的不断增大,电池颗粒达到破坏时的最大主应力都呈递增趋势,这意味着随着临界能量释放率的提高,电池颗粒在充放电过程中承受的应力也相应增加,进而会导致电池颗粒不容发生破坏。但是,当充放电频率增大时,颗粒表面的通量也会增加,这会加快脱锂的速率。这样的结果会导致电池颗粒表面和内部中心区域之间形成更大的浓度梯度,进而增加颗粒的收缩变形程度,会造成电池颗粒达到破坏时的最大主应力降低。如图26所示为相同充放电频率,不同临界能量释放率下,电池颗粒表面的损伤程度随时间的变化曲线,这也呈现出了相同的趋势,降低临界能量释放率会大幅度降低电池的使用寿命。

如图27所示不同充放电频率,不同临界能量释放率的锂离子电池颗粒损伤图。结果表明,提升充放电频率和降低临界能量释放率都会对电池的使用寿命造成影响。随着锂离子的脱出,电池颗粒内部与表面形成更大的浓度梯度,导致颗粒表面的拉伸应力不断增加,而内部的中心区域形成环向的压缩应力,会在颗粒表面形成严重的径向裂纹。随着充放电频率的升高或者临界能量释放率的降低,原始裂纹进一步向内扩展,裂纹尖端延伸在某些区域产生了向前或分支裂缝。此外,随着裂纹向颗粒中心前进,尖端遇到更高的锂浓度,断裂韧性更高,而驱动裂纹生长的拉伸应力降低或变为压缩应力。因此,裂缝增长较慢,直到应变引起的应变能不足以推动裂纹扩展,裂纹停止生长。从图27所示数值模拟结果可以看出,电池颗粒的损伤断裂是径向和环向的组合裂纹形式。这与Tian等[64]所做的研究是一致的。

在锂离子电池颗粒的相场模型研究中,通常使用的假设是将活性材料颗粒建模为球形、椭圆形或片状体,并假定裂纹按照模式Ⅰ(即法向张开型)的方式扩展。裂纹的扩展主要是由力-化耦合效应驱动的,这是因为锂离子的扩散速率会受到应力场的影响,以及锂扩散或由具有不同锂浓度的相之间的晶格失配引起的拉伸环向应力。力-化耦合断裂相场模型能够更准确地描述裂纹扩展的行为,并考虑到力学、电化学和扩散效应之间的相互作用。

根据通量边界条件可知,活性材料颗粒的断裂行为主要受到充放电频率Crate和活性材料颗粒的尺寸和形状以及微缺陷的初始长度的影响。实际上,这些因素都是影响活性材料颗粒内的应力分布,进而影响断裂行为。Ai等[71]基于力化耦合疲劳断裂相场模型研究了圆柱形锂离子电池电极颗粒在不同充放频率Crate、初始裂纹长度a0和电池颗粒半径R0对裂纹扩展的影响,结果如图28所示。

如图28(a)所示为相同循环加载次数下不同充放电频率对电池颗粒表面的损伤程度。可以看出电池颗粒表面的损伤程度随着充放电频率的增大而增大,这是因为随着充放电频率的增大,边界通量增大,锂离子扩散的速度会显著加快,所造成的浓度梯度越高,引起的表面拉伸应力越高,且充放电频率会大幅度缩减电池颗粒损伤程度的加载次数。对于高频率充放电达到损伤最大仅仅只需要10个循环,而低频率下需要数百个循环才能达到最大的损伤程度,且损伤程度远远低于高频率加载下。在低频率0.5Crate和0.2Crate加载下,电池颗粒并没有引起损伤,这也意味着初始裂纹前方的断裂能将显著小于新生成的裂纹表面所需要的能量。吴建营[1]针对此问题做了针对性的研究,并且给出了最大主应力随着充放电时间的变化曲线,在低频率下,随着充放电频率的加载,直到锂离子浓度为0时,所造成的最大主应力还远远低于材料的强度极限,此时认为材料处于线弹性。

Cao等[9]也针对锂离子电池的充电过程在不同充放电频率下进行了模拟。充放电频率对断裂行为有两种不同的促进作用,在充放电频率较低时,会形成单裂纹,当增大频率时,裂纹会进一步扩大,裂纹尖端会向前延伸并出来分叉,当进一步增大,新的裂纹形成并扩展通至边界或与其他裂缝相连,最后会贯穿电池颗粒。Zhang等[63]对充电/放电期间沿着NMC电池颗粒的裂纹成核和萌生进行了建模。发现裂纹起始于颗粒中心,然后向颗粒表面生长并在两个方向上分支,在主裂纹的两个初始分支附近出现额外的小裂纹。这与Sun等[72]的实验研究结果是一致的。

如图28(b)所示为相同循环加载次数下,由于初始裂纹的存在对电池颗粒表面的损伤程度。可以看出,较小的初始裂纹在低频率下,对电池颗粒的损伤程度是可忽略不计的,这也说明了在低频率下,预制裂纹长度有一个阈值,对于低于该阈值的预制裂纹长度对裂纹的扩展造成的影响是十分小的。此外,在许多研究中[67,70-71,73-74]发现,当裂纹长度突然增加时,裂纹会发生不稳定扩展。

图28(c)说明了电池颗粒的不均匀膨胀和收缩引起的变形不仅受到充放电频率和预制裂纹长度的影响,而且还受电池颗粒半径大小的控制。结果表明,电池颗粒尺寸越大,颗粒体积变化的影响越显著,造成的疲劳损伤程度越高。同样,与充放电频率和预制裂纹长度一样,存在一个临界颗粒尺寸,低于该临界尺寸时,增大循环加载的次数也不会导致裂纹扩展。

在锂离子电池疲劳充放电方面,力-化耦合断裂相场模型能够准确地描述锂离子的脱出和嵌入过程,数值模拟结果显示预测不仅能够描述具有初始裂纹的电极颗粒的损伤过程,而且同样适用于无初始裂缝的完整电池颗粒。重点概述了锂离子电池颗粒充放电频率、预制裂纹长度、临界能量释放率、电池颗粒尺寸大小对电池颗粒的损伤,结论如下。

(1)电池颗粒的疲劳损伤程度会随着充放电频率和电池颗粒半径的增大而增大,增大充放电频率会加剧裂纹的开裂和不稳定的裂纹生长。

(2)增大预制裂纹的长度会造成更大的应力集中,但裂纹能否继续扩展基于能量准则,当边界处的

锂离子浓度降低为0时,所产生的应变能仍然不能推动裂纹扩展时,电池颗粒则处于安全状态。

(3)不同电池颗粒的临界能量释放率也会对电池颗粒的损伤程度造成影响,较小的临界能量释放率在同等加载次数下对电池颗粒的损伤程度越大。

(4)电池颗粒充放电频率、预制裂纹长度、临界能量释放率和电池颗粒尺寸大小对材料造成损伤都会有一个预制,低于该阈值对电池颗粒的损伤可忽略不计。

5 结论与展望

力-化耦合断裂相场模型基于断裂力学的变分原理,结合相场理论和扩散过程,通过引入扩散耦合项将材料的扩散行为考虑到模型中。因此,该模型在模拟由浓度变化引起的裂纹扩展问题上具有重要的应用前景,特别是在仿真锂离子电池颗粒的充放电过程以及金属材料的氢脆方面。该模型能够准确地描述锂离子的脱出和嵌入过程,数值模拟结果显示预测不仅能够描述具有初始裂纹的电极颗粒的损伤过程和浓度变化的分布情况,而且同样适用于无初始裂缝的完整电池颗粒。在金属材料的氢脆研究中,该模型能够很好地预测金属材料在氢浓度影响下的断裂行为。由于氢增强局部塑性效应和氢增强脱黏效应会降低材料的力学性能,裂纹的扩展速度随着氢浓度的升高而加快,且静水应力的存在会促进氢浓度的扩散。虽然该模型去得了一定的进展,但力-化耦合断裂相场模型在许多方面仍然有值得进一步研究的价值。

(1)化学反应与裂纹扩展的耦合:在氢辅助开裂过程中,氢与材料中的化学物质发生反应,形成氢化物或氢化物相。这些化学反应会影响材料的力学性能和氢的扩散行为。可以研究化学反应与裂纹扩展的耦合效应,以及不同化学反应条件下的裂纹扩展行为。

(2)多尺度建模:氢辅助开裂问题涉及多个尺度,从原子尺度到宏观尺度。可以通过多尺度建模方法,将不同尺度的信息进行耦合,来更准确地描述氢辅助开裂的机理和行为。

(3)求解效率:力化耦合断裂相场模型需要进行耦合的有限元模拟计算,需要较长的计算时间。需要进一步改进和优化数值求解算法,提高模型的计算效率。

参考文献

[1]

吴建营. 固体结构损伤破坏统一相场理论、算法和应用[J]. 力学学报, 2021, 53(2): 301-329.

[2]

Wu Jianying. On the unified phase-field theory for damage and failure in solids and structures: theoretical and numerical aspects[J]. Chinese Journal of Theoretical and Applied Mechanics, 2021, 53(2): 301-329.

[3]

Du X, Jin L, Zhang R. Modeling the cracking of cover concrete due to non-uniform corrosion of reinforcement[J]. Corrosion Science, 2014, 89: 189-202.

[4]

Golahmar A, Kristensen P K, Niordson C F, et al. A phase field model for hydrogen-assisted fatigue[J]. International Journal of Fatigue, 2022, 154: 106521.

[5]

Mandal T K, Nguyen V P, Wu J Y. Comparative study of phase-field damage models for hydrogen assisted cracking[J]. Theoretical and Applied Fracture Mechanics, 2021, 111: 102840.

[6]

Isfandbod M, Martínez-Pañeda E. A mechanism-based multi-trap phase field model for hydrogen assisted fracture[J]. International Journal of Plasticity, 2021, 144: 103044.

[7]

Singh A, Pal S. Coupled chemo-mechanical modeling of fracture in polycrystalline cathode for lithium-ion battery[J]. International Journal of Plasticity, 2020, 127: 102636.

[8]

Xue R, Li X, Zhao H, et al. Phase field model coupling with strain gradient plasticity for fracture in lithium-ion battery electrodes[J]. Engineering Fracture Mechanics, 2022, 269: 108518.

[9]

Pistorio F, Clerici D, Mocera F, et al. Review on the numerical modeling of fracture in active materials for lithium ion batteries[J]. Journal of Power Sources, 2023, 566: 232875.

[10]

Cao X, Lu Y, Chen Z, et al. A chemo-mechanical phase-field framework for dynamic fracture with viscoplastic flow for large-deformed electrode in lithium-ion batteries[J]. Journal of Alloys and Compounds, 2023, 965: 171387.

[11]

Mohtarami E, Baghbanan A, Eftekhari M, et al. Investigating of chemical effects on rock fracturing using extended finite element method[J]. Theoretical and Applied Fracture Mechanics, 2017, 89: 110-126.

[12]

Luo Z, Zhang N, Zhao L, et al. An extended finite element method for the prediction of acid-etched fracture propagation behavior in fractured-vuggy carbonate reservoirs[J]. Journal of Petroleum Science and Engineering, 2020, 191: 107170.

[13]

Chen Z, Dai Y, Liu Y. Life prediction of corrosion-fatigue based on a new crack growth rate model with damage and the extended finite element method[J]. Engineering Fracture Mechanics, 2023, 289: 109445.

[14]

Dinachandra M, Alankar A. Adaptive finite element modeling of phase-field fracture driven by hydrogen embrittlement[J]. Computer Methods in Applied Mechanics and Engineering, 2022, 391: 114509.

[15]

Zhang X, Jin H, Yu S, et al. Analysis of bending deflection of tunnel segment under load- and corrosion-induced cracks by improved XFEM[J]. Engineering Failure Analysis, 2022, 140: 106576.

[16]

Scheider I, Pfuff M, Dietzel W. Simulation of hydrogen assisted stress corrosion cracking using the cohesive model[J]. Engineering Fracture Mechanics, 2008, 75(15): 4283-4291.

[17]

Zheng T, Chen N Z. A cyclic cohesive zone model for predicting hydrogen assisted fatigue crack growth (FCG) of subsea pipeline steels[J]. International Journal of Fatigue, 2023, 173: 107707.

[18]

Chen Z, Yang D, Bian H. Peridynamic modeling of crack propagation driven by hydrogen embrittlement[J]. Engineering Fracture Mechanics, 2023, 293: 109687.

[19]

Wang X, Tong Q. Peridynamic modeling of delayed fracture in electrodes during lithiation[J]. Computer Methods in Applied Mechanics and Engineering, 2023, 404: 115774.

[20]

袁训锋, 杨燕, 刘飞. 多晶硅凝固界面形态演化过程的相场法模拟分析[J]. 科学技术与工程, 2019, 19(34): 123-128.

[21]

Yuan Xunfeng, Yang Yan, Liu Fei. Research on interface morphology evolution process of solidification for multicrystalline silicon using phase-field method[J]. Science Technology and Engineering, 2019, 19(34): 123-128.

[22]

王辉明, 单禹宁, 王国政. 基于统一相场理论和内聚力模型的钢纤维混凝土损伤破坏细观研究[J]. 科学技术与工程, 2024, 24(32): 14003-14012.

[23]

Wang Huiming, Shan Yuning, Wang Guozheng. Research on mesoscopic damage and failure of steel fiber reinforced concrete based on the unified phase field theory and cohesive zone model[J]. Science Technology and Engineering, 2024, 24(32): 14003-14012.

[24]

李明耀, 李绍金, 左建平, 等. 基于相场法的岩石三点弯曲细观破坏数值模拟[J]. 科学技术与工程, 2022, 22(30): 13441-13449.

[25]

Li Mingyao, Li Shaojin, Zuo Jianping, et al. Numerical simulation of the three-point bending failure behavior of basalt based on the phase field method[J]. Science Technology and Engineering, 2022, 22(30): 13441-13449.

[26]

Zhao J, Cheng Y F. A phase field method for predicting hydrogen-induced cracking on pipelines[J]. International Journal of Mechanical Sciences, 2024, 283: 109651.

[27]

Suvin V S, Ooi E T, Song C, et al. Adaptive scaled boundary finite element method for hydrogen assisted cracking with phase field model[J]. Theoretical and Applied Fracture Mechanics, 2024, 134: 104690.

[28]

Sun B, Zhao H, Dong X, et al. Current challenges in the utilization of hydrogen energy—a focused review on the issue of hydrogen-induced damage and embrittlement[J]. Advances in Applied Energy, 2024, 14: 100168.

[29]

Singh V, Raj A, Charles Y, et al. Hydrogen embrittlement in nickel oligocrystals: experimentation and crystal plasticity-phase field fracture modeling[J]. International Journal of Hydrogen Energy, 2024, 84: 667-681.

[30]

Zhang Y, Xu P, Ding W, et al. Phase-field simulation of dual-crack system hydrogen embrittlement in metallic materials[J]. Theoretical and Applied Fracture Mechanics, 2024, 131: 104332.

[31]

Noii N, Milijasevic D, Khodadadian A, et al. An efficient fenics implementation for coupling lithium-ion battery charge/discharge processes with fatigue phase-field fracture[J]. Engineering Fracture Mechanics, 2024, 306: 110251.

[32]

Chen W X, Peng X L, Wu J Y, et al. A length-scale insensitive cohesive phase-field interface model: application to concurrent bulk and interface fracture simulation in lithium-ion battery materials[J]. Journal of the Mechanics and Physics of Solids, 2025, 196: 106013.

[33]

Wu W, Lin Y, Hu Y, et al. Phase-field modelling for degradation/failure research in lithium battery: progress and prospects[J]. Journal of Energy Chemistry, 2025, 102: 628-650.

[34]

Bourdin B, Francfort G A, Marigo J J. Numerical experiments in revisited brittle fracture[J]. Journal of the Mechanics and Physics of Solids, 2000, 48(4): 797-826.

[35]

Amor H, Marigo J J, Maurini C. Regularized formulation of the variational brittle fracture with unilateral contact: numerical experiments[J]. Journal of the Mechanics and Physics of Solids, 2009, 57(8): 1209-1229.

[36]

Miehe C, Welschinger F, Hofacker M. Thermodynamically consistent phase-field models of fracture: variational principles and multi-field FE implementations[J]. International Journal for Numerical Methods in Engineering, 2010, 83(10): 1273-1311.

[37]

Pham K, Amor H, Marigo J J, et al. Gradient damage models and their use to approximate brittle fracture[J]. International Journal of Damage Mechanics, 2011, 20(4): 618-652.

[38]

Wu J Y. A geometrically regularized gradient-damage model with energetic equivalence[J]. Computer Methods in Applied Mechanics and Engineering, 2018, 328: 612-637.

[39]

Wu J Y. A unified phase-field theory for the mechanics of damage and quasi-brittle failure[J]. Journal of the Mechanics and Physics of Solids, 2017, 103: 72-99.

[40]

Wu J Y. Robust numerical implementation of non-standard phase-field damage models for failure in solids[J]. Computer Methods in Applied Mechanics and Engineering, 2018, 340: 767-797.

[41]

Wu J Y. Unified analysis of enriched finite elements for modeling cohesive cracks[J]. Computer Methods in Applied Mechanics and Engineering, 2011, 200(45/46): 3031-3050.

[42]

Alessi R, Marigo J J, Vidoli S. Gradient damage models coupled with plasticity: variational formulation and main properties[J]. Mechanics of Materials, 2015, 80: 351-367.

[43]

Karma A, Kessler D A, Levine H. Phase-field model of mode III dynamic fracture[J]. Physical Review Letters, 2001, 87(4): 045501.

[44]

Borden M J, Hughes T J R, Landis C M, et al. A phase-field formulation for fracture in ductile materials: finite deformation balance law derivation, plastic degradation, and stress triaxiality effects[J]. Computer Methods in Applied Mechanics and Engineering, 2016, 312: 130-166.

[45]

Kuhn C, Schlüter A, Müller R. On degradation functions in phase field fracture models[J]. Computational Materials Science, 2015, 108: 374-384.

[46]

Lorentz E. A nonlocal damage model for plain concrete consistent with cohesive fracture[J]. International Journal of Fracture, 2017, 207(2): 123-159.

[47]

Lorentz E, Cuvilliez S, Kazymyrenko K. Convergence of a gradient damage model toward a cohesive zone model[J]. Comptes Rendus Mécanique, 2011, 339(1): 20-26.

[48]

Lorentz E, Cuvilliez S, Kazymyrenko K. Modelling large crack propagation: from gradient damage to cohesive zone models[J]. International Journal of Fracture, 2012, 178(1/2): 85-95.

[49]

Ambati M, Gerasimov T, De Lorenzis L. A review on phase-field models of brittle fracture and a new fast hybrid formulation[J]. Computational Mechanics, 2015, 55(2): 383-405.

[50]

Lancioni G, Royer-Carfagni G. The variational approach to fracture mechanics: a practical application to the French Panthéon in Paris[J]. Journal of Elasticity, 2009, 95(1/2): 1-30.

[51]

Zhang F, Huang W, Li X, et al. A study on phase-field models for brittle fracture[J]. International Journal of Numerical Analysis and Modeling, 2018, 19(6): 793-821.

[52]

Wu J Y, Nguyen V P, Zhou H, et al. A variationally consistent phase-field anisotropic damage model for fracture[J]. Computer Methods in Applied Mechanics and Engineering, 2020, 358: 112629.

[53]

Anand L, Mao Y, Talamini B. On modeling fracture of ferritic steels due to hydrogen embrittlement[J]. Journal of the Mechanics and Physics of Solids, 2019, 122: 280-314.

[54]

Alvaro A, Thue J I, Kheradmand N, et al. Hydrogen embrittlement in nickel, visited by first principles modeling, cohesive zone simulation and nanomechanical testing[J]. International Journal of Hydrogen Energy, 2015, 40(47): 16892-16900.

[55]

Scheider I, Pfuff M, Dietzel W. Simulation of hydrogen assisted stress corrosion cracking using the cohesive model[J]. Engineering Fracture Mechanics, 2008, 75(15): 4283-4291.

[56]

Serebrinsky S. A quantum-mechanically informed continuum model of hydrogen embrittlement[J]. Journal of the Mechanics and Physics of Solids, 2004, 52(10): 2403-2430.

[57]

Tehranchi A, Zhang X, Lu G, et al. Hydrogen-vacancy-dislocation interactions in α-Fe[J]. Modelling and Simulation in Materials Science and Engineering, 2017, 25(2): 025001.

[58]

Jiang D E, Carter E A. First principles assessment of ideal fracture energies of materials with mobile impurities: implications for hydrogen embrittlement of metals[J]. Acta Materialia, 2004, 52(16): 4801-4807.

[59]

Martínez-Pañeda E, Golahmar A, Niordson C F. A phase field formulation for hydrogen assisted cracking[J]. Computer Methods in Applied Mechanics and Engineering, 2018, 342: 742-761.

[60]

Hondros E D, Seah M P. The theory of grain boundary segregation in terms of surface adsorption analogues[J]. Metallurgical Transactions A, 1977, 8(9): 1363-1371.

[61]

Dinachandra M, Alankar A. Adaptive finite element modeling of phase-field fracture driven by hydrogen embrittlement[J]. Computer Methods in Applied Mechanics and Engineering, 2022, 391: 114509.

[62]

Si Z, Hirshikesh, Yu T, et al. An adaptive phase-field simulation for hydrogen embrittlement fracture with multi-patch isogeometric method[J]. Computer Methods in Applied Mechanics and Engineering, 2024, 418: 116539.

[63]

Wu J Y, Mandal T K, Nguyen V P. A phase-field regularized cohesive zone model for hydrogen assisted cracking[J]. Computer Methods in Applied Mechanics and Engineering, 2020, 358: 112614.

[64]

Yang G, Yang L, Liu Z, et al. Phase field simulation of hydrogen-assisted cracking with length-scale insensitive degradation function[J]. Computational Materials Science, 2023, 228: 112309.

[65]

Kristensen P K, Niordson C F, Martínez-Pañeda E. A phase field model for elastic-gradient-plastic solids undergoing hydrogen embrittlement[J]. Journal of the Mechanics and Physics of Solids, 2020, 143: 104093.

[66]

Huang C, Gao X. Phase field modeling of hydrogen embrittlement[J]. International Journal of Hydrogen Energy, 2020, 45(38): 20053-20068.

[67]

Zhang Y, Zhao C, Guo Z. Simulation of crack behavior of secondary particles in Li-ion battery electrodes during lithiation/de-lithiation cycles[J]. International Journal of Mechanical Sciences, 2019, 155: 178-186.

[68]

Tian H, Gao L T, Guo Z S. Microstructural adjusting crack evolution of polycrystalline NCM particle during charge/discharge cycle[J]. Journal of the Electrochemical Society, 2022, 169(9): 090513.

[69]

Tian H, Gao L T, Huang P Y, et al. Simulation of intergranular fracture behavior inside randomly aggregated LiNixCoyM${{{{n}_{1-}}_{x}}_{-}}_{y}$O2 polycrystalline particle[J]. Engineering Fracture Mechanics, 2022, 266: 108381.

[70]

Xu R, Zhao K. Corrosive fracture of electrodes in Li-ion batteries[J]. Journal of the Mechanics and Physics of Solids, 2018, 121: 258-280.

[71]

Zhu X, Chen Y, Chen H, et al. The diffusion induced stress and cracking behaviour of primary particle for Li-ion battery electrode[J]. International Journal of Mechanical Sciences, 2020, 178: 105608.

[72]

Alessi R, Vidoli S, De Lorenzis L. A phenomenological approach to fatigue with a variational phase-field model: the one-dimensional case[J]. Engineering Fracture Mechanics, 2018, 190: 53-73.

[73]

Carrara P, Ambati M, Alessi R, et al. A framework to model the fatigue behavior of brittle materials based on a variational phase-field approach[J]. Computer Methods in Applied Mechanics and Engineering, 2020, 361: 112731.

[74]

Klinsmann M, Rosato D, Kamlah M, et al. Modeling crack growth during Li insertion in storage particles using a fracture phase field approach[J]. Journal of the Mechanics and Physics of Solids, 2016, 92: 313-344.

[75]

Ai W, Wu B, Martínez-Pañeda E. A coupled phase field formulation for modelling fatigue cracking in lithium-ion battery electrode particles[J]. Journal of Power Sources, 2022, 544: 231805.

[76]

Sun H H, Manthiram A. Impact of microcrack generation and surface degradation on a nickel-rich layered Li(Ni0. 9Co0. 05Mn0. 05)O2 cathode for lithium-ion batteries[J]. Chemistry of Materials, 2017, 29(19): 8486-8493.

[77]

Mesgarnejad A, Karma A. Phase field modeling of chemomechanical fracture of intercalation electrodes: role of charging rate and dimensionality[J]. Journal of the Mechanics and Physics of Solids, 2019, 132: 103696.

[78]

Klinsmann M, Rosato D, Kamlah M, et al. Modeling crack growth during Li extraction in storage particles using a fracture phase field approach[J]. Journal of the Electrochemical Society, 2016, 163(2): 102-118.

基金资助

国家自然科学基金(12002230)

PDF (19082KB)

39

访问

0

被引

详细

导航
相关文章

/

〈 〉