過渡型艇數(shù)值模擬的校核與驗(yàn)證

余澤爽 毛筱菲

(武漢理工大學(xué)交通學(xué)院1) 武漢 430063) (武漢理工大學(xué)教育部高性能船重點(diǎn)實(shí)驗(yàn)室2) 武漢 430063)

過渡型艇數(shù)值模擬的校核與驗(yàn)證

余澤爽1)毛筱菲1,2)

(武漢理工大學(xué)交通學(xué)院1)武漢 430063) (武漢理工大學(xué)教育部高性能船重點(diǎn)實(shí)驗(yàn)室2)武漢 430063)

基于粘性CFD軟件STAR-CCM+的結(jié)構(gòu)化笛卡爾網(wǎng)格和重疊網(wǎng)格技術(shù)，通過改變船體表面和加密區(qū)網(wǎng)格大小，形成四種不同密度的網(wǎng)格形式，對過渡型艇在靜水中從排水到過渡航行階段(Fr▽=0.664,0.860,1.934)進(jìn)行數(shù)值模擬，對其阻力和縱傾角在數(shù)值模擬中對網(wǎng)格的不確定度進(jìn)行分析，采用兩種方式對結(jié)果隨網(wǎng)格的收斂情況進(jìn)行分析，并采用試驗(yàn)結(jié)果對其進(jìn)行驗(yàn)證.由校核與驗(yàn)證的結(jié)果看，在排水航行階段，笛卡爾網(wǎng)格在一定程度上能夠滿足模擬需求；但在過渡航行階段，重疊網(wǎng)格具有十分明顯的優(yōu)勢.

STAR-CCM+；重疊網(wǎng)格；不確定度分析；校核與驗(yàn)證

0 引 言

計算流體力學(xué)(CFD)和實(shí)驗(yàn)?zāi)M(EFD)是研究水動力問題的主要手段，其中高效的CFD技術(shù)具有無可比擬的優(yōu)勢，但對數(shù)值模擬精確性與不確定度的分析與評估，顯得比數(shù)值模擬本身同等重要[1].

Stern等[2]在數(shù)學(xué)框架下提出不確定度與誤差分析的概念，并提出校核與驗(yàn)證的具體內(nèi)容.為ITTC提出了數(shù)值模擬校核與驗(yàn)證的草案，并提出了校核與驗(yàn)證的定義與理論基礎(chǔ)[3].他們將該理論對一艘貨船和集裝箱船RANS模擬進(jìn)行了校核和驗(yàn)證，提出了針對實(shí)例的精確性和不確定度的分析與評估[4].Jin等[5]應(yīng)用ReFRESCO和PARNASSOS兩種求解器對下臺階流動問題進(jìn)行模擬，基于該草案的方法對結(jié)果進(jìn)行校核和驗(yàn)證，采用了最小二乘法對數(shù)值模擬的結(jié)果進(jìn)行擬合，以校核其數(shù)值模擬的收斂性與可靠性.Zhu等[6-8]基于URANS求解器的Star-CCM+對FLNG船在迎浪中的運(yùn)動的數(shù)值模擬，并采用該方法對所受載荷、波浪等進(jìn)行校核，并通過模型試驗(yàn)，確定試驗(yàn)不確定度，對數(shù)值模擬進(jìn)行驗(yàn)證.

文中從ITTC關(guān)于數(shù)值模擬的校核與驗(yàn)證的角度，以純數(shù)學(xué)的角度，證明與詮釋了對于過渡型艇的中高速數(shù)值模擬，重疊網(wǎng)格比普通笛卡爾網(wǎng)格的優(yōu)勢.

1 理論背景

1.1 定義與概念

在對某物理問題的求解和模擬時，出于可行性的要求，會根據(jù)實(shí)際問題提出概念模型，一方面通過制作船模，在水池中進(jìn)行模型實(shí)驗(yàn)(EFD)；另一方面通過建模軟件制作船模，通過代碼建立數(shù)值水池，通過求解N-S方程進(jìn)行數(shù)值模擬(CFD)，見圖1[9-12].出于成本的考慮，數(shù)值模擬成為重要的研究手段.

圖1 校核與驗(yàn)證(V&V)

V&V的目的是要確定數(shù)值模擬的不確定度Uφ，對真實(shí)值φexact有95%置信區(qū)間的數(shù)值解φi，即

φi-Uφ≤φexcat≤φi+Uφ

(1)

數(shù)值模擬誤差δs為數(shù)值模擬S與真實(shí)值T之差，建模誤差δsM與數(shù)值誤差δsN之和，不確定度為

