LiAlH<sub>4</sub>晶体结构及稳定性的第一性原理研究

张艺龙; 崔慢爱; 刘艳辉

doi:10.11858/gywlxb.20170561

Abstract: The structural stability of LiAlH₄, a promising hydrogen storage material, under high pressure was researched using the ab initio pseudopotential plane wave method.It is found that the phase transition occurs at 1.6 GPa from the α-LiAlH₄ phase to the β-LiAlH₄ (space group I2/b) phase.This phase transition is identified as first-order in nature with volume contractions of 18%.Moreover, the analysis of the phonon dispersion curves suggests that phase transition is related to the phonon softening.Mulliken population analyses indicated that the ambient phase (α-LiAlH₄) is expected to be the most promising candidate for hydrogen storage.

Key words:

crystal structure /
hydrogen storage materials /
density functional theory /
pressure-induced structural transition /
electronic structure

摘要: 基于密度泛函理论的第一性原理赝势平面波方法，研究高压下三元碱金属氢化物LiAlH₄的相变行为，分析了LiAlH₄高压相变的物理机制。研究表明，在1.6 GPa时LiAlH₄发生了相变，从α-LiAlH₄转变为空间群为I2/b的β-LiAlH₄，相变时伴随18%的体积坍塌，即一级相变。通过分析声子色散曲线得出，相变与声子软化有关。Millikan布局分析表明，常压相(α-LiAlH₄)是很有潜力的储氢材料。

关键词:

随着煤炭开采深度逐渐增加，井下工程面临诸多问题，如爆破、采掘活动扰动诱发矿柱失稳、煤层顶板垮落、坍塌等矿井灾害，这些工程问题不仅是单一煤体或岩体破坏失稳的结果，而且与煤岩层复合整体结构的力学行为密切相关^[1-3]。此外，煤层并非单一赋存于地层中，而是与顶底板岩体共存。采矿工程中常见煤岩层复合结构，如井下半煤岩巷道的掘进、采煤工作面推进过程中顶底板与煤层所形成的受力系统、煤柱与顶板岩层形成新的围岩承载结构等。煤层和岩层在开采过程中共同承受载荷，使得煤岩层复合结构体的力学性质与单一煤体或岩体有很大区别。在这种情况下，将地层视为不同黏聚力和不同黏结方式组成的复合体更合理，即煤岩组合体。图1为常见煤岩层复合结构工程体示意图。受高地应力及采掘扰动的影响，煤岩组合体难免受到不同应力波（动荷载）的作用，组合体结构的复杂性导致应力波传播规律、破坏模式与煤、岩单体存在差异，因此研究煤岩组合体的动态响应特征对其工程结构体系统的稳定性十分重要。

图 1 常见煤岩层复合结构工程体示意图

Figure 1. Schematic diagram of common coal-rock composite structure

下载: 全尺寸图片幻灯片

国内外学者通过室内试验、数值计算、理论分析等方法，对煤岩层状复合结构体的力学特性开展了深入而全面的研究。在试验方面，左建平等^[4]发现单轴压缩荷载下煤岩组合体的强度略高于单一煤体，但低于单一岩体的强度，组合体以煤体破坏为主，整体呈拉伸劈裂破坏模式。苗磊刚等^[5]通过霍普金森压杆（split Hopkinson pressure bar，SHPB）试验分析了不同冲击速度下煤岩组合体的动态力学行为，结果表明，组合体的变形具有明显的弹性后效特性，且试样的分形维数与冲击速度近似呈线性关系。此外，煤岩组合体的力学行为还与其组合方式^[6-7]、组合角度^[8-10]、煤厚比^[11-13]、强度比^[14-15]等其他因素相关。理论分析方面，Liu等^[6]以不同的煤岩组合方式为研究对象，分析了静载下组合体中煤的力学性质，并提出了计算组合试样中煤体应力和应变的方法，建立的“串联损伤-牛顿体”本构模型可以很好地反映煤岩组合体的静态力学行为。解北京等^[13]对不同煤岩厚度比的组合试样进行了冲击动力学试验，测得了其动态应力-应变曲线，结合层叠模型理论，构建了多参数煤岩组合体的动态力学本构模型，该模型可以很好地反映煤岩组合体的动态力学行为及破坏特性。在数值模拟方面，杨磊等^[15]通过颗粒流离散元PFC^2D程序研究了不同强度比的组合煤岩的能量演化规律，并从细观力学角度分析了组合体的应力-应变关系及微观破坏过程。Zhang等^[16]借助UDEC离散元数值分析软件，定量分析了煤岩组合体破坏过程中裂纹的数量、长度及宏观破坏面积随时间变化的关系，捕捉到了组合试样裂纹的萌生、扩展、聚结、成核，直至界面贯穿破坏的全过程。

