潛艇受撞數(shù)值仿真的附加質(zhì)量模型

葉劍平，劉俊杰，萬正權(quán)，趙 南

（1海軍駐七一九所軍事代表室，武漢430064；2中國船舶科學研究中心，江蘇 無錫214082）

1 引 言

潛艇在航行或作業(yè)時，有可能受到來自其它物體的撞擊威脅，例如受到其它航行器、大型海洋生物以及水面艦船墜落重物的撞擊等，在潛艇上浮過程中，還有可能受到水面浮冰的撞擊。撞擊事故輕者會造成潛艇結(jié)構(gòu)的損傷變形，從而影響結(jié)構(gòu)的結(jié)構(gòu)強度和疲勞壽命，重者會造成潛艇殼體破裂，從而引起艙室進水，造成電子儀器的損壞以及人員傷亡，其平衡系統(tǒng)也有可能因此而破壞，導致潛艇的沉沒，后果不堪設(shè)想。

潛艇通常采用金屬材料建造，撞擊事故發(fā)生時，在巨大沖擊載荷作用下，受撞區(qū)的金屬構(gòu)件往往會迅速跨越彈性變形階段而進入塑性流動狀態(tài)，產(chǎn)生塑性大變形。撞擊是在很短時間內(nèi)完成的，伴隨整個撞擊過程會出現(xiàn)撞擊物和受撞潛艇的運動非線性、材料塑性變形非線性以及接觸非線性等問題，因此，潛艇的受撞過程是一個非常復雜的非線性瞬態(tài)響應(yīng)過程。區(qū)別于陸上結(jié)構(gòu)物撞擊的一個顯著特征是：潛艇受撞不能忽視周圍流體介質(zhì)的影響。撞擊發(fā)生時，潛艇在撞擊載荷下，將會發(fā)生局部變形和剛體運動，這種變形或運動會改變周圍流體介質(zhì)的流場，流場的變化又會反過來影響到潛艇的變形或運動，這種相互作用，形成了撞擊過程中的流—固耦合現(xiàn)象，參見圖1。

從廣義上講，潛艇受撞也屬于艦艇碰撞研究領(lǐng)域。大量的研究工作表明[1-2]，在眾多艦艇碰撞問題的研究方法中，有限元數(shù)值仿真方法是最有效的。艦艇碰撞有限元數(shù)值仿真方法的計算精度主要取決于對問題描述（即建立的計算模型）的準確程度，文獻[3]給出了水面船舶碰撞數(shù)值仿真的三種建模方法，即：流—固耦合法、附加質(zhì)量法和等效船體梁法，通過對三種方法進行比較表明，計算精度和計算時間上三種方法都是依次降低的，合理的建模方法和仿真模型能夠在保證計算精度的前提下，有效地降低計算時間，從而提升計算效率。本文建立了兩種附加質(zhì)量模型分別對潛艇受物體撞擊的損傷過程進行了有限元數(shù)值仿真，通過與流—固耦合模型仿真結(jié)果的對比，給出了合理的潛艇受撞有限元數(shù)值仿真計算的附加質(zhì)量建模方法。

圖1 潛艇受撞的流—固耦合系統(tǒng)Fig.1 Fluid-structure coupling system of submarine collision

2 計算對象及有限元模型

本文選取某型潛艇作為受撞物體，簡化為由環(huán)肋圓柱殼和環(huán)肋圓錐殼連接而成的單殼體結(jié)構(gòu)，環(huán)肋骨為T型材，潛艇的橫剖面是軸對稱的。本文假設(shè)該潛艇的重心在中心軸線上，潛艇的水下排水量為14 335 t，全長138 m，全寬12.5 m；撞擊物為一圓柱體，撞擊端為半圓球形，質(zhì)量為56 t，全長12 m，直徑2.3 m。潛艇處于水面以下且靠近水面的位置，忽略靜水壓力對其的影響，其軸線平行于水平面，處于靜止狀態(tài)，撞擊物從距離潛艇受撞點30 m的高空在重力作用下自由下落并撞擊潛艇，其軸線通過潛艇的重心并與潛艇的軸線垂直。