(2)

對于修正后模擬值Sc，其不確定度為

(3)

1.2 校核

在數(shù)值模擬中，數(shù)值誤差是最主要的部分.

δSN=δI+δG+δT+δP

(4)

式中：δI為迭代誤差，采用雙精度一般小于0.2%δsN，可以忽略；δG為網(wǎng)格收斂誤差；δT為時間步長收斂誤差；δP為其他參數(shù)誤差.不確定度為

(5)

1.2.1收斂性研究

改變所研究的參數(shù)K，得到不同的算例因子NKn以及模擬結(jié)果SKn，則收斂因子RK為

(6)

當(dāng)參數(shù)K不滿足固定變化率時，

(7)

當(dāng)模擬結(jié)果SKn為向量時，

(8)

對于收斂因子RK：當(dāng)01為發(fā)散.為了更好的判斷收斂因子，因此至少需要三個算例作為評估.

1.2.2單調(diào)收斂

對于單調(diào)收斂的參數(shù)K，對其不確定度進(jìn)行求解和分析.收斂階數(shù)PK用于校核收斂速度：

(9)

修正系數(shù)CK用于校核計算結(jié)果與漸進(jìn)值的趨近程度：

(10)

式中：Pest=2為參考收斂精度.

(11)

(12)

則經(jīng)過修正后的結(jié)果SC為

(13)

當(dāng)修正系數(shù)CK遠(yuǎn)大于或者小于1時，則對參數(shù)K數(shù)值模擬的不確定度UK為

(14)

當(dāng)修正系數(shù)CK接近1時，則對參數(shù)K數(shù)值模擬修正的不確定度UCK為

(15)

1.2.3收斂評估

當(dāng)參數(shù)增長率不固定時的收斂評估：采用最小二乘法將結(jié)果與增長率進(jìn)行擬合，通過比較收斂階數(shù)，以評估結(jié)果的收斂情況.

1.3 驗(yàn)證

模型實(shí)驗(yàn)的比較誤差E為實(shí)驗(yàn)結(jié)果D與數(shù)值結(jié)果S之差，包括實(shí)驗(yàn)誤差δD，歷史數(shù)據(jù)誤差δsPD、建模假定誤差δsMA和數(shù)值模擬誤差δsN，即有

E=|D-S|=δD-(δSPD+δSMA+δSN) (16)

驗(yàn)證不確定度Uval為

(17)

式中：USMA一般難以確定；USPD一般不予考慮，故有保守估計

(18)

當(dāng)|E|UV時，需要比較USMA與E量級的關(guān)系，確定是否需要考慮USMA，從而進(jìn)行驗(yàn)證.

2 過渡型艇數(shù)值模擬

2.1 模型參數(shù)

過渡型艇主尺度見表1.

表1 過渡型艇主尺度

2.2 計算說明

計算域取為船前2倍船長、船后3.5倍船長、船側(cè)1倍船長、自由面下2倍船長、自由面上1.5倍船長，選取半模進(jìn)行計算.湍流模型為k-ε；時間步長為0.05 s.

圖2a)為普通笛卡爾網(wǎng)格，在自由面出采用各向異性的加密網(wǎng)格，船體周圍采用各向同性的加密處理；圖2b)為重疊網(wǎng)格，在自由面和重疊網(wǎng)格區(qū)域外進(jìn)行加密處理之外，為了捕捉船后的興波，在船后的梯形區(qū)域內(nèi)采用各向異性的加密網(wǎng)格.

圖2 網(wǎng)格圖

對普通結(jié)構(gòu)化笛卡爾網(wǎng)格，保持基礎(chǔ)尺寸為船體表面尺寸的8倍，船體表面網(wǎng)格尺寸分別為0.010 6，0.015，0.021 2和0.03 m.對重疊網(wǎng)格，保持背景網(wǎng)格基礎(chǔ)尺寸為重疊網(wǎng)格的2倍，為船體表面網(wǎng)格的8倍，重疊網(wǎng)格基礎(chǔ)尺寸分別為0.03，0.05，0.075，0.1 m.網(wǎng)格數(shù)量見表2.

表2 網(wǎng)格數(shù)量

2.3 數(shù)值模擬

圖3為笛卡爾網(wǎng)格在V=1.26,1.68,3.78 m/s時，以及重疊網(wǎng)格在V=3.78 m/s時的航行模擬圖.

圖3 航線模擬圖