综上所述，尽管许多学者已经对煤岩组合体的静、动态力学行为特征开展了大量的研究，但是对煤岩组合体力学性质的研究主要集中在压缩荷载方面，而对煤岩组合体抗拉性能的研究甚少，对煤岩工程结构体动态破坏的关键影响因素认识不够全面。同时，抗拉性能对于煤岩层复合工程体稳定性控制、冲击地压灾害机理研究等均具有重要意义。此外，在地下工程中，煤岩体会受到不同强度的应力波扰动作用，且受到破坏冲击的方向与煤岩体的赋存状况和工程需要相关。因此，考虑煤岩组合体的冲击速度、冲击方向等因素更加符合工程实际。本研究以室内力学试验测得的煤、岩基本力学参数为依据，确定单一煤、岩体材料的Holmquist-Johnson-Cook（HJC）模型参数，并验证材料模型的有效性，然后通过显式动力学分析软件LS-DYNA研究不同的冲击荷载、冲击方向和加载角度下煤岩组合体的破坏过程及破坏特征。

1. SHPB数值模型及验证

1.1 SHPB数值模型

根据SHPB实验装置的实际尺寸，基于Hypermesh前处理软件，建立SHPB有限元模型，如图2所示。其中，单元类型选择规则六面体8节点Solid164实体单元。数值模型中入射杆和透射杆均选用理想线弹性材料，其主要参数为：密度7850 kg/m³，弹性模量210 GPa，泊松比0.25^[17]。数值试样由煤、岩两部分半圆盘组合而成，材料选用HJC本构模型。杆件与试件的接触类型选择自动面面接触，不计接触面间的摩擦。选取入射杆和透射杆中点处的两个单元作为监测点记录冲击过程中的应力波数据。

图 2 SHPB数值模型及试样分组示意图

Figure 2. SHPB numerical model and sample grouping diagram

下载: 全尺寸图片幻灯片

SHPB数值模型中入射杆的长度为2 000 mm，径向网格等分为80份，轴向网格等分为200份；透射杆的长度为1500 mm，径向网格等分为80份，轴向网格等分为150份。煤岩组合体试样为直径50 mm，厚度25 mm 的圆盘，为了更精确地分析试样的冲击破坏过程，试样的网格划分更精细。其中试样沿径向划分200份网格，轴向划分50份网格。通过旋转数值试样，以达到不同冲击方向和不同加载角度的目的。为了简化建模过程、提高计算效率，以入射杆端部直接施加入射波的方式替代子弹冲击，以约束透射杆末端的方式等效替代吸能杆的作用，并约束入射杆和透射杆的转动自由度，以达到理想状态下的一维应力条件。采用后处理软件LS-PrePost提取应力波曲线等数据，获取不同冲击荷载、不同冲击方向和加载角度条件下煤岩组合体动态劈裂过程中的应力波传播特征、试样破坏过程及破坏特征。

1.2 HJC材料模型

HJC本构方程既考虑了材料的应变率效应，又考虑了损伤效应、孔隙压实和块体压碎效应等因素的影响，能够很好地反映中高应变率、大变形和高静水压力条件下岩石/混凝土类材料的力学行为，因此被LS-DYNA有限元软件引入，广泛应用于冲击、爆炸和侵彻等强动载工程问题，符合煤岩组合体SHPB冲击劈裂数值分析条件^[18-20]。此外，HJC材料模型的形式简单，相关参数的物理意义明确，主要包括基本力学参数、状态方程相关参数和损伤因子等21个参数。材料模型的压力p与体积应变µ用3个不同阶段（线弹性阶段、塑性变形阶段和完全压密阶段）的方程式描述，如图3(a)所示。

图 3 HJC材料模型本构方程描述

Figure 3. Constitutive equation of HJC model

下载: 全尺寸图片幻灯片

(1) 线弹性加卸载阶段（p<p_c）

$p = {K\rm_e}\mu ,\;\;\;\;\;\;{K\rm_e} = {p\rm_c}/{\mu \rm_c}$

(1)

式中：K_e为体积模量，p_c和μ_c分别为材料压碎时的应力和应变。

(2) 塑性变形阶段（p_c≤ p ≤ p_l）

$\text{加载阶段}\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad p = {p\rm_l} + \left[ {\left( {{p\rm_l} - {p\rm_c}} \right) + \left( {\mu - {\mu \rm_c}} \right)} \right]/\left( {{\mu \rm_l} - {\mu \rm_c}} \right) \quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad$

