王小波， 韓端鋒

(1. 航運(yùn)技術(shù)與安全國家重點(diǎn)實(shí)驗(yàn)室， 上海 200135； 2. 哈爾濱工程大學(xué), 哈爾濱 150001)

0 引 言

回轉(zhuǎn)橢球體與潛艇主體相似，雖然外形簡單，但其繞流場卻是一個(gè)典型的三維流動(dòng)實(shí)例，對(duì)這一流動(dòng)預(yù)報(bào)的好壞可以間接地檢驗(yàn)數(shù)值計(jì)算方法對(duì)潛艇主體繞流流動(dòng)的預(yù)報(bào)能力。因此，作為數(shù)學(xué)船型的典例，橢球體的基礎(chǔ)研究有一定的必要性和重要意義。近年來，橢球體粘性分離流研究受到普遍重視，如國外學(xué)者對(duì)長細(xì)比為6∶1橢球體的粘性繞流進(jìn)行了實(shí)驗(yàn)[1-3]和數(shù)值計(jì)算研究[4-5]，主要集中于湍流模式，對(duì)模型前期網(wǎng)格的研究較少；國內(nèi)對(duì)橢球體模型的研究也致力于湍流模型等方向，如李佳等采用不同的湍流模型對(duì)某一潛器進(jìn)行了計(jì)算對(duì)比，并在入口湍流強(qiáng)度方面做了深入研究，相對(duì)而言對(duì)網(wǎng)格劃分的相關(guān)研究較少。

參考相關(guān)文獻(xiàn)及著作，對(duì)回轉(zhuǎn)橢球體粘性三維繞流場進(jìn)行了數(shù)值模擬，集中于繞流阻力的計(jì)算；以試驗(yàn)值為參考標(biāo)準(zhǔn)，研究了網(wǎng)格劃分對(duì)繞流場模擬及計(jì)算結(jié)果的影響規(guī)律。

1 橢球體幾何建模及控制域生成

根據(jù)試驗(yàn)結(jié)果，建立長細(xì)比為6∶1的回轉(zhuǎn)橢球體模型，其長半軸a=750 mm，短半軸b=c=125 mm，橢球面通過定軸旋轉(zhuǎn)而成，幾何模型見圖1。

模擬無界流場中橢球體三維繞流。為消除壁面影響，建立長方體橢球體流場控制域(見圖2)，采用笛卡爾直角坐標(biāo)系，橢球體幾何中心位于坐標(biāo)系原點(diǎn)?？刂朴蛟O(shè)定：X向：-5.5L～2.5L，Y向：-2L～2L，Z向與Y向同。其中L為橢球體特征長度，L=1 500 mm。

圖1 橢球體幾何建模

圖2 橢球體幾何控制域

2 控制域網(wǎng)格劃分

從CFD模擬現(xiàn)狀分析，網(wǎng)格質(zhì)量對(duì)CFD計(jì)算精度和計(jì)算效率有重要影響。以網(wǎng)格形式入手，逐步深入探討網(wǎng)格數(shù)目以及邊界層等相關(guān)因素對(duì)計(jì)算結(jié)果的影響。

目前，網(wǎng)格形式主要集中為結(jié)構(gòu)化網(wǎng)格和非結(jié)構(gòu)化網(wǎng)格兩種。首先應(yīng)用結(jié)構(gòu)化網(wǎng)格，在橢球體兩端采用O型結(jié)構(gòu)化網(wǎng)格形式，這樣可極大地提高網(wǎng)格質(zhì)量，使控制域內(nèi)網(wǎng)格尤其是近壁面網(wǎng)格不至于出現(xiàn)過小尖角等極限情況，邊界層與體網(wǎng)格均為六面體；之后使用非結(jié)構(gòu)化網(wǎng)格，則模型的體網(wǎng)格為四面體單元，邊界層為棱錐形的五面體單元。兩種網(wǎng)格劃分形式分別見圖3、圖4，非結(jié)構(gòu)化網(wǎng)格模型的表面網(wǎng)格見圖5。

