基于動網格的兩棲車航行姿態數值模擬

趙彬，張敏弟，劇冬梅

（1.北京理工大學機械與車輛學院，北京100081；2.中國兵器科學研究院，北京100089）

基于動網格的兩棲車航行姿態數值模擬

趙彬1，張敏弟1，劇冬梅2

（1.北京理工大學機械與車輛學院，北京100081；2.中國兵器科學研究院，北京100089）

針對兩棲車水上航行姿態,構建描述兩棲車水上運動的動力學模型,采用混合耦合算法和動網格技術研究靜水直航狀態下,兩棲車航行姿態的變化規律,模擬結果與實驗值吻合較好。研究結果表明：車體達到穩定航行姿態會經歷大幅振蕩調整階段和平穩運動階段,車體平衡后,升沉和縱傾仍會有小幅變化；兩棲車由排水航行狀態逐步增速到滑行狀態的過程中,車體重心逐漸升高,動升力在支撐車重成分中所占的比例越來越大,靜浮力則越來越小。

兵器科學與技術；兩棲車；航行姿態；動網格技術；數值模擬

0　引言

兩棲車在水上實際航行過程中伴隨著升沉、搖擺等復雜運動，其水上航行特性對其安全性能的保證具有重大作用［1-4］。對于兩棲車水上動態性能的數值計算，主要有兩種計算方法：第一種是靜態網格技術，第二種是動態網格技術。

關于靜態網格技術，近年來學者已經做了大量的工作。韓占忠等［5］對兩棲車的黏性阻力進行數值計算，并對車首防浪板的設計進行了計算分析；王濤等［6］分析兩棲車各阻力成分隨速度變化的規律，為分析兩棲車水上快速性提供了依據；吳珂等［7］運用計算流體力學軟件Fluent，對兩棲車車輪收起前后的阻力進行了仿真分析。靜態網格技術的航行姿態調整方法，主要是手動調整方法。對兩棲車航行姿態的數值模擬主要是采用車體相對于流場靜止不動，賦予來流一定的速度，從而得出車體的繞流場情況。根據車體受力調整航行姿態，直到航行姿態滿足要求。這種模擬方法的不足是：工作量大，車體靜止不動，沒有考慮浮態的變化。

由于靜態網格技術的局限性，動態網格技術應運而生。兩棲車在水中的運動是一動態過程，采用動網格模型可以有效解決這一問題。相比于靜態網格技術，兩棲車動態網格技術的研究國內還較少。徐國英等［8］通過網格分域方法，對兩棲車進行了阻力和浮態的數值模擬；孫偉等［9］通過應用動網格技術和VOF模型，計算了兩棲車在一定速度的航態和水動力特性；萬曉偉等［10］提出可以通過任意拉格朗日-歐拉方法來解決流-固動力耦合問題，從而獲取車輛最終航態。動網格技術主要是對車體實施主動或被動運動，從而實現車體的動態運動過程，并在過程中監測車體姿態的變化，該方法可充分利用計算機的資源，根據每一個時間段的受力狀態，基于軟件對車體的姿態進行自動調整，但網格的精度不易于保證。

兩棲車在靜水中航行時，由于受到外力的作用，會產生6個自由度的運動，這幾種運動之間存在非線性耦合作用。針對兩棲車水上航行姿態，本文構建了兩棲車水上航行數學模型，采用動網格技術研究了靜水直航狀態下，兩棲車航行姿態的變化規律。

1　數學模型和數值模擬方法

1.1研究對象

本研究采用簡化的兩棲車車體為研究對象，車質量為6.5 t，在水中行駛時車輪處于提升狀態，簡化三維模型如圖1所示。

圖1　簡化兩棲車模型Fig.1 Modified model of amphibious vehicle

1.2湍流模型

兩棲車水上航行屬于典型的流-固耦合運動，不僅伴隨水和空氣兩種黏性流體的相互作用，還伴隨著水和車體、空氣和車體之間的相互作用。車體在水中運動時，由于水和空氣兩種流體是有黏性的，所以車體受力滿足黏性流體力學的基本方程。

