基于二階雙漸近法的雙層圓柱殼在水下爆炸作用下的鞭狀運(yùn)動(dòng)＊

劉云龍，汪 玉，張阿漫

（1.武漢第二船舶設(shè)計(jì)研究所，湖北 武漢 430205；2.哈爾濱工程大學(xué)船舶工程學(xué)院，黑龍江 哈爾濱150001；3.中國人民解放軍92857部隊(duì)，北京100161）

圓柱殼作為潛艇典型的結(jié)構(gòu)形式，常被用作潛艇結(jié)構(gòu)強(qiáng)度機(jī)理性研究的對象。在很多情況下對圓柱殼結(jié)構(gòu)抗沖擊性能深入研究的成果可以推廣到實(shí)際潛艇上。因此，本文中以雙層圓柱殼結(jié)構(gòu)為研究對象，對水下爆炸載荷下圓柱殼結(jié)構(gòu)的鞭狀運(yùn)動(dòng)特性進(jìn)行分析，旨在為潛艇抗沖擊設(shè)計(jì)和研究提供參考。水下爆炸問題具有瞬態(tài)強(qiáng)非線性特點(diǎn)，其涉及的流固耦合問題一直是艦船抗沖擊研究的熱點(diǎn)和難點(diǎn)。早在20世紀(jì)60年代，研究人員就嘗試通過解析的方法求解一些簡單的沖擊問題。比較有代表性的是H.Huang等早期的工作，其中得到了一些球殼以及圓柱殼等簡單結(jié)構(gòu)在平面階躍波［1－2］及球面波［3］作用下的彈性響應(yīng)。盡管此類方法對于實(shí)際復(fù)雜結(jié)構(gòu)并不適用，但在后來的研究中通常用來驗(yàn)證理論和數(shù)值模型的精度和有效性。隨著計(jì)算機(jī)性能的提高和邊界元的發(fā)展，流固耦合問題的數(shù)值求解技術(shù)得到了很大的提高。對現(xiàn)代水下沖擊流固耦合研究影響較大的是T.L.Geers在1971年提出的雙漸近（double asymptotic approximation，DAA）法［4］。該方法分別采用平面波假設(shè)和勢流假設(shè)對水下沖擊問題的前期和后期響應(yīng)進(jìn)行近似，而中間頻段則通對兩者進(jìn)行匹配得到。隨后的研究中，T.L.Geers［5］對DAA法［4］進(jìn)行了改進(jìn)，考慮了結(jié)構(gòu)曲率的影響以及前后期的非線性匹配問題，形成了二階雙漸近（second－order doubly－asymptotic approximation，DAA2）法。DAA2法已被廣泛應(yīng)用于水下沖擊防護(hù)研究。針對潛艇結(jié)構(gòu)聲學(xué)覆蓋層問題，劉建湖［6］對DAA2法進(jìn)行改進(jìn)，考慮波在聲學(xué)覆蓋層中的反射和透射特性，將該過程分為聲學(xué)反射階段和DAA2階段分別求解，發(fā)現(xiàn)用傳統(tǒng)DAA2法處理聲學(xué)覆蓋層問題時(shí)低估了沖擊波載荷的作用。隨后，孫士麗［7］、王詩平等［8］通過計(jì)入非線性的二階速度項(xiàng)，對DAA2法進(jìn)行進(jìn)一步改進(jìn)，得到了非線性雙漸近（non－linear second－order doubly－asymptotic approximation，NDAA2）法。宗智等［9］結(jié)合DAA2法和有限元法，計(jì)算了氣泡脈動(dòng)載荷下簡單船體結(jié)構(gòu)的動(dòng)態(tài)響應(yīng)和應(yīng)力集中現(xiàn)象。董海等［10］針對潛艇的鞭狀運(yùn)動(dòng)問題，將DAA2法和船體梁模型結(jié)合在一起，采用DAA2法計(jì)算流體載荷，用變截面梁模型計(jì)算結(jié)構(gòu)響應(yīng)，對結(jié)構(gòu)的時(shí)域和頻域響應(yīng)均進(jìn)行了深入分析。然而，由于DAA2方程是基于無限域三維波動(dòng)方程推導(dǎo)的，因此對內(nèi)流場問題不適用。盡管M.A.Sprague等、T.L.Geers等對經(jīng)典DAA2法進(jìn)行修改并推廣到內(nèi)域問題11－13，但在實(shí)際計(jì)算中經(jīng)典DAA2法仍有較多限制。因此，本文中對雙層圓柱殼結(jié)構(gòu)舷間水對應(yīng)的內(nèi)流場采用聲固耦合求解，而對外流場采用DAA2法求解，兩者在非耐壓殼處耦合，形成適用于雙層圓柱殼結(jié)構(gòu)的流固耦合分析方法，并據(jù)此分析水下爆炸載荷下雙層圓柱殼結(jié)構(gòu)的鞭狀運(yùn)動(dòng)特性。

