船體三維變形響應的數值預報

任慧龍，于鵬垚，李 輝，馮國慶

（哈爾濱工程大學船舶工程學院，黑龍江哈爾濱150001）

船體三維變形響應的數值預報

任慧龍，于鵬垚，李 輝，馮國慶

（哈爾濱工程大學船舶工程學院，黑龍江哈爾濱150001）

針對航行中艦船變形問題，采用三維水彈性理論實現了波浪中艦船的動態變形預報。通過擬合不同模態下船體有限元模型的節點位移，實現了將三維模態信息應用于廣義流固耦合界面條件。建立了艦船在規則波中的廣義水彈性運動方程，求解出各階振動模態的主坐標。將各階彈性模態主坐標與預報位置的角變形矩陣相乘，得到角變形響應。針對某實船，分別采用三維水彈性方法和三維剛體方法進行計算，通過對船體結構角變形響應的對比，驗證了三維水彈性方法的正確性，同時也證明了該方法在計算效率上的優勢。

艦船變形；流固耦合；三維水彈性方法；模態分析；變形響應；數值預報

隨著船舶的大型化，船體在波浪中航行時，由周圍波動壓力作用引起的彈性變形愈發明顯。對于某些特殊用途的船舶（如科考船、軍艦等），較大的動態變形會影響設備的正常使用，進而影響其正常作業。因此，國內外相關學者提出一些方法對船舶營運時的波浪載荷誘導的動態變形響應進行預報。陳曙梅等［1?4］分2步進行載荷預報和結構響應計算，這些學者先將艦船看作三維剛體結構來計算不同工況下波浪載荷，再將載荷加載到三維有限元結構上來獲取船體變形，這種方法在后續論述中簡稱為剛體方法。陳敏康等［5?6］等采用二維水彈性理論對船體的變形進行預報。

剛體方法是在準靜態假定下進行計算，無法反映真實船體結構動力響應特性，而且在計算船體角變形幅頻響應時，需要針對較多頻率進行壓力加載和靜力分析，工作量較大。二維水彈性理論是將船體結構簡化為一根變截面的船體梁，不能夠考慮船寬方向上船體結構的變形。本文利用三維有限元模型的模態分析結果，采用的三維水彈性預報方法實現了船體三維變形響應的數值預報。

1 船體角變形矩陣的建立

對角變形響應而言，關心的是一個位置相對另一位置的角變形，如船體上某位置相對于羅經平臺的變形。本文引入角變形矩陣X3n×m來表達不同模態下不同位置間的相對角變形振型。為了便于進行描述，引入如圖1所示坐標系，其中原點位于重心在靜水面上的投影，x軸有船舵指向船艏，y由右舷指向左舷，z軸數值向上。

圖1 坐標系Fig.1 Coordinate system

本文采用三維有限元方法對船體結構進行干模態分析，并以艉垂線與上甲板交點的線位移進行歸一化，得到各階彈性模態對應的振型，廣義質量矩陣a，廣義剛度矩陣c。

2 流固耦合界面條件的數值處理

基于勢流理論求解艦船周圍的流場時，需要求解流場的繞射勢φd和輻射勢φk。通常可采用源匯分布法［7?9］，結合自由面條件、遠方條件和廣義的流固耦合的界面條件［10］進行求取。

廣義的流固耦合界面條件：

式中：φ0為已知的單位波幅入射勢；ω為波浪自然頻率；n為面元單位法向量；uk為第k模態下的面元中心點位移；W表示船體周圍流場的定常速度，忽略定常速度勢時，W＝-Ui，i為x坐標軸的單位方向向量，U為航速。

從式（2）可以看出，船體結構的振型信息是通過面元中心點的位移和該位移對x、y、z方向的偏導數體現在彈性結構的水動力求解中。因而，問題的關鍵就是如何通過三維模型在不同模態下的節點的位移信息來獲得水動力求解中的面元中心點的位移和位移對x、y、z方向的偏導數。本文采用擬合位移形函數的方法來解決這一問題。

船體結構三維有限元模型的節點線位移矢量u＝（ux，uy，uz）包含3個方向的位移分量，本文以ux方向位移為例，推導出如何利用三維有限元模型外殼的振型信息確定對應模態下面元中心點的位移及其導數，其他兩個方向位移的推導方法與其相同，不再贅述。

為便于后續描述，作如下定義：將面元中心點到面元節點的最遠距離定義為面元控制域半徑，以該半徑作球的區域定義為面元控制域。則面元中心點的位移和位移的導數由面元控制域內節點的振型信息來確定。

