内爆炸下管道裂纹渐进式动态扩展机制研究

杜 洋，郑津洋，马 利

(1.中国石油大学(华东) 机电工程学院，山东青岛 266580；2.浙江大学 化工过程机械研究所，杭州 310027；3.浙江工业大学 机械工程学院，杭州 310014)

0 引言

近年来发生的重/特大爆炸事故(青岛输油、天津港、江苏响水等地)造成了极为严重的安全事故。爆炸产生的冲击波往往导致周围结构的变形与动态破坏，且两者之间相互耦合、相互影响[1-2]。掌握爆炸载荷下结构动态断裂行为与特性是进行定量风险评估与事故分析推演的重要条件，相应的研究对提升我国生产安全保障与事故防控水平具有重要意义[3-5]。

内爆炸下管道失效研究可追溯到MOTT[6]在二战时期有关动态碎裂的研究工作。21世纪以来，更多学者开始关注管道裂纹的全周期扩展特性。代表性试验研究如CHAO等[7-8]开展的气体爆轰加载铝合金管道动态断裂试验，该研究系统获得了裂纹的典型动态扩展与分叉行为。2017年，MIRZAEI等[9]在采用线状炸药加载金属管道试验中首次清晰捕捉了Ⅰ型裂纹的渐进/周期式扩展现象(见图1[9])，并认为这与管道应力波相关，但尚未探明具体的机制与内在规律。内爆炸下管道动态断裂是一个高度瞬态的过程，裂纹扩展速率高达200～300 m/s且呈现高度不稳定的特性[10]。目前，通过试验方法对裂纹扩展及应力波-裂纹相互作用过程进行准确且详尽地记录存在一定困难，因此构建可靠的数值分析模型至关重要。

图1 内爆炸下管道裂纹渐进式扩展“阶梯”状痕迹Fig.1 Staircase markings on the fracture surfaceof pipe due to internal explosion

对于数值模拟研究方面，LIU等[11-15]将内部爆炸/爆轰载荷作为压力边界条件施加于内壁面，以分析管道的动态断裂行为，由于忽略了爆轰波-管道间的流固耦合作用，模拟的管道断裂形貌往往与试验存在较明显差别[16]。在流固耦合模拟研究方面，由于研究对象的流-固-断裂-化学反应强耦合特性，相应数值模型构建难度较大，目前主要是国际顶尖研究机构开展了一些研究工作，如CIRAK等[17]提出的Level-set-based方法和WANG等[18]提出的FIVER-based 计算框架，均实现了气体爆轰下管道裂纹扩展的流固耦合模拟。无论是解耦还是耦合模拟，目前相关研究主要侧重于算法实现，模型中缺乏对高应变率下材料失效机制的细致考究，大多采用临界应力、应变、内聚能等失效准则判断裂纹的起始与扩展，这使得模拟的裂纹扩展行为与试验结果仍存在一定差别。

已有研究[19-21]表明，爆炸载荷下金属材料的失效模式往往是绝热剪切，而细观剪切带的演化同时依赖于应变与应变率。这表明目前相关研究中广泛采用的单一临界值失效准则可能无法较好地表征高应变率下材料的失效行为与特征。

本文在前期内爆炸下管道失效模式研究基础上，考虑高应变率下材料失效应变的应力三轴度与应变率依赖特性，建立相应的损伤演化模型，结合简化气体爆轰模拟算法与改进的流固耦合算法，构建内爆炸下管道裂纹动态扩展分析模型，深入解析管道应力波传播、管道径向振动与裂纹动态扩展间的互作关系，研究阐明管道裂纹渐进式动态扩展的物理机制，进而揭示内爆载荷速率与裂纹渐进扩展长度间的定量关系。

1 典型试验概述

为便于模型验证，本文以加州理工学院开展的典型试验为研究对象构建数值分析模型。试验装置如图2所示。含外表面初始缺陷的6061-T6铝合金试验管道通过法兰与钢制爆轰管道相连接，初始时刻，两管道内充满按化学计量比混合的乙烯和氧气，试验开始时，通过爆轰管道左端安装的火花塞将可燃混合气体点燃，火焰在爆轰管道内加速并经历爆燃转爆轰(DDT)后形成爆轰波进入试验管道，在爆轰波瞬态高压作用下，初始缺陷处形成贯穿裂纹并扩展。铝合金试验管道右端采用聚酯膜片密封，以消除爆轰波反射的影响。铝合金试验管道壁厚0.89 mm，初始缺陷宽度0.30 mm，深度0.56 mm，长度12.7～76.2 mm。试验中通过控制内部混合气体初始压力和初始缺陷长度，以获得不同的裂纹动态扩展行为与管道断裂形貌。

