潛艇旋臂回轉試驗數值模擬

肖昌潤，劉瑞杰，許可，劉洋，徐亞運

(1.海軍工程大學艦船工程系，湖北武漢430033)(2.中國人民解放軍91039部隊，北京102401)

(3.中國人民解放軍92339部隊，廣東湛江524000)

潛艇是海軍重要的作戰(zhàn)武器，對于維護國家海權、協同水面艦船作戰(zhàn)、進行戰(zhàn)略威懾具有重要意義.獲取潛艇水動力系數是評價其操縱性的重要途徑.以節(jié)省經費和縮短研究周期為目的，計算流體力學(computational fluid dynamics，CFD)預報潛艇操縱性越來越受到重視，且已經成為重要學科.隨著應用和研究的深入，用計算流體力學方法獲得的水動力系數越來越準確，計算精度越來越受到認同.

文獻［1-5］中對船舶定常計算進行了研究，分析了網格和湍流模型對計算精度的影響，為模擬回轉運動起到了指導作用.國內通過數值計算方法模擬回轉運動的文獻較少，其中文獻［6］中通過MRF方法對回轉運動進行了數值模擬，分析了回轉運動的原理，并進行了初步研究.國外對于回轉運動的研究已經比較成熟，文獻［7］中通過變形方法預報流動分離情況，為研究回轉運動提供了新的思路.文獻［8-10］中通過自編程序對回轉運動進行了研究，分析了湍流模式、網格結構、邊界層等對計算精度的影響.文獻［11］中通過穩(wěn)流NS方程和旋轉坐標系的方法計算了最大漂角達到18°的回轉運動，計算結果誤差在11%以內.

基于6自由度的潛艇操縱性方程涉及到固定于地球的固定坐標系和固聯于潛艇的運動坐標系.通過CFD對潛艇非定常運動進行計算時需采用動網格.動網格參考坐標系可以是固定坐標系也可以是固聯于潛艇的運動坐標系.文中的Mesh Motion方法采用地球坐標系，而基于源項的方法采用位于潛艇的運動坐標系.文中以潛艇SUBOFF模型為研究對象，通過流體計算軟件FLUENT，對潛艇回轉運動進行了研究，分析了兩種方法對潛艇水動力計算精度的影響，結果表明計算精度滿足工程預報的要求.

1 基本理論與控制方程

1.1 控制方程

基于雷諾時均方法，用張量形式表述的控制方程［12］為:

式中:ui為速度;t為時間;ρ ui′uj′為雷諾應力項;ρ為水的密度.

1.2 Mesh Motion方法模擬旋臂回轉運動

當采用地球坐標系時，網格整體隨潛艇做回轉運動.Fluent軟件的Mesh Motion功能可以幫助實現這個過程.回轉運動是非定常的，通過這種方式進行迭代時，網格不會發(fā)生任何變形，避免了計算時網格發(fā)生拉伸、重構，使網格質量變差，甚至產生負網格，影響計算精度.

基于6自由度潛艇操縱性方程，關于無因次化的橫向力Y′和偏航力矩N′，有方程組:

式中:r為角速度;Y′r，N′r，N′r|r|為旋轉水動力導數.

1.3 基于源項法的潛艇旋臂回轉試驗

由于潛艇做回轉運動，在運動坐標系下的速度可以分解為線速度V和角速度Ω.此時，對于任意流體單元，其在地球坐標系下的速度可以分解為U=Ur+Ue，Ur為相對速度，Ue為牽連速度.Ue可以表示為Ue=V+Ω×r′，r′為相對于運動坐標系的坐標向量.則任意流體單元在固定坐標系下的絕對加速度為a=ar+ae+ac.

又Ue是關于時間的函數，故2U≡0，則上式可以簡化為:

將 -ρae-ρac定義為源項，用MS表示.

通過引入MS源項，即可將運動坐標系下非定常的旋臂回轉運動等效為地球坐標系下的定常運動.

