郜冶,譚大力,李海旭,王金玲,劉長(zhǎng)猛
(1.哈爾濱工程大學(xué)航天與建筑工程學(xué)院,黑龍江哈爾濱150001;2.海軍裝備研究院,北京100192;3.中國(guó)船舶工業(yè)系統(tǒng)工程研究院,北京100094)
大型水面艦船的機(jī)動(dòng)性對(duì)軍事戰(zhàn)術(shù)的實(shí)施尤其對(duì)艦載機(jī)的起飛降落過(guò)程起著決定性作用,這就要求艦船在航行過(guò)程中必須確定相對(duì)于優(yōu)勢(shì)風(fēng)的正確航向和特定的前進(jìn)速度,確保艦載機(jī)安全著艦[1]。目前大型艦船主要依靠安裝在船橋中央桅桿上的測(cè)風(fēng)儀器測(cè)量風(fēng)速,風(fēng)速儀本身的測(cè)量精度符合要求,但由于受艦船甲板及其艦上島型建筑周圍渦流的影響安裝在渦流區(qū)的風(fēng)速儀測(cè)得的風(fēng)場(chǎng)與實(shí)際真風(fēng)存在很大的誤差。在國(guó)外,M.L.Thiebaux利用風(fēng)洞試驗(yàn)研究了艦載風(fēng)速儀的誤差校正問(wèn)題[2];P.K.Taylor和 E.C.Kent等針對(duì)商用船,利用實(shí)驗(yàn)和CFD方法研究了受船體影響其周圍風(fēng)場(chǎng)產(chǎn)生的扭曲和變形對(duì)風(fēng)速儀測(cè)量精度的影響,并與衛(wèi)星測(cè)報(bào)數(shù)據(jù)進(jìn)行了對(duì)比[3-4];Philippe L.Nacass采用 CFD計(jì)算船載風(fēng)速儀周圍的氣流失真現(xiàn)象[5]。在國(guó)內(nèi),劉連吉等用實(shí)驗(yàn)方法研究了船上風(fēng)速儀測(cè)量風(fēng)速與浮標(biāo)測(cè)量風(fēng)速的誤差,并提出了測(cè)量風(fēng)速的改進(jìn)意見(jiàn)[6-7];劉慧等將船測(cè)風(fēng)場(chǎng)資料與QuikSCAT衛(wèi)星遙感風(fēng)場(chǎng)資料進(jìn)行對(duì)比分析研究了船測(cè)風(fēng)場(chǎng)偏差統(tǒng)計(jì)特征[8]。近年來(lái)郜冶等致力于艦船流場(chǎng)的CFD計(jì)算,文獻(xiàn)[9]中運(yùn)用CFD計(jì)算空氣正向流過(guò)滑躍式艦船產(chǎn)生的甲板氣流場(chǎng),討論分析不同網(wǎng)格劃分形式對(duì)計(jì)算結(jié)果的影響,通過(guò)對(duì)全尺寸的滑躍式甲板艦船的數(shù)值計(jì)算分析得到艦船甲板氣流場(chǎng)渦旋結(jié)構(gòu)特征,并將不同湍流模型應(yīng)用于流場(chǎng)計(jì)算,分析了艦船流場(chǎng)的特點(diǎn)。文獻(xiàn)[10]利用FLUENT軟件UDF接口將LK和MMK模型引入FLUENT進(jìn)行流場(chǎng)計(jì)算改進(jìn)了計(jì)算結(jié)果。
本文采用數(shù)值模擬方法對(duì)不同工況下風(fēng)速儀周圍的氣流場(chǎng)進(jìn)行計(jì)算分析,得出風(fēng)速儀更加合適的安裝位置。如計(jì)算條件允許,理論上可以研究流體在任何條件下的運(yùn)動(dòng),再現(xiàn)風(fēng)速儀周圍的風(fēng)場(chǎng)。
為便于計(jì)算,將風(fēng)速儀安裝位置處物理模型(包括梯形柱臺(tái)、雷達(dá)、天線和其他附屬結(jié)構(gòu))進(jìn)行合理簡(jiǎn)化,忽略掉不影響流場(chǎng)特征的細(xì)微結(jié)構(gòu)如毫米級(jí)的天線、柱臺(tái)上的小孔等,簡(jiǎn)化后的模型如圖1所示。
圖1 計(jì)算模型Fig.1 The calculation model
計(jì)算區(qū)域X軸表示船長(zhǎng)方向,由船頭指向船尾為正方向,Y軸表示船寬方向,由左舷指向右舷為正方向,Z軸表示豎直方向,由下方指向上方為正方向,模型關(guān)于中心位置左右對(duì)稱。坐標(biāo)原點(diǎn)位于梯形柱體底部中心處,現(xiàn)有方案風(fēng)速儀安裝在如圖2中所示監(jiān)測(cè)點(diǎn)位置,其坐標(biāo)分別為D(0.0 m,-3.68 m,13.23 m)、D'(0.0 m,3.68 m,13.23 m)。
圖2 風(fēng)速儀監(jiān)測(cè)點(diǎn)位置俯視圖Fig.2 The location of the anemometers in the top viewport
艦船及其風(fēng)速儀周圍的流場(chǎng)屬于三維復(fù)雜湍流場(chǎng),為簡(jiǎn)化計(jì)算設(shè)空氣為不可壓縮流體且符合Boussinesq假設(shè),流動(dòng)為穩(wěn)態(tài)紊流,忽略固體壁面間的熱輻射。采用標(biāo)準(zhǔn)k-ε兩方程模型對(duì)氣流場(chǎng)進(jìn)行數(shù)值模擬。建立流場(chǎng)穩(wěn)定后的湍流控制通用方程[11-12]:
式中:div(ρVφ)表示對(duì)流項(xiàng),div(Γgrad φ)表示擴(kuò)散頂,S表示源項(xiàng)??刂品匠贪ㄟB續(xù)方程如式(2)、動(dòng)量方程如式(3)、能量方程如式(4)、k方程如式(5)和ε方程如式(6)。
式中:Gk是由平均速度梯度引起的湍動(dòng)能k的產(chǎn)生項(xiàng),μ是湍動(dòng)粘度,可表示t成湍動(dòng)能k和耗散率 ε 的函數(shù),即 μt=ρCμk2/ε;C1ε、C2ε為經(jīng)驗(yàn)常數(shù);σk、σε分別是湍動(dòng)能k和耗散率ε對(duì)應(yīng)的Prandtl數(shù)。在數(shù)值模擬過(guò)程中,各系數(shù)取為Cμ=0.09,C1ε=1.44,C2ε=1.92,σk=1.0,σε=1.3,σT=0.9。
流場(chǎng)計(jì)算區(qū)域取風(fēng)速儀前為10倍柱臺(tái)底面邊長(zhǎng)(柱臺(tái)底面邊長(zhǎng)為5.75 m),風(fēng)速儀后為20倍柱臺(tái)底面邊長(zhǎng),左右均為5倍柱臺(tái)底面邊長(zhǎng),垂向?yàn)?倍柱臺(tái)底面邊長(zhǎng),基于柱臺(tái)底面邊長(zhǎng)為特征尺度的雷諾數(shù)量級(jí)為1.0×107。由于計(jì)算模型比較復(fù)雜采用非結(jié)構(gòu)網(wǎng)格,對(duì)梯形柱臺(tái)、雷達(dá)、天線及附屬結(jié)構(gòu)區(qū)域進(jìn)行網(wǎng)格加密。全部網(wǎng)格由ANSYS ICEM生成,網(wǎng)格總數(shù)為8.77×106,計(jì)算區(qū)域及全流場(chǎng)網(wǎng)格如圖3所示,局部網(wǎng)格加密如圖4所示。
圖3 計(jì)算區(qū)域及全流場(chǎng)整體網(wǎng)格Fig.3 The calculation region and grids for the whole flow field
圖4 局部網(wǎng)格加密Fig.4 The refinement of local grid
計(jì)算中入口設(shè)為速度入口(具體速度方向大小隨不同計(jì)算工況改變),梯形柱體、雷達(dá)、天線及附屬結(jié)構(gòu)所有表面均設(shè)定為無(wú)滑移壁面,其余邊界均設(shè)為自由滑移壁面,出口設(shè)為壓力出口。
采用文獻(xiàn)[11]中的方法驗(yàn)證CFD計(jì)算方法的準(zhǔn)確性,將相同計(jì)算方法及網(wǎng)格劃分方法的美國(guó)LHA級(jí)艦船的CFD計(jì)算結(jié)果與文獻(xiàn)[15]中的風(fēng)洞實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比分析。LHA模型為1/48縮比模型,全長(zhǎng)為 5.2 m,甲板寬為 0.75 m,網(wǎng)格數(shù)量約為 2.0×106,U0=6.858 m/s,使用標(biāo)準(zhǔn)k-ε 模型進(jìn)行計(jì)算。驗(yàn)證計(jì)算中,著艦點(diǎn)位置如圖5中所示。
圖5LHA著艦點(diǎn)的位置(m)Fig.5 The touchdown points of LHA(m)
圖6為CFD計(jì)算結(jié)果與文獻(xiàn)[15]中風(fēng)洞試驗(yàn)V-22旋翼通過(guò)艦點(diǎn)2的傾斜軌道上垂向速度對(duì)比。橫坐標(biāo)為坐標(biāo)值與甲板寬度d的比值,縱坐標(biāo)為Z向速度W與入口速度U0的比值。圖6中數(shù)據(jù)對(duì)比顯示CFD計(jì)算結(jié)果與風(fēng)洞試驗(yàn)數(shù)據(jù)吻合較好,誤差不超過(guò)7%。證明論文中運(yùn)用的網(wǎng)格劃分方法及標(biāo)準(zhǔn)k-ε模型可以用于預(yù)測(cè)船體繞流場(chǎng)的特征。
圖6 著艦點(diǎn)2傾斜軌道上垂向速度對(duì)比Fig.6 The vertical velocity of inclined orbit comparison at touchdown point 2
為討論不同船速、風(fēng)速及不同風(fēng)向角對(duì)風(fēng)速儀測(cè)量位置處流場(chǎng)的影響,對(duì)現(xiàn)有安裝方案設(shè)置3種工況進(jìn)行對(duì)比,當(dāng)船速3 m/s,風(fēng)速3 m/s(低航速低風(fēng)速)時(shí)記為工況A;當(dāng)船速14 m/s,風(fēng)速3 m/s(高航速低風(fēng)速)時(shí)記為工況 B;當(dāng)船速14 m/s,風(fēng)速14 m/s(高航速高風(fēng)速)時(shí)記為工況C。并對(duì)每種工況設(shè)置4種不同的風(fēng)向角進(jìn)行計(jì)算,圖3中給出了風(fēng)向角示意表示,從船體正前方吹風(fēng)時(shí)設(shè)為是0°風(fēng),從船體左舷吹風(fēng)時(shí)設(shè)是左舷風(fēng),由于模型關(guān)于中心位置左右對(duì)稱,因此只計(jì)算左舷風(fēng)。計(jì)算時(shí)風(fēng)向角分別選取 0°、10°、20°、30°。將每種工況左右監(jiān)測(cè)點(diǎn)D、D'處數(shù)據(jù)進(jìn)行提取,分離出此處風(fēng)向角α'及風(fēng)速υ',并與無(wú)窮遠(yuǎn)處真風(fēng)風(fēng)向角α及風(fēng)速υ進(jìn)行對(duì)比,得到風(fēng)向角差值 Δα和風(fēng)速差值的模,其計(jì)算方法如圖7所示?,F(xiàn)定義風(fēng)向角差值Δα方向:與實(shí)際真風(fēng)相比,監(jiān)測(cè)點(diǎn)真風(fēng)逆時(shí)針偏轉(zhuǎn)為正,順時(shí)針偏轉(zhuǎn)為負(fù)。
圖7 風(fēng)向角差值和風(fēng)速差值計(jì)算示意圖Fig.7 The calculation schematic diagram for wind angle and vector difference
圖8和圖9分別為工況A、B、C在不同風(fēng)向角α?xí)rα'及e(υ')計(jì)算結(jié)果,其中e(υ')表示風(fēng)速值誤差。風(fēng)速儀在現(xiàn)有安裝位置得到的α'及e(υ')均較大,其中工況A與工況C誤差相差不大,而工況B中α'及e(υ')遠(yuǎn)遠(yuǎn)高于工況A和工況C,即低風(fēng)速高船速時(shí)風(fēng)速儀的測(cè)量誤差急劇升高;隨α值增大,各工況近風(fēng)側(cè)D測(cè)點(diǎn)(左舷測(cè)點(diǎn))α'誤差逐漸減小,遠(yuǎn)風(fēng)側(cè)D'測(cè)點(diǎn)(右舷測(cè)點(diǎn))α'誤差逐漸增大,且遠(yuǎn)風(fēng)側(cè)誤差明顯小于近風(fēng)側(cè)誤差;但隨風(fēng)向角α增大,工況B及工況C中風(fēng)速值誤差e(υ')減小,工況A在α=30°時(shí)風(fēng)速值誤差e(υ')明顯增大,且遠(yuǎn)風(fēng)側(cè)風(fēng)速值誤差始終略大于近風(fēng)側(cè)誤差。
圖8 3種工況不同位置處α'計(jì)算結(jié)果Fig.8 The calculation result α'for the 3 conditions
圖9 3種工況不同位置風(fēng)速值誤差計(jì)算結(jié)果Fig.9 The wind speed error for the 3 conditions
鑒于風(fēng)速儀在現(xiàn)有安裝位置風(fēng)向角及風(fēng)速值測(cè)量誤差均較大,特設(shè)置風(fēng)速儀不同安裝位置進(jìn)行計(jì)算,試圖找到更合適的安裝位置以減小測(cè)量誤差。
不改變風(fēng)速儀安裝高度,分別令x為0 m、-4.6 m 時(shí),y為±3、±4、±5、±6、±7,±8 m,即將風(fēng)速儀安裝位置在x=0 m與z=13.23 m及x=-4.6 m與z=13.23 m交線上分別由內(nèi)向外進(jìn)行平移。風(fēng)速儀監(jiān)測(cè)點(diǎn)分布線如圖2中所示。
3種工況都分別對(duì)α為0°和30°時(shí),風(fēng)速儀安裝在不同測(cè)點(diǎn)時(shí)的流場(chǎng)風(fēng)速進(jìn)行計(jì)算。圖10和11分別為3種工況在不同風(fēng)向角α值時(shí)不同測(cè)點(diǎn)的α'計(jì)算結(jié)果,圖12和13分別為3種工況在不同風(fēng)向角α值時(shí)不同測(cè)點(diǎn)的風(fēng)速值誤差e(υ')。
圖10 α為0°時(shí)3種工況不同位置處α'計(jì)算結(jié)果Fig.10 α'calculation result for the 3 conditions in 0°
圖11 α為30°時(shí)3種工況不同位置處α'計(jì)算結(jié)果Fig.11 α'calculation result for the 3 conditions in 30°
圖12 α為0°時(shí)3種工況不同位置風(fēng)速值誤差計(jì)算結(jié)果Fig.12 Wind speed error for the 3 conditions in 0°
圖13 α為30°時(shí)3種工況不同位置風(fēng)速值誤差計(jì)算結(jié)果Fig.13 Wind speed error for the 3 conditions in 30°
計(jì)算結(jié)果顯示隨著風(fēng)速儀安裝位置由內(nèi)向外移動(dòng),風(fēng)向角α'及風(fēng)速值誤差e(υ')均隨y值絕對(duì)值(風(fēng)速儀監(jiān)測(cè)點(diǎn)距中心的距離)增大而減小,即將風(fēng)速儀安裝位置向外移動(dòng)有助于減小測(cè)量誤差。在遠(yuǎn)風(fēng)側(cè),3種工況x=-4.6 m處的α'誤差普遍略大于x=0處,且隨α增大2位置處的α'差距加大。在近風(fēng)側(cè),工況A和C在x=-4.6 m處的α'誤差小于x=0處,且隨風(fēng)向角α增大兩位置處的α'差距加大;而工況B在x=-4.6 m處的風(fēng)向角α'誤差略大于x=0 m處,且隨風(fēng)向角α增大兩位置處的α'差距減小。對(duì)于風(fēng)速值,3種工況在x=-4.6m處誤差均明顯小于x=0處,且隨風(fēng)向角增大兩位置處風(fēng)速值誤差差距加大。即將風(fēng)速儀安裝位置向前移到天線位置處可以明顯減小風(fēng)速值測(cè)量誤差且風(fēng)向角測(cè)量誤差不會(huì)顯著增大。
以上計(jì)算中僅考慮風(fēng)速儀安裝位置的物理模型,忽略了船體結(jié)構(gòu)的影響。為討論船體結(jié)構(gòu)對(duì)風(fēng)速儀測(cè)量位置處流場(chǎng)結(jié)構(gòu)的影響,將船體模型和風(fēng)速儀模型合并,如圖14所示。當(dāng)船速為0,風(fēng)速為5 m/s時(shí)為工況D,當(dāng)船速為0,風(fēng)速為20 m/s時(shí)為工況E,且分別對(duì)風(fēng)速儀模型及合并模型進(jìn)行計(jì)算。工況D對(duì)風(fēng)向角為0°,工況E對(duì)風(fēng)向角為0°和30°時(shí),風(fēng)速儀不同測(cè)點(diǎn)的流場(chǎng)進(jìn)行計(jì)算。
圖14 船體和風(fēng)速儀合并模型Fig.14 The merge model of ship and anemoscope
圖15~20分別為工況D、E風(fēng)向角及風(fēng)速值誤差計(jì)算結(jié)果。計(jì)算結(jié)果顯示2種工況風(fēng)速儀單獨(dú)模型與合并模型測(cè)量誤差隨y值的變化趨勢(shì)是相同的,隨風(fēng)速儀安裝位置向外延伸兩模型的風(fēng)向角及風(fēng)速值測(cè)量誤差差距均減小。
圖15 α為0°時(shí)工況D風(fēng)向角差值α'計(jì)算結(jié)果Fig.15 α'calculation results for condition D in 0°
圖16 α為0°時(shí)工況E風(fēng)向角差值α'計(jì)算結(jié)果Fig.16 α'calculation results for condition E in 0°
圖17 α為30°時(shí)工況E風(fēng)向角差值α'計(jì)算結(jié)果Fig.17 α'calculation results for condition E in 30°
圖18 α為0°時(shí)工況D風(fēng)速值誤差計(jì)算結(jié)果Fig.18 e(υ')calculation results for condition D in 0°
圖19 α為0°時(shí)工況E風(fēng)速值誤差計(jì)算結(jié)果Fig.19 e(υ')calculation results for condition E in 0°
圖20 α為30°時(shí)工況E風(fēng)速值誤差計(jì)算結(jié)果Fig.20 e(υ')calculation results for condition E in 30°
α=0°時(shí),距離中心位置3 m以外,在左右舷對(duì)稱位置2種模型測(cè)得的風(fēng)向角α'數(shù)值基本相等但方向相反,且2種模型測(cè)量結(jié)果的差值基本維持在2°以內(nèi);對(duì)于風(fēng)速值測(cè)量誤差,在距離中心位置1~4 m,風(fēng)速儀單獨(dú)模型測(cè)量結(jié)果明顯大于合并模型的測(cè)量結(jié)果,在距離中心位置4 m以外,合并模型測(cè)量的風(fēng)速值誤差略微大于風(fēng)速儀單獨(dú)模型的測(cè)量結(jié)果但兩模型的測(cè)量誤差在3%以內(nèi)。
α=30°時(shí),遠(yuǎn)風(fēng)側(cè)的風(fēng)向角及風(fēng)速值誤差測(cè)量結(jié)果均大于近風(fēng)側(cè)的測(cè)量結(jié)果,近風(fēng)側(cè)在4 m以外兩模型風(fēng)向角測(cè)量值的差值在5°以內(nèi),風(fēng)速值誤差的差值在2%以內(nèi),而在遠(yuǎn)風(fēng)側(cè),2種模型的測(cè)量誤差差距相對(duì)較大。
1)風(fēng)速儀在現(xiàn)有安裝位置測(cè)量誤差較大,其中低風(fēng)速高船速時(shí)風(fēng)速儀的測(cè)量誤差急劇升高。
2)隨吹風(fēng)角度增大,近風(fēng)側(cè)風(fēng)向角測(cè)量誤差逐漸減小,遠(yuǎn)風(fēng)側(cè)風(fēng)向角測(cè)量誤差逐漸增大,且遠(yuǎn)風(fēng)側(cè)風(fēng)向角測(cè)量誤差明顯小于近風(fēng)側(cè)風(fēng)向角測(cè)量誤差。
3)隨吹風(fēng)角度增大,風(fēng)速值測(cè)量誤差變化不大,且遠(yuǎn)風(fēng)側(cè)風(fēng)速值測(cè)量誤差始終略大于近風(fēng)側(cè)風(fēng)速值測(cè)量誤差。
4)將風(fēng)速儀安裝位置向外平移有利于減小風(fēng)向角和風(fēng)速值測(cè)量誤差;將風(fēng)速儀安裝位置向前移到天線位置處可以明顯減小風(fēng)速值測(cè)量誤差且風(fēng)向角測(cè)量誤差不會(huì)顯著增大。
5)船體本身會(huì)對(duì)風(fēng)速儀附近的流場(chǎng)產(chǎn)生較大的影響。在此區(qū)域內(nèi),合并模型風(fēng)速儀測(cè)量精度明顯較風(fēng)速儀單獨(dú)模型低,隨風(fēng)速儀安裝位置外移,2種模型的測(cè)量差距減小,尤其在近風(fēng)側(cè)2種模型的測(cè)結(jié)果基本相同。此時(shí)可以使用風(fēng)速儀單獨(dú)模型代替合并模型進(jìn)行定性分析,但是定量分析時(shí)應(yīng)該使用合并模型。
[1]宋劍,何建忠.航空母艦在隨機(jī)海況下的運(yùn)動(dòng)統(tǒng)計(jì)特性[J].艦船電子工程,2011,31(3):67-69.SONG Jian,HE Jianzhong.Statistical properties of the movement of aircraft carrier in the random sea conditions[J].Ship Electronic Engineering,2011,31(3):67-69.
[2]THIEBAUX M L.Wind tunnel experiments to determine correction functions for shipborne anemometers[J].Hydrog and Ocean Sci,1990,36(2):57-61.
[3]TAYLOR P K,KENT E C,YELLAND M J,et al.The accuracy of wind observations from ships[J].Advances in the Applications of Marine Climatology,1995,2(2):132-155.
[4]TAYLOR P K,KENT E C.The accuracy of meteorological observations from voluntary observing ships-present status and future requirements[C]//WMO,F(xiàn)irst Session of the CMM Subgroup on Voluntary Observing Ships.Athen,Geneva,1999:12.
[5]PHILIPPE L N.Shipborne wind measurements corrected for airflow distortion by computational fluid dynamics[R]//Brétigny-sur-Orge:Centred'Aviation Météorologique,2001:228-233.
[6]劉連吉,趙永平.船上氣溫、濕度和風(fēng)速觀測(cè)的代表性[J].山東海洋學(xué)院學(xué)報(bào),1981,11(3):24-31.LIU Lianji,ZHAO Yongping.The reliability of the wind speed air temperature and humidity observed on ship[J].Journal of Shandong College of Oceanology,1981,11(3):24-31.
[7]劉連吉.對(duì)現(xiàn)用海洋水文氣象要素調(diào)查規(guī)范的改進(jìn)意見(jiàn)[J].海洋技術(shù),1982,1(2):70-73.LIU Lianji.Suggestions on the instruction manual for obtaining oceanographic and meteorological data[J].Acta Oceanologica Sinica,1982,1(2):70-73.
[8]劉慧,胡松,鄒曉榮.船測(cè)資料與智利外海QuikSCAT風(fēng)場(chǎng)比較分析[J].遙感技術(shù)與應(yīng)用,2012,27(5):763-769.LIU Hui,HU Song,ZOU Xiaorong.Comparison between the fishing vessels and QuikSCAT scatterometer wind data of the offshore Chile[J].Remote Sensing Technology and Application,2012,27(5):763-769.
[9]劉長(zhǎng)猛,郜冶.滑躍式甲板氣流場(chǎng)數(shù)值模擬[J].華中科技大學(xué)學(xué)報(bào),2012,40(10):68-71.LIU Changmeng,GAO Ye.Numerical simulation of airwake on ski jump deck[J].Journal of Huazhong University of Science and Technology,2012,40(10):68-71.
[10]郜冶,劉長(zhǎng)猛.護(hù)衛(wèi)艦氣流場(chǎng)數(shù)值計(jì)算研究[J].哈爾濱工程大學(xué)學(xué)報(bào),2013,34(1):1-5.GAO Ye,LIU Changmeng.Numerical simulation of frigate ship airwake[J].Journal of Harbin Engineering University,2013,34(1):1-5.
[11]LAUNDER B E,SPALDING D B.The numerical computation of turbulent flows[J].Comput Methods Appl Mech Eng,1974,80(3):269-289.
[12]COMINI G,GIUDICE S A.(k-ε )model of turbulent flow[J].Numer Heat transfer,1985,90(8):299-316.
[13]RAJAGOPALAN G,SCHALLER D,WADCOSK A J,et al.Experimental and computational simulation of a model ship in a wind tunnel[C]//43rd AIAA Aerospace Sciences Meeting.[S.l.],2005.