黃曉明，朱錫，牟金磊，李海濤

（海軍工程大學船舶與動力學院，湖北武漢 430033）

0 引言

水中爆炸載荷是海戰(zhàn)中各類艦艇面臨的主要威脅，艦艇在水中爆炸作用下的響應是一個強烈的動態(tài)過程，艦船自身在水中的總體振動特性，特別是其低階模態(tài)，對響應結果有著顯著的影響。因此，對艦艇水中振動特性的分析是進行艦艇抗爆研究的基礎。

艦艇水中振動模態(tài)也被稱為艦艇濕模態(tài)，船舶工程領域對船舶濕模態(tài)問題開展了長期的研究，形成了一系列行之有效的計算方法。早期，劉易斯、托德等人在試驗研究基礎上提出了附加質量計算公式，將航道、船體等對周圍流體三維流動的影響，繪制成系列的圖譜進行修正［1］;此后，科文－克勞夫斯基進一步發(fā)展了附連水計算的切片理論［2］;錢勤等利用勢流理論，考慮結構變形對附連水質量的影響，推導出無限流場中無限長圓柱殼的附連水質量和圓柱殼振動波數的關系［3］。J．A．Deruntz和T．L．Geers利用邊界積分法計算附連水質量，用結構浸水的邊界元代替三維流場，從而將三維問題降為二維問題，減少了計算量［4］;邊界元法僅限于規(guī)則結構在簡單水域中的求解，而試驗研究的成本較高，且受自然條件的影響，使用有限元方法對結構濕模態(tài)特性進行預報，即可考慮復雜的結構和水域等因素，同時可節(jié)省大量高昂的試驗費用，是當前的研究熱點［5－7］。濕模態(tài)有限元分析的關鍵在于水域的建模與水體網格劃分，為模擬無限流場需要建立較大尺寸的水域，為得到精確的結果則要求網格細分，然而對于三維水體上述2個方面的要求會直接導致計算規(guī)模的驟然增大。濕模態(tài)仿真計算的難點就在于在計算精度與計算資源中尋求平衡。

本文針對船體整體結構低階垂向濕模態(tài)特性，使用Abaqus軟件對有限元濕模態(tài)仿真的參數控制問題進行了研究，具體比較了網格尺度，水域形狀，網格劃分方式對計算結果的影響，本文計算方法適用于水中結構的低階濕模態(tài)分析，對工程實際問題有較大的參考價值。

1 有限元聲固耦合理論基礎

Abaqus把流體當作聲學介質，在滿足壓力的小范圍變化假設條件下，能模擬濕模態(tài)以及水下爆炸等多種流固耦合問題。其計算理想流體中結構模態(tài)的方法是，將流體壓力作為流場中的基本未知量來求解結構的水下模態(tài)。其基本原理是求解方程（1）［8］:

式中:M，K和P分別為質量矩陣、剛度矩陣和載荷矩陣;s為結構;f為流體;C為阻尼矩陣;S為流固耦合矩陣，代表流固耦合界面上每個節(jié)點處的有效面積，它將界面上的流體壓力轉換成結構所受載荷;U和P分別為節(jié)點位移矩陣和流體聲學壓力矩陣。

2 基準模型試驗

仿真計算的準確性需以試驗結果作為基準，本文試驗模型為一細長加筋圓柱殼，如圖1所示。截面為圓形，直徑為85 mm，模型長度3 860 mm，殼板厚度1 mm，均勻設置10根外肋骨用來模擬艙壁，肋骨間距為322 mm，肋骨尺寸為2 mm×10 mm，采用Q235鋼制造。

按照船體振動學中方法估算其一階干、濕模態(tài)頻率，如式（2）［1］:

圖1 試驗模型照片Fig.1Photo of model

式中:fi為第i階的固有頻率，Hz;l為梁的長度;E為彈性模量;I為截面慣性矩;μi為第i階的頻率參數; m為單位長度的質量。對于Q235鋼，彈性模量E= 2.05e11Pa，材料密度ρ=7 850 kg/m3。當計算干模態(tài)頻率時，m=ρA，A為梁的截面積，當計算濕模態(tài)頻率時，需計入附連水質量。全浸沒圓形剖面柱體垂向運動時，單位長度上的附加水質量為