由圖3可知，隨著網(wǎng)格數(shù)量的增加，數(shù)值模擬對流場信息的捕捉能力隨之提高.由圖3a)～b)可知，在排水航行階段，隨著網(wǎng)格數(shù)的提高，對船后興波的橫波系與散波系的捕捉與表達(dá)也更為清晰；由圖3c)～d)可知，在過渡航行階段，船后興波的橫波系消失，由于重疊網(wǎng)格在船后興波的扇形區(qū)域內(nèi)進(jìn)行了網(wǎng)格加密，其散波系的捕捉更為精確.

圖4a)為重疊網(wǎng)格在過渡航行階段的模擬.過渡型艇在高速航行時產(chǎn)生的雞尾流現(xiàn)象以及明顯的縱傾航態(tài)變化，因此，該艇在高速航行時會產(chǎn)生與重力相當(dāng)?shù)乃畡由Γ兄谔w，從而減少濕表面積降低摩擦阻力，見圖4b).隨著航速增加時，雞尾流與尾板之間形成的“空穴”，該“空穴”的長度Δl為“虛長度”，相當(dāng)于增加船長，使船體的修長度系數(shù)提高，從而降低航行的興波阻力.

圖4 過渡航行階段的現(xiàn)象

圖5為過渡型艇在排水和過渡航行階段，船體表面Y+值的分布情況.總體來看，當(dāng)使用k-ε湍流模型時，Y+值的分布基本滿足30≤Y+≤200.此外可見在排水航行階段，Y+值主要為100～200；當(dāng)航速增加時，在過渡航行階段，Y+值急劇下降，保持在30左右.這與文獻(xiàn)[13]關(guān)于滑行艇Y+值分布規(guī)律的結(jié)論基本一致.說明，數(shù)值模擬的邊界層處理，基本滿足湍流模型的需求.

圖5 Y+值分布

數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果見表3，阻力系數(shù)CT為

(19)

式中：RT為總阻力；V為航速.

表3 數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果

3 校核與驗(yàn)證

3.1 校核

3.1.1網(wǎng)格不確定度分析

表4為使用笛卡爾網(wǎng)格時，在V=1.26,1.68,3.78 m/s時，以及重疊網(wǎng)格在V=3.78 m/s時阻力系數(shù)和縱傾角的不確定度.

表4 不確定度計算結(jié)果

注:①V=1.26 m/s,Fr▽=0.644(笛卡爾網(wǎng)格)；②V=1.68 m/s,Fr▽=0.860(笛卡爾網(wǎng)格)；③V=3.78 m/s,Fr▽=1.934(笛卡爾網(wǎng)格)；④V=3.78 m/s,Fr▽=1.934(重疊網(wǎng)格).

由表4可知，在排水航行階段，笛卡爾網(wǎng)格計算結(jié)果的不確定度均小于15%SG，經(jīng)過修正后的不確定度則降低至不超過6%SG，而在過渡航行階段，其不確定度則高達(dá)90%SG(阻力系數(shù))，修正后的不確定度則依然高達(dá)77%SG，說明在該階段，笛卡爾網(wǎng)格的計算結(jié)果存在極大地不確定性；而在該階段對于重疊網(wǎng)格而言，其不確定度均小于6%SG，經(jīng)過修正的不確定度則低至3%SG以下.同時經(jīng)過修正的不確定度均得到有效降低，說明該數(shù)值模擬的結(jié)果存在漸進(jìn)值，并且沒有出現(xiàn)橫穿漸進(jìn)值的現(xiàn)象，那么可以排除震蕩收斂的可能性，如果結(jié)果收斂，則基本可以斷定該計算是單調(diào)收斂的.

因此從數(shù)值模擬的校核的角度而言，在排水航行階段，笛卡爾網(wǎng)格計算的不確定度基本屬于可以接受的范圍；在過渡航行階段，笛卡爾網(wǎng)格則不可接受，而重疊網(wǎng)格則是可以接受的.

3.1.2網(wǎng)格收斂性分析

圖6為笛卡爾網(wǎng)格在V=1.26,1.68,3.78 m/s時，以及重疊網(wǎng)格在V=3.78 m/s時阻力系數(shù)和縱傾角的模擬結(jié)果隨網(wǎng)格因子R的收斂情況.網(wǎng)格因子為特征網(wǎng)格尺寸之比，即網(wǎng)格數(shù)量之比的三次根.

(20)

圖6 模擬結(jié)果圖

式中：Ni為第i種網(wǎng)格數(shù)量.