图2 气体爆轰下管道裂纹动态扩展试验装置结构示意Fig.2 Schematic diagram of experimental setup for dynamiccrack propagation of pipe due to gaseous detonation

2 数值分析模型

2.1 本构与损伤模型

爆炸载荷下结构应变率高达102～107s-1，相应材料力学响应由应变硬化、应变率硬化和温度软化3种效应共同表征[22]。本文采用Johnson-Cook型热黏塑性本构模型[23]描述材料应力-应变关系：

(1)

研究[22-23]表明，高速冲击/爆炸载荷下金属材料往往以绝热剪切模式发生失效，此时塑性功转变的热量会在局部积聚导致变形高度局域化，形成宽10～100 μm的绝热剪切带(ASB)，最终裂纹沿剪切带扩展而发生破坏。无论是细观上对剪切带演化过程的研究；还是宏观上热黏塑性本构的失稳分析结果，均表明绝热剪切失效由应变与应变率共同表征[24-25]。更广泛的研究[26]表明，高应变率下材料失效应变pf由应力三轴度、应变率和温度共同决定：

(2)

式中，D1～D5为材料失效参数；σ*为应力三轴度。

学者们从微观、细观层面对材料损伤演化机制开展了大量研究，为便于应用，往往需要回到宏观层面上，根据连续损伤力学理论，在式(1)基础上，耦合宏观损伤的热黏塑性本构方程为：

(3)

式中，D为宏观连续损伤变量；r为损伤等效塑性应变。

为便于损伤演化参数的确定，文中假定损伤阈值pd与临界损伤DC均为材料常数，并采用以下方程描述损伤演化过程：

(4)

其中损伤阈值pd取值0.04[27]，根据模拟的管道断裂形貌与试验的对比，确定DC=0.45。其余材料参数由文献[28]中获得，列于表1中。

表1 管道6061-T6铝合金材料参数Tab.1 Material parameters for aluminum alloy 6061-T6 pipe

2.2 气相爆轰模拟

爆轰波由冲击波及其后方紧邻的化学反应区薄层组成，在管道内往往以超声速稳定传播。根据CJ爆轰理论[29-30]，若不关注爆轰波的细观胞格结构，可将其简化为一包含化学反应薄层的强间断面，化学反应的效果可用一外加能量项反映到间断面关系中。基于上述思想，本文采用一种自定义燃烧模型模拟气体爆轰过程[31]，忽略气体黏性并考虑其可压缩性，基本控制方程为质量、动量以及能量守恒方程：

(5)

式中，ρ为密度;u为速度矢量;E为比总能量。

为使式(5)闭合，还需要燃烧产物的状态方程，对于气相爆轰，可采用理想气体状态方程p=(γ-1)ρe，(其中，e为比内能，e=E-u2/2)。为将化学反应与流体运动耦合，引入表征化学反应能量释放的变量λ到状态方程中，即：

p=λ(γ-1)ρe

(6)

λ由以下简化化学反应率函数计算：

(7)

式中，t为时间;tl为每个单元的点燃时间，由单元中心距离起爆中心的长度除以爆轰速度DCJ得到；Lmin为沿爆轰波传播方向的最小单元长度。

式(7)中爆轰波厚度为1.5Lmin，当爆轰波传播经过某处单元时，认为该单元内可燃混合气体被点燃，相应化学反应能量加入到爆轰产物状态方程(见式(6))中，使得该处压力升高至爆轰压力pCJ。

针对本文模拟试验的工况(p0=180 kPa)，相应的气体爆轰压力pCJ=6.1 MPa，爆轰速度DCJ=2 400 m/s，可燃混合气体的比化学反应能量为8.93×106J/kg。

2.3 流固耦合实现

采用浸入边界法(IBM)实现流体与固体的耦合，即固体附着于Lagrange单元并为Euler单元内流体提供几何边界条件，反之，Euler单元内流体则为固体提供压力边界条件。流体与固体之间接触力采用罚函数法[32]计算，即接触力正比于接触刚度与接触深度。