圖3 結(jié)構(gòu)化網(wǎng)格及近壁面邊界層

圖4 非結(jié)構(gòu)化網(wǎng)格及邊界層

對(duì)邊界層網(wǎng)格的設(shè)定亦進(jìn)行了定量研究。不同的邊界層數(shù)以及底層邊界層尺寸的變化都會(huì)給計(jì)算結(jié)果帶來一定的波動(dòng)。根據(jù)現(xiàn)有的一些研究成果可知，針對(duì)同一網(wǎng)格劃分策略，改變邊界條件對(duì)計(jì)算結(jié)果精確性的影響相對(duì)較大，而改變模型引起的影響則相對(duì)較小[8]。在邊界層的設(shè)置上，近壁面的網(wǎng)格處理采用壁面函數(shù)法，根據(jù)壁面函數(shù)法的要求設(shè)定底層網(wǎng)格的尺寸。圖6為局部空間網(wǎng)格及邊界層網(wǎng)格。

圖5 模型表面的非結(jié)構(gòu)網(wǎng)格

圖6 結(jié)構(gòu)網(wǎng)格邊界層橫截面

3 邊界條件及求解器參數(shù)

欲尋求可靠、高效的數(shù)值模擬方法，湍流模型的選取也尤為重要。根據(jù)李佳的結(jié)論，不同的湍流模型對(duì)計(jì)算結(jié)果會(huì)有影響，但是差別不大甚至很小。通過初期計(jì)算發(fā)現(xiàn)，選取標(biāo)準(zhǔn)k-e湍流模型的阻力計(jì)算結(jié)果與試驗(yàn)值是非常吻合的。除去網(wǎng)格和湍流模型，入口邊界條件以及求解器參數(shù)也是影響求解的重要因素。目前針對(duì)k-e模型，主要使用的湍流參數(shù)有兩種：湍動(dòng)能K和湍動(dòng)耗散率e，入口湍流強(qiáng)度I和湍動(dòng)粘度比μt/μ[6]。對(duì)于不可壓縮外流而言，入口湍流強(qiáng)度的取值并不像內(nèi)流一樣有經(jīng)驗(yàn)公式可循，只能根據(jù)試驗(yàn)值找出一定的規(guī)律[9]。入口處湍動(dòng)能和湍動(dòng)能耗散率根據(jù)經(jīng)驗(yàn)公式設(shè)定：Kin=3.75×10-3U02，εin=K3/2/0.03[10]；入口處湍流強(qiáng)度和湍動(dòng)粘度比根據(jù)文獻(xiàn)[11]可取I=0.5%，μt/μ=5，變化范圍為I=0.1%～0.5%；李佳總結(jié)出自由來流湍流強(qiáng)度取值的函數(shù)為I=2×107×Re-1.229 8(%)[9]。

經(jīng)過初步的計(jì)算比較，確定的數(shù)值方法為：求解器采用分離式求解器，基于有限體積法來離散控制方程；求解穩(wěn)態(tài)三維RANS方程、單相流模型、粘性模型選擇標(biāo)準(zhǔn)k-e雙方程湍流模型。邊界條件：采用速度進(jìn)口入口條件，出口條件為outflow，物體壁面滿足無滑移條件。入口邊界條件的湍流參數(shù)項(xiàng)為K and Epsilon；壓力速度耦合方式為SIMPIEC，壓力插值項(xiàng)為PRESTO，對(duì)流項(xiàng)及擴(kuò)散項(xiàng)均采用二階迎風(fēng)格式離散。

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

4.1 網(wǎng)格形式

首先，固定在某一雷諾數(shù)下計(jì)算，與試驗(yàn)值比較，初步確定適宜的網(wǎng)格形式；在此基礎(chǔ)上，通過考察網(wǎng)格數(shù)量、邊界層網(wǎng)格數(shù)以及底層邊界層尺寸等參數(shù)進(jìn)一步尋求高精度高效率的網(wǎng)格劃分策略；最后，通過改變計(jì)算工況驗(yàn)證其準(zhǔn)確度。

