水翼三維非穩(wěn)態(tài)空化流動(dòng)特性數(shù)值研究

鄧 飛， 張 松， 梁 興， 李偉平

（1.南昌工程學(xué)院南昌市低品位能源轉(zhuǎn)化與系統(tǒng)節(jié)能技術(shù)重點(diǎn)實(shí)驗(yàn)室，南昌 330099；2.渝能（集團(tuán)）有限責(zé)任公司，重慶 401121）

0 引言

水翼空化具有強(qiáng)烈的非定常特性，易引起流場(chǎng)空泡團(tuán)產(chǎn)生，形成周期性脈動(dòng)造成水翼水動(dòng)力載荷變化，從而出現(xiàn)水翼振動(dòng)和流場(chǎng)噪聲等問(wèn)題［1］，影響流體機(jī)械的安全可靠運(yùn)行。一直以來(lái)，研究空化發(fā)生的機(jī)理、特性和抑制技術(shù)都是流體機(jī)械研究的重點(diǎn)和難點(diǎn)［2］。研究表明，水翼安放角度［3］、雷諾數(shù)、繞流速度［4］等水力因素決定了水翼的繞流特征。謝春梅等［5］針對(duì)高雷諾數(shù)下三維水翼的梢渦流場(chǎng)，采用基于浸沒(méi)邊界法的壁模型大渦模擬方法，對(duì)NACA16020 水翼進(jìn)行了數(shù)值研究，系統(tǒng)觀察了流場(chǎng)的各類統(tǒng)計(jì)特性。Wei 等［6］對(duì)翼型在6°攻角下的水動(dòng)力及空化性能進(jìn)行數(shù)值研究，當(dāng)空化發(fā)生時(shí)水動(dòng)力特性將劇烈惡化。劉胡濤等［7］對(duì)高雷諾數(shù)的NACA0012 水翼繞流開(kāi)展了數(shù)值計(jì)算，探討了初始攻角條件、繞流速度以及剛心位置三個(gè)因素對(duì)水翼振動(dòng)的影響，水翼振動(dòng)狀態(tài)依賴于系統(tǒng)的初始條件，當(dāng)繞流速度不斷增大時(shí)水翼將發(fā)生顫振現(xiàn)象。田北晨等［8］基于建立的多尺度算法對(duì)繞NACA66 水翼空化流動(dòng)進(jìn)行模擬，在回射流發(fā)展階段，離散空泡逐漸增加并分布在回射流擾動(dòng)區(qū)。Yang等［9］采用MRT-LBM方法對(duì)空泡潰滅進(jìn)行熱力學(xué)分析并驗(yàn)證了MRT-LBM空泡潰滅的有效性。

目前，針對(duì)水翼空化現(xiàn)象的機(jī)理研究主要集中在空泡初生階段［10-11］，對(duì)不同繞流速度下的水翼空化特性及其準(zhǔn)周期性演變過(guò)程研究較少。本文針對(duì)不同繞流速度下的NACA6408 水翼非穩(wěn)態(tài)繞流過(guò)程，數(shù)值研究各流動(dòng)因素的變化規(guī)律以及其對(duì)水翼空化過(guò)程的影響，旨在獲取水翼發(fā)生空化的臨界繞流速度和壓力等重要參數(shù)，以及空化過(guò)程的準(zhǔn)周期性特征等。

1 數(shù)值模型與網(wǎng)格劃分

本文采用計(jì)算流體動(dòng)力學(xué)軟件Fluent 16.0 進(jìn)行數(shù)值模擬研究，如圖1 所示為NACA6408 三維水翼計(jì)算域模型，依據(jù)該水翼的相對(duì)彎度和相對(duì)厚度，對(duì)其及流體域進(jìn)行建模。整個(gè)流體域長(zhǎng)度為10b，其中b為水翼弦長(zhǎng)，進(jìn)口距水翼前緣距離為3b，出口距水翼尾緣為6b，計(jì)算域的進(jìn)口面、出口面尺寸為150 mm×150 mm。

圖1 NACA6408三維水翼計(jì)算域模型

本文利用ICEM 網(wǎng)格劃分軟件對(duì)三維NACA6408水翼計(jì)算域進(jìn)行結(jié)構(gòu)化網(wǎng)格劃分，針對(duì)水翼的頭部、尾部及上下壁面等部位進(jìn)行局部網(wǎng)格加密。經(jīng)過(guò)網(wǎng)格無(wú)關(guān)性計(jì)算驗(yàn)證，當(dāng)網(wǎng)格尺寸小于0.07 mm 時(shí)，水翼的升力系數(shù)趨于穩(wěn)定［12］，為合理利用計(jì)算機(jī)資源，提高數(shù)值模擬的準(zhǔn)確性、經(jīng)濟(jì)性，采用邊界層為0.07 mm的網(wǎng)格尺寸，網(wǎng)格總數(shù)為5501430。采用k-ε 湍流模型［13］和Zwart-Gerber-Belarnri（ZGB）空化模型［14］進(jìn)行模擬計(jì)算，與Schnerr-Sauer（SS）空化模型比較，本文采用的Zwart-Gerber-Belarnri（ZGB）空化模型使用空泡密度表示氣相的體積分?jǐn)?shù)［15］，能夠更好地模擬出非定常空化的準(zhǔn)周期性，捕捉到更強(qiáng)烈的云狀空化。