由上述最小二乘法擬合的結(jié)果，可以對計算結(jié)果的收斂速度與網(wǎng)格收斂速度進(jìn)行比較，當(dāng)收斂階數(shù)大于2時，說明該結(jié)果的收斂速度大于網(wǎng)格的收斂速度，則可以考慮進(jìn)一步增加網(wǎng)格數(shù)量，以提高數(shù)值模擬結(jié)果的可靠性；反之，則說明其收斂速度小于網(wǎng)格的收斂速度，如果結(jié)果的不確定度已經(jīng)在可接受范圍內(nèi)，則沒有必要再通過增加網(wǎng)格數(shù)量來優(yōu)化結(jié)果的不確定度.由上述可知，除了在V=1.26 m/s時的縱傾角外，數(shù)值模擬結(jié)果的收斂速度均小于2，鑒于不確定度的分析結(jié)果，即可判斷該數(shù)值模擬的網(wǎng)格數(shù)量是可以接受的.綜上所述，在數(shù)值模擬階段的校核情況比較理想，需要對結(jié)果作進(jìn)一步的驗(yàn)證.

3.2 驗(yàn)證

表5為采用笛卡爾網(wǎng)格時，在V=1.26，1.68，3.78 m/s時，以及重疊網(wǎng)格在V=3.78 m/s時阻力系數(shù)與縱傾角的實(shí)驗(yàn)驗(yàn)證.其中，驗(yàn)證不確定度Uval中基于保守估計，忽略了USPD和USMA，而實(shí)驗(yàn)誤差不確定度以試驗(yàn)結(jié)果的百分比給出為

UD=2.5%D

(21)

式中：D為實(shí)驗(yàn)數(shù)據(jù).

已知驗(yàn)證條件為比較誤差E應(yīng)小于驗(yàn)證不確定度Uv，由上述計算，可知笛卡爾網(wǎng)格在排水航行階段時，阻力系數(shù)和縱傾角均得到驗(yàn)證(achieved validation)，而在過渡階段，雖然比較誤差小于驗(yàn)證不確定度，但是不難發(fā)現(xiàn)，此時驗(yàn)證不確定度過大，說明網(wǎng)格需要進(jìn)一步優(yōu)化；重疊網(wǎng)格的計算可見，對于網(wǎng)格(1-3)，經(jīng)過修正后結(jié)果均得到驗(yàn)證，而網(wǎng)格(2-4)，有未得到驗(yàn)證的結(jié)果，說明網(wǎng)格(1)比(2)更具優(yōu)勢，在過渡航行階段，推薦使用重疊網(wǎng)格(1)進(jìn)行模擬.

表5 試驗(yàn)驗(yàn)證誤差計算結(jié)果 %

注：①V=1.26 m/s,Fr▽=0.644(笛卡爾網(wǎng)格1-3)；②V=1.26 m/s,Fr▽=0.644(笛卡爾網(wǎng)格2-4)；③V=1.68 m/s,Fr▽=0.860(笛卡爾網(wǎng)格1-3)；④V=1.68 m/s,Fr▽=0.860(笛卡爾網(wǎng)格2-4)；⑤V=3.78 m/s,Fr▽=1.934(笛卡爾網(wǎng)格1-3)；⑥V=3.78 m/s,Fr▽=1.934(笛卡爾網(wǎng)格2-4)；⑦V=3.78 m/s,Fr▽=1.934(重疊網(wǎng)絡(luò)1-3)；⑧V=3.78 m/s,Fr▽=1.934(重疊網(wǎng)絡(luò)2-4).

4 結(jié) 論

1) 在排水航行階段，采用笛卡爾網(wǎng)格時，對該艇的阻力、縱傾和流場信息的表達(dá)和模擬，基本滿足需求，并隨著網(wǎng)格數(shù)量增加，流場細(xì)節(jié)的表現(xiàn)越清晰；在過渡航行階段，由于該艇的船型特點(diǎn)，其航態(tài)發(fā)生較大變化，此時笛卡爾網(wǎng)格已不能準(zhǔn)確模擬，而重疊網(wǎng)格對變航態(tài)問題的模擬和流場細(xì)節(jié)的捕捉，均可得到較好的結(jié)果.

2) 對計算結(jié)果的不確定度分析可見，兩種網(wǎng)格均能得到收斂的結(jié)果.但是在過渡航行階段，笛卡爾網(wǎng)格不確定度過大，不具實(shí)用價值；此外，經(jīng)過修正的不確定度均得到優(yōu)化，進(jìn)一步說明該計算滿足單調(diào)收斂；采用最小二乘法對收斂速度進(jìn)行了分析，可知不再需要增加網(wǎng)格來優(yōu)化不確定度了，說明網(wǎng)格數(shù)量已是較優(yōu)方案.