湍流模型選取標準k-ε模型［11］，k和ε的輸運方程為

式中：t為時間；uj和xj分別為速度分量和坐標分量；ρm和μ分別為密度和分子黏性系數；μt為湍流黏性系數，μt=Cμρmk2/ε，其中Cμ為經驗常數，Cμ=0.09；σk和σε分別為k和ε的湍流普朗特數，σk=1.0，σε=1.3；Cε1和Cε2為ε方程常數，Cε1=1.44，Cε2=1.92；Pt為平均速度梯度引起的湍動能k產生項，

1.3兩棲車動力學模型

滿足兩棲車升沉和縱傾兩個自由度的三維模型見圖2.圖中：L表示車體y軸所受升力，y軸正向為正；M表示車體z軸所受轉動力矩，正負方向遵守右手法則。

圖2　滿足h和θ兩個自由度的兩棲車三維模型Fig.2 A 3-DOF model representing h and θ degrees of freedom of an amphibious vehicle

升沉和縱搖的運動幅度通過重心升沉h和縱傾角θ來定義。整個系統的控制方程［12］為

式中：ms為車體的質量矩陣；Cs為車體的阻尼矩陣；Ks為車體的剛度矩陣；｛Mex｝為外部激勵矩陣；為在結構求解器計算出網格變形后通過使用計算流體力學（CFD）求解器計算出來的非線性黏性流體表面力；｛x｝為車體的變形量矩陣。

1.4流-固耦合混合算法

流-固耦合問題研究的是流體和固體結構之間的相互作用，包括固體結構在流體載荷作用下產生的變形或運動及其對流場的影響，因此，求解流-固耦合問題，需要同時考慮流場和結構場的求解及其耦合［13-15］。一般用到4種數值解法［16-17］，分別為完全耦合（FC）算法、松式（LC）算法、緊式（TC）算法、混合耦合（HC）算法。FC算法是完全耦合算法，對流體和結構建立統一的耦合方程，在一個時間步內對流體域和固體域中所有的未知量進行求解，求解精度較高，但在工程應用上求解范圍較窄；LC算法和TC算法將模型的結構和流體部分獨立進行求解，通過交互程序實現結構求解器和流體求解器之間的數據交換，但求解誤差較大；HC算法在（3）式左右兩邊增加附加質量力項，適當預測了流體對運動的慣性、阻尼以及剛度，減小了對物體運動位移的過預測，從而提高了數值算法穩定性。相比之下LC算法和TC算法需要更小的時間步長或者更長的計算時間，甚至計算發散。鑒于此，本次計算采用的是HC數值算法。HC算法離散后的方程為

式中：n為時間步增量；i為每一個時間步的子迭代步數；｛Mfluid｝是附加質量力項。｛Mfluid｝［16］假設為

式中：2b為兩棲車的特征長度，2b=5.04 m；ρf為流體密度。

（5）式的求解采用2階精度的半隱式Crank-Nicholson方法。

本次計算中，忽略車體本身的阻尼項，同時給出車體的質量矩陣和剛度矩陣，根據每一步計算得到的車體升力L和作用力矩M，得到車體升沉和縱傾的變化規律。

混合耦合算法的流程圖如圖3所示。其中：δ是一個介于0～1之間的調節系數，δ=0.5；ε是對于子迭代的容差值，ε=10-5；上標i表示相同時間步的子迭代步數；imax是子迭代的最大步數；下標n表示不同的時間步；x0和分別是初始變形量和變形速度。計算時，將imax值設為1在第一個子迭代步（i=1），附加質量力是從上一個時間步的最后一個子迭代步的變形算起的：

混合耦合算法具體步驟如下：

首先，從t=0，子迭代步為1開始計算，假設初始變形量x0，初始變形速度x·0，通過CFD軟件計算出，并將其代入（5）式中進行迭代求解，解得下一子迭代步的變形量，變形速度以及變形加速度；然后，將相鄰子迭代步求解的變形量進行容差比較，若滿足要求，則輸出變形結果，進行下一時間步的迭代計算，若不滿足，則重新進行（5）式的子迭代計算得到，直到滿足容差要求為止，得到，并輸出結果；最后，當計算達到預設定總時間時停止計算，輸出模型變形結果。

