基于ALE 方法的航行體高速入水緩沖降載性能數值研究*

魏海鵬，史崇鑌，孫鐵志,3，鮑文春，張桂勇,3,4

（1. 北京宇航系統工程研究所，北京 100076；2. 大連理工大學船舶工程學院，遼寧 大連 116024；3. 大連理工大學工業裝備與結構分析國家重點實驗室，遼寧 大連 116024；4. 高新船舶與深海開發裝備協同創新中心，上海 200240）

結構物以一定的速度穿過自由液面的過程被稱為入水過程，入水是一個短暫、變化劇烈而又復雜的力學過程，特別是在入水沖擊瞬間通常會形成很大的沖擊載荷[1]。空投魚雷入水即為一個典型常見的入水問題，入水沖擊過程中由于介質密度的突變形成的軸向力作用可能導致航行體頭部變形、儀器設備失效等，并可能使結構彎曲破壞[2]。因此，降低入水沖擊載荷受到廣泛關注，且具有重要的工程意義。

目前用于入水緩沖降載的措施主要有利用結構物外形降載、主動噴氣降載和利用緩沖組件降載等。利用結構物外形降載指的是通過改變入水結構物外形，從而降低入水過程中的砰擊載荷峰值，主要的措施有改變結構物頭型、在結構物表面設計凹槽等。石漢成等[3]研究了半球頭、錐形頭結構入水載荷，發現錐形頭相比于半球頭結構能有效降低入水瞬時載荷。Sharker 等[4]通過3D 打印入水模型開展入水沖擊動力學實驗研究，發現基于跳水類海鳥的頭型可有效降低入水加速度。Shi 等[5]采用數值模擬方法，通過對比阻力系數研究了不同頭型參數的自主水下航行器(autonomous underwater vehicles, AUV)入水抨擊載荷，發現入水載荷隨著頭部半球角的增大而增大。此外，入水時向外通出非凝結氣體可改變結構物與周圍的接觸環境，從而改善入水沖擊載荷和空泡載荷的影響。張學廣等[6]研究了具有一定深度的凹槽圓柱體入水載荷，相比于平底圓柱體，凹槽內存在的空氣可降低入水沖擊載荷，延長了沖擊作用時間。潘龍等[7]在對頭部噴氣平頭圓柱體入水緩沖機制的研究中也得到了類似的結論。Sun 等[8]對比了通氣空泡和自然空泡對結構表面載荷的影響，發現通氣后可改變自然空化過程與結構表面的劇烈作用，從而降低空泡對結構表面的脈動壓力。利用緩沖組件進行降載是最常見的方法，主要的措施為安裝緩沖罩殼吸收入水產生的能量。緩沖罩殼降載最早由Howard[9]提出，通過在魚雷等入水結構物頭部固定一個流線型頭罩，在入水過程中罩殼破碎，從而降低入水載荷。宣建明等[10]對魚雷頭部模型開展了有緩沖罩殼入水實驗，發現緩沖罩殼可以有效地降低入水沖擊載荷并延長作用時間。錢立新等[11]研究了魚雷罩殼破壞機制和理論設計方法，并進行了入水動力學數值模擬和實驗驗證，提出“撐進破壞”是罩殼的有效破壞模式。Shi 等[12]設計了一種由緩沖罩殼、緩沖墊、定位結構以及連接器組成的緩沖組件，通過數值計算得到安裝該緩沖組件的航行體入水沖擊降載率可以達到42.2%。Li 等[13]通過對使用不同緩沖泡沫材料的緩沖組件進行高速入水數值計算研究，對比不同材料下緩沖罩殼的動態破壞過程，為緩沖泡沫材料的選擇提供參考。

本文中，首先，基于有限元分析軟件LS-DYNA 中的任意拉格朗日-歐拉(arbitrary Lagrangian-Eulerian, ALE)模型，建立高速入水流固耦合數值計算模型。ALE 方法由于結合了拉格朗日方法和歐拉方法的優點，被廣泛應用于帶自由液面的液體晃動問題和固體材料的大變形問題。然后，設計緩沖罩殼模型，對內部結構進行改進并設計罩殼內部緩沖泡沫模型。最后，利用已建立的數值模型對加裝緩沖罩殼后航行體入水過程進行數值計算研究，分析不同緩沖方案下緩沖罩殼和緩沖泡沫的破壞過程以及整體流場演化過程，并對航行體的運動參數進行分析，從而對比不同方案的降載效果，為入水緩沖罩殼設計提供參考。

