接触爆炸作用下特大型LNG储罐的动力响应分析

翁大根 葛庆子

同济大学土木工程防灾国家重点实验室

特大型（16×104m3及以上）液化天然气储罐是国家重要的生命线工程，现代局部战争或恐怖活动中，其极有可能遭到精确制导武器的打击或汽车炸弹等恐怖袭击，故评估其安全性至关重要［1－3］。特大型LNG储罐储液量为16×104～20×104m3，其内罐由厚度仅为12～26mm的9%镍合金钢构成，对储罐整体抗爆性能贡献不大。储罐的外罐由预应力钢筋混凝土构成，厚度一般为800～1 000mm，其作用类似于核反应堆安全壳，是LNG储罐抗爆和防止储液泄漏的主要构件［4］。精确制导武器打击或汽车炸弹袭击一般多为接触爆炸，这与外部爆炸产生的冲击波作用于罐体是不同的。国内外一些学者已经对爆炸荷载作用下壳体结构的动力响应进行了理论和试验研究，取得了一些有价值的成果，但是尚未有关于特大型LNG储罐抗爆的详细研究［5－7］。国外现行规范中规定的特大型LNG储罐抗冲击波荷载设计方法为一种拟静力分析方法，将动力荷载简化为等效静力荷载施加于罐体进行分析，且其规定验算的荷载较小、荷载作用类型单一，不能完全满足生命线工程安全评估的要求［8－9］。为此，本文采用流固耦合算法，对接触爆炸作用下特大型LNG储罐的动力响应和破坏形态进行了数值模拟分析，并对储罐整体抗爆安全性能做出评估，得到了一些对工程有实际意义的结论。

1 分析方法与计算模型

接触爆炸作用下特大型LNG储罐的动力响应分析涉及储罐内外空气、罐体和储液之间的流固耦合计算，是典型的多物质流体弹塑性力学问题。

1.1 分析方法

目前，多物质流固耦合算法主要包括Lagrange方法、Euler方法、ALE（Arbitrary Lagrangian Eulerian）方法 和 CLE（Coupled Lagrangian Eulerian）方 法等［10］。为了准确模拟空气中爆炸冲击波的传播过程，以及空气、罐体和储液之间的流固耦合相互作用，笔者采用LS－DYNA有限元程序中的ALE方法建立多物质分析模型进行数值分析。该方法既可以克服Lagrange单元严重畸变引起的数值计算困难，又可以通过多物质单元来实现流固耦合的动态分析［11］。

1.2 计算模型

笔者以我国沿海地区某全容式特大型LNG储罐为计算原型，设计工作容积160 000m3，工作温度－161℃，设计液位高度34.76m，正常操作最高液位34.26m，设计液体密度480kg／m3。

储罐外罐为预应力钢筋混凝土结构，混凝土强度等级为C50，内径80m，罐体混凝土厚度0.8m，罐高（至穹顶）53.023m，穹顶高度11.117m，穹顶为0.4～0.795m均匀变化的变截面，顶部厚度0.4m，加腋处0.795m；预应力筋为7股钢绞线，公称直径15.7mm，极限强度标准值1770N／mm2；普通钢筋为HRB400级。

储罐内罐由9%镍合金钢组成，内径78m，高度37.896m，钢板厚度从下至上从26～12mm均匀变化。内罐与外罐之间密实填充珍珠岩形成保温层，填充振实密度为70kg／m3。由于本文爆炸点位于罐顶，该处没有内罐和保温层，故本文计算模型中未考虑两者作用，储液直接与外罐进行流固耦合计算。

有限元计算模型中，混凝土采用三维实体8节点单元模拟；预应力筋和普通钢筋采用3节点梁单元；空气、储液和炸药采用ALE单元网格，三者材料在整体网格中自由流动，图1为计算模型示意图。考虑到常规导弹爆炸威力一般为500～1 000kg TNT当量，故本文选取的TNT等效爆炸当量为900kg，考虑最不利荷载情况，故接触爆炸点选取在外罐厚度最小处，即穹顶顶部中心。分析时，安全壳上选择6个特征点，用以考察储罐整体抗爆性能（如图1所示）。