圖3　混合耦合迭代算法具體流程Fig.3 Flow chart of hybrid algorithm

2　網格劃分和邊界條件設置

2.1網格劃分及網格無關性驗證

網格劃分采用全結構化網格，網格總數約為100萬。在車體周圍進行網格加密，并設置邊界層，網格的劃分見圖4.

圖4　車體附近網格劃分Fig.4 Mesh near amphibious vehicle

網格無關性從兩棲車輛阻力系數Cd進行分析，Cd定義為

式中：R為車身所受阻力；ρ為水密度；v為車體航行速度；V為排水體積。通過計算，獲得網格數量n與阻力特性關系如圖5所示。可以看出，使用阻力系數作為網格無關性評判標準時，網格數量在100萬以后，阻力特性變化不大，鑒于此，文中所使用的數值模擬網格約100萬。

圖5　網格無關性驗證Fig.5 Grid-independent validation

2.2邊界條件設置

如圖6所示，計算區域以車體為中心向前延伸了1.5個車長，向后延伸了3個車長，向下延伸了10倍于車體靜吃水深。邊界條件：入口邊界給定為速度入口；出口和周邊外域給定為開放邊界；車體給定為移動固壁邊界。

圖6　邊界條件設置Fig.6 Boundary conditions

3　數值計算結果及分析

3.1相同弗勞德數Fr車體航行姿態分析

以車體Fr=0.93的計算過程為例進行分析。如圖2，車體初始姿態為0.2 m的重心升沉，8°的縱傾角。不同時刻，升沉和縱傾的變化曲線見圖7和圖8.可知車體在運動響應計算過程中，經過一段時間的大幅振蕩調整，找到了受力平衡點，達到動平衡狀態。車體結束大幅振蕩調整大概需要3.5 s的時間，初始大幅振蕩調整階段過后逐步微調進入穩定航行姿態，經過初期短暫振蕩調整，車體后續航行狀態基本穩定。

圖7　車體升沉變化圖Fig.7 Time-history of the predicted center-of-gravity

圖8　車體縱傾變化圖Fig.8 Time-history of the predicted pitch angle

為了說明車體航行姿態變化的原因，與此同時，對車體的單位排水量總阻力CD、單位排水量總升力CL和縱傾力矩系數CN也進行監測，從而確定車體達到平衡的時間和平衡時的具體狀態。

3個無量綱系數的定義為

式中：D為靜水時車體排水重力；N為車體所受實際繞z軸轉矩。

3個系數的變化情況見圖9.為了下文表述方便，將圖中單位排水量總阻力CD、單位排水量總升力CL和縱傾力矩系數CN采用6個時間點進行劃分，分別為t1～t6.

對比圖7、圖8和圖9曲線可知，車體的運動分為兩個階段：第1個階段（t1～t4時刻）為大幅振蕩調整階段，在此階段中車身的運動幅度較大，車體的升沉變化量最大達到了0.25 m，縱傾變化量最大達到了0.13 rad，動態響應特性明顯；第2個階段（t4～t6時刻）為車體的平穩運動階段，在此階段中車身的運動幅度較小，升沉變化量最大只有0.08 m，縱傾變化量最大只有0.05 rad，動態響應特性減弱。以初始姿態運動時，在運動初期，航行中車體所受浮力大于靜止時所受浮力，故車體重心會提升直到浮力和重力再度平衡；設計姿態航行時，縱傾力矩并非為0，故車體縱傾角會減小直到達成新的力矩平衡。經過這兩個階段之后，車體已進入平穩運動階段。在平穩運動階段中，車體的航行姿態也會略有改變，繼而總升力和縱傾力矩也會改變，但變化幅度已不大。由此可以得出車體平衡實際是動態平衡的概念，從而彌補了靜態網格技術車體姿態靜止不動的不足。