在气体爆轰加载铝合金管道过程中，流体与固体之间接触力经历了高幅、瞬态的极端变化过程，且两种材料间密度、刚度差异较大，耦合模拟中往往出现穿透、泄漏、计算不稳定的问题。为解决上述问题，本文改进了接触刚度k的计算方法，使其综合考虑计算稳定性、总体时间步长以及节点质量，得到下式[31]：

(8)

式中，K为材料的体积模量；s为接触段面积；v为单元体积；m，Δt分别为节点质量和时间步长。

2.4 有限元分析模型

爆燃转爆轰不是本文研究重点，故仅对铝合金试验管道部分按照1∶1尺寸比例进行建模，如图3所示。模型由含缺陷试验管道、内部可燃混合气体和外部空气3部分组成。其中管道采用Lagrange单元离散，内、外气体采用Euler单元离散。在裂纹可能经过的区域对单元进行加密，单元特征尺寸为0.3～0.9 mm，模型共划分单元约120万个。

图3 气体爆轰下管道动态断裂流固耦合分析模型Fig.3 Fluid-structure coupling analysis model for dynamicfracture of pipe due to gaseous detonation

针对试验条件，试验管道左右两端施加固定边界条件，左侧出口施加壁面边界条件，右侧出口及模型空气域外表面施加无反射边界条件。对应试验获得的3种典型裂纹扩展行为，内部可燃混合气体初始压力为180 kPa，起爆点位于模型左端中心，管道初始缺陷长度分别为25.4，50.8，76.2 mm。

3 结果与分析

3.1 模型验证

数值模拟借助LS-DYNA求解器完成，相关模型及算法通过编写自定义材料模型与编辑K文件(Keyword file)实现。

图4[7]示出了模拟得到的3种典型裂纹扩展行为及其与试验的对比。定义扩展方向与爆轰波传播方向(自左向右)一致的裂纹为前向裂纹;反之为后向裂纹。管道初始缺陷长度为25.4 mm时，前后向裂纹均沿管道轴线直向扩展；初始缺陷长度为50.8 mm时，前向裂纹扩展一段距离后发生分叉，分叉裂纹沿管道环向扩展；初始缺陷长度增加至76.2 mm时，前后向裂纹均发生了分叉，4条环向扩展的分叉裂纹将管道切分为3部分。总体来说，数值模拟较好地复现了裂纹的典型动态扩展行为与管道最终断裂形貌，主要区别在于试验中，后向裂纹最终发生了弯折(初始缺陷长度为25.4 mm和50.8 mm时)，这可能是试验中管材力学性能并不完全均一引起的。

图4 气体爆轰下管道典型裂纹扩展行为模拟与试验结果对比Fig.4 Comparison between simulated and experimentalresults of propagations of typical cracks of pipe due togaseous detonation

3.2 应力波传播及管道动态响应

爆轰波在管道内形成后，首先以球面波形式向外传播，作用至管道左端后，使该位置产生瞬态径向鼓胀变形，该变形同时包含了弯曲与剪切成分，分别形成弯曲波和剪切波自左向右传播，其中剪切波强度较低，在图5绘制的应变标尺下不可见，而弯曲波强度较高，其传播速率与气体爆轰速度DCJ一致，为2 400 m/s。图5(c)示出了距左端25 cm处管道的环向应变响应，由于弯曲波的连续形成与传播，管道环向应变历程呈现高度震荡的特性。对环向应变数据进行快速傅里叶变换(FFT)，得到其震荡频率为38.8 kHz，这一频率理应与弯曲波的形成频率(管道径向振动频率)相一致。管道轴对称径向振动频率可通过理论公式[33](见式(9))计算得到，为40.7 kHz，两者相差4.7%，具有较好吻合度，也说明了本文构建模型的准确性。

图5 爆轰波传播及其作用下管道动态响应(管道几何变形放大比例50∶1)Fig.5 Propagation of detonation wave and the resultingpipe responses (Pipe displacement scale factor 50∶1)

(9)

3.3 裂纹渐进式动态扩展

管道裂纹扩展经历了加速-高速扩展-减速3个阶段，其中高速扩展阶段贡献了裂纹最终长度的绝大部分，也是本文分析的重点。

图6示出了3种初始缺陷长度下管道裂纹高速扩展阶段的裂纹长度-时间历程曲线。可以看出，裂纹扩展平均速率为200～300 m/s，约为管道弯曲波传播速率(2 400 m/s)的1/10，并且前向裂纹速率高于后向裂纹，这是由于前向裂纹扩展方向与弯曲波传播方向一致，使得弯曲波作用于前向裂纹裂尖的时间长于后向裂纹。此外，图6(b)(c)还表明，主裂纹的扩展速率总体上高于分叉裂纹。