(2)

$\text{卸载阶段} \quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad p - {p_{\max }} = \left[ {\left( {1 - F} \right){K\rm_e} + F{K\rm_l}} \right]\left( {\mu - {\mu _{\max }}} \right)\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad$

(3)

$F = \left( {{\mu _{\max }} - {\mu \rm_c}} \right)/\left( {{\mu \rm_l} - {\mu \rm_c}} \right)$

(4)

式中：p_l为压实压力，μ_l为压实应变，K_l为塑性体积模量，p_max和μ_max分别为岩石类材料在应力卸载前达到的最大应力和应变，这个阶段材料内部的微小孔隙被逐渐压密压实，试样开始起裂破坏。

(3) 完全密实段（p > p_l）

$\text{加载阶段} \quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad p = {K_1}\overline \mu + {K_2}{\overline \mu ^2} + {K_3}{\overline \mu ^3} \quad\quad\quad\quad\quad\quad\quad\;\quad\quad\quad\quad\quad\quad\quad\quad\quad$

(5)

$\text{卸载阶段} \quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad p - {p_{\max }} = {K_1}\left( {\overline \mu - {{\overline \mu }_{\max }}} \right)\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad$

(6)

$\overline \mu {{ = }}\left( {\mu {{ - }}{\mu _1}} \right)/\left( {1 + {\mu _1}} \right)$

(7)

式中：K₁、K₂和K₃为常数， $\overline {\mu}$ 为修正的体积应变。该阶段材料完全破裂，完全压密。其本构方程的等效强度充分考虑了材料应变率、损伤等对屈服应力的影响。

模型屈服面方程如图3(b)所示，其曲线通过无量纲等效应力描述，表达式为

${\sigma ^ * }{\text{ = }}\left[ {A\left( {1 - D} \right) + B{p^{ * N}}} \right]\left( {1 + C\ln {{\dot \varepsilon }^ * }} \right)$

(8)

式中： $\sigma ^ *$ 、 $p ^ *$ 分别为无量纲的等效应力和等效静水压力； ${\dot \varepsilon ^*} = \dot \varepsilon /{\dot \varepsilon _0}$ 为等效应变率， $\dot \varepsilon$ 为应变率， $\dot \varepsilon _0 = 1.0\;{{\rm{s}}^{ - 1}}$ 为参考应变率； $A$ 为特征化黏聚力参数；B为特征化压力硬化系数；C为应变率效应系数；N为特征化压力硬化指数；D为损伤因子。

材料损伤效应以塑性应变的累积过程来表示，其损伤模型见图3(c)，损伤效应演化方程为

$D = \sum {({{\Delta {\varepsilon \rm_p} + \Delta {\mu \rm_p}}})/({{\varepsilon \rm_p^f + \mu \rm_p^f}})}$

(9)

${\varepsilon \rm_p^f + \mu \rm_p^f }= {D_1}{\left( {{p^ * } + {T^ * }} \right)^{{D_2}}} \geqslant {\varepsilon _{\rm f,\min }}$

(10)

式中：Δε_p和Δμ_p分别为一个积分循环内材料的等效塑性应变和塑性体积应变，T^*为材料的最大特征化等效拉应力，D₁ 和D₂为模型的损伤因子，ε_f,min为材料断裂时的最小塑性应变， $\varepsilon _{\rm{p}}^{\rm{f}} + \mu _{\rm{p}}^{\rm{f}}$ 为损伤塑性应变。

通过室内静态和动态基本力学试验，对单一煤体和岩体的材料参数进行标定，其中HJC材料模型的大部分参数可根据试验及理论直接获得，其余参数取经验值，并通过试错法验证模型，以确定最终参数，如表1和表2^[21-23]所示，其中：ρ为密度，G为剪切模量，f_c为材料静态单轴抗压强度，T为最大拉应力，S_max为特征化等效应力所能达到的最大值。

表 1 煤的HJC模型参数

Table 1. Parameters of HJC model for coal

ρ/(kg·m^–3)	G/GPa	f_c/MPa	A	B	C	N	S_max	D₁	D₂
1594	0.58	9.0	0.400	0.700	0.005	0.500	7.000	0.031	1.00
T/MPa	p_c/MPa	μ_c/10^–3	p_l/GPa	μ_l	K₁/GPa	K₂/GPa	K₃/GPa	$\dot \varepsilon$ ₀/s⁻¹	${\varepsilon \rm_{f,\min } }$
1.86	3.0	0.8	1.0	0.121	81	−170	208	1.00	0.005