1 數值方法與計算模型

1.1 基本控制方程與耦合方法

LS-DYNA 中的ALE 算法兼具拉格朗日方法和歐拉方法的優點并改進二者的缺點，在結構邊界運動的處理上引進了拉格朗日方法的特點，能夠有效跟蹤物質結構邊界的運動，在內部網格的劃分上結合了歐拉方法的優點，在對本文中的航行體高速入水并伴隨罩殼和緩沖泡沫破壞的過程模擬時非常有利。流體區域涉及氣相和液相流動，使用歐拉網格進行描述，航行體及緩沖組件由拉格朗日網格進行描述。ALE 算法下的控制方程包括質量、動量和能量守恒方程。

在建立幾何模型及進行網格劃分時，固體結構與流體的幾何形狀以及網格可以重合。通過罰函數耦合約束方法將結構與流體耦合在一起，實現力學參數的傳遞。罰函數約束方式是通過追蹤主從節點間的相對位移從而判斷是否引入界面力，如果出現相對位移即表面的貫穿，則界面力F就會分布到歐拉流體的節點上，界面力的大小與發生貫穿的數量成正比[14]，即：

式中：ki為基于主、從節點質量模型特性的剛度系數；d為結點對主物質表面的貫穿量。由于在每一個時間積分步上都要對等式中的界面力進行求解，可以認為F是等式中一個外力，因此每一時間積分上都可以對總節點力進行求解，由總節點力所引起結構加速度、速度和位移等的變化。

1.2 流體材料模型

在LS-DYNA 中，通過關鍵字*MAT_NULL 來定義流體黏性應力的本構關系，流體介質的壓力由狀態方程進行描述。水介質選擇Grüneisen 狀態方程，并通過關鍵字*EOS_GRüNEISEN 進行設置：

式中：c為聲音在水中的傳播速度；α 為對Grüneisen系數γ0的一階修正；S1、S2、S3為us-up曲線斜率無量綱系數，us為沖擊波速度，up為流體質點的速度；Ew為水的體積內能；μw為水的體積變化率。表1 中給出了水狀態方程的參數，其中Ew0為水的初始體積內能，Vw0為水的初始相對體積。

表1 水狀態方程參數Table 1 Equation-of-state parameters for water

對于空氣介質使用LINER-POLY-NOMIAL 線性狀態方程，通過*EOS_LINEAR_POLYNOMIAL 關鍵字施加，狀態方程壓力公式用下式表示：

式中：μa為氣體體積變化率；Ea為氣體的體積內能；C0～C6為線性多項式狀態方程系數，對于空氣而言，C0=C1=C2=C5=C6=0，C3=C4=0.4。

1.3 計算模型與材料參數

本文的數值模擬與緩沖組件設計方案如圖1所示，航行體最大直徑為324 mm，長2 720 mm。由于本文中只考慮緩沖組件的破壞變形和航行體的運動參數，不考慮自身變形，因此將其視為剛體模型。

圖1 航行體模型示意圖Fig. 1 Schematic diagram of the vehicle

罩殼模型如圖2 所示。罩殼壁厚為5 mm，在罩殼內部預設一定數量的深度為2 mm 的溝槽，以便于控制罩殼的破壞特性并使罩殼及時脫離航行體，同時罩殼中有凸臺用于固定泡沫的位置。

圖2 緩沖罩殼模型Fig. 2 A buffer cover model

材料選擇方面，緩沖罩殼的材料選擇了以不飽和聚氨酯樹脂為基體、短切玻璃纖維為增強基體的復合材料，在瞬時強沖擊載荷下，該復合材料易發生脆性破壞。這種復合材料的參數取自文獻[13]：密度為1 200 kg/m3，楊氏模量為8.5 GPa，泊松比為0.33，屈服應力為45 MPa，切線模量為9 MPa，在LS-DYNA 中利用PLASTIC KINEMATIC 關鍵字定義罩殼材料參數。