潛艇的殼體、艙壁以及環(huán)肋骨均采用殼單元模擬，由于結(jié)構(gòu)的損傷變形主要集中在撞擊接觸區(qū)域，因此對該區(qū)域構(gòu)件的網(wǎng)格進行了細化，遠離接觸區(qū)的結(jié)構(gòu)在保證外形特征以不影響流—固耦合計算精度的前提下，采用較粗的網(wǎng)格劃分。考慮到潛艇金屬材料的應(yīng)變率敏感性，采用Cowper-Symonds本構(gòu)方程來計算材料的動屈服應(yīng)力

式中：σk是材料在塑性應(yīng)變率時的動屈服應(yīng)力，σ0是材料的靜屈服應(yīng)力，D和q是與材料相關(guān)的常數(shù)，船用鋼D取40.4，q取5[4]。為了真實地模擬該潛艇殼體結(jié)構(gòu)的損傷變形以及整體的剛體運動，本文通過布置質(zhì)量點并調(diào)節(jié)其分布和單元屬性，使得模型的重心、質(zhì)量以及轉(zhuǎn)動慣量與實際潛艇對應(yīng)的參數(shù)相同；考慮到本文的主要研究對象是潛艇結(jié)構(gòu)的損傷，并且撞擊物撞擊部位的結(jié)構(gòu)剛度要比潛艇受撞部位的結(jié)構(gòu)剛度大，因此，撞擊物采用剛體材料模擬，采用MSC.Patran[5]建立的有限元模型如圖2所示。

圖2 潛艇受撞有限元模型Fig.2 The FE model of submarine collision

潛艇受物體撞擊屬于接觸力學研究范疇，接觸面的位置和區(qū)域、接觸模式、接觸力的方向和大小都隨著撞擊過程的進行而在發(fā)生變化。本文采用MSC.Dytran的主-從面接觸算法來模擬撞擊物與受撞潛艇之間的接觸、分離，這種算法允許兩個面之間連續(xù)不斷的接觸和滑動。由于該接觸算法是不對稱的，程序只檢查從面是否穿透主面，而不檢查主面是否穿透從面，要求從面的網(wǎng)格要比主面網(wǎng)格細，因此，從面定義在撞擊物外表面，而將主面定義在可能與撞擊物發(fā)生接觸的潛艇構(gòu)件上，參見圖2（b）。對撞擊過程中接觸面之間差生的摩擦力，采用經(jīng)典的庫倫摩擦定律來近似計算，摩擦系數(shù)計算公式為：

式中：μs為靜摩擦系數(shù)，μk為動摩擦系數(shù)，β為指數(shù)衰減系數(shù)，v為主、從接觸面之間的相對滑動速度。

3 流—固耦合模型

本文采用MSC.Dytran的一般耦合算法對潛艇受撞過程中的流-固耦合現(xiàn)象進行模擬[6]，將潛艇結(jié)構(gòu)及其周圍的流體介質(zhì)作為一個整體來考慮，潛艇結(jié)構(gòu)采用Lagrange模型來模擬，而流體介質(zhì)則采用Euler模型模擬，通過定義在Lagrange模型上的耦合面建立聯(lián)系，參見圖3，耦合面是Lagrange網(wǎng)格與Euler網(wǎng)格之間的相互作用力的傳遞者。

圖3 一般耦合示意圖Fig.3 The general coupling

從本文所研究的潛艇的受撞模式可以看出，撞擊主要引起潛艇在中縱剖面內(nèi)運動（垂蕩和縱搖），其它方向的運動很小，因此，長度取潛艇長度的1.3倍，寬度取潛艇寬度的2倍，高度取潛艇高度（寬度）的3倍，建立Euler域，近似地作為受撞潛艇影響的流場范圍；受撞潛艇處在Euler域的中心，并用Lagrange單元模擬，耦合面定義在潛艇的外殼體封閉表面上，建立的流—固耦合模型如圖4所示。

圖4 流—固耦合模型Fig.4 The FE model of fluid-structure coupling