下载: 导出CSV

| 显示表格

表 2 岩的HJC模型参数

Table 2. Parameters of HJC model for rock

ρ/(kg·m^–3)	G/GPa	f_c/MPa	A	B	C	N	S_max	D₁	D₂
2670	10.8	135	0.762	1.65	0.009	0.74	4.00	0.045	1.00
T/MPa	p_c/MPa	μ_c/10^–3	p_l/GPa	μ_l	K₁/GPa	K₂/GPa	K₃/GPa	$\dot \varepsilon$ ₀/s⁻¹	${\varepsilon \rm_{f,min } }$
7.48	45.27	2.9	1.0	0.101	85	−150	208	1.00	0.005

下载: 导出CSV

| 显示表格

1.3 数值模型验证

在SHPB数值模拟中冲击荷载的施加有两种方式：一种是按照实验装置的原比例建立子弹模型并赋予其一定的速度冲击入射杆，另一种是将实验中的入射波曲线直接施加在数值模型中的入射杆端部。相较于第1种方法，后者的建模简单，计算效率高，且与实验条件的吻合度最高，因此本研究采用第2种方式施加冲击荷载。首先选中入射杆起始端部的所有节点，建立Segment组，然后将实验中测得的入射波荷载曲线施加到该组上，则入射杆即可获得初速度。考虑到实际实验中纺锤形冲头冲击入射杆时产生的波为类半正弦形状，为此在模拟时将正弦函数拟合出的波形作为冲击波形施加在入射杆起始端。其中实验测得的一组入射波数据与正弦函数拟合的入射波如图4(a)所示，可见，拟合的半正弦波形与实验测得的入射波形的吻合度较高，此方法适用于冲击荷载的施加。此外，为了验证上述煤、岩的材料参数在煤岩组合体动态劈裂冲击模拟中的有效性，对煤、岩单体试样进行SHPB冲击压缩实验和数值计算，将数值模拟得到的应力-应变曲线与实验结果进行比较，如图4(b)和图4(c)所示，数值计算结果与实验结果比较接近，说明本材料模型参数可以用于煤岩组合体冲击劈裂数值计算及分析。

图 4 冲击荷载下应力波加载方式及材料模型验证

Figure 4. Verification of stress loading mode and material model under impact loading

下载: 全尺寸图片幻灯片

1.4 应力平衡验证

在煤岩组合体冲击劈裂数值计算中，试件两端达到动态应力平衡是一维应力波假设与应力、应变均匀性假设的必要前提，同时也决定着数值计算结果的可靠性。典型煤岩组合体巴西劈裂动态冲击模拟试样两端的应力变化如图5所示。由图5可知，在煤岩组合体数值试样达到峰值应力之前，入射应力与反射应力的叠加之和与透射应力近似相等，即组合体试样两端的动态应力变化基本一致，表明冲击加载过程中组合体试样始终处于动态应力平衡状态。至此，在其余条件下组合体动态拉伸过程中的动态受力平衡均得到了验证。

图 5 煤岩组合体的动态应力平衡验证

Figure 5. Dynamic stress equilibrium verification of coal-rock combination bodies

下载: 全尺寸图片幻灯片

2. 结果及分析

2.1 应力波传播特征

对比不同冲击荷载下C-R与R-C组合体的波形发现，C-R与R-C组合体的应力波表现出相同的性质：不同冲击荷载下，入射波、反射波和透射波的波形基本一致，且入射波与透射波的波形类似于正弦波形；组合体的入射波、反射波和透射波幅值均随着冲击荷载的增大而增大，但透射波幅值的变化不明显。值得注意的是，在相同的冲击荷载条件下，C-R与R-C组合体的入射波幅值基本一致，且R-C组合体的反射波幅值相较于C-R组合体整体偏大，但透射波幅值相较于C-R组合体整体偏小。当冲击荷载较小时，这种差异较明显；随着冲击荷载的增大，差异逐渐减小，具体变化情况如图6所示。

图 6 不同冲击荷载下C-R与R-C组合体的应力波传播特征

Figure 6. Stress waves of C-R and R-C combination bodies under different impact loads

下载: 全尺寸图片幻灯片