图1 计算模型立面图

2 材料参数

在数值模拟中，炸药采用JWL状态方程［12］，其形式为：

式中p为压力；V为相对体积；E0为初始比内能，A1、B1、R1、R2、ω 为材料常数。

各状态方程参数见表1。

表1 JWL状态方程参数表

空气和储液均采用LS－DYNA中的MAT＿NULL材料模型，空气采用式多项式状态方程描述，储液采用式GRUNEISEN状态方程描述。

混凝土采用Johnson－Holmquit本构模型，模型参数见表2。该模型因综合考虑了损伤、应变率、静水压以及压碎、压实效应，适合大应变、高应变率和高静水压情况，被广泛地应用于混凝土撞击、侵彻和爆炸等强动载问题的数值模拟。钢筋材料模型采用双线性随动硬化模型，用两个斜率（弹性和塑性）来表示钢筋的应力应变特性，具体参数见表3。

表2 Johnson－Holmquit本构模型参数表

表3 钢筋模型参数表

3 数值模拟结果及分析

3.1 冲击波荷载和冲击波流场

储罐各特征点处在爆炸作用下的荷载时程曲线如图2所示，由于E、F两点距离爆心较远，其超压现象不明显，本文并未示出这两点的相应曲线。比较前3个特征点的荷载时程曲线，其极值分别为59.15MPa、0.28MPa和0.18MPa；峰值出现时间分别为0.5 ms、15.1ms和35ms；正压区作用时间分别为2ms、8 ms和10.5ms。本文参考文献［14］中给出了装药密度为1 630kg／m3时TNT炸药所形成空气冲击波的初始值，其初始波后压力为54.35MPa，传播速度为6 450m／s。

图2 爆炸作用下储罐特征点处荷载时程曲线图

从图2－a中可以看出，在A点，超压在爆炸发生后的极短时间内上升到极值59.15MPa，误差9%，满足工程精度要求。其余各点的极值压力时刻随偏离爆心（外罐穹顶顶心）距离的增大而滞后，超压幅值随爆心距的增加而迅速降低，另外超压的正压区作用时间也随着爆心距的增加而延长，各超压变化符合三角形变化规律。在D点，荷载曲线有剧烈的脉动现象，其原因主要是其位于穹顶与罐壁的连接处，结构形式特殊，应力波传播到此处发生了一系列反射作用，在这种反复作用下，混凝土结构发生了疲劳破坏。爆炸初期，冲击波并未以球面波的形式向外传播，这是因为药包形状是立方体装药。然而随着传播距离的不断增大，冲击波逐渐形成典型的球面波形式向外扩张，正压区随着冲击波的传播不断拉宽。

为了进一步验证荷载输入的正确性从而保证本文数值模拟结果的可靠性，现将本文计算的超压结果与4个常用经验公式结果进行对比［14－15］，各经验公式均来自试验结果，如图3所示。

图3 衰减曲线对比图

从图3中可以看出，在近爆区（爆心距小于10 m），由于爆炸空气冲击波流场的复杂性，试验条件的差异性以及研究手段的不同且各公式均未考虑装药位置、装药形状的影响，故给出的衰减关系存在较大的差异性，本文数值模拟结果较各公式稍大，但仍在工程精度范围内（除Henrych公式，其余公式误差均小于10%）；对离爆心较远处的超压衰减规律，本文数值模拟与各公式的结果十分接近。由以上分析结果可以看出本文的荷载输入较各经验公式略有偏大，满足工程要求。

3.2 外罐破坏形态

