劉云龍，汪 玉，張阿漫

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

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

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

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

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

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

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

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

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

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

1.2 水下爆炸載荷

水下爆炸氣泡脈動(dòng)是一個(gè)強(qiáng)非線(xiàn)性的物理化學(xué)過(guò)程，針對(duì)不同的問(wèn)題已提出了很多氣泡模型。其中T.L.Geers等基于DAA2法提出的單個(gè)水下爆炸氣泡脈動(dòng)模型［15］考慮了流場(chǎng)的可壓縮性和氣泡內(nèi)部氣體的聲學(xué)特性，因此在氣泡各次脈動(dòng)過(guò)程中，能夠真實(shí)反映氣泡最大半徑和氣泡脈動(dòng)壓力峰值的衰減。另外考慮了氣泡中心的上浮，并通過(guò)以往水下爆炸氣泡實(shí)驗(yàn)數(shù)據(jù)對(duì)上浮中的阻力進(jìn)行了匹配，得到結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好。在隨后的研究中，T.L.Geers等對(duì)Geers－Hunter模型進(jìn)行了一些優(yōu)化調(diào)整［16］。采用Geers－Hunter模型計(jì)算得到氣泡的脈動(dòng)和遷移運(yùn)動(dòng)后，則可通過(guò)以下方法計(jì)算流體中任意一點(diǎn)的流體速度勢(shì)［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)的速度勢(shì)，第2項(xiàng)為氣泡遷移等效的偶極子所誘導(dǎo)的速度勢(shì)，后面2項(xiàng)為自由面鏡像的氣泡所誘導(dǎo)的速度勢(shì)。其中e1和e2分別為源強(qiáng)和偶極子強(qiáng)度；kf為自由面的反射因數(shù)，由于本文中分析深水中結(jié)構(gòu)對(duì)水下爆炸的動(dòng)態(tài)響應(yīng)，因此不考慮自由面影響，kf取值為0；r1和r2分別為氣泡中心和鏡像氣泡中心的距離；θ1和θ2分別為2連線(xiàn)與鉛垂線(xiàn)的夾角，見(jiàn)圖2。

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

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

1.3 BEM與FEM耦合接口

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

圖3中的數(shù)字表示時(shí)間增量步，在同一個(gè)時(shí)間增量步中，F(xiàn)EM計(jì)算得到的結(jié)構(gòu)變形、速度以及加速度通過(guò)接口程序傳遞到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用戶(hù)子程序數(shù)據(jù)流示意圖Fig.3 Data flow in the user－subroutine of ABAQUS

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

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

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

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

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

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

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

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

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

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

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

對(duì)艦艇在水下爆炸作用下的總體響應(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)周期可通過(guò)庫(kù)爾估算公式獲得［1］：

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

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

2.2.1 周期比的影響

首先考慮周期參數(shù)對(duì)鞭狀運(yùn)動(dòng)的影響，工況設(shè)置見(jiàn)表1，其中TNT炸藥質(zhì)量均為1t，爆距均為20m，來(lái)分析不同水深d對(duì)鞭狀運(yùn)動(dòng)的影響。事實(shí)上，方位角對(duì)結(jié)構(gòu)鞭狀運(yùn)動(dòng)的影響：一方面，體現(xiàn)在氣泡向上遷移時(shí)不同方位角下的爆距改變趨勢(shì)不同；另一方面，則是由于氣泡遷移引起的流場(chǎng)非對(duì)稱(chēng)性導(dǎo)致的。方位角為90°時(shí)氣泡向上遷移對(duì)結(jié)構(gòu)鞭狀運(yùn)動(dòng)響應(yīng)影響最小。因此，為減少變量數(shù)量，本文的分析均是針對(duì)方位角為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)撓度來(lái)評(píng)估圓柱殼鞭狀運(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)過(guò)程中在240和450ms出現(xiàn)了2個(gè)相鄰的波峰，間隔時(shí)間較短，3個(gè)鞭狀運(yùn)動(dòng)幅值沒(méi)有明顯增大。

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

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

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

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

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

從圖10可看出：與工況1相同，在從第1個(gè)波峰回落的過(guò)程中，由于同時(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；在接下來(lái)的第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)幅值開(kāi)始回落。

在工況5中，周期參數(shù)約為1.1，所對(duì)應(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è)周期的波峰同樣來(lái)自沖擊波載荷，氣泡前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)載荷基本對(duì)系統(tǒng)不做功。同時(shí)由于系統(tǒng)能量的耗散，使第4次的運(yùn)動(dòng)幅值開(kāi)始減小。隨后的氣泡脈動(dòng)壓力峰值開(kāi)始出現(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í)歷曲線(xiàn)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)幅值下降階段，對(duì)系統(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)的形式和幅值。為分析爆距對(duì)圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)的影響，進(jìn)行一系列的計(jì)算，工況設(shè)置見(jiàn)表2，TNT炸藥的質(zhì)量均為1t，水深均為80m，氣泡脈動(dòng)周期為490ms。不同工況下圓柱殼鞭狀運(yùn)動(dòng)的幅值曲線(xiàn)見(jiàn)圖14～15。從圖14～15可知：隨著爆距的增大，鞭狀運(yùn)動(dòng)幅值迅速減小。其次，爆距越小的工況，鞭狀運(yùn)動(dòng)第1個(gè)周期越長(zhǎng)，如R＝10m的工況第1個(gè)鞭狀運(yùn)動(dòng)周期約為550ms，超過(guò)圓柱殼結(jié)構(gòu)固有周期380ms約45%，因此使后續(xù)周期相對(duì)于其他工況后延了約1／4個(gè)周期。而對(duì)于R≥20m的工況，鞭狀運(yùn)動(dòng)每個(gè)周期的長(zhǎng)度基本都是均勻的。圖14給出了80m水深下各工況最大鞭狀運(yùn)動(dòng)幅值隨爆距的變化趨勢(shì)。從圖14可以看出，隨著爆距的增大，圓柱殼結(jié)構(gòu)鞭狀運(yùn)動(dòng)幅值迅速減小。這種減小一方面來(lái)自于水下爆炸載荷隨爆距的倒數(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 爆距對(duì)鞭狀運(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)流場(chǎng)問(wèn)題，采用ABAQUS接口耦合實(shí)現(xiàn)了圓柱殼水下爆炸沖擊響應(yīng)的計(jì)算分析。通過(guò)簡(jiǎn)單算例的驗(yàn)證表明，本文中建立的數(shù)值模型可用于圓柱殼水下沖擊的工程和理論研究。據(jù)此分析了雙層圓柱殼結(jié)構(gòu)在水下爆炸載荷下的總體響應(yīng)特性，發(fā)現(xiàn)對(duì)于本研究的參數(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］.無(wú)錫：中國(guó)船舶科學(xué)技術(shù)研究所，2002：89－112.

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

［8］王詩(shī)平，孫士麗，張阿漫，等.沖擊波和氣泡作用下艦船結(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ì)長(zhǎng)加筋圓柱殼的鞭狀響應(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］王詩(shī)平，孫士麗，張阿漫，等.可壓縮流場(chǎng)中氣泡脈動(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.