圖1 DAA2法后期時(shí)域近似示意圖Fig.1 Scheme of later approximation in time domain of DAA2

1 理論和數(shù)值模型

1.1 雙漸近法理論模型

關(guān)于DAA2法，對前期近似采用在平面波假設(shè)下的流固耦合模型。在該假設(shè)下，任意時(shí)刻時(shí)邊界上的流體動(dòng)壓力不受附近其他點(diǎn)的影響，推導(dǎo)過程相對簡單，本文不再贅述。對后期近似，有不同的推導(dǎo)方法。王詩平等［14］通過將延遲積分方程進(jìn)行拉普拉斯變換，并在拉氏空間內(nèi)對方程進(jìn)行簡化，最后通過拉氏逆變換變回時(shí)域空間，為DAA2法提供了新的切入點(diǎn)。采用常規(guī)格林函數(shù)法，重新推導(dǎo)DAA2法后期近似公式。該推導(dǎo)方法更簡便，物理意義更明確，旨在為DAA2法后續(xù)研究提供新的視角。假設(shè)流場速度勢φ滿足波動(dòng)方程：

式中：cf為流體中的聲速。若取格林函數(shù)，則有：

式中：R＝p－q為流場中從源點(diǎn)q指向場點(diǎn)p的矢量。于是，令式（1）×G－式（2）×φ，得：

對式（3）在整個(gè)流域V內(nèi)積分并應(yīng)用格林第2公式可得：

式中：S為流場V的邊界，n為其單位法向量。式（4）即為求解波動(dòng)方程的DBEM（domain boundary element method）積分方程。對式（4）左端積分，可采用常規(guī)邊界元法離散實(shí)現(xiàn)，而右側(cè)體積分項(xiàng)，由于積分域是自由場中的無界空間，只能進(jìn)行近似。根據(jù)DAA2法的思想，考慮t cf?L時(shí)對方程進(jìn)行近似，其中L為結(jié)構(gòu)的特征尺度，此時(shí)被結(jié)構(gòu)擾動(dòng)的流場可認(rèn)為是一個(gè)半徑為t cf的球體，結(jié)構(gòu)的影響可看作球體中心的點(diǎn)源和偶極子的貢獻(xiàn)。由于點(diǎn)源和偶極子的誘導(dǎo)速度勢分別按照進(jìn)行衰減，因此略去偶極子影響所產(chǎn)生的二階項(xiàng)，結(jié)構(gòu)對流場的貢獻(xiàn)相當(dāng)于球體中心一個(gè)點(diǎn)源的影響。于是根據(jù)波動(dòng)方程三維空間的時(shí)域基本解，流場中任意點(diǎn)的速度勢可以近似表示為：

式中：m （t）為源強(qiáng)的任意函數(shù)。于是式（4）右端的體積分可化為：

式中：C為待定常數(shù)。為確定C的取值，將其化為積分形式代入式（4），當(dāng)t＝0時(shí)，顯然C＝0。記u·n＝將式（6）代入式（4），有：

可見式（7）即為DAA2后期近似積分方程。通過平面波假設(shè)，沖擊問題的前期近似可表示為［5，14］：

式中：χ為邊界的局部曲率。將前期近似和后期近似在頻域內(nèi)匹配，整理可得DAA2邊界積分方程［5，14］：

1.2 水下爆炸載荷

水下爆炸氣泡脈動(dòng)是一個(gè)強(qiáng)非線性的物理化學(xué)過程，針對不同的問題已提出了很多氣泡模型。其中T.L.Geers等基于DAA2法提出的單個(gè)水下爆炸氣泡脈動(dòng)模型［15］考慮了流場的可壓縮性和氣泡內(nèi)部氣體的聲學(xué)特性，因此在氣泡各次脈動(dòng)過程中，能夠真實(shí)反映氣泡最大半徑和氣泡脈動(dòng)壓力峰值的衰減。另外考慮了氣泡中心的上浮，并通過以往水下爆炸氣泡實(shí)驗(yàn)數(shù)據(jù)對上浮中的阻力進(jìn)行了匹配，得到結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好。在隨后的研究中，T.L.Geers等對Geers－Hunter模型進(jìn)行了一些優(yōu)化調(diào)整［16］。采用Geers－Hunter模型計(jì)算得到氣泡的脈動(dòng)和遷移運(yùn)動(dòng)后，則可通過以下方法計(jì)算流體中任意一點(diǎn)的流體速度勢［17］：