為了進一步分析車體運動規律變化原因，表1中列出了典型時刻車體仰視圖的壓力變化云圖和縱剖面速度矢量圖，結合圖9縱傾力矩系數的變化，選定這幾個時刻分別為力矩最大值時刻t2，力矩最小值時刻t3以及平穩運動階段中的時刻t5.

圖9　不同時刻車體單位排水量總阻力、單位排水量總升力和力矩系數變化趨勢圖Fig.9 Resistance，lift and torque coefficients at different moments

表1　典型時刻車體壓力云圖及速度矢量圖Tab.1 Contours of static pressures and velocity vectors at different moments

大幅振蕩調整時期，縱傾力矩系數的兩個峰值點，最大峰值t2及最小峰值t3時刻，對比其壓力分布，可發現各處的壓力分布與其水上興波密切相關。在t2時刻的仰視圖中，后底部的高壓區域對車體產生了一個繞z軸正的轉矩，仰視圖中其余位置的壓力分布幾乎沿縱剖面對稱；而t3時刻的仰視圖中，首前部產生壓力較高的區域，將產生一個繞z軸負的轉矩。但是在平穩運動時期的時刻t5，車體表面的壓力基本趨于穩定，不會導致大的波動出現。

觀察縱剖面的速度矢量圖，在t2和t3時刻，車體首部、尾部附近有大量旋渦存在，在車尾處水流有明顯向上移動的趨勢。旋渦的出現會增加渦流損失，增大兩棲車在水中航行時的阻力，并且在旋渦區域內部會形成低壓區，增大車體的黏壓阻力，進而改變車體表面的壓力分布情況，從而影響車體航行姿態的調整；而當車處于平穩運動t5時刻，車體首部、尾部附近只有少量興波存在，不再有大尺度的旋渦出現，車體姿態趨于穩定。

3.2不同Fr數車體航行姿態對比

不同Fr下車體穩定后最終航行繞流形態圖見表2。由表2車體繞流形態圖得知，隨著Fr數的增大車體重心位置逐漸升高，車體逐漸抬起。Fr= 0.93時，車體處于排水航行狀態，靜浮力是支持車重的主要成分，車體大部分在水面以下，車首興波較小，尾部流場沒有明顯的“雞尾流”形態［18］。Fr= 1.40時車體穩定后重心位置高于Fr=0.93工況，但仍有部分車體在水面以下，即車的重量由靜浮力和動升力共同支持。這表明隨著航速的提高，動升力的比重將越來越大，尾部流場開始出現“雞尾流”形態。Fr=1.87時，大部分車體從水中抬起，車體基本處于在水面滑行的狀態，車體首部興波上卷較大，兩側興波向后噴濺，尾部出現典型的“雞尾流”形態，動升力為支持車重的主要成分。

3.3數值模擬結果與實驗結果對比

實驗水池長160 m，寬7 m，水深3.7 m；拖車為空腹梁結構，全數字式直流調速系統，速度范圍0.01～8 m/s，速度精度為±1 mm/s；船模采用鋼板制成，內外涂防腐漆。車模拖曳方式示意圖如圖10所示。

表2　不同Fr下車體穩定航行繞流形態圖Tab.2 Water volume fractions round vehicle at different speeds

圖10　車輛拖曳方式圖Fig.10 Amphibious vehicle model in the towing tank

如圖10所示，拖點位置為推力軸線和重心縱向位置線的交點；阻力儀用于記錄車體航行阻力；首尾升沉記錄裝置用于記錄車體航行時首尾升沉位移，并可由此計算出車體航行縱傾角；首尾導航板用于防止車體偏航。具體實驗方法執行我國船舶行業標準CB/Z 244—1988滑行艇船模阻力測試方法

將自動調整方法模擬值與拖模實驗值進行對比（車體最終穩定航行繞流形態和航行姿態），見表3和圖11.由表3和圖11可知，模擬結果與實驗結果吻合較好，但是受網格尺度影響，無法描述車首周圍向外拋出的噴濺水體與空氣混合后的流動結構。

表3　數值模擬與實驗的繞流形態圖對比Tab.3 Comparison of simulation and experimental results for sailing

