考虑初始应力的混凝土重力坝水下爆炸毁伤特性研究

王高辉， 高 政， 卢文波， 陈叶青， 吕林梅

(1.武汉大学 水资源与水电工程科学国家重点实验室，武汉 430072； 2.军事科学院 国防工程研究院，北京 100850)

近些年，国内外恐怖袭击和意外爆炸事件频繁发生，重要结构的抗爆设防给工程界和学术界带来巨大挑战。这些恐怖袭击的重要目标往往处于复杂的初始应力状态下，如深埋洞室存在较大的初始地应力，典型构件结构承受上部结构的重力荷载，大坝结构承受重力、静水压力、泥沙压力和扬压力等初始荷载。当结构遭受爆炸冲击荷载作用时，由于初始静应力的存在，结构的毁伤发展过程及毁伤模式将发生一定的改变。因此，开展初始应力对结构抗爆性能的影响将有利于进一步认识结构的爆炸毁伤机理。

目前有关初始应力对结构抗爆性能的影响研究主要集中在典型构件、地下洞室以及采矿工程等方面。如Xie等[1]对地应力和爆炸荷载耦合作用下的巷道损伤进行了研究，研究表明地应力的大小会影响爆破荷载下岩体的破坏特征；Yang等[2]的模型试验研究表明单轴压力的大小和方向会对预裂爆破产生显著的影响；Dong等[3]的研究得出柱面应力波在完整岩石中的物理衰减随着围压的增加先减小后增大；陈明等[4]研究了初始应力对隧洞开挖爆生裂隙区的影响，发现在压剪破坏破坏模式下，爆生裂隙区比例半径随地应力的增大而减小；方秦等[5]认为，爆炸作用下，RC柱轴力越大越易发生弯曲破坏；瞿海雁等[6]的研究表明，在抵抗侧向爆炸荷载时，轴向压力可以提高立柱的抗剪刚度和抗压刚度；杨立云等[7]对初始应力场对爆生裂纹行为演化效应进行了试验研究，结果表明垂直于裂纹扩展方向的压应力会降低裂纹尖端的应力集中程度，从而阻碍裂纹的扩展；白羽等[8]对不同地应力条件下双孔爆破进行了数值模拟研究，得出随着埋深增加和初始地应力加大，裂纹扩展半径和裂纹区面积减小，地应力对爆破至裂的抑制作用明显。以上研究表明，初始应力状态对构件的抗爆性能、爆炸毁伤特性有一定影响。

高坝由于其显著的政治经济效益，无疑成为局部战争和恐怖袭击的重点打击目标，大坝一旦失事，不仅使工程本身遭受巨大损失，而且将给下游生命财产带来巨大灾难。因此，近年来有关大坝抗爆防护的研究越来越受到重视。如张社荣等[9]对不同爆炸方式作用下混凝土重力坝的抗爆性能和可能破坏模式进行了研究；王高辉等[10]通过建立混凝土重力坝水下爆炸全耦合模型，分析了水下爆炸荷载下混凝土重力坝的破坏效应；Linsbauer等[11]的研究给出了库底爆炸荷载作用下混凝土重力坝挡水坝段的动力响应和稳定性；李麒等[12]研究了库前水位对混凝土重力坝抗爆安全性能的影响。然而，由于水下爆炸冲击荷载持续时间短、峰值大，其对大坝结构的冲击毁伤破坏往往在短时间内完成，上述研究中均忽略了初始应力的影响。

本文基于显式动力学有限元分析程序，采用显隐转换技术实现了初始应力从静力学分析到动力学分析的平稳过渡，通过对比分析不考虑初始应力状态、只考虑重力以及综合考虑重力和水压力三种不同初始应力状态下大坝的水下爆炸毁伤效应，研究了不同初始应力状态下大坝的水下爆炸毁伤特性以及重力、水压力对混凝土重力坝抗爆性能的影响，同时探讨了初始应力对混凝土重力坝接触和非接触爆炸毁伤特性的影响。研究成果将有助于进一步认识混凝土重力坝的毁伤机理和破坏模式。