采用非線性瞬態(tài)動力學分析軟件MSC.Dytran作為求解器，定義撞擊物與潛艇接觸時刻為仿真計算開始時間t=0 s，終止時間為t=0.25 s，計算得到耦合面在三個方向上的耦合力如圖5所示，其中耦合力正值表示潛艇殼體受到流體介質(zhì)的壓力，負值則表示受到流體介質(zhì)的吸力（負壓）。從圖中可以看出，潛艇在長度方向和寬度方向上的耦合力較小，而在受撞方向（垂向）上的耦合力較大，這是因為撞擊主要造成潛艇在垂向上的變形和運動，而在其余兩個方向上的變形和運動卻很小；耦合力曲線呈現(xiàn)一定的波動性，這是由撞擊引起的潛艇在垂向上的震蕩運動產(chǎn)生的，隨著潛艇運動幅度的減弱，耦合力的幅值也隨之降低，最終會趨于零。

圖5 耦合力曲線Fig.5 The curves of coupling force

圖6 流場壓力分布Fig.6 The distribution of fluid pressure

潛艇受撞過程中，不同時刻流場中縱剖面內(nèi)的流體壓力分布如圖6所示，從圖中可以看出，潛艇受撞產(chǎn)生的變形和剛體運動改變了周圍流體介質(zhì)的流場，壓力分布發(fā)生了變化，應(yīng)力波以潛艇為“源”向周圍擴散，由潛艇中部向兩端擴散；潛艇的垂蕩運動使得流體壓力循環(huán)作用在其上下表面，隨著其震蕩運動的減弱，流場也逐漸趨于平穩(wěn)，并恢復到初始狀態(tài)。

4 附加質(zhì)量模型

流—固耦合數(shù)值仿真模型考慮了潛艇受撞過程中周圍流體介質(zhì)的影響，能夠反映真實的受撞場景，因此計算精度較高，但耦合計算需要花費巨大的計算時間，以本文計算為例，用于Euler單元計算的CPU時間占到了總計算時間的將近1/3。開展大規(guī)模潛艇受撞數(shù)值仿真分析時，在對受撞區(qū)復雜、細小構(gòu)件的損傷變形情況需要進行詳細、精確描述的前提下，采用流—固耦合模型必然會使仿真計算時間急劇增加。王自力等人[7]在開展船舶碰撞研究時，提出了一種附加質(zhì)量模型，采用船體附加質(zhì)量的形式來考慮周圍流體介質(zhì)的動力影響，避開了耗時的流—固耦合計算，顯著提高了計算效率。本文采用兩種附加質(zhì)量模型對潛艇受撞進行了數(shù)值仿真計算。

附加質(zhì)量主要用來反映受撞潛艇與周圍流體介質(zhì)的相互作用，其大小與受撞潛艇的結(jié)構(gòu)特征以及受撞后的運動特征等相關(guān)，目前，附加質(zhì)量可以采用切片法計算，也可以通過數(shù)值方法或試驗獲得。由流—固耦合模型計算結(jié)果可以看出，本文研究的受撞潛艇主要產(chǎn)生垂蕩運動，因此只考慮該方向上的附加質(zhì)量，結(jié)合該潛艇結(jié)構(gòu)特征，采用文獻[8]給出的值

式中：mz為垂蕩運動附加質(zhì)量，m為潛艇的質(zhì)量。

4.1 附加質(zhì)量模型一

本文建立的潛艇有限元模型的質(zhì)量是由兩部分組成，一部分是由殼體、艙壁以及環(huán)肋骨等構(gòu)件根據(jù)自身材料密度和體積得到的質(zhì)量；另一部分是由分布在殼體上的、用來代表潛艇其余構(gòu)件以及儀器、儀表設(shè)備等的質(zhì)量點單元的質(zhì)量。附加質(zhì)量模型一是通過增加潛艇有限元模型質(zhì)量的第一部分質(zhì)量，即增加材料的密度使得附加質(zhì)量添加到艇體上，密度增加值的計算公式為

式中：ρz為垂蕩運動的附加密度，ρ為潛艇實際材料密度，mP為質(zhì)量點單元的總質(zhì)量。

4.2 附加質(zhì)量模型二