對于空間運動的剛體，當給定剛體質心處的6自由度運動，便可以得到任意位置的位移：

式中：un為待求任意一點的線位移；uc、ωc分別為質心處的線位移和角位移；r為待求位置與質心之間的向徑；xn、yn、zn為待求位置的坐標；xc、yc、zc為質心坐標；i、j、k分別為x、y、z坐標軸的單位方向向量。

可見，任意節點線位移可由節點坐標與剛體質心處運動相乘的形式來表示。借鑒這種思路，當面元不與坐標軸垂直時，假定船體外殼上，面元對應區域內不同位置的x方向位移形函數具有如下形式：

式中：a0、a1、a2、a3為待定系數。

當面元與坐標軸垂直時，如與z軸垂直，控制域內的各個節點可能為同一個值，若同時為零，則a3為任意值，可導致后續求解方程組時系數矩陣奇異。故當面元與坐標軸垂直時，假定型函數的形式如下：

可以看出想要確定位移型函數，面元控制域內節點不得少于4，當少于4時可以適當增加控制域半徑來增加節點。

若控制域內包含n個節點，可根據振型信息和位移型函數確定n個方程。方程的個數往往大于待定系數，即n＞4，此時可以基于最小二乘理論確定方程的系數，進而確定型函數，代入面元中心點坐標得到面元中心點位移。型函數分別對x、y、z方向求導，代入面元中心點坐標便得到面元中心點位移的導數。

對于與坐標軸垂直的面元，這里定性的確定該面元中心點位移對所垂直坐標軸方向的導數。如當面元垂直x軸時，該面元應位于結構的端部，通常對應于船體艏艉兩端，該處對應的剖面內力為零，則各方向位移對x的導數應為零。同理對于垂直于y軸或z軸的面元，各方向位移對y軸或z軸的導數應為零。

3 船體角變形響應的求解

根據水彈性力學建立起適用于規則波激勵下的廣義線性水彈性運動方程［11］為

式中：a為船體廣義質量矩陣，b為船體廣義剛度矩陣，a和c可由前面的三維有限元結構模態分析獲得；b為船體廣義阻尼矩陣，可通過試驗或經驗公式獲得；pa為主坐標幅值列向量；ωe為遭遇頻率；A、B、C為廣義水動力系數，具體含義分別為廣義流體附加質量、廣義流體附加阻尼和廣義流體恢復力系數矩陣；F為廣義波浪激勵力。

對于廣義水動力系數A、B、C和廣義波浪激勵力F，具體可依據下式求得：

式中：ρ為流體的密度；ρb為船體結構材料的密度；ds為船體濕表面的面積微元；dV為船體結構的體積微元；n為面元中心點法向量；ur為第r階模態下的面元中心點位移；φk為第k階模態對應的輻射勢；wk為第k模態下面元中心點的垂向位移；gs表示重力加速度矢量，gs＝（0，0，-g）；θk表示第k階船體振型產生的結構轉角變形，θk＝（?×uk）／2；ζa為入射波波幅；φ0為已知的單位波幅入射勢；φd為單位波幅繞射勢；U為航速。

綜上，在確定方程（7）的各項系數后，可求得在波浪激勵下不同模態的主坐標響應。主坐標響應與角變形矩陣X3n×m相乘便可得到船體的任意位置間的角變形響應。

4 應用實例

利用提出的基于三維水彈性理論預報變形的方法，對某艦船在波浪作用下的角變形響應進行預報。首先，將三維水彈性得到的運動響應與傳統的三維剛體載荷方法進行比較，然后再將剛體方法得到的波浪載荷加載到三維有限元模型上進行準靜態分析，并提出取船體上任意位置間的相對角變形響應，最后將兩種方法得到的角變形響應進行比較。

考慮到對于船體上任意位置間的相對角變形響應預報方法相同，本文僅以船體甲板上某一位置相對于船舯附近的羅經平臺的變形為例進行角變形響應預報。為了便于設備的使用，預報點位置和羅經平臺都在船體的強力構件附近，并采取了局部加強。

4.1 艦船主尺度數據及關注位置

選取某船舶進行船體三維變形響應的數值預報的相關研究，船舶的主尺度數據和關注位置坐標分別如表1、2所示。

表1 艦船主尺度數據Table1 Principal parameters of the ship

表2 關注位置坐標Table2 Coordinates of concerned positions

4.2 干模態分析