1 混凝土重力坝初始应力的显隐转换分析

为了探究初始应力对大坝水下抗爆性能的影响，选取某混凝土重力坝典型坝段为分析对象，最大坝高为100 m，坝顶宽度为16 m，坝头高度为20 m，坝底宽度为73.6 m，假定大坝上游水深90 m，单坝段宽度为15 m。为了获得大坝的初始应力，建立了混凝土重力坝的有限元模型如图1所示。值得注意的是：结构静力计算模型和动力计算模型的有限元网格必须完全一致。考虑到后续会进行爆炸动力学分析计算，故对距离爆源位置较近的网格进行加密，将网格尺寸控制在0.4 m左右，单坝段网格总数449 040个；由于坝基岩体不作为主要研究对象，故基岩网格划分较粗，以减少计算时间；为使大坝与基岩之间传力正常，将坝体网格与基岩网格共节点。

据已有的研究成果，在解决与时间无关的静力学问题时常使用隐式求解方法，而与时间相关的动力学问题则会选用显式求解方法进行求解。大坝初始应力的求解可以看做是与时间无关的静力学问题，而大坝受到水下爆炸打击则是与时间相关的动力学问题。若要求解考虑初始应力的大坝遭水下爆炸打击的问题，则必须使大坝的初始应力状态由隐式分析平稳过渡到显式分析中，因此须使用显隐转换技术。本文利用ANSYS Workbench良好的交互性来实现初始应力状态由隐式到显式的转换，进而通过显式动力学分析程序对水下爆炸问题进行求解。

如图2所示，建立从静力学模块到显示动力学模块的计算分析流程，以达到显隐转换的目的；依据图1所示的静力计算模型，在静力学模块中对坝体施加重力及水压力，在基岩的底部施加全约束以限制其位移和转角，坝体和基岩的截断面处施加法向支撑约束以限制法向位移；对坝体及基岩进行结构静力计算，图3给出了大坝的竖向应力云图和顺水流向位移云图，可以看出，在重力和水压力的作用下，大坝上游面坝踵部位受到的竖向压应力最大约3.2 MPa，大坝顶部的顺水流向位移达到最大约2 mm；通过将静力结果转入到显式动力有限元中，得到了显式动力初始应力结果如图4所示，可以看出显示动力分析中的初始应力与静力分析结果完全一致，初始应力由隐式分析平稳过渡到了显式分析中。至此实现了大坝初始应力的显隐转换。

图2 初始静力分析及显隐转换Fig.2 Initialstatic analysis and explicit-to-implicit conversion

(a) 大坝竖向正应力云图(kPa)

(b) 大坝顺水流向位移云图(mm)图3 大坝静力分析结果Fig.3 Static analysis results of the dam

(a) 大坝竖向正应力云图(kPa)

(b) 大坝顺水流向位移云图(mm)图4 显式动力分析中大坝的初始应力及位移Fig.4 The initial stress and displacement of dam in explicit dynamic analysis

2 考虑初始应力的混凝土重力坝水下爆炸全耦合模型的建立及验证

2.1 水下爆炸耦合模型