圖2 水下爆炸氣泡脈動(dòng)載荷計(jì)算示意圖Fig.2 Scheme of the calculation of underwater explosion bubble pulsation load

式中：右側(cè)第1項(xiàng)為氣泡脈動(dòng)等效的點(diǎn)源所誘導(dǎo)的速度勢，第2項(xiàng)為氣泡遷移等效的偶極子所誘導(dǎo)的速度勢，后面2項(xiàng)為自由面鏡像的氣泡所誘導(dǎo)的速度勢。其中e1和e2分別為源強(qiáng)和偶極子強(qiáng)度；kf為自由面的反射因數(shù)，由于本文中分析深水中結(jié)構(gòu)對水下爆炸的動(dòng)態(tài)響應(yīng)，因此不考慮自由面影響，kf取值為0；r1和r2分別為氣泡中心和鏡像氣泡中心的距離；θ1和θ2分別為2連線與鉛垂線的夾角，見圖2。

然后通過非定常伯努利方程得到流場任意位置壓力。氣泡等效的變強(qiáng)度點(diǎn)源和偶極子強(qiáng)度需要滿足氣泡表面的連續(xù)性條件，可分別表示為：

式中：a為氣泡等效半徑，vm為氣泡垂向遷移速度，d為水深。

1.3 BEM與FEM耦合接口

對FEM（finite element method）求解器選擇ABAQUS／Explicit，可通過時(shí)間顯式積分的方法計(jì)算結(jié)構(gòu)在外載荷作用下的響應(yīng)。由于不需要迭代過程，計(jì)算一般不存在收斂性問題，而時(shí)間積分的誤差對于爆炸沖擊這種瞬態(tài)分析來講通常為小量，不會(huì)對結(jié)果產(chǎn)生明顯影響。

圖3中的數(shù)字表示時(shí)間增量步，在同一個(gè)時(shí)間增量步中，F(xiàn)EM計(jì)算得到的結(jié)構(gòu)變形、速度以及加速度通過接口程序傳遞到BEM程序中，然后根據(jù)DAA2法計(jì)算得到結(jié)構(gòu)濕表面的流體動(dòng)壓力，并在下一個(gè)時(shí)間增量步中傳遞給ABAQUS的FEM求解器，更新FEM中的載荷。重復(fù)以上流程，實(shí)現(xiàn)BEM和FEM之間的數(shù)據(jù)傳遞和沿時(shí)間的顯式推進(jìn)。

圖3 ABAQUS用戶子程序數(shù)據(jù)流示意圖Fig.3 Data flow in the user－subroutine of ABAQUS

FEM同BEM的數(shù)據(jù)傳遞過程中存在一個(gè)延遲問題，即FEM求解器在計(jì)算時(shí)采用的是上一個(gè)時(shí)間增量步的壓力數(shù)據(jù)，對此可采用前2個(gè)增量步的載荷數(shù)據(jù)外插的方法得到當(dāng)前增量步的載荷：

式中：下標(biāo)表示增量步數(shù)。對于BEM和FEM間的數(shù)據(jù)傳遞過程，本文中采用非共節(jié)點(diǎn)的網(wǎng)格映射方法，一方面減小BEM的計(jì)算規(guī)模，另一方面盡可能地保持較高的精度。

2 計(jì)算結(jié)果分析與討論

2.1 數(shù)值模型精度驗(yàn)證

通過和P.Zhang等［18］針對充水球殼經(jīng)典算例推導(dǎo)的解析解進(jìn)行對比，驗(yàn)證本文計(jì)算模型的有效性。該算例旨在考察厚度為1mm、內(nèi)部充滿水的彈性鋼質(zhì)球殼在1MPa平面階躍波作用下的響應(yīng)。