2 計(jì)算結(jié)果與討論

計(jì)算模型的b＝150 mm，計(jì)算域長(zhǎng)度則為1.5 m，在繞流速度υin＝1、5、10、15 和20 m／s 5 種工況下分別進(jìn)行數(shù)值模擬計(jì)算。在沿翼型主流方向等間距截取11 個(gè)截面i（i＝1，2，…，11）。由于翼型結(jié)構(gòu)特點(diǎn)，同一個(gè)截面將被翼型分成上下兩部分（見(jiàn)圖1），令水翼吸力面以上的截面為B 部分，即Bi，壓力面以下的截面為A部分，即Ai。

通過(guò)對(duì)繞流攻角為4°時(shí)的5 種不同繞流速度工況依次計(jì)算，數(shù)值結(jié)果表明翼型在υin＝1、5 和10 m／s工況下不會(huì)發(fā)生空化，在υin＝15 和20 m／s 時(shí)產(chǎn)生空化現(xiàn)象。在前3 種工況下，截面A 部分的壓力在不同位置處變化較小，主要是由于NACA6408 翼型結(jié)構(gòu)的壓力面曲率很小；截面B部分的平均壓力pav先減小后增加，即pav從截面B1減小至截面B4后增大，直至截面B11，原因是截面面積的變化使速度呈現(xiàn)出先增大后減小的規(guī)律，pav則是先減小后增大。此時(shí)同一截面的A與B部分的pav和湍動(dòng)能κ 變化平穩(wěn)，不存在空化現(xiàn)象。

如表1 所列為不同υin工況下以截面1、4、6、8、11為代表的B部分截面pav值。從表中可以看出，在同一υin下，不同時(shí)間（t＝0.05 和0.1 s）與不同截面間的pav略有不同，但是維持在同一個(gè)數(shù)量級(jí)。υin＜10 m／s時(shí)沒(méi)有發(fā)生空化，截面B1～B4間的pav隨著速度的增加而逐漸減小，相同速度下截面B6～B11間的pav逐漸增大；υin＝5 m／s時(shí)，pav≥89.528 kPa；當(dāng)υin＝10 m／s時(shí)，pav≥52. 192 kPa；而當(dāng)υin＝15 m／s 時(shí)，pav≥11.890 kPa，當(dāng)υin＝20 m／s 時(shí)，pav≥9.763 kPa。因此發(fā)生空化時(shí)，pav≤52.192 kPa。

表1 不同工況下截面B部分的pav變化kPa

如圖2 所示為υin＝15 m／s下翼展中間截面體積分?jǐn)?shù)及壓力云圖。由圖2（a）可知，在此工況下翼型的空化呈現(xiàn)一定的周期性，且周期T＝75 ms。由于此工況為非定常流動(dòng)，無(wú)法對(duì)每個(gè)時(shí)刻進(jìn)行描述，故將一個(gè)周期均分為8 等分以方便分析。在T／8 時(shí)，空化開(kāi)始出現(xiàn)于沿翼型吸力面主流方向b／4 ～b／2 處且附著于吸力面，隨著時(shí)間推移空化過(guò)程越來(lái)越明顯；在2T／8時(shí)，翼型吸力面b／4 ～b／2 處附著型空泡開(kāi)始向鉤狀空泡轉(zhuǎn)變；在3T／8 時(shí)，附著型空泡長(zhǎng)度達(dá)到最大，其在鉤狀空泡向后延伸至翼型吸力面4b／5 處開(kāi)始脫落；在4T／8 ～5T／8 時(shí)，空泡在翼型吸力面4b／5 處部分脫落，脫落后的空泡流向下游；在6T／8 時(shí)，剩余空泡繼續(xù)附著于翼型吸力面b／3 ～5b／6 處且空泡長(zhǎng)度隨著時(shí)間推移而逐漸變小；在7T／8 時(shí)，空泡長(zhǎng)度縮短至翼型吸力面b／3 ～2b／3 處，最終在8T／8 時(shí)逐漸減小為零并進(jìn)入下一個(gè)周期，空化過(guò)程呈現(xiàn)出準(zhǔn)周期性規(guī)律。由圖2（b）中可知，在6T／8 時(shí)，低壓區(qū)的覆蓋范圍最大，隨后低壓區(qū)的覆蓋范圍開(kāi)始減小；在8T／8 時(shí)，低壓區(qū)的覆蓋范圍達(dá)到最小，且與初始時(shí)刻的壓力p相同，從而形成一個(gè)完整的周期。結(jié)合圖2（a）可知，此時(shí)翼型吸力面低壓區(qū)的p低于空化發(fā)生的臨界壓力。