上述显隐转换过程将大坝初始应力作为动力计算的初始计算条件，为了探讨水下爆炸冲击荷载作用下大坝的非线性动态响应及毁伤破坏过程，在含初始应力条件的模型基础上，建立库水、炸药和空气单元，由于静水压力已通过荷载方式施加，故显式动力计算模型中只需要考虑库水与大坝-地基系统的动态耦合作用。如图5所示，空气、库水和TNT均采用Euler网格进行建模；在计算中采用自动时间步长，以便捕获爆炸冲击波的强间断性，同时可提高数值计算效率，该程序在自动选择时间步长时，与有限元模型的网格尺寸大小相关，已有研究表明[13]，在水下爆炸冲击数值分析中，网格尺寸采用炸药半径尺寸的1/3时，其计算精度基本可以满足工程要求，因此将炸药附近区域的网格加密，网格尺寸控制在0.2 m之内；库水与大坝和地基之间采用流固耦合算法，需要说明的是，在流固耦合算法中模型中的水体不会对大坝和基岩产生水压力的作用；为了缩短计算时间，本文在建模过程中只截取了部分的库水、空气和基岩，为了模拟库水、空气和基岩的半无限区域，对其所有的人工截断面施加无反射边界条件；同时为使冲击波在坝段两侧不发生反射，亦对坝段两侧施加无反射边界；经显隐转换后，基岩底部的固定约束以及基岩四周和坝体两侧的法向支撑约束会由隐式分析导入到显式分析中，因此不需要对坝体、基岩重复施加固定约束和法向支撑约束。

图5 考虑初始应力的混凝土重力坝水下爆炸全耦合模型Fig.5 Fully-coupledcalculation model of concrete gravity dam subjected to underwater explosion considering initial stress

为了充分研究考虑初始应力的混凝土重力坝水下爆炸毁伤特性，本文分为水下非接触打击和水下触坝打击两种打击工况进行分析研究。其中，水下非接触打击时，炸药选用TNT炸药1 000 kg，爆心距为5 m，起爆水深选为水下30 m；水下接触打击时，选用TNT炸药1 000 kg，起爆水深选为水下30 m。

2.2 混凝土动态损伤本构模型

本文中，坝体混凝土和基岩均选用RHT模型。RHT本构模型是由Riedel等[14]基于HJC本构模型[15]提出的。该模型综合考虑了混凝土材料的大应变、高应变率、高压效应、应变硬化、应变软化和应力偏量第三不变量的影响，能很好地模拟混凝土从弹性到失效的过程。如图6所示，RHT模型引入了弹性极限面、失效面和残余失效面3个极限面，分别对混凝土的初始屈服强度、失效强度和残余强度的变化规律进行了描述。

图6 RHT模型的3个失效面Fig.6 Three failure surfaces of RHT constitutive model

弹性极限面Ye用于考虑材料的应变硬化效应，其通过失效面来确定，可以表示为：Ye=Yf·ω·Fcap(p)，其中，Yf是RHT模型的破坏面，ω是弹性拉伸强度或弹性压缩强度与相应的极限强度之比，Fcap(p)是用于限制静水压力作用下弹性偏应力的盖帽函数，其取值在0～1之间。

α是压缩应变率因子，δ是拉伸应变率因子。

残余失效面Yr是用来描述完全粉碎材料强度的，可表示为：Yr=B·(p*)M，其中，B为残余失效面常数；M为残余失效面指数。

如表1所示为坝体混凝土的材料参数[16]；基岩材料密度为2.66 g/cm3，抗压强度为40 MPa，剪切模量为2.19×107kPa，体积模量为2.57×107kPa，其余参数与坝体混凝土一致。

表1 坝体混凝土材料参数Tab.1 Concrete material parameters of dam

2.3 状态方程

本文中，炸药使用JWL状态方程[17]进行描述：p=A1(1-ω/R1V)e-R1V+B1(1-ω/R2V)e-R2V+ωE0/V。其中，p为爆轰压力，V表示爆轰产物相对体积,初始比内能E0=6.0 GPa，特征参数A1=373.77 GPa、B1=3.75 GPa、R1=4.15、R2=0.90、材料常数ω=0.35,e表示底数。

水体采用Polynomial状态方程[17]，其在不同压缩状态下具有不同的形式。当水压缩时(μ>0)，状态方程为：p=A1μ+A2μ2+A3μ3+(B0+B1μ)ρ0e；当水膨胀时(μ<0)，状态方程为：p=T1μ+T2μ2+B0ρ0e；当水既不压缩也不膨胀时，可简化为：p=A1B0ρ0e。其中，p为水中的压力；μ为压缩比，μ=ρ/ρ0-1；e为水的内能；水密度ρ0=1 g/cm3；A1=T1=2.2×106kPa，A2=9.54×106kPa，A3=1.46×106kPa，B0=B1=0.28，T2=0。

