相场断裂方法发展概况

张豪; 于继东; 裴晓阳; 彭辉; 李平; 蔡灵仓; 汤铁钢

doi:10.11858/gywlxb.20190777

相场断裂方法发展概况

doi: 10.11858/gywlxb.20190777

中国工程物理研究院流体物理研究所，四川绵阳　621999

基金项目: 国家自然科学基金面上项目（11772067）；国家自然科学基金青年科学基金（11702277）；国家自然科学基金重点项目（11532012）

详细信息

作者简介:
张　豪（1988－），男，博士研究生，助理研究员，主要从事高压物理与力学研究. E-mail：zhanghao17@gscaep.ac.cn

通讯作者:
蔡灵仓（1964－），男，博士生导师，研究员，主要从事凝聚态物理研究. E-mail：cai_lingcang@aliyun.com

中图分类号: O347.3
计量
- 文章访问数: 9041
- HTML全文浏览量: 3706
- PDF下载量: 188
出版历程
- 收稿日期: 2019-05-15
- 修回日期: 2019-05-22
- 刊出日期: 2019-05-25

An Overview of Phase Field Approach to Fracture

Institute of Fluid Physics, CAEP, Mianyang 621999, China

摘要

摘要: 相场断裂方法自21世纪初开始发展以来，一直备受关注，在裂纹扩展模拟方面表现出了一定的特点，取得了一定的研究成果。本文分析了相场断裂方法较其他断裂模拟方法的优势，简单介绍了相场断裂方法的发展现状和发展趋势：目前脆性断裂相场方法已较为成熟，能够模拟诸多脆性断裂中的经典问题，在此基础上正在朝着解决多场耦合情况下的断裂问题发展，且也取得了一定的研究成果。最后，简单介绍了延性断裂相场方法的发展现状，提出在该方向进行深入研究的展望。
- 相场断裂方法 /
- 数值模拟 /
- 脆性断裂 /
- 延性断裂
Abstract: Phase field modeling to fracture has received much attention since the beginning of this century, which exhibits an advantage in fracture propagation simulation. In this work, we compare the phase field approach to fracture with other simulation methods, and show an overview and development of phase field approach to fracture. Up to now, the phase field method has been successfully applied to the brittle fracture and could simulate some classical crack problems. Based on this, the multi-fields problem coupled with the fracture is currently pursued. Furthermore, we introduce the study situation of the phase field simulation to the ductile fracture and put forward its development in the future.
- phase field approach to fracture /
- simulation /
- brittle fracture /
- ductile fracture

HTML全文

随着人们对安全要求的不断提高和人为爆炸事故的频繁发生，发生爆炸时工程的安全可靠性已逐渐成为人们关注的重要问题，深入了解爆炸载荷下结构的动态响应具有重要意义。与圆管相比，填充可压缩泡沫芯层的夹芯管能够通过冲击衰减和塑性功耗散来吸收大量爆炸能量，具有更好的抗爆性能。因此，泡沫夹芯结构作为抗爆吸能装置，引起了研究者们的广泛关注^[1-3]。

Fan等^[4]研究了短空心圆管和夹芯管在横向冲击载荷下的动态响应，发现结构的能量耗散主要通过塑性铰的形成来实现。Hall等^[5]研究了横向冲击载荷下泡沫铝、泡沫黄铜和泡沫钛填充管的能量吸收性能，结果表明，泡沫铝填充管的能量吸收特性最佳。Song等^[6]提出了横向爆炸载荷作用下钢圆管塑性变形的解析模型，预测了钢管的挠度和变形角，发现理论预测结果与实验结果吻合较好。Yuen等^[7]对横向爆炸载荷作用下的圆管进行了实验和数值模拟分析，并与Wierzbicki等^[8]提出的理论模型进行了对比，发现实验结果与数值模拟结果吻合良好，但理论模型预测的最终挠度偏低。Jing等^[9-11]研究了爆炸载荷下泡沫铝夹芯壳结构的动态响应和吸能特性，结果表明，载荷强度和几何构型影响结构的变形/破坏模式和吸能特性，增加面板厚度或泡沫铝芯层的密度可提高其抗爆性能，并建立了流固耦合、芯层压缩和整体弯曲拉伸的三阶段框架理论模型预测结构的动态响应。Hoo Fatt等^[12]考虑了泡沫芯层的弹塑性行为，对外爆炸载荷作用下的复合材料夹芯壳进行了应力波分析，结果表明，高峰值、低持续时间的冲击波作用下结构能量吸收最高，低峰值、长持续时间的冲击波反之。Liu等^[13]通过实验和数值模拟研究了两端固支的泡沫铝填充管在横向外爆炸载荷作用下的动力响应和能量吸收性能，结果表明，泡沫铝芯层的引入减小了结构跨中挠度，提高了其能量吸收能力。上述研究表明，填充泡沫的管状结构具有良好的抗冲击性能，是一种理想的吸能装置。

多目标优化可以平衡各个目标，使目标之间作出协调和折中处理，最终总体的目标尽可能达到最优。Baroutaji等^[14-15]通过多目标优化设计方法对准静态横向压缩下的圆管和泡沫铝夹芯管进行了优化，发现圆管在最大壁厚和最小直径下耐撞性最好，而夹芯管在最大芯层厚度和最小内管直径下耐撞性最好。Qi等^[16]对爆炸载荷下的泡沫铝夹芯板进行了多目标优化，通过非支配排序遗传算法（non-dominated sorting genetic algorithm, NSGA-Ⅱ）对结构的几何尺寸和芯层的相对密度进行优化，得到了最佳参数，平衡了结构的吸能能力和结构变形量。

综上所述，已有的研究主要集中于横向冲击载荷下薄壁圆管、填充管和夹芯壳结构的动态响应和能量吸收性能的研究。但是，关于泡沫金属夹芯管结构在横向爆炸载荷作用下的变形失效模式、吸能机制及抗冲击性能等方面的研究比较有限，已有的研究也局限于带曲率的夹芯壳结构，针对爆炸载荷下夹芯管结构的多目标优化的研究尚不多见。本研究基于动力显式有限元法，探究泡沫铝夹芯管在横向爆炸载荷下的动态响应，并对其抗爆性能进行多目标优化。建立横向爆炸作用下泡沫铝夹芯管的有限元模型，研究泡沫铝夹芯管的几何参数、芯层的相对密度和爆炸加载条件对内外管跨中挠度和吸能性能的影响；以内管、外管壁厚和芯层相对密度为设计变量，以芯层比吸能和外管跨中挠度为优化指标，采用响应面法（response surface method, RSM）和NSGA-Ⅱ算法实现泡沫铝夹芯管抗爆性能多目标优化设计，为抗爆结构设计提供技术指导。

1. 有限元模拟

1.1 有限元模型