緩沖泡沫的材料選擇了聚甲基丙烯酰亞胺(polymethacrylimide, PMI)泡沫，該材料能在保證撞水后吸收能量的同時也較為容易破壞[15]。并在對緩沖泡沫進行分層設計時，選擇了不同力學參數的PMI 泡沫，將單層泡沫模型按長度三等分進行不同參數的設置。PMI 泡沫材料參數見表2。

表2 聚甲基丙烯酰亞胺泡沫(PMI)材料參數Table 2 Material parameters of polymethacrylimide(PMI) foam

網格劃分時，空氣域與水域共節點，航行體、罩殼、緩沖泡沫與流域均采用實體單元進行劃分，并對航行體入水路徑與緩沖組件進行網格加密，提高計算效率，網格劃分完成后共2 104 138個單元，2 945 706 個節點，截取的罩殼與水域部分網格如圖3 所示。

圖3 局部網格劃分示意圖Fig. 3 Part of the finite element mesh

進行數值模擬時，航行體均以150 m/s 的初速度垂直入水，為了對比安裝緩沖組件后的緩沖性能并探討對緩沖泡沫進行分層設計帶來的改變，本文中選取了3 種工況進行數值模擬，具體工況如表3 所示，其中內層、中層和外層泡沫距離航行體頭部由近及遠排布。

表3 工況Table 3 Simulation cases

2 計算結果及分析

2.1 計算模型有效性驗證

為保證數值模型的準確性，下面分別對材料建模、材料破壞模型及ALE 模型參數等進行驗證。首先驗證泡沫材料建模和破壞模型，基于楊洋等[16]開展的密度為110 kg/m3的聚甲基丙烯酰亞胺泡沫（110WF）平面斷裂實驗研究，建立了與實驗相同的泡沫模型并通過關鍵字進行材料參數的定義，對泡沫模型施加與實驗相同的載荷，得到泡沫斷裂形態如圖4 所示。從圖4 可以看到所建立的泡沫模型和破壞模型能夠較好地模擬材料的破壞過程。

圖4 材料破壞模型驗證Fig. 4 Validation of the material failure model

驗證材料建模和材料破壞模型后，下面將驗證文中使用的ALE 數值計算方法的合理性。基于Chen 等[17]開展的不同頭型彈體高速入水實驗，選擇扁平頭型彈體，入水初速度為106.8 m/s的垂直入水工況進行對比，計算模型尺寸和初始條件均與實驗相同。實驗結果、數值計算結果對比如圖5 所示，可以看到本文中所使用的ALE數值計算模型能夠較好地模擬結構物高速入水過程。綜上所述，本文中所使用的材料建模、材料破壞模型和數值算法具有較高的有效性。

圖5 數值計算模型驗證Fig. 5 Validation of the numerical calculation model

2.2 航行體入水流場演化特性與緩沖組件整體破壞過程

首先從航行體入水空泡演化過程與緩沖組件整體破壞特性進行分析。圖6 給出了使用單層71WF泡沫作為緩沖材料的工況空泡演化與緩沖組件破壞的過程，將罩殼與水接觸的時刻定義為時間零點。可以看到入水時罩殼頭部首先接觸自由面，在強沖擊載荷的作用下發生破壞，隨后緩沖泡沫進一步吸收撞擊產生的能量也迅速破壞，在t= 1.5 ms 時看到罩殼在航行體、水與泡沫的共同作用下發生了變形彎曲，罩殼此時沒有發生開裂；而在t= 2.0 ms 時，罩殼已經沿預設溝槽開裂，發生撐進破壞；由于71WF 泡沫最易破壞，隨著航行體的下落，泡沫也同時完全失效，因此起到的吸能緩沖效果有限，隨著緩沖過程的結束，罩殼殘片向周圍散開雖然沒有完全破壞，但對空泡影響較小。

圖6 工況2 入水流場演化與破壞過程Fig. 6 Flow field evolution and destroyed process in case 2

圖7 為將緩沖泡沫分為3 層作為緩沖材料的入水空泡演化與緩沖組件破壞過程，將泡沫模型三等分并分別賦予3 種泡沫材料參數，靠近航行體部分選擇最容易破壞的71WF 泡沫，靠近罩殼部分選擇最難破壞的200WF 泡沫，中間段為110WF 泡沫。在罩殼頭部破壞后，由于泡沫硬度的增加與其流線型的結構，在t= 1.0 ms 可以看到泡沫仍保持較為完整的形態不被破壞；t= 1.5 ms 時中后段泡沫在水與航行體的共同作用下擠壓變形，導致罩殼沿溝槽開裂；t= 2.0 ms 時中后段泡沫完全破壞，到t= 3.0 ms 時頭部200WF 泡沫仍保持一定的形態，由此可見分層后的泡沫可以延長緩沖作用的時間，出現二次緩沖的效果。