空气采用Ideal Gas状态方程：p=(γ-1)ρe。空气密度ρ=1.225 kg/m3，空气初始内能e=2.068×105kJ/kg，材料常数γ=1.4。

3 考虑初始应力的大坝水下接触爆炸毁伤特性

大坝在正常运行时主要会受到重力和水压力的作用。为了研究初始应力对大坝水下爆炸毁伤特性的影响，本文通过数值计算研究了三种不同初始应力状态下大坝水下爆炸的毁伤特性，分别为不考虑初始应力、只考虑重力以及综合考虑重力和水压力作用三种情况。

图7给出了不考虑初始应力的大坝水下接触爆炸的毁伤过程。可以看出，触坝打击时巨大的爆炸冲击荷载直接作用在坝体上游面，对坝体上游面造成了爆炸成坑破坏；随着爆炸冲击波传入坝内，坝身处出现冲切破坏，其不断由坝体上游表面斜向坝体内部发展，最终几乎贯穿坝体上部；爆炸冲击波传入大坝内部后会衰减为压缩应力波，压缩波传至坝体下游面后又被反射为拉伸波，因此坝体下游折坡处产生了拉伸破坏，并逐渐向上游发展，直至贯穿坝体；最终由于坝体的整体响应，坝踵位置出现不断向下游发展的整体拉伸破坏；由于大坝整体受到巨大的冲击荷载，大坝上游面局部受弯，因此产生了弯曲破坏。

(a) t=10 ms

(b) t=20 ms

(c) t=30 ms

(d) t=50 ms图7 不考虑初始应力的大坝水下接触爆炸毁伤过程Fig.7 Damage evolution process of the dam without initial stress subjected to underwater explosion

图8给出了考虑重力和水压力作用下大坝的水下接触爆炸毁伤过程。可以看出，考虑初始应力条件时大坝的毁伤破坏模式和部位与不考虑初始应力工况下基本一致，但毁伤发展过程及破坏形态有所不同。

(a) t=10 ms

(b) t=20 ms

(c) t=30 ms

(d) t=50 ms图8 重力和水压力作用下大坝的水下接触爆炸毁伤过程Fig.8 Damage evolution process of the dam with initial stress subjected to underwater explosion

图9分别给出了三种不同应力状态下大坝的水下接触爆炸毁伤形态。可以看出，在触坝打击情况下，大坝的主要毁伤模式有上游面的爆炸成坑破坏、坝身处的爆炸冲切破坏、下游折坡处的拉伸破坏以及坝踵位置的整体拉伸破坏。

(a) 不考虑初始应力时

(b) 综合考虑重力和水压力时

(c) 只考虑重力时图9 不同初始应力状态下大坝水下接触爆炸的毁伤形态Fig.9 Damage of dams subjected to underwater contact explosion in different initial stress states

如图9(a)和(b)所示，不考虑初始应力和综合考虑重力和水压力情况下大坝的毁伤形态和范围有所不同。其中，当不考虑初始应力状态时，大坝上游面爆炸成坑破坏的爆坑直径约为7.9 m，爆坑深度约3.8 m；坝身爆炸冲切破坏范围约22.2 m，且冲切破坏基本贯穿坝体上部；由于坝体的整体响应，坝踵位置的整体拉伸破坏达29.5 m；由于边界效应，坝段的两侧出现了一定的毁伤。当综合考虑重力和水压力作用后，大坝上游面混凝土处于初始压实状态，故上游面的爆炸成坑破坏较不考虑初始应力时有所减小，其爆炸成坑破坏的爆坑直径约为6.1 m，减小约22.8%，爆坑深度约2.8 m，减小约26.3%；同时由于考虑大坝重力和水压力后，坝身受到合力方向为斜向下，与大坝冲切破坏的方向基本一致，故爆炸冲切破坏相较于不考虑初始应力时有所增大，增大了约8.1%，其毁伤范围约24.0 m；坝踵位置整体受压，故整体拉伸破坏相对于不考虑初始应力情况减小了8.5%，毁伤范围约27.0 m；相较于不考虑初始应力的情况，此种情况下坝体两侧的毁伤有较明显加重。由此可以看出，在水下触坝打击情况下，大坝的初始应力主要会对坝体上游面爆炸成坑破坏以及坝身爆炸冲切破坏的范围及形态产生一定影响。