进一步在微秒尺度上解析裂纹动态扩展过程，发现裂纹扩展是渐进式的。为阐明其特性与规律，对图6数据进行微分处理，得到裂纹扩展瞬态速率(见图7)，对各裂纹速率历程进行FFT处理，表明裂纹的渐进式扩展存在某一主导频率(尤其对于前后向裂纹)，范围为33～45 kHz，恰恰将管道弯曲波的形成频率(38.8 kHz)包含在内。因此，认为裂纹渐进式扩展与管道弯曲波传播高度相关。

图6 管道裂纹扩展长度-时间历程与平均速度

图7 管道裂纹扩展速率与震荡频率Fig.7 Crack growth rates and oscillation frequency of pipe

结合数值模拟结果解析弯曲波与裂纹相互作用过程(见图8)，可以看出，当管道弯曲波传播经过裂尖时，即对该裂纹的扩展产生一次加速/促进作用，而在裂纹整个扩展历程中，弯曲波是不断形成并沿管道轴向传播的，所以导致了裂纹的渐进式动态扩展。基于以上认识，可提出以下理论公式计算裂纹的渐进式扩展频率：

图8 管道弯曲波-裂纹相互作用过程(以初始缺陷长度25.4 mm为例)Fig.8 Interaction between flexural waves and crack of pipe(Take initial flaw length of 25.4 mm as an example)

(10)

表2列出由式(10)计算的裂纹渐进扩展频率及其与数值模拟结果的对比，结果显示两者具有较好的吻合度，最大差值在6%以内，进一步说明了上述认识的正确性。

以上研究表明,内爆炸下管道裂纹实际上是在某一主导频率下有规律的渐进式扩展过程，这一特性是管道轴向传播弯曲波与裂纹尖端相互作用的结果，相应的渐进扩展频率可由式(10)计算。

3.4 内爆载荷速率与裂纹渐进扩展长度的关系

绘制管道裂纹扩展过程的损伤云图(见图9)，可见裂纹扩展路径上损伤是相间分布的，不难理解这与裂纹的渐进式扩展相对应，同样也对应试验[9]中管道断面的“阶梯”状痕迹。

图9 管道裂纹扩展路径损伤相间分布云图Fig.9 Distribution nephogram of spaced damages alongthe crack growth path of pipe

(11)

式中，变量下标1,2分别代表前、后向裂纹。

表3列出了由式(11)计算的裂纹平均渐进扩展长度与数值模拟结果的对比，其中由于初始缺陷为76.2 mm时的裂纹前后向渐进扩展次数较少，数据参考价值较低，表中未予列出。数据表明，式(11)的计算值与数值模拟结果的相对差值在9%以内，吻合度较好。

表3 裂纹平均渐进扩展长度理论公式计算值与数值模拟结果对比Tab.3 Comparison between calculated and simulatedresults of average incremental growth length of crack

(12)

式(12) 实际上给出了由管道前后向裂纹渐进扩展长度计算内爆载荷速率的新方法。由于目前爆炸事故分析中，爆炸强度(等效TNT当量)的推断主要依据事故中建筑物或其他结构的毁坏程度，往往难以给出较准确的结果，特别对于气体发生了爆燃(亚声速)还是爆轰(超声速)的分辨也较为困难，因此，式(12)的应用价值在于：在管道爆炸事故分析中，可通过测量管道断面裂纹渐进扩展长度计算得到爆炸发生时的载荷速率，根据载荷速率及其与声速大小的关系，从而对爆炸强度与类型(爆燃、爆轰)进行较准确地推断。

4 结论

(1)构建了计及损伤演化与流固耦合的内爆炸下管道动态断裂分析模型，较好地复现了裂纹典型扩展行为与管道最终断裂形貌，为在微秒尺度解析裂纹渐进式扩展提供了模型基础。

(2)阐明了内爆炸下管道裂纹是在某一主导频率下的渐进式动态扩展，其实际上是管道轴向传播弯曲波与裂纹尖端相互作用的结果，并给出了裂纹渐进扩展频率的理论计算公式。

(3)对于管道这一典型结构，揭示了内爆载荷速率与裂纹渐进扩展长度间的内禀关系，建立了相应的理论关系式，为管道爆炸事故分析中爆炸强度与类型的推演提供了新方法。