圖2 在υin ＝15 m／s時(shí)，翼展中間截面的體積分?jǐn)?shù)及壓力云圖

2.1 水翼各截面壓力p變化分析

針對(duì)υin＝15 m／s工況水翼表面發(fā)生空化時(shí)A、B兩部分各截面處的壓力p、湍動(dòng)能κ、速度分量u、v和w等參數(shù)變化進(jìn)行分析。圖3 所示υin＝15 m／s時(shí)，翼型主流方向i位置的A和B 部分在一個(gè)周期內(nèi)最大、最小和平均壓力折線圖。由圖可知，A 部分截面的最大壓力pmax和最小壓力pmin隨著截面呈現(xiàn)逐漸減小的趨勢(shì)，即在A2～A4基本趨于平穩(wěn)；A5～A11呈現(xiàn)逐漸下降。結(jié)合NACA6408 翼型結(jié)構(gòu)以及在υin＝1、5 和10 m／s工況下無(wú)空化現(xiàn)象，可以判定A1～A4壓力下降不是截面面積變化造成的，而是受空化發(fā)展的影響。B部分截面的pmax和pmin在B1～B4間呈斷崖式下降，原因是翼型結(jié)構(gòu)的影響。B4～B11呈現(xiàn)上升趨勢(shì)，上升過(guò)程波動(dòng)幅度存在較大的差異。

圖3 在υin ＝15 m／s時(shí)，11個(gè)不同截面位置的壓力變化折線圖

在0 ～2T／8 段，pav呈現(xiàn)大幅度上升，整體平均增加了12 kPa 左右，原因是因?yàn)榭张莞街谝硇臀γ妫铱张菡急群苄。瑢?duì)流場(chǎng)影響較小，從而使pav逐漸上升；在3T／8 ～6T／8 段，B4～B11平均壓力pav在降至最低后幾乎不再變化，中間位置的截面壓力維持在110 kPa左右，這是由于此時(shí)附著型空泡開(kāi)始向鉤狀空泡轉(zhuǎn)變，并且鉤狀空泡脫離翼型吸力面形成獨(dú)立的空泡，此時(shí)為空化的核心階段，壓力整體處于低壓狀態(tài)；在7T／8 ～8T／8 時(shí)刻B1～B11的pav開(kāi)始逐漸升高，并恢復(fù)至1T／8 的壓力水平，這是由于脫落的空泡隨主流流向下游并逐漸潰滅的原因。

綜合分析圖3 可知，A 部分截面的pmax、pmin和pav的波動(dòng)幅值都明顯大于B 部分，因此水翼翼型結(jié)構(gòu)對(duì)壓力脈動(dòng)的幅值有著重要影響，空化顯著的區(qū)域其壓力脈動(dòng)更加明顯，可以利用壓力脈動(dòng)變化來(lái)探究空化現(xiàn)象特性。

2.2 各截面κ分布特性

圖4 所示υin＝15 m／s下A、B部分各截面κ折線圖。由圖可知，在一個(gè)周期內(nèi)，A部分截面的κ幾乎不存在變化，。此時(shí)流場(chǎng)A 部分截面相對(duì)穩(wěn)定，可以排除翼型結(jié)構(gòu)對(duì)κ的影響。B部分截面中，在B1～B5區(qū)段內(nèi)湍動(dòng)能隨時(shí)間變化較小，流場(chǎng)相對(duì)穩(wěn)定；而在B6～B11內(nèi)湍動(dòng)能變化較為明顯，在B11處湍動(dòng)能最大值增加至3.6 J／kg，說(shuō)明B6～B11內(nèi)部空泡的產(chǎn)生使得流場(chǎng)紊亂。空化發(fā)生區(qū)域在0 時(shí)刻湍動(dòng)能變化明顯，為空化初始階段，此時(shí)流場(chǎng)相對(duì)紊亂，能量耗散大；在2T／8 和3T／8 時(shí)產(chǎn)生附著型空泡，流場(chǎng)相對(duì)穩(wěn)定；4T／8至8T／8 時(shí)空化由附著型空泡向鉤狀空泡轉(zhuǎn)化，在脫離翼型吸力面后最終潰滅，流場(chǎng)相對(duì)紊亂，能量耗散大，如圖4 的B部分截面所示。