針對(duì)結(jié)構(gòu)化網(wǎng)格和非結(jié)構(gòu)化網(wǎng)格做了大量計(jì)算，從中選取了有代表性的五組，同組的網(wǎng)格數(shù)相同，而網(wǎng)格形式不同，其余參數(shù)均一樣，計(jì)算結(jié)果見圖7。

圖7中：縱向看，非結(jié)構(gòu)化網(wǎng)格模型的計(jì)算值基本上隨著網(wǎng)格數(shù)量的增加呈減小的趨勢(shì)，而結(jié)構(gòu)化網(wǎng)格模型的計(jì)算值在網(wǎng)格數(shù)量多于20萬以后變化非常小。橫向?qū)Ρ?，結(jié)構(gòu)化網(wǎng)格模型的結(jié)果普遍低于同組的非結(jié)構(gòu)化網(wǎng)格模型，且更接近于試驗(yàn)值(0.63)；在所有的非結(jié)構(gòu)化網(wǎng)格模型中，最好的結(jié)果是0.684 5，誤差為8.65%。可以說，結(jié)構(gòu)化網(wǎng)格形式整體上更有利于計(jì)算結(jié)果的收斂和計(jì)算精度的提高，因此下面的工作基于結(jié)構(gòu)化網(wǎng)格進(jìn)行。

4.2 網(wǎng)格數(shù)量

由于數(shù)值模擬的精度受到網(wǎng)格數(shù)量的影響，從常規(guī)的角度思考，研究人員往往以增加網(wǎng)格數(shù)量的方式保證計(jì)算精度。此時(shí)，雖然數(shù)值模擬不失真，但會(huì)增大計(jì)算工作量和CPU消耗。在實(shí)際操作中，事前并不知道流場的確切情況，因此很難根據(jù)流場來建造網(wǎng)格，勢(shì)必造成網(wǎng)格數(shù)量的急劇上升[7]。針對(duì)同一流場域，劃分了不同網(wǎng)格數(shù)量，經(jīng)過逐一計(jì)算，得出了網(wǎng)格數(shù)量對(duì)計(jì)算結(jié)果的影響關(guān)系。

從實(shí)際計(jì)算中抽取5個(gè)網(wǎng)格數(shù)量不同且增量較為均勻的計(jì)算模型，對(duì)比計(jì)算結(jié)果，見圖8。

圖7 兩種網(wǎng)格形式的計(jì)算結(jié)果對(duì)比

圖8 網(wǎng)格數(shù)對(duì)結(jié)果的影響

由圖8可知，各組計(jì)算值與試驗(yàn)值較為接近，誤差<6%，滿足工程精度要求。隨著網(wǎng)格數(shù)目的增加，計(jì)算值小幅度減小，并逐漸靠近試驗(yàn)值。按照此趨勢(shì)，更多的網(wǎng)格數(shù)目將取得更為理想的結(jié)果。然而，龐大的網(wǎng)格數(shù)量必然耗費(fèi)大量計(jì)算時(shí)間，繼而影響工作進(jìn)度，因此，對(duì)邊界層展開了研究。

4.3 邊界層及底層網(wǎng)格尺寸

圖9為同一模型在邊界層設(shè)置不同網(wǎng)格層數(shù)時(shí)的阻力計(jì)算結(jié)果，5組網(wǎng)格層數(shù)分別為0、6、8、13、25。為考察邊界層的必要性，特設(shè)定一個(gè)無邊界層的模型作為比較，其余4組的網(wǎng)格層數(shù)逐漸增加，對(duì)于回轉(zhuǎn)橢球體，25層的邊界層已經(jīng)足夠。

