应变梯度长度的试验标定：从经典断裂韧性到简化SGE参数的映射与精化

1  引言

经典线弹性断裂力学（LEFM）以应力强度因子KIc为判据，在裂纹长度远大于材料微观特征尺寸时具有很高的预测精度。然而当裂纹尺度接近材料细观结构尺度时，材料表现出明显的尺寸效应：名义强度不再满足σ∝a⁻¹/²的经典幂律规律，而出现由短裂纹向长裂纹的渐变过渡。这类非经典尺寸效应在准脆性材料（如陶瓷、混凝土、纤维增强复合材料、多孔材料）中尤为突出，传统LEFM框架难以统一描述。

为弥补这一不足，Toupin与Mindlin最早发展了应变梯度弹性理论（SGE），其应变能密度同时依赖于应变本身及其梯度，使裂纹应力场正则化，避免了经典弹性理论中的r⁻¹/²奇异性。后续的广义及简化梯度理论将附加材料常数从最初的5个逐步减至3、2乃至1个。其中Aifantis形式的简化SGE仅引入单一长度尺度参数l即可在保留正则化效应的同时极大简化计算，因而在断裂力学问题中获得了广泛应用。

然而，SGE长度尺度参数l的实验辨识一直是制约该理论工程应用的核心难题：l既无法通过常规单轴拉伸试验获得，也不便直接由压痕或疲劳试验反推。Askes和Susmel曾基于隐式梯度法给出l与LEFM特征长度LLEFM≈KIc²/(πσult²)间的比例关系，但所得到的比例系数(1/8)依赖于简化假设。针对以上问题，Yury Solyaev、 Kirill Shelkov、Pavel Polyakov在“Materials Science"期刊上发表了题为“Identification of the length scale parameter of simplified strain gradient elasticity from standard Mode I fracture tests"的论文，该文基于精确满足C¹-连续性要求的富集有限元法，对三种标准I型断裂试样进行系统数值模拟，显著细化了该比例系数，并将该方法推广到短裂纹过渡区。本文按其原文结构，依次梳理简化SGE理论框架、富集C¹-FEM数值实现、收敛性验证、两类参数辨识方法及多类准脆性材料的应用结果。

图1  (a)裂纹富集有限元网格（蓝色为富集单元）；(b)CCT、SENT与SENB三种试样几何与边界条件

2  简化应变梯度弹性理论

简化SGE的应变能密度w定义为应变εij与其一阶梯度κijk=εij,k的函数，w=½Cijklεijεkl+½Aijklmnκijkκlmn。其中Cijkl为经典四阶弹性张量，Aijklmn为含附加材料常数的六阶张量。在Aifantis简化形式中，Aijklmn=l²Cijlmδkn，使理论中仅含一个具有长度量纲的附加常数l，l=0时自动退化到经典线弹性。

由w对εij和κijk分别求导，得到经典应力τij=λδijεkk+2μεij与偶应力µijk=l²τij,k，后者物理上等价于经典应力梯度乘以l²。由此建立的边值问题除满足经典平衡方程与外力-位移边界条件外，还需补充双应力型边界条件mi、ui,jnj及棱线载荷si。正是这些附加项在裂尖附近自然地产生类Barenblatt-Dugdale型内聚牵引效应，使裂尖应力场有界。

作者指出，在大多数工程边界条件下只需规定面力ti和零双应力mi=0；对称面则需额外约束法向位移梯度（剪力分量为零），以保证C¹-连续性。简化SGE形式的简洁性使其可被嵌入到标准有限元框架中，并成为本文数值方案的理论基础。

3  富集C¹-有限元数值实现

作者将变分原理作用于SGE的弱形式，得到以节点位移u、其一阶与二阶梯度为自由度的离散方程Kbu=f。域内采用Bell三角形上的五阶多项式作为形函数，使C¹-连续性在各节点上得到严格保证；每个非富集单元具有36个自由度。围绕裂尖的若干单元则采用基于SGE平面应变裂尖渐近解的富集策略，引入两个幅值因子K1、K2作为附加自由度（共38个），其几何与边界条件按图1(b)分别对应CCT、SENT与SENB三种试样。

J积分由所得幅值因子按下式计算：J=(1+η)/(8µπl²)·[(3K1+K2)²+8K2²(η+2)/(η+2)]，其中η=3−4ν为Kolosov常数。该方法在Abaqus中通过UEL子程序实现，沿对称轴和y轴的位移、位移一阶和二阶切向导数均被显式约束。

4  收敛性分析与J积分验证