圖4為不同時(shí)刻球殼內(nèi)流場的聲壓云圖。從圖4可知：在相互作用初期，內(nèi)流場透射沖擊波幅值與入射波基本一致，且波陣面以平面向前推進(jìn)；隨著時(shí)間的增加，由于結(jié)構(gòu)中應(yīng)力波的波速大于流場波速，因此遠(yuǎn)離殼體的球殼中部的沖擊波明顯慢于外側(cè)靠近球殼部分流體，從而使波陣面不再保持平面，見圖4（b）；當(dāng)透射波到達(dá)球殼下部時(shí)，又產(chǎn)生了一定的反射，透射壓力疊加，形成高壓區(qū)，見圖4（c）和（d）。

圖4 不同時(shí)刻內(nèi)流場聲壓云圖Fig.4 Acoustic pressure of the inner field at different times

圖5分別是不同方法計(jì)算結(jié)果中球殼上迎爆點(diǎn)和背爆點(diǎn)在沖擊波方向的速度對比。可以看到，在圖5（a）中內(nèi)外域均為ABAQUS／聲固耦合計(jì)算結(jié)果的高頻部分衰減過快，因此在后期的速度曲線較光滑。在圖5（b）中，背爆點(diǎn)產(chǎn)生了2次明顯的加速過程：第1次發(fā)生在t＝0.6ms左右，為應(yīng)力波通過殼體傳遞到背爆面所致；第2次發(fā)生在t＝1.3ms左右，為內(nèi)流場透射壓力到達(dá)背爆面所致。3種方法的計(jì)算結(jié)果基本一致，表明本文中所建立的計(jì)算模型可用于水下沖擊問題的工程和理論研究。

圖5 數(shù)值模型對比驗(yàn)證Fig.5 Comparison among different methods

2.2 圓柱殼結(jié)構(gòu)總體響應(yīng)特性研究

計(jì)算模型如圖6～7所示。計(jì)算模型長約75m，非耐壓殼直徑約8m，耐壓殼直徑約6m，肋距0.6m，一階濕模態(tài)頻率約為2.6Hz。對外流場采用DAA2法模擬，對內(nèi)流場采用聲固耦合法模擬，兩者在非耐壓殼處滿足連續(xù)性條件。由沖擊動(dòng)力學(xué)和結(jié)構(gòu)動(dòng)力學(xué)的相關(guān)理論可知，舷間流場的存在并未改變細(xì)長圓柱殼的縱向體質(zhì)量和剛度的分布情況，因此對總體低階振動(dòng)特性不會(huì)產(chǎn)生明顯影響。但舷間流場對高頻沖擊波的透射作用會(huì)顯著改變結(jié)構(gòu)的局部高頻響應(yīng)特征，尤其是前期沖擊波作用階段更明顯。

圖6 雙層圓柱殼總體響應(yīng)計(jì)算示意圖Fig.6 Configuration of the total response simulation of double cylindrical structure

圖7 雙層圓柱殼計(jì)算模型Fig.7 Calculation models for double cylindrical shell structures

對艦艇在水下爆炸作用下的總體響應(yīng)已進(jìn)行了大量的理論和實(shí)驗(yàn)研究，通常認(rèn)為當(dāng)氣泡脈動(dòng)頻率接近或略小于艦艇一階模態(tài)頻率，即氣泡脈動(dòng)周期略大于結(jié)構(gòu)固有周期時(shí)，鞭狀總體響應(yīng)最明顯。氣泡脈動(dòng)周期可通過庫爾估算公式獲得［1］：

式中：K 為材料常數(shù)，對于 TNT，K≈2.11s·kg－1／3·m5／6；W 為裝藥質(zhì)量，kg；d 為水深，m；10.3m 為大氣壓對應(yīng)的壓力水頭。

為了深入研究不同工況參數(shù)對鞭狀運(yùn)動(dòng)的影響，定義周期參數(shù)γT＝Tb／Ts和爆距參數(shù)γR＝R／L，作為衡量水下爆炸氣泡誘導(dǎo)鞭狀運(yùn)動(dòng)的量綱一參量。Tb為氣泡脈動(dòng)周期，可通過式（14）計(jì)算；Ts為結(jié)構(gòu)一階濕模態(tài)對應(yīng)的固有頻率；R為爆距；L為結(jié)構(gòu)特征長度，對加筋圓柱殼，本文中L取殼體長度。

2.2.1 周期比的影響