由圖9可知，不設(shè)置邊界層的模型，計(jì)算結(jié)果與試驗(yàn)值相差較多，這與相關(guān)文獻(xiàn)的邊界層的重要性結(jié)論一致。其余4個(gè)模型的計(jì)算結(jié)果較為理想，最大誤差為5.2%，最小為2%。

根據(jù)覃文潔的研究結(jié)論，選取合適的壁面層網(wǎng)格的y+值能夠有效提高解的精度；在處理近壁面時(shí)采用壁面函數(shù)法，壁面層網(wǎng)格的y+值選取的合理性將會(huì)影響計(jì)算結(jié)果的合理性[8]。當(dāng)然，根據(jù)壁面函數(shù)法確定的y+有一個(gè)取值范圍，為尋求更為適宜的近壁面網(wǎng)格尺寸，進(jìn)行了對(duì)比計(jì)算，結(jié)果見圖10。

圖9 邊界層數(shù)對(duì)計(jì)算結(jié)果的影響

圖10 近壁面底層網(wǎng)格尺寸對(duì)計(jì)算結(jié)果的影響

經(jīng)壁面函數(shù)法確定的底層網(wǎng)格尺寸范圍是2.8～14.2 mm。圖中給出了五組不同網(wǎng)格尺寸的計(jì)算結(jié)果，誤差最小的是徑向尺寸為8.2 mm的模型。另外，超出范圍的三組模型計(jì)算結(jié)果較差，可知底層網(wǎng)格尺寸過小或過大都不利于計(jì)算。

經(jīng)過以上研究，最有利的網(wǎng)格劃分策略為:結(jié)構(gòu)網(wǎng)格形式，邊界層的網(wǎng)格層數(shù)為6層，底層徑向尺寸為8.2 mm，滿足壁面函數(shù)法的要求，總網(wǎng)格數(shù)約為7.5萬。求解器參數(shù)方面，壓力速度耦合采用SIMPLEC格式，各方程均采用二階迎風(fēng)格式，入口邊界的湍流參數(shù)為K和E。

4.4 邊界條件

入口邊界條件的湍流參數(shù)尚有多種結(jié)論，應(yīng)用上述研究結(jié)論，再考察李佳設(shè)定的入口邊界條件是否有利于計(jì)算。分別賦予湍流強(qiáng)度I以0.1、0.5、1.271 37、5、10，而湍動(dòng)粘度比恒定為5，計(jì)算結(jié)果見圖11。

由圖11可知，隨著I的增加，阻力計(jì)算值呈上升趨勢(shì)，并與試驗(yàn)值不斷接近，且在I取值為10時(shí)，達(dá)到較為理想的結(jié)果，誤差僅為0.59%。而李佳給出的I的經(jīng)驗(yàn)公式并不完全適用于此計(jì)算模型，也反映了CFD的不確定性。

為驗(yàn)證研究成果，對(duì)雷諾數(shù)為1.42×106(速度為1.0 m/s)及2.13×106(速度為1.5 m/s)都做了數(shù)值計(jì)算。應(yīng)用以上結(jié)論，使用相同的計(jì)算模型、保持統(tǒng)一網(wǎng)格劃分策略、微調(diào)網(wǎng)格數(shù)量以及邊界層設(shè)置，均取得了較為良好的計(jì)算結(jié)果，見圖12。

圖11 湍流強(qiáng)度對(duì)計(jì)算結(jié)果的影響

圖12 不同雷諾數(shù)下計(jì)算值與試驗(yàn)值的比較

圖12中：連接正方塊的曲線為試驗(yàn)值，連接菱形的曲線為數(shù)值模擬的計(jì)算結(jié)果。計(jì)算值與試驗(yàn)值的吻合程度較好，誤差很小，滿足工程應(yīng)用的要求，驗(yàn)證了研究成果的準(zhǔn)確性及價(jià)值。同時(shí),在這兩工況計(jì)算中，為了達(dá)到更好的計(jì)算結(jié)果，改變了入口邊界條件的湍流參數(shù)設(shè)定，與速度為0.5 m/s時(shí)略有不同。