圖4 在υin ＝15 m／s時(shí)，11個(gè)不同截面位置的κ變化折線圖

因此，可看出湍動(dòng)能代表流場(chǎng)的紊亂程度，空化的發(fā)生與湍動(dòng)能大小有關(guān)，同時(shí)與翼型結(jié)構(gòu)和回射流等有關(guān)。湍動(dòng)能又隨速度增加而增大，因此流場(chǎng)速度也是判斷流場(chǎng)紊亂的因素之一。

2.3 水翼各截面速度分布

為了分析各截面的速度對(duì)流場(chǎng)脈動(dòng)特性的影響規(guī)律，如圖5 所示υin＝15 m／s 下A、B 部分截面在同1個(gè)周期內(nèi)X、Y、Z3 個(gè)方向的平均分速度u、v和w折線圖。

圖5 在υin ＝15 m／s時(shí)，X、Y、Z方向11個(gè)不同截面位置的平均分速度變化

（1）A部分截面。X方向上A1～A3的主流方向速度u基本平穩(wěn)；A4～A11間速度u有小幅度上升，原因是翼型截面積變化造成的影響。Y方向上速度v在A1～A4速度為正且變化緩慢，而在A5～A11間速度v方向?yàn)樨?fù)且逐漸增大，原因是水流以較大的υin沖擊翼型頭部結(jié)構(gòu)從而產(chǎn)生反方向的速度。Z方向的速度w基本維持在0 附近，整體比較穩(wěn)定。

（2）B部分截面。X方向上B1～B4區(qū)間主流方向速度u逐漸增大；在B5～B11區(qū)間u逐漸減小且變化幅度較大，原因是翼型截面積變化以及空化現(xiàn)象的產(chǎn)生。Y方向上在B1～B4速度v為正且逐漸減小；在B5～B11速度v方向?yàn)樨?fù)，原因是空化的產(chǎn)生形成復(fù)雜旋渦；隨著時(shí)間推移，Y方向上速度零點(diǎn)的出現(xiàn)時(shí)間沿截面逐漸后移，說(shuō)明產(chǎn)生的旋渦逐漸往后移動(dòng)，證明了空泡的脫離并且向下游移動(dòng)，B9～B11速度v方向?yàn)樨?fù)且逐漸增大，主要是由于空泡隨主流流向下游或者發(fā)生潰滅，導(dǎo)致Y方向速度v最后逐漸減小。Z方向上速度w維持在0 附近且有正有負(fù)，說(shuō)明空化過(guò)程產(chǎn)生的旋渦并不是平面型，而是復(fù)雜的立體型旋渦。

綜上所述，X主流方向速度u主要受υin的影響，Y和Z方向的速度v和w主要受空化過(guò)程的影響，速度方向發(fā)生改變說(shuō)明有空泡的脫離、潰滅和立體型旋渦的產(chǎn)生。針對(duì)水翼空化過(guò)程的X、Y、Z3 個(gè)方向分速度分析是探討空化現(xiàn)象特性的重要研究思路之一。

3 結(jié)語(yǔ)

本文基于計(jì)算流體動(dòng)力學(xué)技術(shù)數(shù)值研究了NACA6408 水翼在繞流攻角4°時(shí)不同繞流速度工況的流場(chǎng)特性，通過(guò)數(shù)值分析υin＝15 m／s下水翼空化發(fā)生過(guò)程，獲得流場(chǎng)p、κ以及各方向分速度等參數(shù)的變化規(guī)律。結(jié)果表明：①隨著進(jìn)口速度的增加，各截面上的平均湍動(dòng)能以及速度的波動(dòng)幅值明顯增加，壓力、湍動(dòng)能以及流場(chǎng)速度的變化與空化過(guò)程密切相關(guān)。②水翼空化發(fā)生區(qū)域壓力、湍動(dòng)能以及流場(chǎng)速度變化明顯，水翼壓力面以下區(qū)域的參數(shù)變化相對(duì)緩慢，水翼吸力面以上區(qū)域變化相對(duì)明顯。壓力面變化曲率較大導(dǎo)致進(jìn)口水流流經(jīng)翼型吸力面時(shí)流場(chǎng)相對(duì)紊亂。③在υin＝1、5 和10 m／s工況時(shí)水翼表面不會(huì)產(chǎn)生空化現(xiàn)象，在υin＝15 和20 m／s時(shí)發(fā)生空化，同時(shí)空化基本在B5～B11發(fā)生。④空化發(fā)生時(shí)非主流方向（Y和Z方向）的速度存在方向變換，形成了復(fù)雜的立體型旋渦。本文的數(shù)值計(jì)算結(jié)果對(duì)深入認(rèn)識(shí)水翼空化演變機(jī)理，提高水翼抗空化性能具有重要的參考意義。