孫現(xiàn)有，劉顯龍

(1.海軍駐昆明七五〇試驗(yàn)場軍事代表室；2.中國船舶重工集團(tuán)公司七五〇試驗(yàn)場，云南 昆明 650051)

0 引言

拖纜作為水面及水下武器的承載設(shè)備發(fā)揮著巨大的作用。常見的拖纜拖曳形式有水面艦拖曳航行體，水下航行體拖曳聲學(xué)線列陣等。為了計(jì)算整個拖纜的阻力性能，需要將其參數(shù)化，計(jì)算各個的阻力性能，并將其整體矢量求和。

以往常采用經(jīng)驗(yàn)公式法來計(jì)算拖纜微元的阻力性能[1]。經(jīng)驗(yàn)公式是將來流速度按拖纜的切向、法向2個方向進(jìn)行分解，再由公式計(jì)算切向力和法向力[2]。采用這種方法由于未考慮切向與法向速度相互干擾因素，結(jié)果不夠準(zhǔn)確，連璉等人[3]通過試驗(yàn)的方法嚴(yán)格說明了此種方法帶來的誤差，本文不再贅述。隨著計(jì)算流體力學(xué)(CFD)技術(shù)的發(fā)展，采用CFD方法能夠快速地計(jì)算流體性能。另外，CFD方法僅能計(jì)算有限個攻角時拖纜的阻力性能，為了得到實(shí)際不同攻角的速度，需要采用三次樣條插值的方法來插值求得。

1 流體力學(xué)基礎(chǔ)理論

1.1 控制方程

不可壓縮粘性流體的連續(xù)性方程和動量方程為

(1)

(2)

1.2 湍流模型

要使上述方程封閉，須引入新的湍流模型方程，把應(yīng)力項(xiàng)中的脈動值與時均值聯(lián)系起來。RNGk-ε方程中湍流動能k和湍流脈動強(qiáng)度ε的輸運(yùn)方程為

(3)

(4)

常數(shù)為：G1ε=1.42；G2ε=1.68；Cμ=0.0845；

σk=1.0；σε=1.3；η0=4.38；β=0.012。

2 CFD模擬

2.1 流體計(jì)算域

拖纜微元與迎流速度的夾角(即攻角)從0°到90°范圍變化。0°及90°情況較為簡單，可以較快地得到結(jié)構(gòu)化網(wǎng)格，較為復(fù)雜的是0°至90°之間的情況。對于此種情況，將速度入口由平面改為弧形面，可以較好地解決這一問題[4]，其流體計(jì)算域及邊界設(shè)置情況見圖1。

需要說明的是，拖纜微元在不同攻角時阻力性能是不同的，而計(jì)算時又不可能由CFD的方法計(jì)算每種攻角情況，為了解決這一問題，采用三次樣條插值的方法[5]。

2.2 網(wǎng)格劃分

為了得到高質(zhì)量的網(wǎng)格，對整個計(jì)算域采取分塊化處理，再劃分為結(jié)構(gòu)化網(wǎng)格。在拖纜周圍設(shè)置邊界層網(wǎng)格，第一層控制近壁面的網(wǎng)格高度和質(zhì)量，近壁面高度按估算公式(5)確定：

(5)

式中：Δy為第一層邊層界厚度；L為長度值；Re為雷諾數(shù)。實(shí)際取y+為30，整體網(wǎng)格數(shù)量為90萬個。拖纜流場網(wǎng)格見圖2。

2.3 邊界條件

邊界條件的設(shè)定是保證CFD實(shí)現(xiàn)的必要條件之一。本文計(jì)算域的邊界條件設(shè)定見圖1，包括：速度入口、壓力出口、壁面和對稱面。

1)速度入口：在進(jìn)流面設(shè)定迎流速度的大小和方向，Vin=V0；

2)壓力出口：在出流面設(shè)定相對于參考壓力點(diǎn)的流體靜壓力，Pout=0；

3)壁面條件：潛艇表面及圓柱表面設(shè)定為不可滑移邊界條件，u=v=w=0；

2.4 數(shù)值計(jì)算

采用有限體積法離散控制方程和湍流模式，并用耦合隱式求解法計(jì)算。在差分格式中，壓力項(xiàng)使用標(biāo)準(zhǔn)格式，動量項(xiàng)、湍流動能項(xiàng)以及耗散率項(xiàng)先使用一階迎風(fēng)格式迭代500次，再采用二階迎風(fēng)格式[6]。

采用RNGk-ε湍流模型，并采用有限體積法對控制方程進(jìn)行離散。壓力-速度耦合迭代求解采用協(xié)調(diào)一致的SIMPLEC算法，壓力項(xiàng)采用標(biāo)準(zhǔn)離散格式外，動量、湍流動能和湍流耗散率采用二階迎風(fēng)格式，亞松弛因子均采用默認(rèn)值。近壁處區(qū)域采用壁面函數(shù)法[7]。

3 三次樣條插值

在實(shí)際計(jì)算過程中，通過CFD方法計(jì)算拖纜微元在攻角分別為0°、10°、20°、…、90°的阻力情況。實(shí)際計(jì)算時若碰到不同屬于上述攻角情況時，采用三次樣條插值的方法進(jìn)行估算[8]。

曲線插值一般有指數(shù)函數(shù)插值、拉格朗日插值、三次樣條插值、B樣條曲線、NURBS曲線等方法。這里采用三次樣條進(jìn)行插值。