附加質(zhì)量模型二是通過增加潛艇有限元模型質(zhì)量的第二部分質(zhì)量，即增加質(zhì)量點單元的質(zhì)量使得附加質(zhì)量添加到艇體上，該種模型不改變構(gòu)件的材料密度，相同變形下局部結(jié)構(gòu)的動能與實際情況比較接近。為了保證潛艇模型的重心位置和轉(zhuǎn)動慣性半徑不發(fā)生變化，單個質(zhì)量點單元的質(zhì)量增加值應(yīng)在滿足以下公式的基礎(chǔ)上重新進行調(diào)整：

式中：mei、mi和ni分別為潛艇第i個區(qū)域內(nèi)單個質(zhì)量點單元的質(zhì)量增加值、該區(qū)域的附加質(zhì)量以及該區(qū)域內(nèi)包含的質(zhì)量點個數(shù)，riz為第i個區(qū)域?qū)軸的轉(zhuǎn)動慣性半徑，rz為流—固耦合模型中潛艇對z軸的轉(zhuǎn)動慣性半徑。

5 仿真結(jié)果對比分析

5.1 潛艇損傷變形

仿真計算得到了撞擊終止時刻（t=0.25 s時）的潛艇受撞變形位移圖，如圖7所示。從圖中可以看出，本文給定的撞擊場景下，撞擊主要引起艇體受撞區(qū)局部結(jié)構(gòu)的變形，潛艇整體結(jié)構(gòu)的變形和剛體位移都很小；附加質(zhì)量模型二的計算結(jié)果要比附加質(zhì)量模型一的計算結(jié)果更接近流—固耦合模型的計算結(jié)果，因此能較為準確地反映潛艇受撞瞬間的運動狀態(tài)。

圖7 潛艇整體變形位移圖Fig.7 The deformation and displacement of whole submarine

圖8為潛艇受撞區(qū)艙段結(jié)構(gòu)的損傷變形應(yīng)力云紋圖，同樣可以看出，潛艇受撞具有明顯的局部特性，撞擊造成了潛艇受撞位置殼體結(jié)構(gòu)的局部凹陷。對比研究發(fā)現(xiàn)，附加質(zhì)量模型二的殼體結(jié)構(gòu)應(yīng)力分布特征與流—固耦合模型的計算結(jié)果較為接近。

圖8 應(yīng)力云紋圖Fig.8 The stress distribution

5.2 撞擊力

仿真計算得到三種模型下撞擊物與潛艇之間的撞擊力曲線，如圖9所示，從圖中可以看出，三種碰撞模型下的撞擊力曲線都具有明顯的非線性特征，撞擊力的最大峰值出現(xiàn)在0.025 s左右的時刻；附加質(zhì)量模型二計算得到的撞擊力曲線與流—固耦合模型的計算結(jié)果較為接近，尤其是在受撞過程的前段時間（0 s-0.05 s之間），兩條曲線基本吻合，在0.05 s之后，兩條曲線均出現(xiàn)了震蕩現(xiàn)象，主要是由于潛艇屬于細長體，類似于一根梁，中部受撞時產(chǎn)生的梁震動導致的，由于流—固耦合模型考慮了流體介質(zhì)的影響，梁震動會減弱，因此撞擊力曲線也在逐漸趨為平坦；附加質(zhì)量模型一與流—固耦合模型的撞擊力計算結(jié)果存在一定的差異，前者撞擊力最大峰值要比后者小，并且沒有震蕩現(xiàn)象出現(xiàn)。

圖9 撞擊力曲線Fig.9 The curves of collision force

5.3 能量耗散

潛艇在靜止狀態(tài)下受到撞擊物的撞擊，隨著撞擊過程的進行，撞擊物的初始動能將會耗散，耗散的動能將會轉(zhuǎn)變成以下幾種能量：受撞潛艇的動能、受撞潛艇的結(jié)構(gòu)變形能、沙漏能、潛艇剛體運動引起的流體介質(zhì)的內(nèi)能及摩擦產(chǎn)生的熱能。圖10為流—固耦合模型仿真計算得到的撞擊過程中撞擊物動能、潛艇動能、潛艇變形能和沙漏能的變化曲線，可以看出，撞擊物的動能在終止時刻（t=0.25 s時）幾乎被耗散完，此時潛艇結(jié)構(gòu)的變形能幾乎占到了撞擊物耗散掉動能的90.5%，表明撞擊過程中，撞擊物耗散的動能主要轉(zhuǎn)變成了潛艇結(jié)構(gòu)的變形能；撞擊引起潛艇的剛體運動很小，因此產(chǎn)生的動能很?。蛔矒暨^程中存在沙漏模態(tài)，會產(chǎn)生沙漏能量損失，計算得到的沙漏能也很小。

