孫銘澤 王永生 張志宏 戴余良

（海軍工程大學(xué)船舶與動力學(xué)院 武漢 430033）

0 引 言

潛艇水動力系數(shù)的數(shù)值計算主要是基于勢流理論.隨著計算流體力學(xué)和計算機(jī)技術(shù)的飛速發(fā)展，目前基于雷諾時均方法（RANS）的數(shù)值模擬成為操縱性試驗(yàn)數(shù)值計算的熱點(diǎn).潛艇操縱性物理試驗(yàn)主要有直線拖曳試驗(yàn)（ORT）、回轉(zhuǎn)臂試驗(yàn)（RAT）以及平面運(yùn)動機(jī)構(gòu)試驗(yàn)（PMM）.以往操縱性的數(shù)值模擬主要集中于前2種試驗(yàn).文獻(xiàn)［1］系統(tǒng)闡述了基于粘性流的CFD方法在潛體操縱性預(yù)報中的重要作用，并以SUBOFF幾何為研究對象，計算了在不同攻角下的流場，并將不同湍流模型下的計算結(jié)果與試驗(yàn)數(shù)據(jù)以及基于遲滯分離渦（DDES）模型的計算結(jié)果進(jìn)行了比較；文獻(xiàn)［2］從網(wǎng)格第一層節(jié)點(diǎn)距離、網(wǎng)格數(shù)量以及湍流模型等方面論述了操縱性計算的數(shù)值方法，并與試驗(yàn)數(shù)據(jù)進(jìn)行了對比，該文獻(xiàn)認(rèn)為，第一層網(wǎng)格節(jié)點(diǎn)距離的取值應(yīng)使無因次距離y＋介于40～100之間，而縱向力的精確預(yù)報相對于側(cè)向力或力矩需要較多的網(wǎng)格數(shù)量，同時指出了大漂角姿態(tài)下的操縱性精確預(yù)報的缺陷；文獻(xiàn)［3］應(yīng)用CFD 軟件對潛艇的直線拖曳試驗(yàn)和回轉(zhuǎn)臂試驗(yàn)進(jìn)行了大量的數(shù)值計算，得到了在線性或弱非線性前提下的粘性類水動力系數(shù)和舵角系數(shù)；文獻(xiàn)［4］采用6 種不同的湍流模式計算了SUBOFF 2種模型在變漂角和變攻角條件下的操縱性水動力，并與試驗(yàn)值進(jìn)行比較，計算誤差在20%以內(nèi)，該文獻(xiàn)認(rèn)為κ－ω SST 模型較為適合進(jìn)行潛艇操縱性數(shù)值計算.以上所述的潛艇操縱性數(shù)值計算只能獲取潛艇粘性類水動力系數(shù)，而不能計算慣性類水動力系數(shù)，所需要進(jìn)行的計算量較為龐大，限制了潛艇操縱性的快速預(yù)報.而PMM 數(shù)值試驗(yàn)卻能較好地解決這個問題.

文獻(xiàn)［5－6］系統(tǒng)地闡述了利用基于動網(wǎng)格技術(shù)的FLUENT 軟件實(shí)現(xiàn)平面運(yùn)動的方法，該方法通過編寫UDF函數(shù)控制四面體網(wǎng)格節(jié)點(diǎn)位移來實(shí)現(xiàn)網(wǎng)格運(yùn)動，并在運(yùn)動中求解瞬態(tài)雷諾時均方程從而對操縱性進(jìn)行快速預(yù)報.由于結(jié)構(gòu)化網(wǎng)格較非結(jié)構(gòu)化網(wǎng)格在網(wǎng)格正交性、計算收斂性以及對數(shù)值耗散的抑制等方面具有優(yōu)勢，基于結(jié)構(gòu)化網(wǎng)格的數(shù)值計算結(jié)果將更為可信，于是利用六面體網(wǎng)格進(jìn)行PMM 試驗(yàn)的數(shù)值模擬勢必成為主流，然而這方面的發(fā)展較少見于報道.