为了进一步研究重力和水压力分别对坝体水下接触爆炸毁伤特性的影响，图9(c)给出了只考虑重力时大坝的水下接触爆炸毁伤形态。可以得到，此种情况下，大坝上游面的爆坑直径约为6.1 m，爆坑深度约2.3 m；爆炸冲切破坏长度达24.6 m，基本贯穿坝体上部；坝踵处的整体拉伸破坏范围约19.1 m。通过对比图9(a)～(c)可得，在水下接触爆炸打击情况下，重力会抑制大坝上游面爆炸成坑破坏的形成，使爆坑直径和爆坑深度均有一定程度的减小，但其会加重坝身位置的爆炸冲切破坏；而水压力则使大坝上游面爆炸成坑破坏有所加重，其中水压力对爆炸成坑破坏的影响主要体现在使爆坑深度加深；可以看出，在水下接触爆炸条件下，重力对大坝毁伤特性的影响起主导作用。

为探究大坝的动态响应，研究中选取了五个监测点如图10(a)所示，图10(b)～(f)给出了触坝打击下不同监测点的顺水流向速度-时间曲线。可以看出，触坝打击情况下，重力和水压力对坝体动态响应的影响主要体现在坝踵和坝体上部，坝体中部受影响较小，这是由于爆源位置接近坝体中部，触坝打击情况下坝中位置受到的爆炸冲击荷载远远大于重力和水压力的作用。因此，重力和水压力对坝体动态响应的影响主要在坝顶位置(监测点5)和坝踵位置(监测点1)。

(a) 测点布置

(b) 测点1

(c) 测点2

(d) 测点3

(e) 测点4

(f) 测点5图10 水下接触爆炸下大坝各监测点的顺水流向速度时程曲线Fig.10 Dynamic time history response curves of dam subjected to underwater contact explosion

首先分析坝踵位置的监测点1，如图10(b)所示，从峰值速度来看，不考虑初始应力时的峰值速度约0.25 m/s，只考虑重力时的峰值速度约0.29 m/s，而综合考虑重力和水压力时监测点位置的峰值速度达到0.32 m/s，由此可以得出重力和水压力均使坝踵位置响应的峰值速度略有增大，比不考虑初始应力时增大了约28.0%。

监测点5位于坝顶位置，由图10(f)可以看出，不考虑初始应力时，坝顶的峰值速度约0.85m/s；当只考虑重力时，坝顶的峰值速度约0.8m/s，且此时响应的速度整体低于不考虑初始应力情况；综合考虑重力和水压力作用时，峰值速度约0.95m/s。可以看出三种情况下坝顶的响应均很剧烈，重力作用减轻了坝顶位置的响应，而水压力加重了坝顶的响应，且考虑重力和静水压力下响应的峰值速度相比于不考虑初始应力状态时增大了约11.8%。

4 考虑初始应力的大坝水下非接触爆炸毁伤特性

图11给出了三种初始应力状态下大坝水下非接触爆炸的毁伤形态(爆心距5 m，TNT当量1 000 kg)。由图11可知，大坝遭水下非接触爆炸打击时，其主要毁伤形式有坝体上游面的爆炸冲切破坏、下游折坡处的拉伸破坏以及坝踵位置的整体拉伸破坏。

(a) 不考虑初始应力时

(b) 综合考虑重力和水压力时