首先考慮周期參數(shù)對鞭狀運(yùn)動(dòng)的影響，工況設(shè)置見表1，其中TNT炸藥質(zhì)量均為1t，爆距均為20m，來分析不同水深d對鞭狀運(yùn)動(dòng)的影響。事實(shí)上，方位角對結(jié)構(gòu)鞭狀運(yùn)動(dòng)的影響：一方面，體現(xiàn)在氣泡向上遷移時(shí)不同方位角下的爆距改變趨勢不同；另一方面，則是由于氣泡遷移引起的流場非對稱性導(dǎo)致的。方位角為90°時(shí)氣泡向上遷移對結(jié)構(gòu)鞭狀運(yùn)動(dòng)響應(yīng)影響最小。因此，為減少變量數(shù)量，本文的分析均是針對方位角為90°的情況。

表1 柱殼總體響應(yīng)工況設(shè)置Table 1 Case parameters for total responses of cylindrical shell structures

由于鞭狀運(yùn)動(dòng)主要是一階模態(tài)形式下的運(yùn)動(dòng)，因此采用船體梁中點(diǎn)撓度來評估圓柱殼鞭狀運(yùn)動(dòng)的幅值?？紤]到圓柱殼在沖擊載荷下產(chǎn)生的剛體位移，將圓柱殼撓度定義為：

式中：Uc為圓柱殼中部的垂向位移；Ul和Ur分別為圓柱殼兩端的垂向位移。

先以工況1為例，圓柱殼到達(dá)鞭狀運(yùn)動(dòng)幅值各波峰波谷時(shí)的變形和應(yīng)力云圖如圖8所示。

圖8 工況1鞭狀運(yùn)動(dòng)應(yīng)力云圖Fig.8 Mises stress contours of the whipping response in case 1

在工況1下，質(zhì)量為1t的TNT在水下60m深處爆炸的氣泡一次脈動(dòng)周期為610ms，約為圓柱殼結(jié)構(gòu)一階固有周期的1.60倍。從圖8可以看出，工況1的圓柱殼在水下爆炸的前期的沖擊波作用下已被激起了一階鞭狀運(yùn)動(dòng)，但幅值不大。而鞭狀運(yùn)動(dòng)過程中在240和450ms出現(xiàn)了2個(gè)相鄰的波峰，間隔時(shí)間較短，3個(gè)鞭狀運(yùn)動(dòng)幅值沒有明顯增大。

圖9 工況1鞭狀運(yùn)動(dòng)幅值和沖擊載荷曲線Fig.9 Curves of whipping response to shock pressure in case 1

圖9為工況1中鞭狀運(yùn)動(dòng)幅值和沖擊載荷曲線。其中，綠虛線為沖擊波和氣泡脈動(dòng)載荷，實(shí)線為計(jì)算得到的三維圓柱殼結(jié)構(gòu)的鞭狀運(yùn)動(dòng)幅值曲線。該工況下結(jié)構(gòu)最大幅值僅有0.1m左右；另外，在第2次鞭狀運(yùn)動(dòng)周期出現(xiàn)2個(gè)緊鄰的波峰。當(dāng)沖擊波載荷過后，圓柱殼開始第1次鞭狀運(yùn)動(dòng)，并很快達(dá)到第1個(gè)波峰。在回落過程中，受氣泡脈動(dòng)載荷負(fù)壓區(qū)的作用，對圓柱殼結(jié)構(gòu)輸入了較大的能量，使第1個(gè)鞭狀運(yùn)動(dòng)波谷的幅值大于第1個(gè)波峰的幅值，約為0.17m。由于工況1的氣泡脈動(dòng)周期較長，當(dāng)圓柱殼結(jié)構(gòu)第2個(gè)鞭狀運(yùn)動(dòng)峰值出現(xiàn)時(shí)，氣泡一次脈動(dòng)壓力尚未到來，因此在這段時(shí)間內(nèi)，圓柱殼的鞭狀運(yùn)動(dòng)主要是由本身的有阻尼震蕩和脈動(dòng)壓力的負(fù)壓區(qū)共同主導(dǎo)的，幅值僅有0.10m，較第1個(gè)波谷有明顯的減小。氣泡一次脈動(dòng)載荷峰值在鞭狀運(yùn)動(dòng)第2個(gè)波峰的回落過程中出現(xiàn)，從而使結(jié)構(gòu)調(diào)轉(zhuǎn)運(yùn)動(dòng)方向，形成了另一個(gè)波峰。在這種水深較小的工況下，氣泡脈動(dòng)周期較長，使氣泡一次脈動(dòng)壓力峰值到來的相對較晚，若恰巧出現(xiàn)在鞭狀運(yùn)動(dòng)的回落階段，則會(huì)產(chǎn)生類似的鞭狀運(yùn)動(dòng)的雙峰現(xiàn)象，即存在2個(gè)間隔很小的波峰。在這種情況下，通常不會(huì)引起劇烈的鞭狀運(yùn)動(dòng)。

