2.2.3  数值稳定性分析

在显式时间积分中，除非时间步长选得足够小，否则会出现数值不稳定性，从而使计算结果失去意义。尽管对于完全非线性问题，不能进行稳定性分析，但是经验表明线性分析也会给出有用的时间步长估计。一般说来，将线性稳定的时间步长减少10%到20%，可以克服非线性带来的不稳定性。

对于中心差分，线性分析表明时间步长

其中，是系统的最高频率，是最高频率下的临界阻尼系数。

ω是下面问题的特征值：

[K]{X}=ω²[M]{X}                    （2-44）

其中[K]是系统的刚度矩阵，[M]是质量矩阵。

在实际计算中，我们并不是直接计算系统的ω_max，而是计算单元的。在所有的单元的中选择最大值作为ω_max。

对于板单元，有

其中，C_l是平面纵波波速，对四边形板单元来说，L是短边的长度。

2.2.4  气体压力模式

1.3.3节提到了解耦算法中气体压力模式的广泛应用。对于含裂纹扩展的管道，其内部管壁受到的气体压力由两种模式组成：裂纹前面气体减压传播和裂纹后部的气体压力衰减。

Battelle的Maxey等人在实验中采用压缩空气、氮气以及甲烷含量高的天然气观察气体的减压行为，发现理想气体规律足以满足这些气体的状态等式。假设气体的膨胀过程是等熵过程，充满气体的管道突然发生开裂，整个状态中气体成分混和均匀，则裂尖附近气体局部压力与裂纹扩展速度有如下关系:

式中p_d为减压后的压力水平（MPa)；p_i为开裂前管内压力（MPa)；v为裂纹扩展速度（m/s)；v_a为初始压力和温度状态下的声速（m/s)；为起始状态下气体的绝热系数，C_P和C_v分别为定压比热和定容比热。

从（2-46）式可以看出，当裂纹在管道上扩展时，裂纹前面的气体减压传播，其数值低于气体初始压力值p₀，气压减少值主要取决于裂纹的扩展速度，对于稳态扩展的长裂纹，管道内部气体减压趋于稳定的数值。

当裂纹沿管壁快速扩展时，裂纹后面的气体逸出，气压降低，但是与液体压力不同的是，该区域内的气体压力并不是立即降为零值，而是稳定地衰减至大气压力，并有一段明显的衰减区。在这一衰减区内，气体压力仍然作用在已经开裂的管壁上，加速管壁的变形。

根据实验给出的实际管道开裂实验的结果，发展了裂纹后面气体压力衰减模式，将气压的衰减表示成裂纹位置和初始压力的函数。常见的衰减模式有指数衰减的Fourier级数和抛物线衰减函数等。

更为简便的是图2-5所示的气体压力衰减模式是线性衰减函数，图中以裂纹扩展方向z为横轴，压力分布p为纵轴，其计算表达式为

P（z）=P₁（1-），z<L                    （2-47）

式中z是计算截面的位置，L是衰减长度。在计算中，L取1.5倍或2倍的管道外径。

2.3  管道稳态裂纹扩展问题的数值模拟

本节着重介绍稳态扩展条件下的开裂管道计算结果，力求从多个角度模拟管道变形的真实状态，为建立止裂判据与后文的减速机理打下基础。

2.3.1  计算模型与边界条件

管道上的裂纹一般从起裂点同时向两侧扩展。因管道几何形状和载荷分布具有对称性，本文取一侧管道的一半来分析和计算，即四分之一模型，如图2-6所示。采用四节点板壳单元，沿管道轴向和环向分别划分网格，在管道内壁给出气体压力分布。

在对西气东输管道进行计算时，采用的管道模型长度L=37.5m，直径D=1.016m，厚度h=0.0147m，沿轴向划分为250个单元，环向划分为16个单元，共4000个单元，4267个节点，每个四边形单元轴向长度0.15m，环向长度0.lm。下文如无特殊说明，则网格均按此划分。

在图2-7中，C点表示裂尖位置：裂纹沿LINE4扩展，在A、B点和边界LINEl、LINE2、LINE3和LINE4上施加不同的边界条件。I、Ⅱ和Ⅲ区代表卸载区、减压区和裂前区：其中，点A、B和边界LINEl、LINE2、LINE3是约束边界条件，C点和边界LINE4为运动边界条件。I、Ⅱ及Ⅲ区随时间而变化。

根据对称性指定边界条件如下：

B、F点是固定点；

C点随时间而变化，AC是自由边；

边界LINE1（包括A点）与LINE3上约束了z向的位移和绕x、y轴的转动；

边界LINE2上约束x、y向的位移和绕x、z轴的转动；

边界LINE4上约束x向的位移和绕y、z轴的转动。

2.3.2  表征裂纹驱动力的G与CTOA的计算

本节采用节点力释放法计算裂纹驱动力G。根据（2-38）式，以钢制管道为例，指定设计参数：直径为0.72m，壁厚0.015m，管长取37.5m，压力为0.30MPa,裂纹稳定扩展速度为550m/s。图2-8给出了程序计算的结果，反映了裂纹驱动力与裂纹扩展距离的变化关系。在这一计算中，裂纹后面压力衰减长度，即（2-47)式中的L，取为2倍的直径长度。

如图2-8，随着裂纹起裂后的迅速扩展，裂纹驱动力上升直至达到并保持在稳定的数值，这一平台值1.6KN/m就是裂纹稳态（steady state）扩展的驱动力G。裂纹尖端的移动过程中，G在该值附近有一些小的波动，这是由于裂纹长度的突然增加和位移约束的突然解除，引起的有限元解的高频振荡。振荡的幅度与有限元单元网格的疏密和单元内节点力释放的指定规律有关。

图2-9表示了CTOA的计算方法。由于在裂纹后方有很大的非线性变形，计算中张开位移用割线近似，采用8个单元的长度来计算CTOA。这一长度经过与细化后的网格结果比较，误差在10%之内。

2.4  管道稳态裂纹扩展算例分析

如2.1.2节所述，可在Kanninen（1980）按一维梁模型推导的输气管线裂纹扩展极限速度（2-6a）式的基础上，对高延性管道进行修正。按照西气东输管道选用的设计参数，可得到大致的延生断裂稳态裂纹扩展速度，估算裂纹稳态扩展达到该速度时的裂纹驱动力G和CTOA。另外，通过选用不同的设计参数，可给出稳态裂纹驱动力随内压、壁厚、管径以及径厚比等参数变化的规律。

2.4.1  西气东输设计参数下的典型算例

在这一小节，我们研究的管道参数限于：管道长度L=25m，直径D=l.016m,内压10.OMPa，厚度h=0.0147m。网格划分比2.3.1节中的网格两个方向各密一倍，即：沿轴向划分为500个单元，环向划分为32个单元，共16000个单元，16533个节点，每个四边形单元轴向长度0.05m，环向长度0.05m。裂纹稳态扩展速度按（2-6b）式，取为216m/s。轴向压力衰减长度取1.5倍直径，即约1.5m。

首先观察裂尖附近的MISES应力分布，如图2-10。按以针状铁素体和块状铁素体组织为主的X70钢级管线钢的平均屈服强度500～600MPa估算，计算给出的裂尖前部塑性区长度大约在3个单元长度以内，即0.15m左右。

图2-11直观地给出了内压10MPa，裂纹速度216m/s时一维管道气体压力模型（L=1.5D）下的管壁压力沿轴向的分布。每条曲线的衰减起始点即为该时刻的裂尖位置。图中裂尖压力约为6.6MPa。

上一页

下一页