本文采用有限元軟件MSC.nastran對船體整船結構干模態進行分析，考慮到關注的位置位于甲板上，在圖2中列出彈性模態對應甲板位置的振型。表3列出各階模態固有頻率、廣義質量和無因次阻尼系數。表4給出關注位置相對羅經平臺的角變形矩陣。實際計算中，所選取的彈性模態的數目將在后續進行進一步的討論。對于船體自身的結構阻尼，理論計算很難實現精確的確定，而且結構阻尼主要影響船體與波浪共振時的響應峰值，對低頻波浪激勵下的結構響應影響不大，因而本文計算中的阻尼系數主要參照相關文獻［12］的經驗公式進行確定。本文僅以前4階彈性模態為例列出模態分析結果，對于更多的模態只需擴展相應的行或列。

圖2 彈性模態甲板振型Fig.2 Deck shape of elastic principal mode

表3 船體干模態固有頻率、廣義質量和無因次阻尼系數Table3 Eigenfrequencies，generalized mass of flexible modes and nondimensional damping coefficients

表4 角變形矩陣Table4 Angular deformation matrixmin

由船體干模態分析結果可以看出，一階彈性模態振動形式以兩節點的垂向振動為主。若采用二維的船體梁模型進行模態分析，對于對稱船體，該模態下任意位置處繞x和z軸的角變形均為0，而由三維有限元方法得到相對角變形矩陣可以看出，在一階彈性模態下，繞x和z軸的對應元素并不0，說明由于船體結構的局部特點影響，在以垂向振動為主模態下，對于繞x和z軸的角變形也有貢獻。

4.3 彈性模態數目的選取

對于實際的船體結構，從理論上說，它具有無限種振動形式，然而在采用模態疊加方法計算某一種外力激勵下的結構響應時，不可能考慮到結構所有模態，只能夠盡可能的選取對結構響應貢獻較大的模態用于計算。對于本文計算的船體結構，結構響應主要是由于波浪的激勵引起的，而實際船體的固有頻率通常大于波浪頻率的分布范圍，因而固有頻率越低的模態越接近波浪頻率，對結構響應貢獻越大。本文為考察究竟多少個彈性模態能夠充分反映船體關注位置的角變形響應，分別采用前4、8、12階彈性模態用于水彈性計算。計算時，艦船航速為20 kn，浪向角為30°。圖3為選取不同彈性模態數目進行計算時的角變形幅頻響應曲線。

通過計算可以看出，選取4個彈性模態數目的計算與選取8和12個彈性模態的計算結果差別較大，選取8個和12個彈性模態的計算結果幾乎相同，說明選取8個彈性模態已經能夠充分反映出該海況下船體關注位置的角變形響應，其中第8階彈性模態對應的干模態固有頻率為34 rad／s，大于波浪能夠激勵的結構干模態頻率范圍。因而后續采用三維水彈性方法進行計算時均采用8個彈性模態。

圖3 不同模態數目下角變形幅頻響應曲線Fig.3 Angular deformation amplitude?frequency re?sponse curves of different mode numbers

4.4 變形響應的對比

在航速20 kn，浪向角30°時，分別采用三維水彈性方法和加載三維剛體方法的波浪載荷的方法預報船體的角變形的幅頻響應，如圖4。

圖4 2種方法的角變形幅頻響應曲線Fig.4 Angular deformation amplitude?frequency re?sponse curves of the two methods

通過對比可以看出，2種計算結果的趨勢基本一致，而且繞y軸的角變形明顯大于繞x、z軸的角變形。對于繞x、z軸的角變形，剛體方法的計算結果大于水彈性方法，對于繞y軸的角變形，2種方法基本一致。分析如下：三維水彈性方法是一種考慮到流固耦合的動力分析方法，而剛體方法是在船體是剛體，結構響應是準靜態的假定下進行的計算，因而2種方法的計算結果必然存在偏差。總體的看來，對于船體角變形比較大的方向（繞y軸角變形），2種方法預報結果幾乎相同，一定程度上驗證了水彈性方法正確性。

5 結論

對于預報艦船在波浪中角變形這一流固耦合問題，本文采用擬合位移形函數的處理方法，實現將三維有限元的模態分析結果應于與三維水彈性理論來預報結構的角變形響應。通過計算表明：

1）本文所選用的前8階彈性模態已足夠可以反映出關注位置處的角變形響應。

2）通過剛體方法和水彈性方法的變形響應對比，驗證了采用擬合位移形函數的方法處理流固耦合界面條件的正確性。