煤岩组合体中煤体和岩体的力学性质不同，两者的波阻抗存在差异，即岩体的波阻抗大于煤体的波阻抗。波阻抗匹配效应影响波的传递效率。C-R组合体与入射杆的波阻抗匹配效果优于R-C组合体，当入射应力波传到入射杆与煤岩组合体的接触界面时，更优的阻抗匹配组合产生的反射波更少，更多的波通过岩石试样传入透射杆，因此R-C组合体的反射波相较于C-R组合体偏大。但是随着冲击荷载不断增大，入射应力波逐渐增强，波阻抗匹配效应对应力波传播特征的影响逐渐减弱。通过在后处理软件Ls-PrePost中查看煤岩组合体在冲击劈裂过程中的压力云图，可以直观清晰地看到应力波在入射杆、透射杆及组合体试样中的传播过程，详细信息如图7所示。

图 7 煤岩组合体动态劈裂应力波传播过程

Figure 7. Stress wave propagation process of coal-rock mass during dynamic splitting

下载: 全尺寸图片幻灯片

2.2 变形破坏过程

为研究冲击荷载下煤岩组合体动态变形破坏过程及特征，分析煤岩组合体的受力状态十分必要。为简化分析过程，将煤、岩视为一个紧密相连的整体，交界面处有黏聚力，且胶结厚度忽略不计，其力学模型如图8(a)所示。在SHPB冲击中，煤岩组合体两端受轴向压缩动荷载作用，由动态应力平衡验证可知，试样两端动荷载相等，在轴向压应力作用下，煤体和岩体在水平方向均会产生变形，由于煤体和岩体的弹性模量和泊松比不同，由胡克定律可知，两者在水平方向的应变存在差异。这种差异必然导致两者在胶结面处产生相互作用力。在煤岩组合体冲击压缩过程中的线弹性变形阶段，煤、岩的侧向应变表达式为

图 8 煤岩组合体力学模型及受力分析示意图

Figure 8. Mechanical model and stress analysis diagram of coal-rock mass

下载: 全尺寸图片幻灯片

$\begin{cases} {\varepsilon \rm_{2c}} = {\varepsilon \rm_{3c}} = {\nu \rm_c}\dfrac{{{\sigma \rm_{1c}}}}{{{E\rm_c}}} \\ {\varepsilon \rm_{2r}} = {\varepsilon \rm_{3r}} = {\nu \rm_r}\dfrac{{{\sigma \rm_{1r}}}}{{{E\rm_r}}} \end{cases}$

(11)

式中：ν_c、E_1c、σ_1c和ν_r、E_1r、σ_1r分别为同一时刻煤、岩的泊松比、弹性模量和轴向应力，ε_2c、ε_3c和ε_2r、ε_3r分别表示为煤、岩在水平两个垂直方向的应变。

取煤岩组合体胶结面处的一个微单元体进行受力分析，其应力状态如图8(b)所示。由于组合体中岩的弹性模量大于煤，且煤的泊松比大于岩，由式(11)可知，煤的侧向应变大于岩。煤岩组合体中的岩体部分受水平拉应力作用，而煤体部分受大小相等、方向相反的水平方向压应力作用。

鲜学福等^[24]根据岩石经典强度理论，对静态压缩条件下层状复合岩体强度变化规律进行了理论分析，结果表明：三向拉应力状态会使交界面处岩石强度降低，而三向压应力状态会使交界面处的岩石强度提高，且交界面处岩石强度变化的幅度与两种岩石的弹性模量和泊松比有关。煤岩组合体所受的动态轴向压缩荷载作用与静载轴向压缩作用有所区别，但分析方法与结论可类比推理。在动态劈裂冲击数值模拟过程中，煤岩组合体试样两端的轴向应力随时间变化，且其力学行为有明显的应变率效应，因此，煤岩组合体交界面处的轴向应力是一个变量。但组合体中煤、岩的弹性模量、泊松比大小关系始终保持不变，所以交界面处两者的受力方向不会发生改变，只是受力大小发生变化。因此，组合体交界面处的煤体受三向压应力作用，强度提高，大于其他区域煤体强度；交界面处岩体受三向拉应力作用，强度降低，小于其他区域岩体强度。

选取数值模拟中的一组计算结果进行分析，以冲击荷载100 MPa为例，冲击方向由岩入煤，对比分析不同加载角度条件下煤岩组合体的动态变形破坏过程（见图9）。由上述煤岩组合体理论力学模型及交界面处受力分析可知，在煤岩组合体试样受压冲击荷载下，煤体部分相较于交界面的远处区域先发生变形破坏，而交界面处岩体部分相较于其他区域先发生变形破坏，破坏模式以拉伸、剪切为主。从煤岩组合体变形破坏结果可以看出，组合体试样在动态压缩下，先在煤岩交界面远处的煤体部分起裂破坏，然后在交界面近区的岩体部分起裂破坏，随后煤、岩体各部分的裂纹逐渐扩展、贯通、聚合成形，直至最后形成宏观的损伤破裂面。不同加载角度下煤岩组合体的破坏过程如图9所示，其中用红线标出了具有代表性的特征区域。