在區(qū)間[a,b]內(nèi)給定一個分劃Δ：a=x0

1)每個子區(qū)間[xi-1,xi],i=1,2,…,n上為三次多項(xiàng)式；

2)y(x)在整個區(qū)間[a,b]上二階連續(xù)可導(dǎo)。

則稱y(x)為區(qū)間[a,b]上關(guān)于分劃Δ的三次樣條函數(shù)，其中xi(i=0,1,…,n)稱為節(jié)點(diǎn)。插值點(diǎn)如表1所示。

表1 三次樣條曲線插值點(diǎn)

則三次樣條插值函數(shù)為

(6)

式(6)中：

hk=xk+1-xk(k=0,1,…,n-1)

(7)

λkmk-1+2mk+μkmk+1=gk(k=1,2,…,n-1)

(8)

其中：

4 拖纜計(jì)算方法簡介

由于拖纜是柔性的，受力將發(fā)生變形[8]。因此，各處拖纜的來流攻角是不同的。采用有限差分法，將拖纜處理成有限個小段，對各小段分別受力分析，由初始點(diǎn)的力、轉(zhuǎn)角條件求出未知位置的力、轉(zhuǎn)角及坐標(biāo)點(diǎn)。

為了滿足工程中的需要，且基本滿足實(shí)際情況，做出以下幾點(diǎn)假設(shè)：

1)拖纜在張力作用下的伸縮忽略不計(jì)；

2)拖纜為柔性的，不傳遞彎矩；

3)整個流動為二維流動；

4)作用于拖纜微元上的水動力系數(shù)等同于作用于無限長的拖纜。

將一小段拖纜微元進(jìn)行受力分析，見圖3。

按照切向力與法向力方向進(jìn)行分解，由力的平衡可得：

(T+dT)sin(dφ)+pdscos(φ)=Dds

(9)

(T+dT)cos(dφ)=pdssin(φ)+Fds+T

(10)

式 (9)、(10)中：T為拉力；D為單位長度拖纜法線方向的力(由CFD及三次樣條插值計(jì)算)；F為單位長度拖纜切線方向的力(由CFD及三次樣條插值計(jì)算)；p為單位長度拖纜在水中重力。

忽略二階小量，由式(10)得：

(11)

將式(11)代入式(9)得：

(12)

由式(11)、(12)可求出不同位置拖纜拉力、轉(zhuǎn)角等值。

另外：

dx=ds×cosφ

(13)

dy=ds×sinφ

(14)

由式(13)、(14)可求出相應(yīng)位置點(diǎn)的坐標(biāo)值。

5 計(jì)算結(jié)果分析

在實(shí)際計(jì)算時，可能會遇到不同直徑的拖纜。眾所周知，在CFD計(jì)算中，最為復(fù)雜的是采用ICEM來劃分結(jié)構(gòu)化網(wǎng)格，耗費(fèi)了大量的人力、物力。目前，Workbench(15.0版)后已經(jīng)將ICEM集成進(jìn)去，但是當(dāng)修改拖纜直徑后，ICEM的拓?fù)浞謮K處理并不能自動化解決，還是需要大量、重復(fù)性人工干預(yù)。為了較好地解決這一問題，采用ICEM自帶的腳本語言，將劃分網(wǎng)格中的建組、分塊、劃網(wǎng)格、導(dǎo)出網(wǎng)格全過程監(jiān)控、錄制并保存，當(dāng)模型的直徑發(fā)生變化時，重新讀入一次腳本語言，即可一鍵化劃分網(wǎng)格，將原本需要2 h左右的劃分網(wǎng)格時間縮減至1 min以內(nèi)。需要說明的是，建立模型的直徑需要在Workbenchk中的DM模塊下進(jìn)行，建立好模型后，僅修改其直徑等參數(shù)，以便于腳本語言無縫識別，否則在其它建模環(huán)境中，如Pro/E、Rhino等，腳本語言對點(diǎn)、線、面、體的識別極易出錯。

以直徑為10 mm的拖纜(長度為1 m)為例，其在迎流速度為5 kn，其切向力及法向力見表2。

表2 10 mm拖纜在不同攻角時切向力與法向力(v=5 kn)

以攻角為0°、45°和90°為例，其壓力云圖見圖4。需要說明的是，表2中攻角沒有計(jì)算接近0°和90°的角度，因?yàn)榇藭r網(wǎng)格斜角過大，劃分網(wǎng)格質(zhì)量不夠好，影響了計(jì)算的精度，采用三次樣條插值進(jìn)行插值計(jì)算。攻角在5°時，計(jì)算得到切向力0.48 N，法向力0.32 N。

6 結(jié)束語

上述計(jì)算表明，本文采用的方法有以下優(yōu)勢：

1)采用CFD的方法計(jì)算拖纜微元的阻力性能，考慮了切向速度與法向速度間干擾，優(yōu)于以往經(jīng)驗(yàn)公式；

2)采用參數(shù)化的方法建模、劃分結(jié)構(gòu)化網(wǎng)格，大大節(jié)省了前處理時間；

3)采用三次樣條插值的方法，可以快速計(jì)算實(shí)際中任一攻角的阻力性能。

綜合前三者優(yōu)勢，結(jié)合編程的方法，可以將整個拖曳系統(tǒng)(含多節(jié)點(diǎn))在不同速度拖曳時的總阻力及姿態(tài)情況計(jì)算出來。