泡沫铝夹芯管由内、外两个同心的薄壁圆管和泡沫铝芯层3部分组成，内、外薄壁圆管采用钢材，泡沫铝芯层采用相对密度分别为10%、15%和20%的泡沫铝材料。图1为泡沫铝夹芯管横截面的示意图，其中：D为外管外径，H为外管壁厚，d为内管外径，h为内管壁厚。试件长度为L。

图 1 泡沫铝夹芯管

Figure 1. Aluminum foam-filled sandwich tube

下载: 全尺寸图片幻灯片

采用有限元软件ABAQUS/EXPLICIT建立泡沫铝夹芯管的有限元模型，如图2所示。泡沫铝夹芯管的中线与z轴重合，内管和外管采用常规的四节点壳单元（S4R），泡沫铝芯层采用八节点线性实体单元（C3D8R），采用中性轴算法划分网格。内、外管与泡沫芯层之间定义为绑定约束，试件两端约束全部自由度来模拟两端固支的边界条件，分析过程采用动态显式算法。

图 2 有限元模型

Figure 2. Finite element model

下载: 全尺寸图片幻灯片

内外管采用各向同性硬化双线性本构模型，密度ρ=7 830 kg/m³，杨氏模量E=193 GPa，泊松比μ=0.25，屈服应力σ₀=205 MPa，在最大应变ε_UTS=0.4时，最大拉伸强度σ_UTS =520 MPa，应变率敏感性参数C和P分别取100 s⁻¹和10^[17]。泡沫铝芯层采用各向同性硬化塑性模型定义，用可压缩泡沫材料模拟泡沫铝芯层。不同相对密度泡沫铝的材料参数如表1^[18]所示，其中：ρ^*为相对密度，μ_c、μ_p分别为压缩屈服应力比和塑性泊松比。

表 1 泡沫铝材料参数^[18]

Table 1. Material parameters for aluminum foam^[18]

ρ^*/%	ρ/(g·cm⁻³)	E/MPa	μ	μ_c	μ_p
10	0.27	48	0	1.41	0.17
15	0.40	71	0	1.12	0.29
20	0.54	143	0	1.10	0.30

下载: 导出CSV

| 显示表格

1.2 爆炸载荷的假设

球形炸药作用在试件上形成了脉冲载荷，假定作用在泡沫铝夹芯管上的脉冲载荷分布如图3所示，其脉冲近似方程可写为

图 3 爆炸载荷假设^[13]

Figure 3. Blast load assumption^[13]

下载: 全尺寸图片幻灯片