炸弹在储罐穹顶爆炸后，强冲击波首先在顶心及其周围区域形成半径约为3m的大变形区，该区域混凝土与钢筋单元均发生了较大的变形，在2ms时该处竖向极限位移为613mm。从塑性应变分布的变化过程来看，在2ms时，强冲击荷载作用在穹顶，导致其上层混凝土受压而下层则处于受拉状态，由于抗拉强度远小于抗压强度，下层混凝土首先受拉破坏而失去承载能力。多次模拟分析结果表明，当TNT当量超过900kg后，罐顶处混凝土会发生破坏，出现弹坑，故特大型LNG储罐能够抵抗的最大爆炸当量为900kg。

若出现弹坑，空气冲击波会从弹坑中进入罐体内部，其作用范围可从罐顶到达罐体内部约15m处。特大型LNG储罐的设计液位高度为34.760m，距离爆心约18m，空气冲击波的作用范围与此十分接近。特大型LNG储罐穹顶部位结构简图如图4所示。

图4 储罐穹顶部位结构简图

从图4中可以看出，穹顶结构上部仅由外部钢筋混凝土拱顶和厚度为6mm的罐顶衬板组成，其下悬挂吊杆和铝吊顶，铝吊顶上填充珍珠岩或玻璃棉保温层，铝吊顶距离正常操作液位约500mm。如上所述，如果外罐顶出现爆坑，冲击波进入罐内，虽然其冲击荷载已急速衰减，但冲击波会对铝吊顶上部的保温层造成破坏，引起火灾、爆炸等次生灾害。

3.3 混凝土单元应力响应

图5为储罐特征点的混凝土最大主应力时程曲线。在A点，其最大主应力达到85MPa，但由于混凝土动态强度会比静载强度有较大提高，故其并未发生破坏。B点的应力曲线有两个峰值，第一峰值的最大主应力为23.2MPa，出现在2.5ms时，即爆炸初期，与A点超压峰值时间0.5ms十分接近。其原因在于爆炸初期，爆心处产生的初始空气冲击波作用于穹顶中心A点，导致穹顶其余部位受到一个脉冲作用而承受瞬时弯矩，表现为混凝土应力响应出现一个脉冲峰值。在15.5ms时，B点又出现一个峰值为28.7MPa的峰值，这个峰值是由空气冲击波传播到B点所致，与荷载时程曲线中B点在15ms出现峰值相对应，随后混凝土应力出现大幅回落，并在应力水平10MPa左右出现振荡，这一现象是由负压区的作用引起的。

图5 储罐特征点处混凝土最大主应力时程曲线图

C点在4.5ms时也出现了一个10.7MPa的脉冲峰值，原因同上；在24.8ms出现一个11.9MPa的峰值响应，是由于空气冲击波波阵面传播到C点所致，其后应力急剧下降并伴有高频锯齿振荡［16］。从图5中可以看出D点在80ms以前应力持续小幅振荡，80ms以后为大幅低频响应，如3.1节中所述，这一现象与其特殊的位置和构造有关，爆炸冲击波在此处连续作用于罐体，导致罐体疲劳损伤。总体来说，A、B点虽未发生破坏，但均已处于较高的应力水平，C点至F点混凝土的应力响应均在安全范围之内，这也印证了前述900kg当量为罐体抗爆极限值的结论。

3.4 混凝土单元位移和加速度响应

图6 储罐特征点处混凝土节点位移时程曲线图

图6为储罐特征点A～D处混凝土节点位移时程曲线，E、F两处位移响应较小，此处未予列出。从图6中可以看出，A点穹顶处最大位移达到了458mm，出现在25ms左右，说明穹顶处变形并不是在爆炸瞬间达到了最大位移响应，而是有一个逐渐变形的过程，最后其变形维持在100mm左右波动，证明穹顶各处发生了不可恢复变形，即穹顶混凝土已经发生塑性破坏。B、C、D 3处的最大变形分别为164mm、113mm和27 mm，各处均存在不同程度的残余变形，但其值均较小，不影响结构的整体安全。现行国外规范并未对穹顶挠度限制有详细规定，参考我国《混凝土结构设计规范》［17］中的规定，穹顶处挠度限制应小于l0／400，即80 000／400＝200mm。从以上分析中可以看出，穹顶处挠度已超过该值，虽未引起破坏，但结合其他数据可以看出基本已达到其承载力极限。