圖7 工況3 入水流場演化與破壞過程Fig. 7 Flow field evolution and destroyed process in case 3

2.3 緩沖罩殼、泡沫破壞特性與應力云圖

圖8 和圖9 分別給出了工況2、工況3 入水過程中緩沖罩殼破壞過程及應力云圖。從圖8～9 中可以看出在撞水時罩殼頭部均出現了明顯的應力集中現象，同時罩殼尾端邊緣與溝槽處首先產生了應力，應力以應力波的形式在罩殼傳遞。罩殼頭部與預設溝槽處應力更大，頭部破壞后沿預設溝槽發生開裂，并最終完全破壞。使用分層泡沫也會對罩殼的破壞時間產生影響，工況2 在t= 1.5 ms 時罩殼開始出現裂痕，而工況3 在t= 1.1 ms 時罩殼已開始開裂；同時也會影響罩殼的破壞形式，工況2 罩殼完全開裂后，剩下的殘片長度較大，而工況3 罩殼完全開裂后余下殘片被破壞成兩部分。

圖8 工況2 罩殼應力云圖與破壞過程Fig. 8 Stress cloud of the nose cap and destroyed process in case 2

圖9 工況3 罩殼應力云圖與破壞過程Fig. 9 Stress cloud of the nose cap and destroyed process in case 3

圖10 和圖11 分別給出了緩沖泡沫破壞過程與應力云圖，可以看到應力從緩沖泡沫頭部開始沿徑向和軸向傳遞并導致泡沫發生破壞。從時間來看，工況3 泡沫完全破壞的時間遠大于工況2；從破壞形態來看，工況2 的單層71WF 泡沫在罩殼破壞的同時開始破壞，而工況3 對泡沫分層后，頭部200WF 泡沫在罩殼破壞后并不會迅速破壞，中后段泡沫破壞后，在航行體與水的共同作用下才發生破壞，可以說明分層后會出現二次吸能的作用，具有更好的緩沖效果。

圖10 工況2 泡沫應力云圖與破壞過程Fig. 10 Stress cloud of the foam and destroyed process in case 2

圖11 工況3 泡沫應力云圖與破壞過程Fig. 11 Stress cloud of the foam and destroyed process in case 3

圖12 和圖13 為緩沖泡沫中截面的應力云圖，可以清楚地看到應力在不同層泡沫之間的傳遞以及泡沫的破壞規律。應力首先集中在最外層，隨后向內層傳遞。從圖12 可以看到單層泡沫工況應力傳遞比較平滑，且中間應力傳遞速度最快。而圖13中則可以明顯的看到分層后對應力傳遞造成的影響，每層泡沫能承受的最大應力不同，撞水后應力出現分層，當傳遞到內層泡沫時由于其強度最弱，因此最先發生破壞，其次是中間層泡沫，最后為外層泡沫，可以看到分層后的泡沫破壞形式發生了明顯的變化。

圖12 工況2 泡沫內部應力云圖與破壞過程Fig. 12 Stress cloud of the inside foam and destroyed process in case 2

圖13 工況3 泡沫應力云圖與破壞過程Fig. 13 Stress cloud of the inside foam and destroyed process in case 3

為了更直觀地對比不同緩沖方案罩殼與泡沫破壞的相關特征時刻，提取了罩殼頭部開始破壞、罩殼中段開裂、罩殼完全破壞與泡沫頭部開始破壞、泡沫完全破壞的時刻并繪制了散點圖，如圖14 所示。可以看到泡沫分層后罩殼各個時刻均出現提前，總破壞時間縮短；相反泡沫的破壞時刻出現延后，并且總破壞時間延長。

圖14 特征時刻散點圖Fig. 14 Scatter plot of characteristic moments

2.4 航行體運動參數對比分析