圖11　升沉和縱傾角的模擬值與實驗值對比Fig.11 Comparison of simulation and experimental results of the predicted center-of-gravity and the predicted pitch angle

4　結論

本文通過應用動網格技術，對兩棲車靜水直航的航行姿態進行了數值模擬，得出以下結論：

1）結合HC算法的動態網格技術可以有效求解兩棲車的航行姿態，數值計算結果與實驗結果的吻合程度較好。

2）不同Fr數下車體穩定后的航行姿態有顯著差異，具體表現在隨著Fr數的增大車體重心位置逐漸升高，動升力在支持車重成分中所占的比例越來越大，靜浮力則越來越小。

（

）

［1］王濤，徐國英，郭齊勝.兩棲車輛水上動態性能數值模擬方法及其應用［M］.北京：國防工業出版社，2009：1-7. WANG Tao，XU Guo-ying，GUO Qi-sheng.Amphibious vehicle numerical simulation of water［M］.Beijing：National Defense Industry Press，2009：1-7.（in Chinese）

［2］萬曉偉，姚新民，王濤.基于CFD的兩棲車輛水上阻力預測［J］.裝甲兵工程學院學報，2013，27（3）：26-30. WAN Xiao-wei，YAO Xin-min，WANG Tao.Water drag prediction of amphibious vehicle based on CFD［J］.Journal of Academy of Armored Force Engineering，2013，27（3）：26-30.（in Chinese）

［3］王濤，郭齊勝，徐國英，等.基于CFD的兩棲車輛水上耐波性能研究初探［J］.系統仿真學報，2009，21（10）：3142-3145. WANG Tao，GUO Qi-sheng，XU Guo-ying，et al.Preliminary study on seakeeping capacity of sailing amphibious vehicle based on CFD［J］.Journal of System Simulation，2009，21（10）：3142-3145.（in Chinese）

［4］袁益民，郭齊勝.兩棲車輛水上操縱運動建模與仿真［J］.計算機仿真，2011，28（7）：336-339. YUAN Yi-min，GUO Qi-sheng.Modeling and simulation of maneuvering behavior of sailing amphibious vehicle［J］.Computer Simulation，2011，28（7）：336-339.（in Chinese）

［5］韓占忠，王國玉，閆為革.兩棲車輛航行黏性阻力數值分析［J］.車輛與動力技術，2003（2）：6-10. HAN Zhan-zhong，WANG Guo-yu，YAN Wei-ge.Numerical simulation of viscosity resistance around a running amphibian vehicle［J］.Vehicle&Power Technology，2003（2）：6-10.（in Chinese）

［6］王濤，徐國英，姚新民，等.兩棲車兩相繞流場的模擬與水上快速性分析［J］.機械工程學報，2009，44（12）：168-172. WANG Tao，XU Guo-ying，YAO Xin-min，et al.Numerical simulation of two phases flow field around amphibious vehicle and analysis of power performance on water［J］.Journal of Mechanical Engineering，2009，44（12）：168-172.（in Chinese）

［7］吳珂，宋桂霞，趙又群.基于流體仿真的兩棲車車輪收起前后阻力對比分析［J］.系統仿真技術，2007，3（4）：187-191.WU Ke，SONG Gui-xia，ZHAO You-qun.Comparative analysis on resistance of amphibious vehicle before and after retracting wheels based on fluid simulation［J］.System Simulation Technology，2007，3（4）：187-191.（in Chinese）

［8］徐國英，王俊，周景濤.基于CFD的兩棲車輛阻力和浮態數值模擬［J］.艦船科學技術，2006，28（4）：22-25. XU Guo-ying，WANG Jun，ZHOU Jing-tao.Numerical simulation of the amphibious vehicle's drag force and attitude based on CFD［J］.Ship Science and Technology，2006，28（4）：22-25.（in Chinese）

［9］孫偉，韓建禮，劉西俠.基于動網格模型的兩棲車輛數值模擬［J］.艦船科學技術，2009，31（1）：146-150. SUN Wei，HAN Jian-li，LIU Xi-xia.Numerical simulation of amphibious vehicle based on the dynamic-mesh model［J］.Ship Science and Technology，2009，31（1）：146-150.（in Chinese）