令Ur= ［u，v，w］T，V= ［u0，v0，w0］T，Ω = ［p，q，r］T，r′=［x，y，z］T.其中u，v，w分別是縱向速度、橫向速度、垂向速度;p，q，r分別是橫傾角速度、縱傾角速度、偏航角速度;u0，v0，w0表示初始速度.潛艇繞Z軸做回轉運動，則u0=v0=w0=0，u=v=0=0，p=q=0=0.消去這些項可得源項的分量形式:

通過相似的步驟也可以得到繞X，Y軸做回轉運動的源項公式.將式(8)編譯，借助UDF的DEFINE-SOURCE，函數代入Fluent即可對回轉運動進行模擬.

2 計算區(qū)域網格劃分

文中兩種方法使用的網格及邊界條件的設置是相同的.計算區(qū)域是一個旋轉中心距離內壁面1.5L(L=4.356 m是模型的長度)的半環(huán)形.為了降低近壁面網格Y+值，提高網格質量，同時減少網格數量，采用分塊網格劃分技術.使用RNG k-ε湍流模型計算時，網格總數為450萬，Y+值在30左右，并選用標準壁面函數.Mesh Motion方法中，入口及四周壁面均設為速度入口，整個計算區(qū)域設定為旋轉域，潛艇壁面相對旋轉域靜止.旋轉域的旋轉角速度 r分別為 0.08，0.1，0.15，0.2 rad/s.源項法中由于入口處速度是相對速度，且隨半徑變化，Uinlet=-Ω ×r′，速度，故需要用到DEFINE-PROFILE函數，同時將DEFINESOURCE函數通過UDF接口添加到流場中，實現回轉運動模擬.計算區(qū)域網格如圖1.

圖1 計算網格示意Fig.1 Schematic diagram of grids

3 計算結果校驗與分析

3.1 Mesh Motion方法計算結果分析

由于網格繞旋轉中心旋轉，故力矩中心坐標也是旋轉的，通過已知點坐標，換算可得到力矩中心坐標.仿真得到力Y和力矩N后，通過軟件Matlab擬合三次樣條曲線，得到相關水動力系數.力和力矩及其旋轉導數的無因次化表達式分別為:.式中:ρ為流體的密度;U為潛艇旋轉的線速度;L為艇長;Y為回轉試驗中沿Z方向的力;N為繞Y軸的力矩.

圖2 擬合曲線Fig.2 Fitted curves

由圖2的擬合曲線可得旋轉導數計算結果，并與試驗結果進行對比，如表1所示，誤差滿足工程計算需要.圖3顯示不同時刻流場的速度場分布示意圖.

表1力和力矩計算結果與試驗結果對比Table 1 Force and torque calculation results compare with the test results

圖3 基于Mesh Motion方法不同時刻速度場分布及艇體周圍流場細節(jié)Fig.3 Velocity distribution at different time and flow field detail around the submarine based on Mesh Motion

3.2基于源項法的潛艇旋臂回轉試驗計算

Mesh Motion功能計算是一個非定常的過程，其對計算機的性能要求較高，耗費大量時間.基于源項法的方式可以有效擬補Mesh Motion方法的不足.基于源項法的潛艇旋臂回轉模擬是通過Fluent軟件的源項功能，加入UDF自編程序實現的.計算是定常的，極大減少了工作量.

對比RNG k-ε湍流模型和SST k-ω湍流模型的計算精度和計算量.使用SST k-ω湍流模型時，由于SST k-ω要求Y+＜1，故網格數大增，網格數目為930萬.計算結果如表2所示.

表2 計算結果與試驗結果比較Tab.2 Comparison between calculation results and test results

從計算時間上看，基于源項法的回轉運動計算效率明顯高于Mesh Motion方法.由于網格數目大增，SST k-ω湍流模型大大增加了計算量，但是并沒有明顯提高精度.

4 結論

文中通過Mesh Motion方法和基于源項法旋臂試驗對全附體的SUBOFF潛艇回轉運動進行了模擬.對比了RNG k-ε湍流模型和SST k-ω湍流模型的計算精度，并與試驗結果進行對比.結果表明:

1)Mesh Motion方法和基于源項法旋臂試驗的計算精度滿足工程需求.但基于源項法懸臂試驗計算量更小，計算精度更高.

2)SST k-ω湍流模型對于網格Y+值要求較高，造成網格數增加，與RNG k-ε湍流模型對比，計算耗時更長，但并沒有明顯提高計算精度.

References)

［1］ 劉志華，熊鷹，韓寶玉.潛艇流場數值計算網格與湍流模型選?。跩］.華中科技大學:自然科學版，2009，39(7):98-101.Liu Zhihua，Xiong Ying，Han Baoyu.Computational grid and turbulent model for calculating submarine viscous flow field［J］.Journal of Huazhong University of Science and Technology:Natural Science Edition，2009，39(7):98-101.(in Chinese)

［2］ 肖昌潤，劉巨斌，朱建華.DARPA2潛艇模型定常繞流水動力數值計算［J］.華中科技大學:自然科學版，2007，35(8):115-118.Xiao Changrun，Liu Jubin，Zhu Jianhua.Numerical computation of hydrodynamic force of DARPA2 submarine model［J］.Journal of Huazhong University of Science and Technology:Natural Science Edition，2007，35(8):115-118.(in Chinese)

［3］ Xing T，Bhushan S，Stern F.Vortical and turbulent structures for KVLCC2 at drift angle 0，12 and 30 degrees［J］.Ocean Engineering，2012，55:23-43.

［4］ 鄭小龍，黃勝，尚秀敏.基于CFD的船舶阻力預報方法研究［J］.江蘇科技大學學報:自然科學版，2014，28(2):109-113.Zheng Xiaolong，Huang Sheng，Shang Xiumin.Study of ship resistance prediction method based on CFD［J］.Journal of Jiangsu University of Science and Technology:Natural Science Edition，2014，28(2):109-113.(in Chinese)

［5］ 陳淑玲，楊松林，劉智.基于Fluent的五體船靜水中水動力特性數值模擬［J］.江蘇科技大學學報:自然科學版，2012，26(6):541-545.Chen Shuling，Yang Songlin，Liu Zhi.Numerical simulation of hydrodynamic performance of pentamaran in calm water based on Fluent［J］.Journal of Jiangsu University of Science and Technology:Natural Science Edition，2012，26(6):541-545.(in Chinese)

［6］ 劉帥，葛彤，趙敏.基于源項法的潛艇旋臂試驗模擬［J］.大連海事大學學報，2011，37(2):1-4.Liu Shuai，Ge tong，Zhao Min.Simulation for submarine rotating-arm test based on added momentum source method［J］.Journal of Dalian Maritime University:Natural Science Edition，2011，37(2):1-4.(in Chinese)

［7］ Zhang J T，Maxell J A，Gerber A G，et al.Simulation of the flow over axisymmetric submarine hulls in steady turning［J］.Ocean Engineering，2013，57:180-196.

［8］ Toxopeus S，Atsavapranee P，Wolf E，et al.Collaborative CFD exercise for a submarine in a steady turn［C］∥International Conference on Ocean，Offshore and Arctic Engineering.Rio de Janeiro:ASME，2012:761-772.

［9］ Gregory P A，Joubert P N，Chong M S.Flow over a body of revolution in a steady turn［R］.Rockingham:Defence Science and Technology Organisation Victoria Platform Sience Lab，2004.

［10］ Phillips A B，Turnock S R，Furlong M.The use of computational fluid dynamic to aid cost-effective hydrodynamic design of autonomous underwater vehicles［J］.Journal of Engineering for the Maritime Environment，224(4):239-254.

［11］ Marshallsay P G，Eriksson A M.Use of computational fluid dynamics as a tool to assess the hydrodynamic performance of a submarine［C］∥Australasian Fluid Mechanics Conference.Launceston:［s.n.］2012.

［12］ SAS IP.Ansys fluent user′s guide，version 14.5.2013.