$p(x,\theta ,t) = \left\{ \begin{array}{*{20}{l}} {p_0}\left( {1 - t/\tau } \right){{\rm{cos}}^2}(\text{π} x/L){\rm{sin}}\,\theta & \;\;\;t \leqslant \tau \\ 0&\;\;\; t > \tau \\ \end{array} \right.$

(1)

式中：τ为爆炸载荷响应时间，取τ=50 μs^[17]；L为试件总长，L=280 mm；p₀为作用在泡沫铝夹芯管迎爆面处的超压峰值。

根据TNT球形炸药在无限大气中的超压峰值计算公式，结合实验中的实际炸药质量和起爆距离，超压峰值^[19]（单位：10⁵ Pa）可以表示为

${p_0} = \frac{{0.085}}{Z}+\frac{{0.3}}{{{Z^2}}}+\frac{{0.8}}{{{Z^3}}}$

(2)

$Z = R/\sqrt[{{3}}]{{{W_{\rm{TNT}}}}}$

(3)

式中：W_TNT为TNT炸药质量，kg；R为起爆距离，m。实验中所采用的乳化炸药与TNT炸药的换算关系为

${W_{\rm{TNT}}} = 0.7{W_{\rm{s}}}$

(4)

式中：W_s为乳化炸药质量。结合式(1)～式(4)，可得出作用在试件上的总冲量I_T为

${I_{\rm{T}}} = 4\int_0^{L/2} {\int_0^{{\text{π}} /2} {\int_0^\tau {p(x,\theta ,t)\frac{D}{2}{\rm{d}}x{\rm{d}}\theta {\rm{d}}t = } } } \frac{1}{4}{p_0}DL\tau$

(5)

表2列出了横向爆炸载荷下两端固支的两组泡沫铝填充管试件的理论冲量与文献[13]中实验冲量的对比，其中m为泡沫铝填充管的质量，冲量I_E为实验结果，超压峰值p₀由式(2)计算得到，总冲量I_T由式(5)给出。由表2可知，理论计算冲量与实验冲量吻合较好，验证了爆炸载荷假设的合理性。

表 2 理论冲量与实验冲量的对比

Table 2. Comparison of impulse between theoretical solution and experiment

Test piece	W_s/g	R/mm	D/mm	H/mm	m/kg	p₀/MPa	I_E/(N·s)	I_T/(N·s)	$\dfrac{ {I{_{ \rm{T } } }-I{_{ \rm{E} } } } }{ {I{_{ \rm{E} } } } }\Big/$ %
No.1	30	40	76	0.7	96.8	27.7	7.1	7.4	4.2
No.2	50	40	76	0.7	96.7	45.8	12.9	12.2	−5.4

下载: 导出CSV

| 显示表格

1.3 数值模拟验证及泡沫铝夹芯管网格敏感性验证

为了验证本研究中有限元模型的有效性，根据文献[13]中的实验模型，利用ABAQUS建立有限元模型，爆炸载荷以式(1)的形式施加在泡沫铝填充管上。图4给出了试件1的数值模拟和实验的最终变形模式，其中泡沫铝夹芯管截面中的4.1 mm为最终变形区域的径向厚度尺寸。由图4可知，结构的整体变形模式与有限元结果吻合较好。表2中试件1的外管跨中挠度的实验值与模拟值分别为14.0和12.6 mm，试件2的外管跨中挠度的实验值和模拟值分别为22.1和29.0 mm，数值模拟结果与实验测试结果吻合较好，表明所建立的有限元模型是合理的。

图 4 实验^[13]与数值模拟结果的对比

Figure 4. Comparison of experimental^[13] and numerical simulation results

下载: 全尺寸图片幻灯片

为研究横向爆炸载荷下泡沫铝夹芯管的动态响应，对其有限元模型的网格敏感性进行了验证。采用乳化炸药的质量为30 g，起爆位置距离试件跨中40 mm，载荷以式(1)给出的脉冲载荷函数形式施加在泡沫铝夹芯管上。试件的几何参数：外管外径为101 mm，壁厚为0.7 mm；内管外径为54 mm，壁厚为0.7 mm，泡沫铝芯层的相对密度为15%。泡沫铝芯层和内、外管分别采用1.5、2.0、2.5和3.0 mm共4种不同的单元网格尺寸。图5给出了在横向爆炸载荷下4种不同单元网格尺寸泡沫铝夹芯管试件内外管跨中挠度时程曲线，其中虚线和实线分别代表内、外管挠度曲线，外管的挠度大于内管挠度，内外管跨中挠度随时间上下波动，取第一个波谷点后的挠度曲线的平均值作为泡沫铝夹芯管的最终挠度。当网格尺寸为1.5、2.0和2.5 mm时，三者显示了较好的一致性，而网格尺寸为3.0 mm时，数值模拟与其他网格尺寸结果相差较大。综合考虑计算精度和效率，数值模拟中采用2.0 mm的网格尺寸。

图 5 网格敏感性验证

Figure 5. Grid sensitivity verification

下载: 全尺寸图片幻灯片

2. 模拟结果与分析

为研究爆炸载荷下泡沫铝夹芯管结构的动态响应和能量吸收特性，对结构的几何参数、爆炸加载条件进行分析，研究结构的内外管直径与壁厚、泡沫铝芯层的相对密度、炸药起爆距离及炸药质量对泡沫铝夹芯管动态响应的影响。考虑9组不同几何参数的泡沫铝夹芯管试件在横向爆炸载荷下的抗爆性能，其中泡沫芯层的相对密度为15%，试件详细的几何参数如表3所示。

表 3 试件几何参数

Table 3. Geometric parameter of specimens

Test piece	D/mm	d/mm	H/mm	h/mm	L/mm	m/g
NS01	101	76	0.7	0.7	280	1 242.62
NS02	101	65	0.7	0.7	280	1 326.02
NS03	101	54	0.7	0.7	280	1 388.13
NS04	90	54	0.7	0.7	280	1 150.28
NS05	79	54	0.7	0.7	280	933.72
NS06	101	76	0.8	0.7	280	1 312.19
NS07	101	76	0.9	0.7	280	1 381.75
NS08	101	76	0.7	0.8	280	1 294.97
NS09	101	76	0.7	0.9	280	1 347.31

下载: 导出CSV

| 显示表格

泡沫铝夹芯管对于爆炸冲击波能量的吸收主要靠内外管的变形与泡沫铝芯层的坍塌塑性变形来实现，定义比吸能E_a为单位质量泡沫铝夹芯管所吸收的能量

$E_{\mathrm{a}}=\frac{E_{\mathrm{t}}}{m}$

(6)

式中：E_t为泡沫铝夹芯管吸收的总能量，在数值模拟中取2 ms时的塑性能量耗散值为结构整体的吸能，m为泡沫铝夹芯管的质量。

2.1 外管直径和壁厚的影响

为研究泡沫铝夹芯管结构中外管直径和壁厚对泡沫铝夹芯管抗爆性能的影响，对2组试样共计6个试件进行了横向爆炸载荷下的动态响应和能量吸收性能研究，其中：第一组试样NS03、NS04和NS05的内外管壁厚和内管直径保持不变，仅改变外管的直径；另一组试样NS01、NS06和NS07的内外管直径和内管壁厚保持不变，仅改变外管的壁厚。加载条件为：乳化炸药质量30 g，起爆位置距离试件跨中40 mm。表4给出了6个不同几何参数试件吸收的总能量、比吸能和内、外管跨中挠度（w_i和w_o）。由表4可知，随着外管直径的减小，试件的总吸能、比吸能和外管跨中挠度减小，内管跨中挠度增大。原因是随着外管直径的减小，泡沫铝夹芯管迎爆面积减小，结构承受的爆炸载荷降低，外管塑性变形减小，因而结构的总吸能减小；而内管跨中挠度增大是由于泡沫铝芯层厚度减小导致芯层密实化使内管吸收了较多的能量。随着外管壁厚的增大，试件整体结构刚度增加，内外管的塑性变形和芯层的压缩量减小，因此试件的总吸能、比吸能和内、外管跨中挠度均减小。

表 4 外管直径和壁厚的影响

Table 4. Effect of outer tube diameter and wall thickness

Test piece	E_t/J	m/g	E_a/(J·g⁻¹)	w_o/mm	w_i/mm
NS03	539.878	1 388.13	0.389	10.78	1.95
NS04	442.605	1 150.28	0.385	9.06	2.21
NS05	354.622	933.72	0.380	8.88	2.87
NS01	542.923	1 242.62	0.437	10.06	2.45
NS06	456.645	1 312.19	0.348	7.61	2.16
NS07	388.472	1 381.75	0.281	5.81	1.76

下载: 导出CSV

| 显示表格

图6给出了不同外管直径和壁厚的泡沫铝夹芯管试件内外管的跨中挠度时程曲线，其中虚线和实线分别代表内、外管跨中挠度。显然，外管跨中挠度大于内管跨中挠度。外管壁厚对外管跨中挠度的影响程度大于外管直径的影响；随着外管直径的增大及外管壁厚的减小，泡沫铝夹芯管芯层压缩量增大，结构吸收的总能量增加。

图 6 外管直径和壁厚对泡沫铝夹芯管内外管跨中挠度的影响

Figure 6. Effects of outer tube diameter and wall thickness on mid-span deflection of inner and outer tubes of aluminum foam-filled sandwich tubes

下载: 全尺寸图片幻灯片

2.2 内管直径和壁厚的影响

为研究泡沫铝夹芯管结构中内管直径和壁厚对泡沫铝夹芯管动态响应的影响，对2组试样共计6个试件进行了横向爆炸载荷下的力学性能研究，其中：第一组试样NS01、NS02和NS03的外管直径和内、外管壁厚保持不变，仅改变内管的直径；另一组试样NS01、NS08和NS09的内外管直径和外管壁厚保持不变，仅改变内管的壁厚。试件的加载条件为：乳化炸药质量30 g，起爆位置距离试件跨中40 mm。表5给出了试件的总吸能、比吸能及内、外管跨中挠度。由表5可知，随着内管直径的减小，试件的比吸能及内管跨中挠度减小，外管跨中挠度增大。内管直径的减小对泡沫铝夹芯管迎爆面积无影响，因而对结构吸收总能量的影响很小，芯层厚度的增大使内管吸能减少，内管跨中挠度减小。随着内管壁厚的增大，试件整体结构刚度增加，内外管的塑性变形和芯层的压缩量减小，因此试件的总吸能、比吸能和内、外管跨中挠度均减小。

表 5 内管直径和壁厚的影响

Table 5. Effect of inner tube diameter and wall thickness

Test piece	E_t/J	m/g	E_a/(J·g⁻¹)	w_o/mm	w_i/mm
NS01	542.923	1 242.62	0.437	10.06	2.45
NS02	546.692	1 326.02	0.412	10.39	2.32
NS03	539.878	1 388.13	0.389	10.78	1.95
NS08	542.470	1 294.97	0.419	9.92	2.16
NS09	542.161	1 347.31	0.402	9.81	1.92

下载: 导出CSV

| 显示表格

图7给出了不同内管直径和不同内管壁厚的泡沫铝夹芯管试件内外管跨中挠度时程曲线，其中虚线和实线分别代表内管和外管跨中挠度时程曲线，外管挠度大于内管挠度。由图7(b)可知，内管壁厚对外管跨中挠度的影响程度小于对内管跨中挠度的影响；对比图6，外管的直径和壁厚对内、外管跨中挠度的影响程度大于内管直径和壁厚的影响。

图 7 内管直径和壁厚对泡沫铝夹芯管内外管跨中挠度的影响

Figure 7. Effects of inner tube diameter and wall thickness on mid-span deflection of inner and outer tubes of aluminum foam-filled sandwich tubes

下载: 全尺寸图片幻灯片

2.3 芯层相对密度的影响

为研究泡沫铝夹芯管结构中泡沫铝夹芯层的相对密度对泡沫铝夹芯管动态性能的影响，对以NS01试件为基础分别填充相对密度为10%、15%和20%的泡沫铝芯层的3个试件进行横向爆炸载荷下的力学性能研究，试件的加载条件为：乳化炸药质量30 g，起爆位置距离试件跨中40 mm。表6给出了不同相对密度芯层的泡沫铝夹芯管结构的总能量、比吸能和内外管跨中挠度。由表6可知，随着泡沫铝芯层相对密度的增大，结构的整体刚度增大，内、外管的塑性变形和芯层压缩量减小。因此，结构的总吸能、比吸能和内、外管跨中挠度均减小。图8给出了不同相对密度芯层的泡沫铝夹芯管试件的内外管跨中挠度时程曲线，其中虚线和实线分别代表内、外管跨中挠度随时间的变化规律，外管挠度大于内管挠度。由图8可知，相对密度越小，泡沫铝芯层的压缩量越大，结构的总吸能和比吸能越大。

表 6 泡沫铝相对密度的影响

Table 6. Effect of relative density of aluminum foam

ρ^*/%	p₀/MPa	E_t/J	m/g	E_a/(J·g⁻¹)	w_o/mm	w_i/mm
10	36.8	614.301	1 116.12	0.550	14.77	4.99
15	36.8	542.923	1 242.62	0.437	10.06	2.45
20	36.8	379.921	1 378.86	0.276	4.82	2.35

下载: 导出CSV

| 显示表格

图 8 泡沫铝相对密度对泡沫铝夹芯管内外管跨中挠度的影响

Figure 8. Effect of relative density of aluminum foam on mid-span deflection of inner and outer tubes of aluminum foam-filled sandwich tubes

下载: 全尺寸图片幻灯片

图9给出了填充不同相对密度泡沫铝夹芯管内、外管的轴向变形轮廓，其中虚线和实线分别代表内、外管变形轮廓。由图9可知，内、外管的变形可分为3个区域，包括中心部位的塑性大变形区域、围绕塑性铰移动的刚性截面和未变形边界区域^[17]。

图 9 不同相对密度泡沫铝夹芯管的轴向变形轮廓

Figure 9. Axial deformation profiles of aluminum foam-filled sandwich tubes with different relative densities

下载: 全尺寸图片幻灯片

2.4 爆炸加载条件的影响

由2.3节可知，横向爆炸载荷下泡沫铝夹芯管NS01试件在所研究的试件中抗爆性能最佳，因此以NS01试件为研究对象，在距离试件跨中37、40和43 mm处施加炸药质量为30 g的爆炸载荷，研究炸药起爆距离对试件动态响应和吸能特性的影响。在距离跨中40 mm处分别施加乳化炸药质量为30、35和40 g的爆炸载荷，研究炸药质量对试件动态响应和吸能特性的影响。表7给出了不同爆炸加载条件下试件NS01的总能量、比吸能和内、外管跨中挠度。由表7可知，随着炸药起爆距离的减小和炸药质量的增加，泡沫铝夹芯管承受的爆炸载荷增大，内外管的塑性变形区域增大。因此，结构的总吸能、比吸能和内、外管跨中挠度均增大。

表 7 爆炸加载条件的影响

Table 7. Effect of blast loading conditions

R/mm	W_s/g	p₀/MPa	E_t/J	E_a/(J·g⁻¹)	w_o/mm	w_i/mm
37	30	34.9	935.594	0.753	14.96	4.34
40	30	27.7	542.923	0.437	10.06	2.45
40	35	32.3	776.974	0.625	12.60	3.62
40	40	36.8	1058.910	0.852	16.20	4.98
43	30	22.4	312.212	0.251	5.92	1.54

下载: 导出CSV

| 显示表格

图10给出了炸药起爆距离与炸药质量对泡沫铝夹芯管试件内外管跨中挠度的影响，其中虚线和实线分别代表内、外管跨中挠度时程曲线，外管的挠度大于内管的挠度。由图10可知，炸药起爆距离和炸药药量对泡沫铝夹芯管的内、外管跨中挠度影响显著。

图 10 炸药起爆距离和炸药质量对泡沫铝夹芯管内外管跨中挠度的影响

Figure 10. Effect of explosive initiation distance on mid-span deflection of inner and outer tubes of aluminum foam-filled sandwich tubes

下载: 全尺寸图片幻灯片

泡沫铝夹芯管各组成部分的吸能占比如表8所示。由表8可知，外管吸收的能量最多，吸能百分比均超过了50%，其次是泡沫铝芯层，吸能最少的是内管。原因是外管作为迎爆面，首先承受爆炸载荷，之后外管压缩芯层，芯层部分密实化后与外管共同压缩内管变形。随着外管壁厚的减小，外管的吸能百分比减小，泡沫铝芯层的吸能百分比增大；内管壁厚对各部分吸能百分比的影响较小；当夹芯管外管直径或内管直径一定时，随着芯层厚度的增加，内管的吸能百分比减小；芯层的相对密度对泡沫铝夹芯管的吸能百分比影响显著，随着芯层相对密度的减小，芯层和内管的吸能百分比减小，外管的吸能百分比增加。

表 8 泡沫铝夹芯管组成部分的吸能百分比

Table 8. Energy absorption percentage of aluminum foam-filled sandwich tube components

Test piece	$\rho ^*$ /%	R/mm	W_s/g	p₀/MPa	Energy absorption percentage/%
Test piece	$\rho ^*$ /%	R/mm	W_s/g	p₀/MPa	Outer tube	Foam core	Inner tube
NS01	15	40	30	27.7	59.9	32.0	8.1
NS02	15	40	30	27.7	62.6	30.5	6.9
NS03	15	40	30	27.7	61.8	33.7	4.5
NS04	15	40	30	27.7	62.4	30.7	6.9
NS05	15	40	30	27.7	59.7	32.3	8.0
NS06	15	40	30	27.7	64.0	28.1	7.9
NS07	15	40	30	27.7	66.3	25.6	8.2
NS08	15	40	30	27.7	59.5	32.4	8.0
NS09	15	40	30	27.7	59.2	32.9	7.9
NS01	10	40	30	27.7	79.8	17.5	2.7
NS01	20	40	30	27.7	52.0	32.2	15.8
NS01	15	37	30	34.9	61.7	31.0	7.3
NS01	15	43	30	22.4	61.5	29.6	8.9
NS01	15	40	35	32.3	62.0	30.5	7.5
NS01	15	40	40	36.8	61.5	31.2	7.3

下载: 导出CSV

| 显示表格

3. 抗爆性能的多目标优化

3.1 抗爆性能的定义及优化方法

由参数分析可知，内、外管壁厚和泡沫铝芯层的相对密度对泡沫铝夹芯管的抗爆性能具有一定影响，期望芯层比吸能最大化的同时，尽可能将外管跨中挠度控制在合理的范围内。因此，将最大化芯层比吸能E_CSA和最小化外管跨中挠度w_o作为优化目标，通过多目标优化方法研究泡沫铝夹芯管的抗爆性能。以泡沫铝夹芯管NS01内、外管壁厚以及芯层的相对密度作为设计参数（壁厚单位：mm），则多目标优化问题表示为

$\left\{ \begin{gathered} {{f}}\left( x \right){{ = {\rm{min}}}}\left( {{1 \mathord{\left/ {\vphantom {1 {{E_{{\rm{CSA}}}}(H,h,{\rho ^*})}}} \right. } {{E_{{\rm{CSA}}}}(H,h,{\rho ^*})}},{w_{\text{o}}}(H,h,{\rho ^*})} \right) \\ {\rm{s.t.}}\;{\text{ }}10{\text{%}} \leqslant {\rho ^*} \leqslant 20{\text{%}} {\text{ }} \\ \;\;\;\;0.7 \leqslant H \leqslant 1.78 \\ \;\;\;\;0.7 \leqslant h \leqslant 1.78 \\ \end{gathered} \right.$

(7)

多目标优化方案如图11所示。首先，采用正交试验设计方法，在设计变量的范围内选取样本点；然后，采用响应面法（RSM），根据有限元模拟得到样本点数据建立E_CSA和w_o的代理模型，并对代理模型进行评估；最后，基于非支配排序遗传算法（NSGA-Ⅱ）对参数进行优化，得到多目标优化问题的Pareto前沿解集。

图 11 优化方法流程

Figure 11. Flowchart of the optimization method

下载: 全尺寸图片幻灯片

采用正交试验设计的方法选取试验点^[20]，将外管与内管壁厚、泡沫铝芯层的相对密度作为3个影响因子，内、外管壁厚采用4个水平，内、外管壁厚水平分别采用0.70、1.06、1.42、1.78 mm，相对密度采用3个水平，相对密度水平采用10%、15%和20%，按照正交试验设计方法选取48个样本点进行数值计算，乳化炸药的质量为30 g，起爆位置距离试件跨中40 mm，分别得到不同样本点的外管跨中挠度和泡沫铝芯层的比吸能，并以此建立响应面模型。

响应面法具有简单易行、效果良好的特点，在建立高度非线性的代理模型中得到了广泛应用，使用三次多项式构建代理模型，得到48个样本点对应的E_CSA和w_o，通过数据拟合和最小二乘法确定响应面代理模型的系数。由设计变量H、h、 $\,\rho^*$ 构造的E_CSA和w_o三次多项式响应面函数为

$\begin{split} 1/{E_{\rm{CSA}}} =& - 6.477\,3+49.519\,7H - 23.364\,1h+122.008\,0{\rho^ {*} } - 8.864\,1{H^2} - 4.141\,7{h^2} - 481.066\,6{\rho ^ {*2} }+18.808\,4Hh - \\ &739.643\,0H{\rho ^ {*} }+301.309\,9h{\rho ^ {*} }+0.389\,1{H^3} - 0.228\,9{h^3} - 224.429\,9{\rho ^ {*3} }- 4.478\,3{H^2}h+116.427\,8{H^2}{\rho ^ {*} } + \\ & 2.412\,4H{h^2}+21.482\,4{h^2}{\rho^ {*} }+2\,533.473\,6H{\rho ^ {*2} }- 847.899\,4h{\rho ^ {*2} } - 127.966\,0Hh{\rho ^ * } \\[-15pt]\end{split}$

(8)

$\begin{split} {w_{\text{o}}} = & 74.447\,9 - 95.263\,9H - 19.720\,6h - 162.469\,6{\rho ^ * }+44.208\,1{H^2}+3.002\,8{h^2} - 324.652\,3{\rho ^ {*2} }+11.999\,2Hh+ \\ & 151.621\,9H{\rho ^ * }+97.555\,6h{\rho ^ * } - 6.028\,2{H^3} - 0.044\,7{h^3} - 135.309\,7{\rho ^ {*3} }- 2.002\,2{H^2}h - 75.762\,0{H^2}{\rho ^ * } - \\ &1.305\,7H{h^2} - 4.677\,9{h^2}{\rho ^ * }+430.972\,2H{\rho ^ {*2} } - 201.805\,6h{\rho ^ {*2} } - 18.248\,5Hh{\rho ^ * } \end{split}$

(9)

为了评估代理模型的准确性，采用拟合优度R²和相对误差的绝对值 $\left| \delta \right|$ 两种评价指标对多项式响应面进行误差分析，这两种指标的计算方法为

${R^2} = 1 - \frac{{\displaystyle\sum\limits_{i=1}^n {{{({y_i} - {{\hat y}_i})}^2}} }}{{\displaystyle\sum\limits_{i=1}^n {{{({y_i} - \bar y)}^2}} }}$

(10)

$\left| \delta \right| = \frac{{\displaystyle\sum\limits_{i=1}^n {\left| {{y_i} - {{\hat y}_i}} \right|} }}{{\displaystyle\sum\limits_{i=1}^n {\left| {{y_i} - \bar y} \right|} }}$

(11)

式中： ${y_i}$ 和 ${\hat y_i}$ 分别为第i个样本点的有限元模拟和多项式响应面的响应值； $\bar y$ 为 ${y_i}$ 的平均值；n为样本点的个数；R²为评价响应面函数与数值计算数据拟合准确度的指标， $\left| \delta \right|$ 用于评估模型的整体精度， ${R^2} > 0.9$ 、 $\left| \delta \right|$ <0.2表示代理模型有足够高的精确度。结果显示：E_CSA和w_o的R²分别为0.992 6和0.993 8，均大于0.99，接近1，E_CSA和w_o的 $\left| \delta \right|$ 分别为0.096 9和0.084 7，均小于0.1，表明建立的多项式代理模型精度满足要求。

3.2 优化结果分析

最终得到的Pareto前沿如图12所示，空心点为NSGA-Ⅱ得到的数据。当E_CSA增大时，w_o也增大，两个优化目标相互冲突。由图12可知，所有的解集中w_o都小于初始设计点。表9给出了设计点的结构参数，并将有限元模拟结果与响应面函数得到的结果进行比较，误差均小于9.5%，验证了代理模型的准确性。当仅要求芯层比吸能最大化时，最佳参数设计点为Op1，与初始设计点相比，外管跨中挠度减小了28.9%，芯层比吸能提高了11.7%。当仅要求外管跨中挠度最小化时，最佳参数设计点为Op2，与初始设计点相比，外管跨中挠度减小了95.7%，芯层比吸能降低了87.5%。为使E_CSA和w_o达到平衡，采用最小距离选择法设计Op3，最小距离选择法的公式为

图 12 Pareto前沿

Figure 12. Pareto front

下载: 全尺寸图片幻灯片

表 9 最佳设计点

Table 9. Optimum design point

Design point	H/mm	h/mm	$\,\rho ^ *$ /%	w_o			E_CSA
Design point	H/mm	h/mm	$\,\rho ^ *$ /%	FEM/mm	RSM/mm	Error/%	FEM/(J·g⁻¹)	RSM/(J·g⁻¹)	Error/%
Op1	0.70	0.76	13.7	0.457	0.505	−9.5	10.50	11.14	−5.7
Op2	1.78	1.75	18.8	0.051	0.056	−8.9	0.64	0.59	8.5
Op3	0.84	1.76	20.0	0.194	0.205	−5.4	2.85	2.83	0.7

下载: 导出CSV

| 显示表格

$\min S = \sqrt {{{\left[ {\frac{{\left( {1/{E_{\rm{CSA}}}\left( {H,h,{\rho ^*}} \right)} \right) - {{\left( {1/{E_{\rm{CSA}}}\left( {H,h,{\rho ^*}} \right)} \right)}_{\min }}}}{{{{\left( {1/{E_{\rm{CSA}}}\left( {H,h,{\rho ^*}} \right)} \right)}_{\max }} - {{\left( {1/{E_{\rm{CSA}}}\left( {H,h,{\rho ^*}} \right)} \right)}_{\min }}}}} \right]}^2}+{{\left[ {\frac{{{w_{\rm{o}}}\left( {H,h,{\rho ^*}} \right) - {w_{\rm{o}}}{{\left( {H,h,{\rho ^*}} \right)}_{\min }}}}{{{w_{\rm{o}}}{{\left( {H,h,{\rho ^*}} \right)}_{\max }} - {w_{\rm{o}}}{{\left( {H,h,{\rho ^*}} \right)}_{\min }}}}} \right]}^2}}$

(12)

式中：S为Pareto解集中每个非劣解到理想最优解（ideal point）的距离， $\left( {{1 \mathord{\left/ {\vphantom {1 {{E_{\rm{CSA}}}\left( {H,h,{\rho ^*}} \right)}}} \right. } {{E_{{\rm{CSA}}}}\left( {H,h,{\rho ^*}} \right)}}} \right)$ 和 ${w_{\rm{o}}}\left( {H,h,{\rho ^*}} \right)$ 为Pareto解集中每个非劣解对应芯层比吸能的倒数和外管跨中挠度的值， ${\left( {{1 \mathord{\left/ {\vphantom {1 {{E_{{\rm{CSA}}}}\left( {H,h,{\rho ^*}} \right)}}} \right. } {{E_{{\rm{CSA}}}}\left( {H,h,{\rho ^*}} \right)}}} \right)_{\max }}$ 和 ${\left( {{1 \mathord{\left/ {\vphantom {1 {{E_{{\rm{CSA}}}}\left( {H,h,{\rho ^*}} \right)}}} \right. } {{E_{{\rm{CSA}}}}\left( {H,h,{\rho ^*}} \right)}}} \right)_{\min }}$ 分别为所有Pareto解集中芯层比吸能的倒数的最大值和最小值， ${w_{\rm{o}}}{\left( {H,h,{\rho ^*}} \right)_{\max }}$ 和 ${w_{\rm{o}}}{\left( {H,h,{\rho ^*}} \right)_{\min }}$ 分别为所有Pareto解集中外管跨中挠度的最大值和最小值。最佳参数设计点为Op3，与初始设计点相比，外管跨中挠度减小了80.7%，芯层比吸能降低了52.6%。

4. 结　论

采用数值模拟方法研究了泡沫铝夹芯管在横向爆炸载荷下的动态响应，对比分析了泡沫铝夹芯管的几何参数、芯层的相对密度和爆炸加载条件对泡沫铝夹芯管的总吸能、比吸能和内、外管跨中挠度的影响。以泡沫铝夹芯管内管壁厚、外管壁厚和芯层相对密度为设计变量进行多目标优化，以提高泡沫铝夹芯管结构的抗爆性能，得到以下主要结论。

(1) 泡沫铝夹芯管结构的几何参数和芯层的相对密度对其抗爆性能存在显著影响，相对于内管的直径和壁厚，外管直径和壁厚对内外管跨中挠度的影响更显著。随着外管直径的增大和内管直径的减小，泡沫铝夹芯管内管跨中挠度减小，外管跨中挠度增大，比吸能随内外管直径的减小而减小；随着内管和外管壁厚的增大，泡沫铝夹芯管的比吸能和内外管跨中挠度减小。泡沫铝芯层的相对密度越大，泡沫铝夹芯管的比吸能和内、外管跨中挠度越小。

(2) 随着炸药质量的增加和起爆距离的减小，泡沫铝夹芯管的变形区域增大，比吸能和内、外管跨中挠度增大。分析了不同工况下泡沫铝夹芯管各组成部分的吸能情况，其中外管吸收的能量最多，均超过了50%，其次是泡沫铝芯层，吸能最少的是内管；内管壁厚对各部分吸能百分比的影响较小；芯层的厚度和相对密度对泡沫铝夹芯管的吸能百分比影响显著。

(3) 以内管壁厚、外管壁厚和芯层相对密度为设计变量，以芯层的比吸能和外管跨中挠度为优化指标，采用响应面法和非支配排序遗传算法可以实现泡沫铝夹芯管抗爆性多目标优化设计，得到了Pareto前沿图，发现外管跨中挠度取最小值时，对应的芯层比吸能最小。因此，当确定了外管跨中挠度和芯层比吸能的偏好时，可确定泡沫铝夹芯管结构的内外管壁厚和芯层的相对密度。

图 1 I、III型混合模式下锯齿状裂纹的实验结果（a）和相场模拟结果（b）（c）^[21]

Figure 1. Experimental result (a) and phase field simulation (b), (c) of jag fracture under mixed condition of type I & III^[21]

下载: 全尺寸图片幻灯片

图 2 静态/准静态脆性断裂的相场模拟：（a）I、II型裂纹^[27]，（b）L形薄板中的裂纹^[29]

Figure 2. Phase field simulation of static/quasi-static brittle fracture: (a) type I & II fracture^[27], (b) fracture in L-shape sheet^[29]

下载: 全尺寸图片幻灯片

图 3 动态脆性断裂的相场模拟：（a）（b）动态裂纹经典问题^[12]，（c）三维形式的KW试验模拟结果^[32]

Figure 3. Phase field simulation of dynamic brittle fracture: (a), (b) classical problems in dynamic fracture^[12]; (c) 3D KW test^[32]

下载: 全尺寸图片幻灯片

图 4 相场断裂方法的基本框架

Figure 4. Framework of phase field approach to fracture

下载: 全尺寸图片幻灯片

图 5 带裂纹的电场分布^[37]

Figure 5. Electronic field with fracture^[37]

下载: 全尺寸图片幻灯片

图 6 带裂纹场的极化矢量场^[38]

Figure 6. Polarization field with fracture^[38]

下载: 全尺寸图片幻灯片

图 7 水坝（Koyna dam）的水力压裂^[41]

Figure 7. Phase field simulation of hydraulic fracture in Koyna dam^[41]

下载: 全尺寸图片幻灯片

图 8 Ambati等延性断裂相场模拟结果^[44]

Figure 8. Ductile fracture simulated by Ambati et al.^[44]

下载: 全尺寸图片幻灯片

图 9 脆性-延性断裂模型的能量释放率：（a）应变率决定的临界能量释放率转变^[36]，（b）应力三轴度决定的临界能量释放率转变^[46]

Figure 9. Energy release rate in fracture mode for brittle-ductile fracture transition: (a) strain rate induced energy release rate transition^[36], (b) stress triaxiality induced energy release rate transition^[46]

下载: 全尺寸图片幻灯片

图 10 KW试验相场模拟结果^[36]：（a）冲击速度20 m/s产生的脆性断裂，（b）冲击速度39 m/s产生的绝热剪切失效

Figure 10. Phase field simulation to KW test^[36]: (a) brittle fracture induced by impact velocity of 20 m/s, (b) ductile fracture induced by impact velocity of 39 m/s

下载: 全尺寸图片幻灯片

表 1 裂纹扩展的数值模拟方法对比

Table 1. Comparison of methods in fracture propagation simulation

Methods Type I tensile fracture KW test

Experimental result

Element deletion

Inter element

XFEM

Phase field

下载: 导出CSV

表 2 扩展有限元及相场方法处理三维及多裂纹的对比

Table 2. Comparison of 3D and multi-fracture simulation of XFEM and phase field approach

Crack XFEM Phase field

Level sets Fast marching No need crack tracking

3D crack

Multi-crack

下载: 导出CSV

参考文献(46)

[1]	钱学森. 物理力学讲义[M]. 上海: 上海交通大学出版社, 2007.
[2]	BECKER R. How metals fail [R]. Metal Fracture-Lawrence Livermore Nation Laboratory, 2002: 13–30.
[3]	FAN R, FISH J. The rs-method for material failure simulations [J]. International Journal for Numerical Methods in Engineering, 2008, 73(11): 1607–1623. doi: 10.1002/nme.2134
[4]	HALLQUIST J O. LS-DYNA theoretical manual [M]. USA: Livemore Software Technology Corporation, 1998.
[5]	XU X P, NEEDLEMAN A. Numerical simulations of fast crack growth in brittle solids [J]. Journal of the Mechanics and Physics of Solids, 1994, 42(9): 1397–1434. doi: 10.1016/0022-5096(94)90003-5
[6]	CAMACHO G T, ORTIZ M. Computational modelling of impact damage in brittle materials [J]. International Journal of Solids and Structures, 1996, 33(20/21/22): 2899–2938.
[7]	BELYTSCHKO T, BLACK T. Elastic crack growth in finite elements with minimal remeshing [J]. International Journal for Numerical Methods in Engineering, 1999, 45(5): 601–620. doi: 10.1002/(sici)1097-0207(19990620)45:5<601::aid-nme598>3.0.co;2-s
[8]	MOËS N, DOLBOW J, BELYTSCHKO T. A finite element method for crack growth without remeshing [J]. International Journal for Numerical Methods in Engineering, 1999, 46(1): 131–150. doi: 10.1002/(SICI)1097-0207(19990910)46:13.0.CO;2-J
[9]	庄茁, 成斌斌. 发展基于CB壳单元的扩展有限元模拟三维任意扩展裂纹 [J]. 工程力学, 2012, 29(6): 12–21. doi: 10.6052/j.issn.1000-4750.2010.08.0616 ZHUANG Z, CHENG B B. Development of X-FEM on CB shell element for simulating 3D arbitrary crack growth [J]. Engineering Mechanics, 2012, 29(6): 12–21. doi: 10.6052/j.issn.1000-4750.2010.08.0616
[10]	PROVATAS N, ELDER K. Book-phase-field methods in materials science and engineering [M]. Wiley-VCH Press, 2010.
[11]	SONG J H, WANG H, BELYTSCHKO T. A comparative study on finite element methods for dynamic fracture [J]. Computational Mechanics, 2008, 42(2): 239–250. doi: 10.1007/s00466-007-0210-x
[12]	BORDEN M J, VERHOOSEL C V, SCOTT M A, et al. A phase-field description of dynamic brittle fracture [J]. Computer Methods in Applied Mechanics and Engineering, 2012, 217(220): 77–95.
[13]	RAMULU M, KOBAYASHI A S. Mechanics of crack curving and branching-a dynamic fracture analysis [J]. International Journal of Fracture, 1985, 27(3/4): 187–201.
[14]	KALTHOFF J F, WINKLER S. In failure mode transition at high rates of shear loading [C]//Impact Loading and Dynamic Behavioavior of Materials. FRG, 1987: 185–195.
[15]	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. doi: 10.1016/S0022-5096(99)00028-9
[16]	GIACOMINI A, PONSIGLIONE M. A Γ-convergence approach to stability of unilateral minimality properties in fracture mechanics and applications [J]. Archive for Rational Mechanics and Analysis, 2006, 180(3): 399–447. doi: 10.1007/s00205-005-0392-3
[17]	GRAVOUIL A, MOËS N, BELYTSCHKO T. Non-planar 3D crack growth by the extended finite element and level sets-Part II: level set update [J]. International Journal for Numerical Methods in Engineering, 2002, 53(11): 2569–2586. doi: 10.1002/nme.v53:11
[18]	SUKUMAR N, CHOPP D L, BÉCHET E, et al. Three-dimensional non-planar crack growth by a coupled extended finite element and fast marching method [J]. International Journal for Numerical Methods in Engineering, 2008, 76(5): 727–748. doi: 10.1002/nme.v76:5
[19]	MIEHE C, SCHÄNZEL L M, ULMER H. Phase field modeling of fracture in multi-physics problems. Part I. balance of crack surface and failure criteria for brittle crack propagation in thermo-elastic solids [J]. Computer Methods in Applied Mechanics and Engineering, 2015, 294(1): 449–485.
[20]	HAKIM V, KARMA A. Laws of crack motion and phase-field models of fracture [J]. Journal of the Mechanics and Physics of Solids, 2009, 57(2): 342–368. doi: 10.1016/j.jmps.2008.10.012
[21]	PONS A J, KARMA A. Helical crack-front instability in mixed-mode fracture [J]. Nature, 2010, 464(7285): 85–89. doi: 10.1038/nature08862
[22]	FRANCFORT G A, MARIGO J J. Revisiting brittle fracture as an energy minimization problem [J]. Journal of the Mechanics and Physics of Solids, 1998, 46(8): 1319–1342. doi: 10.1016/S0022-5096(98)00034-9
[23]	SPATSCHEK R, BRENER E, KARMA A. Phase field modeling of crack propagation [J]. Philosophical Magazine, 2011, 91(1): 75–95. doi: 10.1080/14786431003773015
[24]	GRIFFITH A A. The phenomena of rupture and flow in solid [J]. Philosophical Transactions of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 1920, 221: 163–198.
[25]	Bourdin B, Francfort G A, Marigo J J. The variational approach to fracture[M]. Springer Science+Business Media B. V., 2008.
[26]	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. doi: 10.1016/j.jmps.2009.04.011
[27]	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. doi: 10.1002/nme.v83:10
[28]	MIEHE C, HOFACKER M, WELSCHINGER F. A phase field model for rate-independent crack propagation: robust algorithmic implementation based on operator splits [J]. Computer Methods in Applied Mechanics and Engineering, 2010, 199: 2765–2778. doi: 10.1016/j.cma.2010.04.011
[29]	MESGARNEJAD A, BOURDIN B, KHONSARI M M. Validation simulations for the variational approach to fracture [J]. Computer Methods in Applied Mechanics and Engineering, 2015, 290(1): 420–437.
[30]	BOURDIN B, LARSEN C J, RICHARDSON C L. A time-discrete model for dynamic fracture based on crack regularization [J]. International Journal of Fracture, 2011, 168(2): 133–143. doi: 10.1007/s10704-010-9562-x
[31]	LARSEN C J, ORTNER C, SÜLI E. Existence of solutions to a regularized model of dynamic fracture [J]. Mathematical Models and Methods in Applied Sciences, 2010, 20(7): 1021–1048. doi: 10.1142/S0218202510004520
[32]	HOFACKER M, MIEHE C. Continuum phase field modeling of dynamic fracture: variational principles and staggered FE implementation [J]. International Journal of Fracture, 2012, 178(1): 113–129.
[33]	HOFACKER M, MIEHE C. A phase field model of dynamic fracture: robust field updates for the analysis of complex crack patterns [J]. International Journal for Numerical Methods in Engineering, 2013, 93(3): 276–301. doi: 10.1002/nme.v93.3
[34]	SCHLÜTER A, WILLENBÜCHER A, KUHN C, et al. Phase field approximation of dynamic brittle fracture [J]. Computational Mechanics, 2014, 54(5): 1141–1161. doi: 10.1007/s00466-014-1045-x
[35]	MIEHE C, MAUTHE S. Phase field modeling of fracture in multi-physics problems. Part III. crack driving forces in hydro-poro-elasticity and hydraulic fracturing of fluid-saturated porous media [J]. Computer Methods in Applied Mechanics and Engineering, 2016, 304(1): 619–655.
[36]	MIEHE C, HOFACKER M, SCHÄNZEL L M, et al. Phase field modeling of fracture in multi-physics problems. Part II. coupled brittle-to-ductile failure criteria and crack propagation in thermo-elastic-plastic solids [J]. Computer Methods in Applied Mechanics and Engineering, 2015, 294(1): 486–522.
[37]	WILSON Z A, BORDEN M J, LANDIS C M. A phase-field model for fracture in piezoelectric ceramics [J]. International Journal of Fracture, 2013, 183(2): 135–153. doi: 10.1007/s10704-013-9881-9
[38]	ABDOLLAHI A, ARIAS I. Phase-field modeling of fracture in ferroelectric materials [J]. Archives of Computational Methods in Engineering, 2015, 22(2): 153–181. doi: 10.1007/s11831-014-9118-8
[39]	MIEHE C, SCHÄNZEL L M. Phase field modeling of fracture in rubbery polymers. Part I: finite elasticity coupled with brittle failure [J]. Journal of the Mechanics and Physics of Solids, 2014, 64: 93–113.
[40]	ZIAEI-RAD V, SHEN L, JIANG J, et al. Identifying the crack path for the phase field approach to fracture with non-maximum suppression [J]. Computer Methods in Applied Mechanics and Engineering, 2016, 312(1): 304–321.
[41]	SANTILLÁN D, MOSQUERA J C, CUETO-FELGUEROSO L. Phase-field model for brittle fracture. validation with experimental results and extension to dam engineering problems [J]. Engineering Fracture Mechanics, 2017, 178: 109–125. doi: 10.1016/j.engfracmech.2017.04.020
[42]	DUDA F P, CIARBONETTI A, SÁNCHEZ P J, et al. A phase-field/gradient damage model for brittle fracture in elastic-plastic solids [J]. International Journal of Plasticity, 2015, 65: 269–296. doi: 10.1016/j.ijplas.2014.09.005
[43]	AMBATI M, GERASIMOV T, LORENZIS L D. Phase-field modeling of ductile fracture [J]. Computational Mechanics, 2015, 55(5): 1017–1040. doi: 10.1007/s00466-015-1151-4
[44]	AMBATI M, KRUSE R, LORENZIS L D. A phase-field model for ductile fracture at finite strains and its experimental verification [J]. Computational Mechanics, 2016, 57: 149–167. doi: 10.1007/s00466-015-1225-3
[45]	MCAULIFFE C, WAISMAN H. A coupled phase field shear band model for ductile-brittle transition in notched plate impacts [J]. Computer Methods in Applied Mechanics and Engineering, 2016, 305: 173–195. doi: 10.1016/j.cma.2016.02.018
[46]	柳占立, 李想, 初东阳, 等. 多物理场耦合断裂的相场方法模拟及工程应用[C]//第十二届全国爆炸力学会议. 桐乡, 2018.