作者针对a/b=0.01的CCT近似无限大板试样，依次考察了高斯积分点数目、裂尖富集单元尺寸R/l与径向层数N对正则化应力τ22、幅值因子K1、K2和J积分的影响。

如图2所示，正则化应力τ22在13点高斯积分下即收敛，但K1、K2两个幅值因子需42点才达到稳定收敛；这一结论细化了文献[26]的早期结果。图3表明，要使τ22达到可靠精度只需R/l=0.1、N=5–10的网格，但K1、K2对N更敏感，至少需N=15。综合考虑精度与效率，作者推荐R/l=0.01、N=26的网格用于后续所有计算。

在以上参数下，作者进一步将所得J积分与Gourgiotis和Georgiadis的半解析解比较（图4a），二者高度吻合；并在图4b中证实了长裂纹极限a/l≫1下JSGE/JLEFM→1，确认SGE在长裂纹下渐近回归LEFM的J积分判据。

图2  富集C¹-FEM解关于积分点数的收敛性：(a)正则化应力τ₂₂；(b)幅值因子K₁；(c)幅值因子K₂

图3  富集C¹-FEM解关于裂尖网格尺寸R/l及径向层数N的收敛性：(a)τ₂₂；(b)K₁；(c)K₂

图4  J积分验证：(a)数值解（点）与Gourgiotis-Georgiadis半解析解（线）比较；(b)J(l)/J(0)随a/l的变化，验证长裂纹极限下回到LEFM

5  长度尺度参数辨识方法

作者给出了两种独立的辨识方法，分别针对裂纹长度远大于l的脆性材料和需要考虑短裂纹过渡的准脆性材料。

5.1  脆性材料（长裂纹）

在长裂纹极限下，简化SGE内的名义强度σn与a/l的关系在双对数坐标中同样呈−1/2的斜率，故可与LEFM的σn=(a/LLEFM)⁻¹/²显式对齐。对b/a=1000的CCT近似无限大板，作者对多个l值进行数值模拟并拟合，最佳比例系数为k=1.137254，从而得到脆性材料的显式公式：l=k·KIc²/(πσult²)≈0.362(KIc/σult)²。拟合在a/l>10的范围内平均绝对百分误差MAPE仅0.24%。

对SENT和SENB试样只需将经典LEFM的KIc替换为KIc/1.1215即保持一致。由图5可见，a/l≫1时CCT曲线与SENT/SENB曲线严格平行，证实该修正的普适性。图6(b)给出的回归式σ∞(a/l)=[1+k(a/l)²/(a/l+c0)]⁻¹/²（k=1.137254, c0=2.088）在10⁻¹<a/l<10²范围内以1.6%的MAPE刻画了SGE数值结果，并同时具有短裂纹极限σn→1（强度平台）与长裂纹极限回归LEFM两个渐近行为。表1列出了多种脆性材料（退火浮法玻璃、AISI 4340钢、Al₂O₃、SiC、PMMA、环氧及其纳米复合体系、淡水冰、海冰、核级石墨IG-110/IG-430等）的KIc、σult、LLEFM和由(18)式求得的l。l值从数十微米跨越到毫米量级，对纳米复合材料尤为明显地大于其代表体积元(RVE)尺寸，提示按经验假设l≈RVE尺寸并不可靠。

图5  a/b=0.1时不同试样的(a)应力集中因子Kt与(b)名义强度σn随归一化裂纹长度a/l的变化

图6  无限大板中心裂纹的名义强度：(a)不同k值的SGE解与LEFM解析解的拟合过程；(b)回归式(19)与数值SGE解比较

5.2  准脆性材料（短-长裂纹全范围）

对有限尺寸试样，名义强度同时依赖a/l与a/b，作者提出总回归式：σn=σ∞(a/l)·F(a/b)⁻¹·G(a/b, a/l)，其中F(a/b)为经典LEFM几何修正因子（CCT、SENT、SENB的解析表达式见表2），G(a/b, a/l)是新引入的SGE附加修正因子，并须满足a/l→∞时G→1（保证与LEFM一致）和a/l→0时G→F（保证强度平台σn→1）的双极限条件。

对CCT和SENT，作者给出G=[X^(c1+c2)]/[X^(c1+c2/F)]形式，X=(a/l)(1−a/b)；对SENB则增加一个衰减型修正项(1+c3(a/b)(a/l)exp(−c4·a/l))，回归系数见表3。图7、图8展示了三种试样上回归式与SGE数值解的吻合程度和最大相对误差分布：CCT的MAPE=1.37%、SENT为2.10%、SENB为2.59%，足以支持实验数据反演l。