图 9 不同加载角度下煤岩组合体的变形破坏过程（100 MPa冲击荷载）

Figure 9. Deformation and failure process of coal-rock mass under different loading angles (impact loading of 100 MPa)

下载: 全尺寸图片幻灯片

在煤岩组合体的冲击劈裂过程中，煤体的强度比岩体的强度低，因而组合体中以煤体破坏为主。从数值模拟结果可以发现，组合体试样的宏观破裂面及损伤破坏区域多出现在煤体部分，岩体的破裂损伤区相对较小。对比宏观破裂区域大小可知，当加载角度为0°和45°时，组合体的破坏程度最大，且岩体的破坏程度也最大。当加载角度为30°时，组合体及其岩体部分的破坏程度最小。随着加载角度的增大，煤岩组合体中岩体的破坏程度呈先减小后增大的趋势，而煤体部分的破坏程度随加载角度的变化不明显。

2.3 破坏特征分析

在不同的冲击荷载作用下 C-R 和 R-C 组合体的变形破坏特征如图10所示。

图 10 不同冲击荷载下C-R与R-C组合体的变形破坏特征

Figure 10. Deformation and failure characteristics of C-R and R-C combination bodies under different impact loads

下载: 全尺寸图片幻灯片

对比不同冲击方向下煤岩组合体的变形破坏模式，可以发现，对于C-R组合体，由煤入岩冲击时，煤体的变形破碎程度更剧烈，岩体的变形破坏并不明显；对于R-C组合体，由岩入煤冲击时，煤体中裂纹扩展的长度、数量和范围相对更小，破坏程度较小。分析认为：由岩入煤冲击时，岩体吸收并集聚了大量的弹性应变能，能量传递至煤体时，能量损耗较多，导致煤体的破坏相对较小；由煤入岩冲击时，由于煤体的强度相对较小，更容易发生变形破坏，大部分入射能量被煤体吸收，能量传递至岩体部分时，剩余的能量难以达到岩体破坏强度，致使岩体的变形破坏并不明显。

对比分析不同冲击荷载下煤岩组合体的变形破坏模式发现，同一加载角度下煤岩组合体裂纹起裂的初始位置及扩展方向基本一致，随着荷载的增大，组合体整体的破坏程度加剧。50 MPa冲击荷载下，裂纹扩展长度、数量及范围相对较小，C-R与R-C组合体的破坏形态基本相同；随着冲击荷载的不断增大，两者的破坏形态差异更明显；75 MPa冲击荷载下，当加载角度为15°和30°时，C-R与R-C组合体形成的宏观裂隙方向大致与交界面平行，其他加载角度下，组合体的破坏形态基本相同；当冲击荷载增大至100 MPa时，组合体中煤、岩部分的裂纹扩展延伸范围最大，其中对于R-C组合体，当加载角度大于或等于45°时，裂纹在煤、岩两侧萌生和扩展后会继续破裂直至贯穿交界面，其具体的破坏特征如图10所示。

上述对于煤岩组合体破坏特征的分析主要以试样的宏观破坏形态为依据，虽然数值模拟获得的试样破坏形态直观清晰，但也只是定性分析。在上述分析的基础上，通过提取煤岩组合体数值模型的损伤失效单元数，以此定量分析组合体试样的破坏程度及破坏规律，结果如图11所示。从图11可以看出，C-R与R-C组合体在同一加载角度下的损伤失效单元数变化规律基本一致，即组合体损伤失效单元数随着冲击荷载的增大而增加。50 MPa冲击荷载下，组合体试样的损伤失效单元数随角度的变化不明显，75和100 MPa冲击荷载下，组合体试样的损伤失效单元数随加载角度的变化幅度较大，且当加载角度为45°时，试样的整体破坏程度最剧烈。基于试样损伤单元数的分析结果与上述理论分析结果基本一致。

图 11 煤岩组合体损伤破坏程度的变化规律

Figure 11. Damage degree of coal-rock mass

下载: 全尺寸图片幻灯片

3. 结　论

通过SHPB实验验证了煤、岩HJC材料模型参数的有效性，采用显式动力学分析软件LS-DYNA研究了不同冲击荷载、冲击方向及加载角度条件下煤岩组合体动态劈裂过程中应力波的传播规律、组合体的变形破坏过程及破坏特征，得到的主要结论如下：