本文以有精細(xì)試驗(yàn)數(shù)據(jù)的DARPA－2全附體SUBOFF潛艇為分析對象，以全局映射式六面體網(wǎng)格進(jìn)行空間離散，利用基于網(wǎng)格變形技術(shù)的ANSYS CFX 軟件，實(shí)現(xiàn)了潛艇小振幅PMM 試驗(yàn)的虛擬仿真，同時也進(jìn)行了直線拖曳試驗(yàn)和回轉(zhuǎn)臂試驗(yàn)的數(shù)值計算，并將計算結(jié)果與試驗(yàn)數(shù)據(jù)進(jìn)行對比，闡述了該方法在操縱性計算領(lǐng)域的適用性，證明了該方法在計算水動力系數(shù)時的可靠性，提供了潛艇操縱性快速預(yù)報的一種新途徑.

1 基本理論與控制方程

1.1 控制方程

目前廣泛使用雷諾時均法來考察瞬時值的影響，用張量的形式表述不可壓湍流作時間平均處理的控制方程為

1.2 網(wǎng)格變形方程

CFX 軟件所提供的網(wǎng)格變形技術(shù)［7］能夠使指定的邊界或子域內(nèi)的節(jié)點(diǎn)按照CEL 語句所規(guī)定的規(guī)律進(jìn)行變形，這項(xiàng)技術(shù)主要是在求解瞬態(tài)RANS方程之前求解位移擴(kuò)散方程，該方程使邊界或子域附近的節(jié)點(diǎn)適應(yīng)性地運(yùn)動

式中：δ為相對位移；（Γdisp為網(wǎng)格剛度系數(shù)，默認(rèn)為1m2／s.

2 PMM 試驗(yàn)的數(shù)值實(shí)現(xiàn)

潛艇PMM 試驗(yàn)主要包括：純升沉運(yùn)動、純俯仰運(yùn)動、純橫蕩運(yùn)動和純搖艏運(yùn)動.利用CFX 軟件實(shí)現(xiàn)各種運(yùn)動，通過編寫CEL語句控制節(jié)點(diǎn)位移，同時控制體進(jìn)口方向模擬無限遠(yuǎn)來流，設(shè)來流速度為v＝3.05m／s.

2.1 數(shù)值計算的準(zhǔn)備

計算時為了保證網(wǎng)格變形的質(zhì)量，將潛艇壁面附近一定范圍內(nèi)的流場單獨(dú)設(shè)置為子域（subdomain），外部流場與子域通過設(shè)置交界面進(jìn)行數(shù)據(jù)傳遞，整個流場邊界范圍為：來流方向取1 L，下游方向取2 L，徑向取1 L，L 為潛艇全長，見圖1.

圖1 計算域的設(shè)置

本文的計算案例中計算域均采用全局映射式六面體網(wǎng)格進(jìn)行離散，這部分工作在專業(yè)網(wǎng)格劃分軟件ICEM CFD 中完成，見圖2.為了適應(yīng)所采用的SST 湍流模型的需要，需控制網(wǎng)格壁面第一層距離.同時，為了適應(yīng)潛艇不同姿態(tài)下的分離流動，在翼面及主艇體背流區(qū)域的網(wǎng)格進(jìn)行局部加密.另一方面，將潛艇艇體附近的網(wǎng)格與遠(yuǎn)離艇體部分的流場分離開來，并使該子域內(nèi)的所有網(wǎng)格節(jié)點(diǎn)按指定規(guī)律同步運(yùn)動以保證網(wǎng)格質(zhì)量，見圖3.網(wǎng)格節(jié)點(diǎn)總數(shù)約為156萬.

圖2 艇體表面網(wǎng)格

圖3 縱中剖面網(wǎng)格變形前后比較（純俯仰運(yùn)動）

2.2 操縱性計算

2.2.1 速度與加速度系數(shù)的求取 基于線性理論的小振幅平面運(yùn)動機(jī)構(gòu)主要是為了測定線性水動力系數(shù)，而純升沉運(yùn)動和純橫蕩運(yùn)動是為了求取速度與加速度系數(shù)，以純升沉運(yùn)動為例，由水動力系數(shù)的意義可得線性方程

設(shè)定振幅a＝300mm，艇長L＝4.356m，模型的運(yùn)動規(guī)律為

式中：ζ為潛艇的垂向位移；ω 為周期運(yùn)動圓頻率，ω＝2πf；f 為運(yùn)動頻率，為了使迭代時間步長取整數(shù)，將頻率設(shè)置為：0.2，0.25，0.3125，0.4，0.5，0.625Hz；w 和˙w 分別為垂向運(yùn)動速度和加速度.