图7  回归关系(线)与SGE数值解(点)的对比：(a)CCT；(b)SENT；(c)SENB

图8  回归关系与SGE数值解的相对误差分布：(a)CCT；(b)SENT；(c)SENB

6  准脆性材料长度尺度参数辨识实例

作者将上述回归流程应用于四类典型准脆性材料，全部无需再行数值模拟。

6.1  短切碳纤维/环氧复合材料的SENT试验

对Ko等[84]报道的75×12、50×8、25×4 mm三种短切碳纤维胶带铺层试样的SENT数据（a/b=0.2，宽度b=6.3–120 mm）进行反演。结果（图9、表4）显示三种胶带对应的l分别为10.7、8.2、5.0 mm，与胶带宽度的数量级高度一致；最大绝对百分误差MAPE分别为9.6%、8.0%、15.3%，NRMSE为10.6%、8.9%、17.0%。这说明在短切纤维铺层结构中，l主要反映最小细观特征尺寸（胶带厚度方向）。

图9  短切碳纤维/环氧复合材料SENT试样断裂强度对裂纹长度的依赖：(a)75×12 mm；(b)50×8 mm；(c)25×4 mm

6.2  高孔隙率SiO₂陶瓷的SENB试验

作者自行加工并测试了孔隙率92%的多孔石英陶瓷SENB试样（图10a-c，b=16 mm、h=12 mm、L=96 mm，跨高比4）。30个试样跨0–9 mm裂纹长度，Instron 4959机1 mm/min十字头速率下测得峰值载荷Pfrac，结合三点弯公式σc=3Pfrac·s/(hb²)反推临界应力。由图10(d)和表5，拟合得到l=0.76 mm，MAPE=17.1%、NRMSE=20.8%。由于该材料内部为连通的延伸孔道通道，非局部效应显著，使l明显大于其颗粒尺寸。

图10  多孔SiO₂陶瓷SENB试验：(a)试样照片；(b)三点弯曲断裂瞬间；(c)缺口与预裂纹的光学显微图；(d)临界断裂应力随裂纹长度的实验点(点)与SGE回归(线)

6.3  闭孔隙SiO₂和SiC耐火陶瓷的SENB试验

Bazant和Kazemi[85]对McKinney-Rice的滑注成型SiO₂(孔隙率12%, 晶粒10–20 µm)和SiC CN-163(孔隙率15%, 最大晶粒约2 mm)SENB试验数据按相同流程处理(图11)，得到l分别为1.3 mm和1.5 mm，均显著大于其晶粒尺寸，再次印证l与孔隙尺寸/形貌的强相关性。

图11  Bazant-Kazemi多孔SiO₂与SiC耐火陶瓷SENB试验反演：(a)SiO₂；(b)SiC CN-163

6.4  致密结构陶瓷的SENB试验

Usami等[87]对Al₂O₃、SiC、Si₃N₄、Sialon四种致密陶瓷的小裂纹SENB数据进行处理(图12)。回归法所得l从5.8 µm(Sialon)到39 µm(Al₂O₃)不等，与各自晶粒尺寸的比值l/d集中在1.95–3.00，与Usami等"晶粒-断裂模型"中r0=2d的设定相符。利用式(18)直接由KIc估算的l则偏大(图12蓝线)，在短裂纹区会高估断裂载荷。

图12  致密结构陶瓷SENB试验反演：(a)Al₂O₃；(b)SiC；(c)Si₃N₄；(d)Sialon

7  结论与展望

Solyaev等的工作表明，简化SGE的单一长度尺度参数l可以由标准I型断裂试验直接获得：对于脆性材料只需KIc和σult两个经典参数即可显式计算l≈0.362(KIc/σult)²；对于准脆性材料则可使用本文给出的统一回归式，无需额外有限元模拟即可对实验数据反演l。多类材料实例进一步显示，l与晶粒尺寸、孔隙率、纤维铺层最小尺度等内部结构特征存在系统联系。

研究中同时观察到SGE幅值因子K1、K2在长裂纹下与a近似呈线性、短裂纹下指数衰减至零的非平凡行为，提示可建立基于幅值因子的K型断裂判据。作者指出未来工作需在以下几方面推进：①将方法推广到混合模式加载及更一般应力集中构型(尖缺口、小孔等)；②发展基于J积分或幅值因子的能量型/类K型判据；③采用含两个或更多长度尺度的广义SGE本构以更精细地描述过渡区尺寸效应。附录B给出的CCT/SENT/SENB三种试样上J积分和幅值因子的全套数值结果，为上述后续研究提供了可直接使用的基准数据。