仿真計算得到三種模型下被撞潛艇結(jié)構(gòu)的動能和變形能的對比曲線，分別參見圖11和圖12，很明顯，采用附加質(zhì)量模型二得到的計算結(jié)果與流—固耦合模型的計算結(jié)果比較接近，而附加質(zhì)量模型一與流—固耦合模型的計算結(jié)果相差較大。值得注意的是，附加質(zhì)量模型二得到的潛艇變形能要比流—固耦合模型的略高（圖12），這主要是由于前者忽略了流體介質(zhì)的內(nèi)能，而后者則考慮了流體介質(zhì)的內(nèi)能，根據(jù)能量守恒定律，在其它能量相同的前提下，前者必然要大于后者。

圖10 能量變化曲線Fig.10 The curves of different energy

圖11 不同模型的潛艇動能變化曲線Fig.11 The kinetic energy of different model

圖12 不同模型的潛艇變形能變化曲線Fig.12 The strain energy of different model

6 結(jié) 論

本文建立了潛艇受物體撞擊的三種有限元數(shù)值仿真模型：流—固耦合模型、附加質(zhì)量模型一和附加質(zhì)量模型二，分別進行了數(shù)值仿真計算，通過對兩種附加質(zhì)量模型與流—固耦合模型仿真計算結(jié)果的對比，可以得到如下結(jié)論：

（1）流—固耦合模型考慮了撞擊過程中結(jié)構(gòu)與流場之間的相互作用，能夠反映潛艇受撞過程中的一些本質(zhì)特征，仿真計算結(jié)果最精確，但花費的計算時間也最長；附加質(zhì)量模型可以獲得與流—固耦合模型相近的計算結(jié)果，降低了計算時間，提高了計算效率；

（2）通過合理布置質(zhì)量點單元或增加質(zhì)量點單元的質(zhì)量來使附加質(zhì)量添加到受撞潛艇上建立的附加質(zhì)量模型（附加質(zhì)量模型二）要比通過增加結(jié)構(gòu)材料密度而使附加質(zhì)量添加到受撞潛艇上建立的附加質(zhì)量模型（附加質(zhì)量模型一）的數(shù)值仿真結(jié)果更加接近流—固耦合模型的計算結(jié)果，因此，前一種附加質(zhì)量模型更適合替代流—固耦合模型對潛艇受撞損傷進行數(shù)值仿真計算。

[1]秦洪德,紀 肖,申靜.船舶碰撞研究綜述[C].船舶結(jié)構(gòu)力學學術(shù)會議暨中國船舶學術(shù)界進入ISSC30周年紀念會,2009:191-205.

[2]王自力,顧永寧.船舶碰撞動力學過程的數(shù)值仿真研究[J].爆炸與沖擊,2001,21(1):29-34.

[3]王自力.船舶碰撞損傷機理與結(jié)構(gòu)耐撞性研究[D].上海:上海交通大學,2000.

[4]Jones N.Structural Impact[M].Cambridge:Cambridge University Press,1989.

[5]MSC Software Corporation.MSC.Patran User’s Guide[K].

[6]丁沛然,錢 純.非線性瞬態(tài)動力學分析MSC.Dytran理論及應(yīng)用[M].北京:科學出版社,2006.

[7]王自力,蔣志勇,顧永寧.船舶碰撞數(shù)值仿真的附加質(zhì)量模型[J].爆炸與沖擊,2002,22(4):321-326.

[8]林 一.潛器耐壓結(jié)構(gòu)碰撞響應(yīng)與強度評估方法研究[D].哈爾濱:哈爾濱工程大學,2009.