［10］萬曉偉，王濤，姚新民.兩棲車輛水上高速航態的數值仿真方法研究［J］.計算機仿真，2012，29（6）：323-327. WAN Xiao-wei，WANG Tao，YAO Xin-min.Study on numerical simulation method for high speed navigation state of amphibious vehicle［J］.Computer Simulation，2012，29（6）：323-327.（in Chinese）

［11］Launder B E，Spalding D B.The numerical computation of turbulent flows［J］.Computer Methods in Applied Mechanics and Engineering，1974，3（2）：269-289.

［12］葉正寅，張偉偉，史愛明，等.流固耦合力學基礎及其應用［M］.哈爾濱：哈爾濱工業大學出版社，2010：37-51. YE Zheng-yin，ZHANG Wei-wei，SHI Ai-ming，et al.Fundamentals of fluid-structure coupling and its application［M］.Harbin：Harbin Institute of Technology Press，2010：37-51.（in Chinese）

［13］Dowell E H，Hall K C.Modeling of fluid-structure interaction［J］.Annual Review of Fluid Mechanics，2001，33：445-490.

［14］Munch C，Ausoni P，Braun O，et al.Fluid-structure coupling for an oscillating hydrofoil［J］.Journal of Fluids and Structures，2010，26（6）：1018-1033.

［15］Akcabay D T，Young Y L.Hydroelastic response and energy harvesting potential of flexible piezoelectric beams in viscous flow［J］.Physics of Fluids（1994-present），2012，24（5）：054106.

［16］Young Y L，Chae E J，Akcabay D T.Hybrid algorithm for modeling of fluid-structure interaction in incompressible，viscous flows［J］.Acta Mechanica Sinica，2012，28（4）：1030-1041.

［17］Young Y L，Chae E J，Akcabay D T.Transient hydroelasticreponse of a flexible hydrofoil in subcavitating and cavitaing flows［C］//29th Symposium on Naval Hydrodynamics.Gothenburg，Sweden：the US Office of Naval Research，2012：26-31.

［18］楊爽，王志東，吳賀賀，等.滑行艇縱向運動響應與靜水航行數值模擬分析［J］.船舶工程，2012，34（4）：31-35. YANG Shuang，WANG Zhi-dong，WU He-he，et al.Numerical simulation and analysis on heave-pitch coupled motion response and motion in calm water of planing ship［J］.Ship Engineering，2012，34（4）：31-35.（in Chinese）

Numerical Simulation of Navigating Pose for Amphibious Vehicle based on Dynamic-mesh Model

ZHAO Bin1，ZHANG Min-di1，JU Dong-mei2
（1.School of Mechanical Engineering，Beijing Institute of Technology，Beijing 100081，China；2.Ordnance Science and Research Academy of China，Beijing 100089，China）

A dynamics model for the navigating pose of amphibious vehicle in water is established.The dynamic mesh technology and the hybrid coupled algorithm are used to study the change rule of navigating pose of amphibious vehicle in sailing，and the navigating pose of amphibious vehicles in still water are discussed.The numerical results are consistent with the experimental results.The result shows that the vehicle undergoes a sharp motion and a smooth motion during movement before achieving the ultimate navigating pose.After achieving the ultimate navigating pose，the center-of-gravity and the pitch angle slightly change.During the process of sailing from drainage state to coasting state，the body center of gravity gradually increases；the dynamic buoyancy to support the vehicle weight increases，and the static buoyancy decreases.

ordnance science and technology；amphibious vehicle；navigating pose；dynamic mesh technology；numerical simulation

U674.78

A

1000-1093（2015）03-0412-09

10.3969/j.issn.1000-1093.2015.03.005

2014-06-04

國家自然科學基金項目（51106009）

趙彬（1989—），男，碩士研究生。E-mail：qzb1989q@163.com；張敏弟（1971—），女，副教授，碩士生導師。E-mail：zmd1971@bit.edu.cn