式中:maV為單位長度上的附加水質量，kg/m;ρ0為水的密度，kg/m3;R為計算剖面處的剖面半徑，m。按照式（2）和式（3）計算得到的模型干模態(tài)一階頻率為35.3 Hz，濕模態(tài)一階頻率為18.9 Hz。采用力錘法測定了模型的低階垂向固有頻率，其干模態(tài)一階頻率為34.2 Hz，濕模態(tài)一階頻率為17.2 Hz。模型干濕模態(tài)的理論計算結果與試驗結果的比較見表1?？梢妼τ趫A柱殼模型，理論估算結果略大于模型實際頻率，但仍能滿足工程預報的精度要求。

3 仿真建模概況與干模態(tài)計算結果

對模型進行無限水域中的模態(tài)仿真計算，結構采用殼單元來模擬，外部水體用聲學單元模擬。結構表面與水體的接觸面用關鍵字TIE，使用該關鍵字可以使水與結構表面始終保持接觸狀態(tài)而不分離。圓柱殼的有限元模型由4 812個4節(jié)點殼單元S4R進行模擬，其有限元模型如圖2（a）所示。模型材料為Q235鋼，E=2.05e11Pa，λ=0.3，ρ=7 850 kg/m3。

首先對結構模型進行干模態(tài)計算，得到模型的干模態(tài)一階頻率為32.3 Hz，振型如圖2（b）所示。對比試驗結果（34 Hz），可見干模態(tài)有限元計算的結果較為準確。

4 仿真參數對計算結果的影響

4.1 網格劃分方式的影響

圖2 圓柱殼有限元模型及干模態(tài)仿真結果Fig.2Finite element model and dry frequency result

對水域中的網格劃分方式進行比較，一種使用結構化的六面體網格劃分，一種使用非結構化的四面體網格劃分。比較其他參數一致時，2種網格劃分方式對計算結果和計算時間的影響，計算結果見表2。由表2可見，結構化六面體網格劃分方式能獲得較高的計算精度同時節(jié)約大量的計算時間。因而對于流固耦合系統(tǒng)的流體部分，應優(yōu)先使用結構化六面體網格劃分方法。

4.2 水域建模方式的影響

流固耦合系統(tǒng)中模擬無限場，流域的邊界距結構的距離至少應為結構尺度的6倍左右［9］。針對模型建立流場，流場半徑取圓柱殼半徑的6倍。以往的流固耦合仿真中［10］流場水域通常劃分為中部圓柱形，兩端均為與中間環(huán)柱體外半徑相同尺寸的半球體的形狀（見圖3（a））。但這種水域模型增加了網格劃分時的復雜度，網格質量較難控制。另外一種簡化的水域形狀如圖3（b）所示。即使用圓柱體代替半球體，該水域模型簡化了端部的網格劃分，但在流體與圓柱殼接觸區(qū)域的網格劃分也較為復雜。對于本文關心的低階垂向振動模態(tài)，嘗試省略端部的水域，即只保留中部的圓柱形水域，而對端部不進行建模，如圖3（c）所示。對3種水域均使用結構化的六面體網格進行劃分，并使用相同的網格尺度。仿真中水的體積模量取值Kf=2.0e9Pa，水的密度ρf=1 000 kg/m3。

3 種水域模型下的圓柱殼一階垂向振動濕模態(tài)振形相近，圖3（d）給出了無封蓋水域中的計算結果，不同水域中濕模態(tài)一階垂向振動頻率結果和計算耗費時間見表3。由表3可知，3種水域模型下計算的結果相近。有封蓋水域，網格劃分復雜，計算結果更加偏離試驗值，同時花費更多的計算時間。而無封蓋水域，網格更加均勻，其精度也更高，因此本文后續(xù)的仿真均使用這種無封蓋的水域建模方式。

圖3 不同水域建模方式及無封蓋水域濕模態(tài)計算結果Fig.3Different water building mode and the wet mode simulate results of water without capping

4.3 水域網格尺度的影響