(c) 只考虑重力时图11 不同初始应力状态下大坝水下非接触爆炸的毁伤形态Fig.11 Damage of dams subjected to underwater non-contact explosion in different initial stress states

当不考虑初始应力状态时，坝体未形成明显的爆炸冲切破坏；坝踵位置的拉伸破坏约26.1 m；且坝体两侧因边界效应会出现一定的毁伤。而当综合考虑重力和水压力作用时，坝身正对爆源位置出现了范围约11.2 m的爆炸冲切破坏，其相较于不考虑初始应力时有明显增大，这是由于考虑重力和水压力后该位置的合力增大，达到了此处混凝土的破坏强度；坝踵整体拉伸破坏约23.3 m，相对于不考虑初始应力时减小了约10.7%；且坝体两侧出现的毁伤明显比不考虑初始应力时更加严重。可以看出，在水下非接触打击下，初始应力主要会影响坝身冲切破坏的形态及范围。由图11(a)～(c)可得，在水下非接触爆炸情况下，重力对坝体毁伤特性的影响依然起主导作用。但对比图9和图11可以看出，水下非接触爆炸下初始应力对坝身冲切破坏的影响大于水下接触爆炸情况，这是由于水下非接触爆炸时大坝受到的爆炸冲击荷载远远小于水下接触爆炸，爆炸冲击荷载越小则初始应力的作用越大，故在水下非接触爆炸时初始应力对坝身冲切破坏的影响更加显著。

图12(b)～(f)给出了水下非接触爆炸时大坝各监测点的顺水流向速度时程曲线。对比各监测点的速度时程曲线可知，水下非接触爆炸下，由于爆源距离上游坝面有一段距离，坝体受到的爆炸冲击荷载相对较小，故重力和水压力对坝体响应的影响相对明显，且主要体现在坝体中部位置(监测点2～4)。在测点3位置，不考虑初始应力时峰值速度为1.30 m/s；只考虑重力时的峰值速度为1.35 m/s；综合考虑重力和水压力情况下的峰值速度约1.70 m/s，相对于不考虑初始应力时增大了30.8%。而测点2和测点4处考虑重力和水压力后的峰值速度相较于不考虑初始应力时分别增大了约27.0%和33.6%。考虑重力和水压力后，坝中部位三测点的峰值速度相较于不考虑初始应力时平均增长了约30%。综上可得，在水下非接触爆炸打击下，初始应力对坝体中部的响应影响较大，且重力和水压力均使坝体中部响应的峰值速度略有增大。

(a) 测点布置

(b) 测点1

(c) 测点2

(d) 测点3

(e) 测点4

(f) 测点5图12 水下非接触爆炸下大坝各监测点x向速度时程曲线Fig.12 Dynamic time history response curves of dam subjected to underwater non-contact explosion

5 结 论

本研究建立了考虑初始应力状态的混凝土重力坝水下接触/非接触爆炸毁伤分析模型，探讨了重力和静水压力对大坝毁伤发展和非线性动态响应的影响。主要结论如下：

(1) 初始应力对混凝土重力坝上游坝面的爆炸成坑破坏、坝身的冲切破坏、坝踵的拉伸破坏等破坏特性和破坏范围均有较明显的影响，在混凝土重力坝抗爆性能评估时应适当考虑初始应力的影响，尤其是重力的影响。

(2) 在水下触坝打击条件下，考虑重力和静水压力时，大坝上游面爆坑范围和坝踵开裂长度有所减小，其中，爆坑深度减小约26.3%，爆坑直径减小约22.8%，坝踵开裂破坏长度减小了8.5%；但加大了坝身的冲切破坏长度，约增大了8.1%，且大坝顶部和坝踵部位的响应峰值速度有所增大。

(3) 在水下非接触打击条件下，考虑初始应力时坝踵的开裂长度有所减小，约10.7%；然而大坝正对爆源部位的冲切破坏长度有所增大，且这种增加较水下接触爆炸情况更加显著；坝身中部质点振动峰值速度平均约增大了30%。