(1) 对比分析由岩入煤和由煤入岩冲击发现，C-R组合体与R-C组合体中入射波的波形及幅值基本一致，但R-C组合体的反射波幅值偏大，透射波幅值偏小，随着冲击荷载增大，差异逐渐减小；

(2) 在不同冲击荷载和不同加载角度条件下，煤岩组合体以煤体破坏为主，岩体的破坏范围较小，并且煤体总是在交界面远处先产生裂隙，而岩体则多在交界面近处先发生起裂破坏；

(3) 在50和75 MPa冲击荷载下，C-R与R-C组合体的破坏形态基本相同，以拉伸、剪切破坏为主，随着荷载的增大，岩体的破坏程度加剧，形成的裂纹逐渐扩展直至穿过交界面与煤体连接贯通；

(4) 以有限元模型损伤失效的单元数为判别依据，定量分析煤岩组合体的破碎程度，当加载角度为45°时，组合体的破坏程度最剧烈，该结果与理论分析吻合。

Figure 1. Volume vs. pressure curves of α-LiAlH₄ and β-LiAlH₄ phases (Enthalpy difference (per formula unit) between α-LiAlH₄ and β-LiAlH₄ as a function of pressure is shown in the insert.)

下载: 全尺寸图片幻灯片

Figure 2. Phonon dispersion relations and density of state in high-symmetry directions for the α-LiAlH₄ structure at 0 and 1.6 GPa

下载: 全尺寸图片幻灯片

Figure 3. Calculated total and partial electronic densities of states for α-LiAlH₄ structure at 0 GPa (a) and β-LiAlH₄ structure at 1.6 GPa (b) respectively

下载: 全尺寸图片幻灯片

Table 1. Optimized structural parameters, atomic position parameters for the α-LiAlH₄ and the β-LiAlH₄ structures

Phase	Unit-cell dimensions	Atom coordinates
α-LiAlH₄ (P2₁/c)	a=0.485 1 nm(0.481 7 nm^) b=0.781 4 nm(0.780 2 nm^) c=0.773 2 nm(0.782 1 nm^*)	Li:(0.585, 0.459, 0.829), (0.560, 0.466, 0.827)^* Al:(0.159, 0.204, 0.938), (0.139, 0.203, 0.930)^* H1:(0.193, 0.102, 0.766), (0.183, 0.096, 0.763)^* H2:(0.377, 0.371, 0.986), (0.352, 0.371, 0.975)^* H3:(0.254, 0.082, 0.119), (0.243, 0.081, 0.115)^* H4:(0.822, 0.268, 0.882), (0.799, 0.247, 0.872)^*
β-LiAlH₄ (I2/b)	a=0.445 2 nm b=0.445 9 nm c=1.010 2 nm β=89.978°	Li:(0, 0.250, 0.125) Al:(0, 0.250, 0.625) H1:(0.259, 0.425, 0.542) H2:(0.324, 0.508, 0.792)
Note:"*" represents experimental data from Ref.[18].

下载: 导出CSV

Table 2. Average net charges, bond length (L) and scaled bond overlap population (P^s) between H, Al, and Li atoms in the α-LiAlH₄ (at 0 GPa) and the β-LiAlH₄ (at 2.0 GPa) structures

Phase	Average net charge			L/nm		P^s
Phase	H	Al	Li	Al─H	Li─H	Al─H	Li─H
α-LiAlH₄	-0.48	0.64	1.29	0.161 5	0.188 6	0.506	0.019
β-LiAlH₄	-0.46	0.55	1.28	0.162 8	0.206 7	0.534	-0.063

下载: 导出CSV