图7为储罐特征点A～D处混凝土节点加速度时程曲线。由图7可见，各点加速度峰值均出现在初始加载阶段，而后反应大幅下降，并伴有锯齿形剧烈振荡，振荡现象随着爆心距的增加而愈加剧烈。这一现象的原因为炸药爆炸后给储罐穹顶一个极大的初始加速度荷载，引起各特征点处加速度出现一个脉冲峰值，而后出现的锯齿形振荡则是由储罐的自振引起的。将图5～7与图2的荷载时程曲线相对比，可以看出，冲击波波阵面的传播对混凝土单元的应力影响很大，但对于其加速度响应影响不大。以C点为例，其荷载峰值出现在35ms处，此时其相应的应力出现了一个脉冲峰值，而其加速度响应未有明显波动。

图7 储罐特征点处混凝土节点加速度时程曲线图

3.5 钢筋应力响应

储罐部分特征点B处钢筋最大主应力时程曲线如图8所示，由于爆心处A点的钢筋瞬间即达到设定的应变极限不具代表性，故此处并未示出。图8显示B点的钢筋一直处于较高的应力水平，最大应力已达到400MPa，但是未达到设置的应变极限0.05，故尚未破坏，初始爆炸冲量并未使其到达应力极值，后续的空气冲击波持续作用使其应力逐渐增加。外罐钢筋除爆心处发生破坏外，其余部分即使有瞬间到达应力限值，但未达到设定的应变极限，总体来说，外罐钢筋应力响应较安全。

图8 储罐特征点处钢筋最大主应力时程曲线图

4 结论

笔者应用流固耦合ALE算法和LS－DYNA有限元程序对接触爆炸作用下特大型LNG储罐的动力响应进行了数值模拟分析，根据前期工作分析，选取900 kg TNT当量炸药作为分析载荷。计算结果表明：

1）本文模拟过程中冲击波荷载的峰值和衰减曲线与4种经验公式基本一致，满足工程精度要求。

2）等效质量900kg TNT炸药在储罐穹顶顶心处接触爆炸后，爆炸瞬间强冲击波即造成外罐穹顶处发生局部塑性变形；在穹顶顶心部位形成一个半径为3 m的大变形区。

3）以爆心为中心，半径为25m以内的混凝土和钢筋在整个作用过程中均处于较高的应力水平，但仍能保证对内罐的保护作用，爆炸对半径25m以外的罐体影响不大。

4）若爆炸荷载当量超过900kg，穹顶部位即会出现弹坑，混凝土剥落，钢筋失效，空气冲击波从此弹坑可进入罐体内部，会造成保温层等关键部位的破坏，同时可能造成火灾和后续爆炸等一系列后果，对场区及生命线工程的安全造成危害。

综上所述，通过本文分析可以确定特大型LNG储罐穹顶处大约可抵抗900kg等效当量TNT炸药的爆炸荷载，基本可抵抗常规导弹的袭击，但考虑到生命线工程的重要性，建议在今后的储罐设计中应适当加大穹顶部位混凝土的截面厚度，以提高穹顶半径25m以内罐体的抗爆性能；其次在铝吊顶上部保温层上可设置一层防护钢板，以期在罐顶部位形成多道抗爆防线，以保证冲击波即使破坏外罐，也不会直接作用于保温层和铝吊顶。

［1］ 张瑞甫，翁大根，倪伟波，等.特大型LNG储罐抗（减）震研究发展综述［J］.结构工程师，2010，26（5）：164－171.ZHANG Ruifu，WENG Dagen，NI Weibo，et al.State of the art of seismic analysis for extra－large LNG storage tanks［J］.Structural Engineers，2010，26（5）：164－171.

［2］ PARFOMAK P W，FLYNN A M.Liquefied natural gas（LNG）import terminals：Siting，safety and regulation［R］.Washington DC：Congressional Research Service，2004.