3）通過水彈性方法，避免了將水動壓力加載到有限元模型上，可以快速的獲得艦船的變形幅頻響應，方便了后續變形響應的長短期分析和可靠性分析。

［1］陳曙梅.航行狀態下艦船結構總體變形分析［D］.哈爾濱：哈爾濱工程大學，2004：7?42.

CHEN Shumei.Deformation analysis of whole vessel in navi?gating state［D］.Harbin：Harbin Engineering University，2004：7?42.

［2］LIU Aili，MA Hongxu，Dai Hongde.Estimation method for ship deformation under wave Loads［C］／／International Con?ference on Mechanic Automation and Control Engineering.Wuhan，China，2010：93?96.

［3］MURAWSKI L.Shaft line alignment analysis taking ship construction flexibility and deformations into consideration［J］.Marine Structures，2005，18：62?84.

［4］YU Pengyao，FENG Guoqing，REN Huilong，et al.Deform?ation analysis and reliability assessment of ship hull in irreg?ular waves［C］／／Proceedings of the International Conference on Offshore Mechanics and Arctic Engineering.Nantes，Frances，2013.

［5］陳敏康，邱強.隨機波浪中船體變形轉角的估算方法及應用［J］.船舶力學，2006，10（2）：87?95.

CHEN Minkang，QIU Qiang.Estimate method and its appli?cations for angular deflection of girder in random wave［J］.Journal of Ship Mechanics，2006，10（2）：87?95.

［6］TAYLOR R E.Reflections on the structural dynamics of floating beams and ship hulls in waves［J］.Journal of Engi?neering for the Maritime Environment，2009，223（3）：257?273.

［7］劉應中，繆國平.船舶在波浪中的運動理論［M］.上海：上海交通大學出版社，1987：56?80.

［8］戴遺山，段文洋.船舶在波浪中的勢流理論［M］.北京：國防工業出版社，2007：63?84.

［9］李輝.船舶波浪載荷的三維水彈性分析方法研究［D］.哈爾濱：哈爾濱工程大學，2009：46?66.

［10］夏錦祝，吳有生.流固耦合問題的一個一般交界面條件［J］.艦船性能研究，1993（2）：15?22.

XIA Jinzhu，WU Yousheng.One common interface condi?tion of fluid solid interaction problem［J］.Ship Perform?ance Research，1993（2）：15?22.

［11］李輝，任慧龍，劉日明，等.FPSO運動與波浪載荷三維水彈性分析方法研究［J］.哈爾濱工程大學報，2008，29（8）：789?793.

LI Hui，REN Huilong，LIU Riming，et al.3?D hydroelas?tic method for analysis of motion and wave loads of an FP?SO ship［J］.Journal of Harbin Engineering University，2008，29（8）：789?793.

［12］董艷秋.船舶波浪外荷和水彈性［M］.天津：天津大學出版社，1991：104?141.

Numerical prediction of three dimensional deformation response of the ship hull

REN Huilong，YU Pengyao，LI Hui，FENG Guoqing
（College of Shipbuilding Engineering，Harbin Engineering University，Harbin 150001，China）

Aiming at the deformation problem of the navigating ship，the three dimensional（3D）hydroelastic meth?od is adopted to predict the dynamic deformation of ships in waves.Through fitting the nodes displacement of differ?ent vibration modes，the modal information of the 3D finite model is applied to the generalized interface boundary condition of fluid?structure interactions.Then，generalized hydroelastic equations of ship motion are established and the principal coordinates of different order vibration modes are solved.Furthermore，the angular deformation re?sponse could be achieved through multiplying the angular deformation matrix and the principal coordinates of differ?ent order elastic vibration modes.Take a real ship as an example，the 3D hydroelastic method and the 3D rigid method are applied to calculate the deformation response of the ship，respectively.Through comparing the results from the two methods，the correctness and high efficiency of the 3D hydroelastic method are proved.

ship deformation；fluid and structural interactions；3D hydroelastic method；modal analysis；fluid?structure interactions；angular deformation matrix；deformation response；numerical prediction

10.3969／j.issn.1006?7043.201311053

U661.3

A

1006?7043（2015）01?0134?05

http：／／www.cnki.net／kcms／doi／10.3969／j.issn.1006?7043.201311053.html

2013?11?15.網絡出版時間：2014?11?07.

國家自然科學基金資助項目（51079034）.

任慧龍（1965?），男，教授，博士生導師.

任慧龍，E?mail：renhuilong＠263.net.