有限元計算中網格尺度對計算結果有較大影響，本文對水域網格尺度對濕模態(tài)計算結果的影響進行分析，對于船體整體結構模型，其長度方向較截面方向的尺度要大得多，因此對長度方向和截面方向的網格尺度分別進行考察。

在水域中沿結構長度方向的網格尺度對模型垂向振動計算結果影響見表4。由表4可見，長度方向網格數在20～80之間，計算結果變化不大，且都與試驗值較為接近。本文后續(xù)計算中水域沿梁長方向劃分為20個網格。

對于圓柱殼模型，定義模型的截面半徑Lmd= 0.042 5 m為結構特征長度。設截面徑向方向的網格邊長為網格特征長度Led，使用無因次量λd=Led/Lmd來描述截面徑向網格劃分的疏密程度。研究網格密度對計算結果的影響。計算結果見表5。由表5可知，當λd＜1/2時，計算結果與試驗結果的誤差約為5%。當λd＜1/4時，計算結果與試驗結果的誤差小于2%。

5 仿真求解策略

綜合以上的計算結果可見，對于無限場水域中圓柱殼結構橫向低階濕模態(tài)仿真計算，使用表6中仿真策略，能在較小的計算資源下獲得合理的計算結果。

6 結語

1）在整體結構低階垂向濕模態(tài)的有限元仿真計算中，水域建?？梢院喕蔀闊o封蓋水域，不但能大為降低水域網格劃分的復雜性，同時由于網格更加均勻，計算精度更高。

2）在水域中沿結構長度方向的網格數取20～80，水域徑向網格尺度小于結構截面特征尺度1/2時，能獲得較為精確的計算結果。

3）結構化六面體網格劃分方式能獲得較高的計算精度。因而對于流固耦合系統(tǒng)的流體部分，應優(yōu)先使用結構化六面體網格劃分方法。

本文所提供的仿真參數控制策略適用于水中結構的低階濕模態(tài)分析，可作為工程實際問題的數值計算方法。

［1］翁長儉，張保玉．船體振動學［M］．大連:大連海運學院出版社，1992．91－95．

［2］董艷秋．船舶波浪外載荷和水彈性［M］．天津:天津大學出版社，1991．36－39．

［3］錢勤，黃玉盈，劉忠族．求附連水質量的一種直接方法［J］．力學實踐，1996，18（5）:19－21．

QIAN Qin，HUANG Yu-ying，LIU Zhong-zu．A direct method for added water mass Calculate．［J］．Mechanics and Engineering，1996，18（5）:19－21．

［4］JOHN，DERUNTZ A，GEERS．Added mass computation by the boundary integral method［J］．International Journal for Numerical Methods in Engineering，1978，12（3）:66－72．

［5］卓頡，李兵．水下復雜聲源結構自振特性分析［J］．計算機輔助工程，2008，17（3）:33－35．

ZHUO Jie，LI Bing．Characteristics of self-vibration analysis on acoustic source with complicated structure under water［J］．Computer Aided Engineering，2008，17（3）:33－35．

［6］劉曉明，張悅，高志龍．集裝箱船振動及響應分析［J］．中國造船，2004，45（4）:23－28．

LIU Xiao-ming，ZHANG Yue，GAO Zhi-long．Structural vibration and response analysis on container ship［J］．Shipbuilding of China，2004，45（4）:23－28．

［7］蘇海東，黃玉盈，陳琴．片狀域流場中三維結構的濕頻率和濕模態(tài)［J］．華中科技大學學報，2007，35（1）:121－124．

SUHai-dong，HUANGYu-ying，CHENQin．Wet frequencies and wet modals of 3-D structures in an infinite liquid with constant depth［J］．Journal of Huazhong University of Science and Technology（Nature Science Edition），2007，35（1）:121－124．

［8］Dassault systems similar Corp．AbaqusTheoryManual［M］．Pawtucker USA:ABAQUS，INC，2008．

［9］張阿漫．水下爆炸氣泡三維動態(tài)特性研究［D］．哈爾濱:哈爾濱工程大學，2006．

［10］董海，劉建湖，吳有生．水下爆炸氣泡脈動作用下細長加筋圓柱殼的鞭狀響應分析［J］．船舶力學，2007，11（2）: 250－258．

DONGHai，LIUJian-hu，WUYou-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．