3) 對數(shù)值模擬的結(jié)果進(jìn)行試驗(yàn)驗(yàn)證可知，笛卡爾網(wǎng)格在排水航行階段的結(jié)果通過了驗(yàn)證；而在過渡航行階段，笛卡爾網(wǎng)格的結(jié)果雖然滿足驗(yàn)證條件，但是不確定度不可接受，故只有重疊網(wǎng)格通過了驗(yàn)證.同時，重疊網(wǎng)格(1)的結(jié)果最為理想，推薦采用.

[1] STERN F, WILSON R V, COLEMAN H W, et al. Verification and validation of cfd simulations [J]. Tools on Martian Rover Software Formal Methods in System Design,1999,57(3):2032-2036.

[2] STERN F, WILSON R V, COLEMAN H W, et al. Comprehensive, approach to verification and validation of CFD simulations - part I: methodology and procedures[J]. Journal of Fluids Engineering,2001,123(4):792-799.

[3] WILSON R V, STERN F, COLEMAN H W, et al. Comprehensive approach to verification and validation of CFD simulations—part 2: application for rans simulation of a cargo/container ship[J]. Journal of Fluids Engineering,2001,123(4):803-810.

[4] ECA L, VAZ G, HOEKSTRA M. A Verification and validation exercise for the flow over a backward facing step[C].Eccomas Cfd，2010.

[5] JIN Y, CHAI S, DUFFY J, et al. Urans prediction of ship hydrodynamics in head sea waves at zero forward speed with model testing validation [C].International Conference on Ocean, Offshore and Arctic Engineering，2016.

[6] ZHU R C，YANG C L，MIAO G P, et al. Computational fluid dynamics uncertainty analysis for simulations of roll motions for a 3D ship[J]. Journal of Shanghai Jiaotong University (Science),2015,20(5):591-599.

[7] CLARK C G, LYONS D G, NEU W L. Comparison of single and overset grid techniques for CFD simulations of a surface effect ship[C].ASME 2014 33rd International Conference on Ocean, Offshore and Arctic Engineering,American Society of Mechanical Engineers,2014.

[8] TAKAI T, KANDASAMY M, STERN F. Verification and validation study of URANS simulations for an axial waterjet propelled large high-speed ship[J]. Journal of Marine Science and Technology,2011,16(4):434-447.

[9] ?MER F S, G?K?E M K, ?MER K K. Overset grid system for large ship motions inside fluid flow[C]. Thirteenth International Conference on Marine Sciences and Technologies, Blacksea,2016.

[10] 湛俊華.基于OpenFOAM數(shù)值波浪水池的過渡型艇水動力性能研究[D].武漢：武漢理工大學(xué)，2015.

[11] 鄒勁,姬朋輝,孫寒冰,等.網(wǎng)格因素對三體滑行艇阻力計算影響探究[J].船舶,2016,27(3):8-14.

[12] 張玲.過渡型艇微氣泡減阻噴氣參數(shù)數(shù)值研究[D].武漢：武漢理工大學(xué),2009.

[13] 張玲,王家楣.離散格式和網(wǎng)格數(shù)量對過渡型艇粘性阻力的影響[J].船海工程,2009(4):9-12.

Verification and Validation for the Numerical Simulation of Semi-displacement Craft

YUZeshuang1)MAOXiaofei1,2)

(SchoolofTransportation,WuhanUniversityofTechnology,Wuhan430063,China)1)(KeyLaboratoryofHighPerformanceShipTechnologyofMinistryofEducation,WuhanUniversityofTechnology,Wuhan430063,China)2)

For an accurate simulation, four different kinds of grids are formed by changing the grid of hull surface and refinement areas based on the structured cartesian grid and overset grid technique in viscous CFD software STAR-CCM+. The numerical simulation from the navigation phases of placement to transition (Fr▽=0.664,0.860,1.934) in the static water is carried out and the uncertainty analysis for mesh of the resistance and trim angle is conducted. Two different approaches to evaluate the convergence of the results with the mesh are used for the verification and the experiment results are adopted to verify the numerical results. From the results of verification and validation, the cartesian grid can meet the requirements to a certain extent in the navigation phase of placement, while in the navigation phases of transition, the overset grid has very obvious advantages.

STAR-CCM+; overset grid; uncertainty analysis; verification and validation

U661.32

10.3963/j.issn.2095-3844.2017.06.022

2017-09-30

余澤爽(1992—)：男，碩士生，主要研究領(lǐng)域?yàn)榇靶阅?/p>