5 結(jié) 語

采用商用CFD分析軟件完成了回轉(zhuǎn)橢球體的阻力計(jì)算，對(duì)其網(wǎng)格劃分及網(wǎng)格質(zhì)量的討論，可以為復(fù)雜模型的網(wǎng)格劃分策略提供借鑒。主要集中探討了幾何模型控制域的網(wǎng)格劃分對(duì)計(jì)算結(jié)果的影響，包括網(wǎng)格形式、網(wǎng)格數(shù)量以及邊界層網(wǎng)格劃分等，并得出以下結(jié)論：

1) 控制域網(wǎng)格劃分對(duì)計(jì)算結(jié)果有明顯的影響，針對(duì)不同的計(jì)算模型有不同的網(wǎng)格劃分策略，沒有通用的一致性原則，但是相近的模型之間可以互相借鑒。

2) 網(wǎng)格數(shù)量的增加一定程度上會(huì)帶來高精度的計(jì)算結(jié)果，但是龐大的網(wǎng)格數(shù)量將給計(jì)算機(jī)帶來過大的負(fù)荷從而影響工作效率，并且網(wǎng)格數(shù)量增加到一定數(shù)量，結(jié)果將趨于穩(wěn)定，變化極小。因此，尋找適宜的網(wǎng)格數(shù)量將有助于提高計(jì)算精度并減少工作量，這在CFD實(shí)際應(yīng)用中有極為重要的意義。

3) 邊界層的存在將極大地提高計(jì)算精度，同時(shí)近壁面底層網(wǎng)格尺寸在提高計(jì)算精度上值得深入探索，目前主要以壁面函數(shù)法為參考。

參考文獻(xiàn)：

[1] Chesnakas C J, Simpson R L. An investigation of the three-dimensional turbulent flow in the cross-flow separation region of a 6:1 prolate spheroid[J]. Experiments in Fluids, 1994(17):68-74.

[2] Chesnakas C J, Simpson R L. Detailed investigation of the three-dimensional separation about a prolate spheroid[J]. AIAA Journal, 1997,35(6):990-999.

[3] Wetzel T G, Simpson R L. Unsteady three-dimensional cross-flow separation measurements on a pro-late spheroid undergoing time-dependent maneuvers[J].21st Symposium on Naval Hydrodynamics,1997.

[4] Rhee S H, Hino T. Computational investigation of 3D turbulent flow separation around a spheroid using an unstructured grid method[J].Journal of Society of Naval Architects of Japan,2000(18):1-9.

[5] Kim S E, Rhee S H, Cokljat D. High-incidence and dynamic pitch-up maneuvering characteristics of a prolate spheroid-CFD validation[J].24th Symposium on Naval Hydrodynamics,2002.

[6] 王福軍.計(jì)算流體動(dòng)力學(xué)分析[M].北京：清華大學(xué)出版社,2004.

[7] 李剛.潛器水動(dòng)力特性的數(shù)值模擬[D].哈爾濱：哈爾濱工程大學(xué),2011.

[8] 覃文潔，胡春光，郭良平,等. 近壁面網(wǎng)格尺寸對(duì)湍流計(jì)算的影響[J]. 北京理工大學(xué)學(xué)報(bào)，2006(5)：388-392.

[9] 李佳,黃德波，鄧銳.載人潛器阻力的數(shù)值計(jì)算方法分析[J].船舶力學(xué)，2010，14(4)：334-339.

[10] 邱遼原.潛艇粘性流場的數(shù)值模擬及其阻力預(yù)報(bào)的方法研究[D].武漢：華中科技大學(xué)，2006.

[11] 韓占忠，王敬，蘭小平．FLUENlL流體工程仿真技術(shù)與實(shí)例應(yīng)用[M]．北京：北京理工大學(xué)出版社，2004．