對于工況3，1t的TNT在水下80m深處爆炸的氣泡一次脈動(dòng)周期為490ms，約為圓柱殼結(jié)構(gòu)一階濕模態(tài)對應(yīng)的周期的1.3倍。工況3鞭狀運(yùn)動(dòng)幅值和沖擊載荷時(shí)歷曲線見圖10。

圖10 工況2鞭狀運(yùn)動(dòng)幅值和沖擊載荷曲線Fig.10 Curves of whipping response to shock pressure in case 2

圖11 工況5鞭狀運(yùn)動(dòng)幅值和沖擊載荷曲線Fig.11 Curves of whipping response to shock pressure in case 5

從圖10可看出：與工況1相同，在從第1個(gè)波峰回落的過程中，由于同時(shí)受到氣泡脈動(dòng)載荷負(fù)壓區(qū)的能量輸入，使第1個(gè)鞭狀運(yùn)動(dòng)波谷的幅值大于第1個(gè)波峰的幅值；而氣泡一次脈動(dòng)載荷又恰好位于鞭狀運(yùn)動(dòng)的第2次上升階段，同樣以做正功的形式向圓柱殼結(jié)構(gòu)輸入能量，使鞭狀運(yùn)動(dòng)幅值繼續(xù)增大，達(dá)到0.28m；在接下來的第3和第4個(gè)周期中，氣泡脈動(dòng)載荷同樣為系統(tǒng)輸送能量，從而使峰值逐個(gè)增大，分別達(dá)到0.35和0.41m；而第5個(gè)周期時(shí)，脈動(dòng)壓力已經(jīng)衰減到很小，所做的功不足以抵消整個(gè)系統(tǒng)的能量耗散，使鞭狀運(yùn)動(dòng)幅值開始回落。

在工況5中，周期參數(shù)約為1.1，所對應(yīng)圓柱殼鞭狀運(yùn)動(dòng)最大幅值僅有0.37m，遠(yuǎn)小于工況3。圖11中共有2個(gè)氣泡脈動(dòng)壓力峰值出現(xiàn)在鞭狀運(yùn)動(dòng)上升階段，鞭狀運(yùn)動(dòng)最大運(yùn)動(dòng)幅值出現(xiàn)在第3個(gè)運(yùn)動(dòng)周期。第1個(gè)周期的波峰同樣來自沖擊波載荷，氣泡前2次脈動(dòng)壓力峰值均在鞭狀運(yùn)動(dòng)上升階段，使運(yùn)動(dòng)幅值逐次增加。而氣泡第3次脈動(dòng)壓力峰值剛好處于圓柱殼鞭狀運(yùn)動(dòng)第3個(gè)波谷，此時(shí)鞭狀運(yùn)動(dòng)速度為零，即此次氣泡脈動(dòng)載荷基本對系統(tǒng)不做功。同時(shí)由于系統(tǒng)能量的耗散，使第4次的運(yùn)動(dòng)幅值開始減小。隨后的氣泡脈動(dòng)壓力峰值開始出現(xiàn)在結(jié)構(gòu)鞭狀運(yùn)動(dòng)的下降階段，從而加快了鞭狀運(yùn)動(dòng)衰減的速度。從圖12可知，隨著水深的增大，圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)周期逐漸減小。綜合所有工況，分析周期比與鞭狀運(yùn)動(dòng)最大幅值Um的關(guān)系，將其繪于圖13。

圖12 不同工況下，鞭狀運(yùn)動(dòng)時(shí)歷曲線Fig.12 Whipping response curves in different cases

圖13 鞭狀運(yùn)動(dòng)最大幅值隨周期比的變化Fig.13 The maximum amplitude of the whipping response varied with period ratio