計算時使子域內(nèi)的網(wǎng)格節(jié)點(diǎn)按正弦規(guī)律作小幅振動，將一個運(yùn)動周期分為40個時間步，通常計算6個周期便可得到穩(wěn)定的垂向力和俯仰力矩關(guān)于時間的變化曲線，見圖4.

圖4 純升沉運(yùn)動垂向力Z 和俯仰力矩M 的迭代曲線（f＝0.3125Hz）

文獻(xiàn)［8］對物理試驗(yàn)的數(shù)據(jù)處理方法做了闡述，依照分量分離的方法，將水動力的同相分量與正交分量分開.可以采用Fourier級數(shù)展開到1階的方法進(jìn)行分量分離，計算曲線見圖5，由該曲線擬合出的速度系數(shù)和加速度系數(shù)為

圖5 純升沉運(yùn)動水動力系數(shù)求取曲線

純橫蕩運(yùn)動與純升沉運(yùn)動非常相似，將潛艇模型置于橫臥狀態(tài)進(jìn)行運(yùn)動即可.在此直接給出水動力系數(shù)計算結(jié)果，見圖6，據(jù)此擬合出的關(guān)于側(cè)向速度v的速度系數(shù)和加速度系數(shù)為

圖6 純橫蕩運(yùn)動水動力系數(shù)求取曲線

2.2.2 角速度和角加速度系數(shù)的求取 純俯仰運(yùn)動和純搖艏運(yùn)動則是為了求取角速度和角加速度系數(shù).同上述分析，在小振幅運(yùn)動的線性前提下可寫出純俯仰運(yùn)動的水動力方程

設(shè)定純俯仰運(yùn)動的振幅θ0＝π／30rad，則純俯仰運(yùn)動的運(yùn)動方程為

式中：θ為潛艇模型的縱傾角.

數(shù)值試驗(yàn)中，通過編寫CEL 語句，一方面使子域內(nèi)的節(jié)點(diǎn)按正弦規(guī)律上下振動，另一方面則繞重心在垂直面內(nèi)按正弦規(guī)律轉(zhuǎn)動，那么兩種運(yùn)動的合成就是純俯仰運(yùn)動.計算時與純升沉運(yùn)動一樣，將一個周期分為40個時間步，計算6～8個周期就可以得到穩(wěn)定的垂向力和俯仰力矩關(guān)于時間的變化曲線，對該曲線按Fourier級數(shù)展開至1階，取不同頻率下的1階系數(shù)擬合出相應(yīng)的角速度和角加速度系數(shù)，圖表數(shù)據(jù)略，直接給出計算得到的水動力系數(shù)如下

2.3 計算結(jié)果校驗(yàn)與分析

美國泰勒研究中心做了SUBOFF 全附體潛艇（AFF－8幾何）的操縱性試驗(yàn)，并公布了部分水動力系數(shù)試驗(yàn)數(shù)據(jù)［9］.將本文計算得到的水動力系數(shù)與試驗(yàn)數(shù)據(jù)進(jìn)行對比［10］，見表1～2.

表1 粘性類水動力系數(shù)對比

表2 慣性類水動力系數(shù)對比

圖7 RANS模擬潛艇渦量場沿流動方向的演變（純橫蕩運(yùn)動，f＝0.2Hz）

圖8 純橫蕩運(yùn)動中潛艇處于平衡位置時渦量向背流面的偏移

文獻(xiàn)［7］指出，由純升沉和純橫蕩試驗(yàn)所測出的水動力系數(shù)中，速度系數(shù)的值較小，所得之測量值往往離散度較大，所以速度系數(shù)一般仍由直拖試驗(yàn)來確定.從表1中也可以看出，由PMM 數(shù)值試驗(yàn)求得的速度系數(shù)與直拖試驗(yàn)或回轉(zhuǎn)試驗(yàn)所得結(jié)果相比誤差稍大，但精度仍稍高于文獻(xiàn)［6］由PMM 數(shù)值計算得到的粘性類系數(shù)的計算精度.另一方面，PMM 試驗(yàn)所采用的數(shù)學(xué)模型是基于小振幅線性假設(shè)的，并沒有考慮非線性項(xiàng)和耦合項(xiàng)的影響，這也是產(chǎn)生誤差的重要原因之一.