[1]	SCHLAPBACH L, ZÜTTEL A.Hydrogen-storage materials for mobile applications[J]. Nature, 2001, 414(6861):353-358. doi: 10.1038/35104634
[2]	SCHÜTH F, BOGDANOVIC' B, FELDERHOFF M.Light metal hydrides and complex hydrides for hydrogen storage[J]. Chemical Communications, 2005, 36(2):2249-2258.
[3]	ZHU C Y, LIU Y H, DUAN D F.Structural transitions of NaAlH₄ under high pressure by first-principles calculations[J]. Physica B:Condensed Matter, 2011, 406(8):1612-1614. doi: 10.1016/j.physb.2011.02.007
[4]	DILTS J A, ASHBY E C.Thermal decomposition of complex metal hydrides[J]. Inorganic Chemistry, 1972, 11(6):1230-1236. doi: 10.1021/ic50112a015
[5]	DYMOVA T N, ALEKSANDROV D P, KONOPLEV V N, et al.Spontaneous and thermal-decomposition of Lithium Tetrahydroaluminate LiAlH₄－the promoting effect of mechanochemical action on the process[J]. Russian Journal of Coordination Chemistry, 1994, 20(4):279-285. http://www.academia.edu/13646109/Hydrogen_production_from_solid_reactions_between_MAlH4_and_NH4Cl
[6]	VAJEESTON P, RAVINDRAN P, VIDYA R, et al.Huge-pressure-induced volume collapse in LiAlH₄ and its implications to hydrogen storage[J]. Physical Review B, 2003, 68:212101. doi: 10.1103/PhysRevB.68.212101
[7]	PITT M P, BLANCHARD D, HAUBACK B C, et al.Pressure-induced phase transitions of the LiAlD₄ system[J]. Physical Review B, 2005, 72:214113. doi: 10.1103/PhysRevB.72.214113
[8]	CHELLAPPA R S, CHANDRA D, GRAMSCH S A, et al.Pressure-induced phase transformations in LiAlH₄[J]. The Journal of Physical Chemistry B, 2006, 110(23):11088-11097. doi: 10.1021/jp060473d
[9]	TALYZIN A V, SUNDQVIST B.Reversible phase transition in LiAlH₄ under high-pressure conditions[J]. Physical Review B, 2004, 70:180101. doi: 10.1103/PhysRevB.70.180101
[10]	HOHENBERG P, KOHN W.Inhomogeneous electron gas[J]. Physical Review B, 1964, 136(3):864-871.
[11]	KOHN W, SHAM L.Self-consistent equations including exchange and correlation effects[J]. Physical Review A, 1965, 140(4):1133-1138. http://tu-freiberg.de/sites/default/files/media/institut-fuer-theoretische-physik-10451/Lehre/Dichtefunktionaltheorie/a9rf1a3.pdf
[12]	PERDEW J P, WANG Y.Accurate and simple analytic representation of the electron-gas correlation energy[J]. Physical Review B, 1992, 45(23):13244-13249. doi: 10.1103/PhysRevB.45.13244
[13]	TROULLIER N, MARTINS J L.Efficient pseudopotentials for plane-wave calculations[J]. Physical Review B, 1991, 43(3):1993-2006. doi: 10.1103/PhysRevB.43.1993
[14]	SEGALL M D, LINDAN P J D, PROBERT M J, et al.First-principles simulation:ideas, illustrations and the CASTEP code[J]. Journal of Physics:Condensed Matter, 2002, 14:2717-2744. doi: 10.1088/0953-8984/14/11/301
[15]	MONKHORST H J, PACK J D.Special points for Brillouin-zone integrations[J]. Physical Review B, 1976, 13(12):5188-5192. doi: 10.1103/PhysRevB.13.5188
[16]	PARLINSKI K, LI Z Q, KAWAZOE Y.First-principles determination of the soft mode in cubic ZrO₂[J]. Physical Review Letters, 1997, 78(21):4063-4066. doi: 10.1103/PhysRevLett.78.4063
[17]	TOGO A, OBA F, TANAKA I.First-principles calculations of the ferroelastic transition between rutile-type and CaCl₂-type SiO₂ at high pressures[J]. Physical Review B, 2008, 78(13):134106. doi: 10.1103/PhysRevB.78.134106
[18]	HAUBACK B C, BRINKS H W, FJELLVG H.Accurate structure of LiAlD₄ studied by combined powder neutron and X-ray diffraction[J]. Journal of Alloys and Compounds, 2002, 346(1):184-189.
[19]	SEGALL M D, PICKARD C J, SHAH R, et al.Population analysis in plane wave electronic structure calculations[J]. Molecular Physics, 1996, 89(2):571-577. doi: 10.1080/002689796173912
[20]	MULLIKEN R S.Electronic population analysis on LCAO-MO molecular wave functions[J]. The Journal of Chemical Physics, 1955, 23(10):1833-1840. doi: 10.1063/1.1740588
[21]	HU C H, CHEN D M, WANG Y M, et al.First-principles investigations of the pressure-induced structural transitions in Mg(AlH₄)₂[J]. Journal of Physics:Condensed Matter, 2007, 19:176205. doi: 10.1088/0953-8984/19/17/176205
[22]	WANG H, LI Q, WANG Y C, et al.High-pressure polymorphs of Li₂BeH₄ predicted by first-principles calculations[J]. Journal of Physics:Condensed Matter, 2009, 21:385405. doi: 10.1088/0953-8984/21/38/385405