［3］ 王冰，陈学东，王国平.大型低温LNG储罐设计与建造技术的新进展［J］.天然气工业，2010，30（5）：108－112.WANG Bing，CHEN Xuedong，WANG Guoping.Design of large low－temperature LNG storage tanks and new progress in its construction technology［J］.Natural Gas Industry，2010，30（5）：108－112.

［4］ 刘勇.液化天然气的危险性与安全防护［J］.天然气工业，2004，24（7）：105－107.LIU Yong.Dangers and safeguards of LNG［J］.Natural Gas Industry，2004，24（7）：105－107.

［5］ 余爱萍，王远功，翁智远.冲击波对核反应堆安全壳的动力响应研究［J］.爆炸与冲击，1992，12（3）：219－227.YU Aiping，WANG Yuangong，WENG Zhiyuan.Dynamic response of the shock wave reactor containment［J］.Journal of Vibration and Shock，1992，12（3）：219－227.

［6］ 张小鹏，葛飞，邢怀念，等.大型爆炸防护结构的动态应力测试与分析［J］.爆炸与冲击，2010，30（1）：101－104.ZHANG Xiaopeng，GE Fei，XING Huainian，et al.Dynamic stress measurement of a large explosion containment vessel［J］.Explosion and Shock Waves，2010，30（1）：101－104.

［7］ HUANG Y.Performance assessment of conventional and base－isolated nuclear power plants for earthquake and blast loadings［D］.Buffalo：State University of New York at Buffalo，2008.

［8］ British Standard Institution.Design and manufacture of site built，vertical，cylindrical，flat－bottomed steel tanks for the storage of refrigerated，liquefied gases with operating temperatures between 0 ℃ and－165 ℃ ［S］.London：British Standard Institution，2006.

［9］ UK Atomic Energy Agency.Guidelines for the design and assessment of concrete structures subjected to impact［S］.Zurich：UK Atomic Energy Agency，1990.

［10］ 岳宝增，李笑天.ALE有限元方法研究及应用［J］.力学与实践，2002，24（2）：7－11.YUE Baozeng，LI Xiaotian.Study of the ALE finite element method and its applications［J］.Mechanics and Engineering，2002，24（2）：7－11.

［11］ 岳宝增，彭武，王照林.ALE迎风有限元法研究进展［J］.力学进展，2005，35（1）：21－29.YUE Baozeng，PENG Wu，WANG Zhaolin.Research progress of ALE upwind finite element methods［J］.Advances in Mechanics，2005，35（1）：21－29.

［12］ 杜修力，廖维张，田志敏，等.炸药爆炸作用下地下结构的动力响应分析［J］.爆炸与冲击，2006，26（5）：474－480.DU Xiuli，LIAO Weizhang，TIAN Zhimin，et al.Dynamic response analysis of underground structures under explosion－induced loads［J］.Explosion and Shock Waves，2006，26（5）：474－480.

［13］ TALASLIDIS D G，MANOLIS G D，PARASKEVOPOULOS E，et al.Risk analysis of industrial structures under extreme transient loads［J］.Soil Dynamics and Earthquake Engineering，2004，24（6）：435－448.

［14］ BRODE H L.Numerical solution of spherical blast waves［J］.Journal of Applied Physics，1955，26（6）：35－59.

［15］ HENRYCH J.The dynamics of explosion and its use［M］.London：Elsevier Scientific Publishing Company，1979.

［16］ 王天运，任辉启，王玉岚.接触爆炸荷载作用下核电站安全壳的动力响应分析［J］.核动力工程，2005，26（2）：187－191.WANG Tianyun，REN Huiqi，WANG Yulan.Analysis of dynamical response of nuclear power plant concrete containment under contact blast load［J］.Nuclear Power Engineering，2005，26（2）：187－191.

［17］ 中国建筑科学研究院，清华大学，天津大学，等.混凝土结构设计规范［S］.北京：中国建筑工业出版社，2010.China Academy of Building Research，Tsinghua University，Tianjin University，et al.Code for design of concrete structures［S］.Beijing：China Building Industrial Press，2010.