從圖13可知，隨著周期比的增大，鞭狀運(yùn)動(dòng)幅值逐漸增大，水深為80～90m，即1.19＜γT＜1.30時(shí)，鞭狀運(yùn)動(dòng)幅值達(dá)到最大，隨后迅速衰減。綜合以上分析可知，當(dāng)1.19＜γT＜1.30時(shí)，可是使連續(xù)多個(gè)氣泡脈動(dòng)峰值出現(xiàn)在圓柱殼鞭狀運(yùn)動(dòng)幅值的上升階段，從而連續(xù)向系統(tǒng)輸送能量，達(dá)到共振。當(dāng)γT從1.19減小時(shí)，氣泡脈動(dòng)峰值處于鞭狀運(yùn)動(dòng)上升階段的個(gè)數(shù)逐漸減小，從而使鞭狀運(yùn)動(dòng)最大幅值隨之減小。當(dāng)γT從1.30增大時(shí)，氣泡第一次脈動(dòng)峰值直接后移至鞭狀運(yùn)動(dòng)幅值下降階段，對系統(tǒng)的能量輸入變?yōu)樨?fù)，在鞭狀運(yùn)動(dòng)的第2個(gè)峰值附近緊接著出現(xiàn)第3個(gè)峰值，從而使鞭狀運(yùn)動(dòng)最大幅值迅速減小。

2.2.2 爆距比的影響

除水深外，爆距也會(huì)明顯影響鞭狀運(yùn)動(dòng)的形式和幅值。為分析爆距對圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)的影響，進(jìn)行一系列的計(jì)算，工況設(shè)置見表2，TNT炸藥的質(zhì)量均為1t，水深均為80m，氣泡脈動(dòng)周期為490ms。不同工況下圓柱殼鞭狀運(yùn)動(dòng)的幅值曲線見圖14～15。從圖14～15可知：隨著爆距的增大，鞭狀運(yùn)動(dòng)幅值迅速減小。其次，爆距越小的工況，鞭狀運(yùn)動(dòng)第1個(gè)周期越長，如R＝10m的工況第1個(gè)鞭狀運(yùn)動(dòng)周期約為550ms，超過圓柱殼結(jié)構(gòu)固有周期380ms約45%，因此使后續(xù)周期相對于其他工況后延了約1／4個(gè)周期。而對于R≥20m的工況，鞭狀運(yùn)動(dòng)每個(gè)周期的長度基本都是均勻的。圖14給出了80m水深下各工況最大鞭狀運(yùn)動(dòng)幅值隨爆距的變化趨勢。從圖14可以看出，隨著爆距的增大，圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)幅值迅速減小。這種減小一方面來自于水下爆炸載荷隨爆距的倒數(shù)關(guān)系，另一方面則是隨著爆距的增大，爆炸載荷在結(jié)構(gòu)縱向的分布逐漸均勻，使結(jié)構(gòu)更多傾向于剛體平動(dòng)，而非鞭狀運(yùn)動(dòng)形式的彎曲響應(yīng)。

表2 柱殼總體響應(yīng)工況設(shè)置Table 2 Case parameters for total response

圖14 爆距對鞭狀運(yùn)動(dòng)幅值的影響Fig.14 Effect of explosion distance on the amplitude of the whipping response

圖15 鞭狀運(yùn)動(dòng)最大幅值隨爆距比的變化Fig.15 The maximum amplitude of the whipping response varied with explosion distance ratio

3 結(jié) 論

從波動(dòng)方程出發(fā)，推導(dǎo)了DAA2法的后期近似，并結(jié)合聲固耦合法初步解決了雙層圓柱殼內(nèi)流場問題，采用ABAQUS接口耦合實(shí)現(xiàn)了圓柱殼水下爆炸沖擊響應(yīng)的計(jì)算分析。通過簡單算例的驗(yàn)證表明，本文中建立的數(shù)值模型可用于圓柱殼水下沖擊的工程和理論研究。據(jù)此分析了雙層圓柱殼結(jié)構(gòu)在水下爆炸載荷下的總體響應(yīng)特性，發(fā)現(xiàn)對于本研究的參數(shù)范圍，周期比為1.2到1.3之間時(shí)，能夠誘發(fā)最嚴(yán)重的鞭狀運(yùn)動(dòng)；當(dāng)周期比從1.2減小時(shí)，鞭狀運(yùn)動(dòng)幅值逐漸減?。欢?dāng)周期比從1.3增大時(shí)，鞭狀運(yùn)動(dòng)幅值迅速減??；圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)幅值隨爆距的增加呈平方關(guān)系減小。

［1］Huang H.Transient interaction of acoustic plane waves with a spherical elastic shell［J］.The Journal of the Acoustical Society of America，1969，45（3）：661－670.