航行體的運動參數能夠反映出加裝緩沖組件后對高速入水沖擊載荷的改善情況，也是判斷緩沖性能的重要指標，因此本部分內容通過提取不同工況運動參數進行對比分析。圖15 為不同工況的航行體位移曲線，可以看到由于入水過程是一個極短暫過程，因此在相同時間內航行體位移大致相同。通過放大t= 10 ms 時刻附近的曲線能看到工況3 具有最大的位移，工況1 位移最小，說明安裝緩沖組件后會增加入水相同時間內的位移，位移增加量與緩沖方案相關。

圖15 不同工況下航行體的位移-時間曲線Fig. 15 Displacement-time curves of the vehicle in different cases

航行體的速度曲線也能體現出緩沖性能的優劣，圖16 給出了不同工況的速度曲線。可以看到航行體速度在入水初期階段衰減最快，沒有緩沖組件時航行體在撞水瞬間受到巨大沖擊載荷，速度以幾乎垂直的斜率下降。工況2 由于緩沖罩殼先吸收了一部分能量，因此在罩殼撞水后航行體速度降低較為平緩，但隨著泡沫與罩殼的完全破壞，速度也出現了快速降低。工況3 由于對泡沫進行了分層處理，各層緩沖能力不同，因此速度曲線始終以較為平緩的趨勢降低，也從側面體現了工況3 的方案具有最好的緩沖效果。

圖16 不同工況下航行體的速度-時間曲線Fig. 16 Velocity-time curves of the vehicle in dfferent cases

最后提取了航行體的加速度曲線，如圖17所示，加速度曲線的峰值直接說明了受到沖擊載荷的大小。在無緩沖作用下入水加速度瞬間上升到峰值，隨后迅速下降，在極短時間內受到巨大的沖擊；工況2 中航行體受到的沖擊明顯減小，說明該方案具有較好的緩沖效果；而工況3 中航行體加速度曲線變化趨勢發生改變，出現了兩個加速度峰值，并且峰值遠小于無緩沖的工況，證明了將緩沖泡沫進行分層設計后出現了二次緩沖的作用，顯著降低了高速入水時受到的沖擊。同時將工況2、工況3 特征時刻的緩沖泡沫進行提取，可以看到載荷最大值出現在泡沫完全破壞時刻，將泡沫分層后，內層泡沫先破壞時航行體也會受到一個較為明顯的沖擊載荷。

圖17 不同工況下航行體的加速度-時間曲線Fig. 17 Acceleration-time curves of the vehicle in different cases

將無緩沖工況1 航行體入水加速度峰值作為基準，計算了工況2 與工況3 緩沖方案的降載率，如圖18所示。可以看到將泡沫分層設計的工況3 降載率為73.2%，緩沖效果遠好于使用單層泡沫的工況2（降載率31.8%）。

圖18 降載率對比Fig. 18 Comparison of load reduction ratios

3 結 論

對航行體高速入水緩沖降載問題進行了研究，設計了緩沖罩殼與相應的緩沖泡沫，利用ALE 方法對入水時罩殼與泡沫的動態破壞過程進行了數值計算研究，并分析了不同緩沖方案的緩沖性能。結果表明：

（1）通過在航行體頭部安裝緩沖罩殼與緩沖泡沫，能夠在撞水時吸收一定的沖擊能量發生破壞，罩殼的破壞形式受泡沫變形的影響，同時材料的選擇也會影響泡沫的破壞形式，將緩沖泡沫分層后，能夠加快罩殼的破壞并延長自身的破壞時間，本文中設計的罩殼與泡沫均能夠及時破壞并脫離航行體。

（2）應力主要以應力波的形式在罩殼與泡沫中傳遞，撞水時應力會明顯集中于罩殼頭部與預設的溝槽處，罩殼會沿預設溝槽進行破壞，可以通過改變預設溝槽的形式調整罩殼破壞特性。泡沫分層后頭部堅硬部分在入水一段時間后在航行體與水的共同作用下破壞，產生了二次緩沖的效果。

（3）安裝緩沖組件后會增加航行體入水后相同時間內的位移，使速度曲線更平緩，通過對比加速度曲線，可知使用分層緩沖泡沫的方案具有最好的緩沖效果，降載率可達到73.2%，而單層泡沫由于其破壞速度較快，吸能效果較差，降載率為31.8%，兩種方案均能有效降低入水沖擊載荷。