從計算效率的角度來看，要想精確獲得全部速度系數(shù)，直拖數(shù)值試驗(yàn)至少需要12個不同的攻角姿態(tài)以及同等數(shù)量的不同的漂角姿態(tài)，計算耗時較長，而PMM 數(shù)值試驗(yàn)只需要在5～6個頻率下進(jìn)行純升沉和純橫蕩試驗(yàn)即可，角速度系數(shù)的計算耗時與此類似.另一方面，直拖數(shù)值試驗(yàn)和回轉(zhuǎn)數(shù)值試驗(yàn)不能在穩(wěn)態(tài)條件下獲得慣性類水動力系數(shù)，這也限制了此二者的應(yīng)用范圍.PMM 數(shù)值試驗(yàn)的最大的優(yōu)勢在于可在一組試驗(yàn)中同時求出慣性類系數(shù)和粘性類系數(shù)，基于大振幅非線性理論的PMM 數(shù)值試驗(yàn)也可獲得非線性系數(shù)和諸多耦合系數(shù)，從而減少了數(shù)值試驗(yàn)次數(shù)，提高了計算效率，這對潛艇操縱性的快速預(yù)報具有重要意義.

3 結(jié) 論

1）潛艇PMM 試驗(yàn)的數(shù)值計算結(jié)果與DTRC所公布的試驗(yàn)數(shù)據(jù)取得較好一致，整體誤差能夠控制在8%以內(nèi)，只有極少數(shù)水動力系數(shù)與試驗(yàn)值偏差超過了20%，從數(shù)值模擬的整體上看計算結(jié)果是真實(shí)可信的.

2）潛艇PMM 數(shù)值模擬可以較快地獲取相關(guān)的水動力系數(shù)，且在一組計算內(nèi)既能夠獲得慣性類水動力系數(shù)也能夠獲得粘性類水動力系數(shù)，這是直拖數(shù)值試驗(yàn)或回轉(zhuǎn)試驗(yàn)數(shù)值計算所不具有的優(yōu)勢.

3）潛艇PMM 試驗(yàn)的數(shù)值模擬所獲得的慣性類水動力系數(shù)具有較高的精度，而粘性類水動力系數(shù)計算精度不如直拖試驗(yàn)或回轉(zhuǎn)試驗(yàn)數(shù)值計算的精度高，若只要求獲得粘性類系數(shù)，本研究建議仍采用常規(guī)試驗(yàn)方法.

［1］VAZ G，TOXOPEUS S，HOLMES S.Calculation of manoeuvring forces on submarines using two viscous－flow solvers［C］∥Proceedings of the ASME 201029thInternational Conference on Ocean，Offshore and Arctic Engineering，OMAE，2010.

［2］ZENG Guanghui，ZHU Jun.Study on key techniques of submarine maneuvering hydrodynamics prediction using numerical method［C］∥2010Second In－ternational Conference on Computer Modeling and Simulation，ICCMS，Taiyuan，China，Oct.21，2010.

［3］詹成勝，劉祖源，程細(xì)得.潛艇水動力系數(shù)數(shù)值計算［J］.船海工程，2008，37（3）：1－4.

［4］柏鐵朝，梁中剛，周軼美.潛艇操縱性水動力數(shù)值計算中湍流模式的比較與應(yīng)用［J］.中國艦船研究，2010，5（2）：22－28.

［5］ZHANG He，YU Yuru，CAI Haopeng.Using CFD software to calculate hydrodynamic coefficients［J］.Marine.Sci，2010（9）：149－155.

［6］張曉頻，龐永杰.多功能潛水器操縱性能與運(yùn)動仿真研究［D］.哈爾濱：哈爾濱工程大學(xué)，2008.

［7］王化明，鄒早建.雙槳多舵船舶操縱性預(yù)報研究［J］.武漢理工大學(xué)學(xué)報：交通科學(xué)與工程版，2006，30（1）：124－127.

［8］陳厚泰.潛艇操縱性［M］.北京：國防工業(yè)出版社，1981.

［9］RODDY R F.Investigation of the stability and control characteristics of several configurations of the darpa suboff model（DTRC model 5470）from captive－model experiments［R］.David Taylor Research Center，Ship Hydrodynamics Department，1990.

［10］ANSYS Inc.Turbulence and near－wall modeling［M］.Canonsburg：ANSTS CFX－Solver Modeling Guide，2009.