［2］Huang H.An exact Analysis of the transient interaction of acoustic plane waves with a cylindrical elastic shell［J］.Journal of Applied Mechanics，1970，37：1091－1099.

［3］Huang H，Lu Y P，Wang Y F.Transient interaction of spherical acoustic waves and a spherical elastic shell［J］.Journal of Applied Mechanics，1971，38：71－74.

［4］Geers T L.Residual potential and approximate methods for three－dimensional fluid－structure interaction problems［J］.Journal of the Acoustical Society of America，1971，49（5B）：1505－1510.

［5］Geers T L.Doubly asymptotic approximations for transient motions of submerged structures［J］.Journal of the A－coustical Society of America，1978，64（5）：1500－1508.

［6］劉建湖.艦船非接觸水下爆炸動(dòng)力學(xué)的理論與應(yīng)用［D］.無錫：中國船舶科學(xué)技術(shù)研究所，2002：89－112.

［7］孫士麗.瞬態(tài)載荷作用下大幅運(yùn)動(dòng)航行體流固耦合方法及應(yīng)用研究［D］.哈爾濱：哈爾濱工程大學(xué)，2009：15－34.

［8］王詩平，孫士麗，張阿漫，等.沖擊波和氣泡作用下艦船結(jié)構(gòu)動(dòng)態(tài)響應(yīng)的數(shù)值模擬［J］.爆炸與沖擊，2011，31（4）：367－372.Wang Shi－ping，Sun Shi－li，Zhang A－man，et al.Numerical simulation of dynamic response of warship structures subjected to underwater explosion shockwaves and bubbles［J］.Explosion and Shock Waves，2011，31（4）：367－372.

［9］張弩，宗智，張文鵬，等.基于雙漸進(jìn)方法的水下爆炸氣泡載荷作用下艦船的動(dòng)態(tài)響應(yīng)分析［J］.振動(dòng)與沖擊，2012，32（23）：50－56.Zhang Nu，Zong Zhi，Zhang Wen－peng，et al.Dynamic response of a ship hull structure subjected to an underwater explosion bubble based on doubly asymptotic approximation method［J］.Journal of Vibration and Shock，2012，32（23）：50－56.

［10］董海，劉建湖，吳有生.水下爆炸氣泡脈動(dòng)作用下細(xì)長加筋圓柱殼的鞭狀響應(yīng)分析［J］.船舶力學(xué)，2007，11（2）：250－258.Dong Hai，Liu Jian－h(huán)u，Wu You－sheng.Whipping response analysis of slender stiffened cylindrical shell subjected to underwater explosion with bubble pulse［J］.Journal of Ship Mechanics，2007，11（2）：250－258.

［11］Geers T L，Hunter K S.An integrated wave－effects model for an underwater explosion bubble［J］.The Journal of the Acoustical Society of America，2002，111（4）：1584－1601.

［12］Geers T L，Park C K.Optimization of the G ＆ H bubble model［J］.Shock and Vibration，2005，12（1）：3－8.

［13］Sprague M A，Geers T L.Response of empty and fluid－filled submerged spherical shells to plane and spherical，step－exponential acoustic waves［J］.Shock and Vibration，1999，6：147－157.

［14］王詩平，孫士麗，張阿漫，等.可壓縮流場中氣泡脈動(dòng)數(shù)值模擬［J］.力學(xué)學(xué)報(bào)，2012，44（3）：513－519.Wang Shi－ping，Sun Shi－li，Zhang A－man，et al.Numerical simulation of bubble dynamics in compressible fluid［J］.Journal of Theoretical and Applied Mechanics，2012，44（3）：513－519.

［15］Geers T L，Zhang P Z.Doubly asymptotic approximations for submerged structures with internal fluid volumes：Formulation［J］.Journal of Applied Mechanics，1994，61（4）：893－899.

［16］Geers T L，Zhang P Z.Doubly asymptotic approximations for submerged structures with internal fluid volumes：Evaluation［J］.Journal of Applied Mechanics，1994，61（4）：900－906.

［17］Wilkerson S A.Elastic whipping response of ships to an underwater explosion loading［D］.Washington：George Washington University，1985：57－72.

［18］Zhang P，Geers T L.Excitation of a fluid－filled，submerged elastic shell by a transient acoustic wave［J］.Journal of the Acoustical Society of America，1993，93（2）：696－705.