DDES方法在復(fù)雜旋翼流場計(jì)算中的應(yīng)用

董軍，葉靚

中國航空工業(yè)空氣動(dòng)力研究院， 沈陽 110034

旋翼流場是一個(gè)典型的復(fù)雜非定常流場。一方面，槳葉的運(yùn)動(dòng)規(guī)律復(fù)雜，在有來流速度時(shí)，需要調(diào)整槳葉的周期變距操縱，實(shí)現(xiàn)旋翼氣動(dòng)力和力矩平衡；同時(shí)在周期性氣動(dòng)力和慣性力作用下，槳葉還要附加揮舞和擺振運(yùn)動(dòng)，槳葉前行側(cè)氣流壓縮、后行時(shí)的動(dòng)態(tài)失速等都是旋翼特有的氣動(dòng)現(xiàn)象。另一方面旋翼的工作狀態(tài)較多，除去常規(guī)的懸停和平飛巡航狀態(tài)，還可能進(jìn)行近地、小速度下降和機(jī)動(dòng)飛行等，在這些工作狀態(tài)下，流場的復(fù)雜度進(jìn)一步增加，預(yù)測的難度進(jìn)一步加大。為準(zhǔn)確計(jì)算旋翼的氣動(dòng)力，旋翼流場中的流動(dòng)細(xì)節(jié)捕捉(激波、槳尖渦和槳葉表面附近空間的氣流分離)則成為關(guān)鍵因素。

在先前的計(jì)算研究中，基于雷諾平均 Navier-Stokes(RANS)方程的流場預(yù)測占據(jù)著主導(dǎo)地位[1]。隨著計(jì)算流體力學(xué)的發(fā)展，一些計(jì)算量更大、計(jì)算精度更高的計(jì)算方法越來越得到重視[2-3]。

近年來，兼顧RANS和大渦模擬(LES)方法的優(yōu)勢，在不同的空間區(qū)域，選擇不同特征長度，能夠更好地預(yù)測分離流動(dòng)的RANS/LES混合方法開始廣泛應(yīng)用于固定翼飛行器流場和氣動(dòng)力的預(yù)測[4-7]。在近壁面區(qū)域，該方法等同于RANS方法，在遠(yuǎn)離物面的區(qū)域，計(jì)算模型為LES方法。該方法基于特征尺度判斷實(shí)現(xiàn)兩種計(jì)算模式的切換。與RANS方法相比，混合方法能夠較好地預(yù)測近壁面附近的分離流動(dòng)，在較大范圍氣流分離發(fā)生時(shí)，預(yù)測效果更好；與LES方法相比，其網(wǎng)格數(shù)量要求大大降低，因而成為一種可接受的計(jì)算手段。目前RANS/LES混合方法發(fā)展出了脫體渦模擬(DES)、延遲脫體渦模擬(DDES)、IDDES(Improved DDES)等多個(gè)版本[8-11]，其中DDES方法設(shè)置了延遲函數(shù)，克服了網(wǎng)格加密時(shí)可能形成的附面層內(nèi)LES方法過早啟動(dòng)的問題。

隨著計(jì)算數(shù)值手段的發(fā)展進(jìn)步，精確預(yù)測流場中的流動(dòng)分離和尾跡發(fā)展成為進(jìn)一步提高旋翼流場和氣動(dòng)力計(jì)算精度的關(guān)鍵制約因素。為更好地模擬在較大的空間范圍內(nèi)螺旋發(fā)展的旋翼尾跡，渦核半徑內(nèi)需要布置多個(gè)網(wǎng)格點(diǎn)，這使得旋翼流場計(jì)算總體網(wǎng)格數(shù)量膨脹速度較快。在計(jì)算機(jī)硬件水平及計(jì)算程序并行能力迅速提升的大背景下，RANS/LES混合方法也開始應(yīng)用在旋翼流場計(jì)算方面。特別是在2014—2017年的AIAA會(huì)議上，陸續(xù)有研究者針對(duì)旋翼的懸停流場預(yù)測問題，采用RANS/LES混合方法進(jìn)行了數(shù)值計(jì)算研究[12-16]。采用混合方法時(shí)，與傳統(tǒng)的基于RANS方法的旋翼流場預(yù)測比較，新方法本身在分離流動(dòng)的預(yù)測能力上有所加強(qiáng)，且預(yù)測工作一般采用了更大的網(wǎng)格數(shù)量，預(yù)測出了更多的流場細(xì)致結(jié)構(gòu)。在國內(nèi)，旋翼流場的預(yù)測一般基于RANS方法[17-20]，較少查閱到基于RANS/LES混合方法的旋翼流場計(jì)算文獻(xiàn)。

鑒于此，本文采用DDES方法進(jìn)行了旋翼流場計(jì)算研究。計(jì)算網(wǎng)格系統(tǒng)構(gòu)建方面，為有效控制計(jì)算規(guī)模，用于尾跡捕捉的背景網(wǎng)格使用了自適應(yīng)直角網(wǎng)格，并根據(jù)特定的策略控制網(wǎng)格分布。為提高計(jì)算速度，采用了分布式并行的計(jì)算方法，無論是在流場求解還是在網(wǎng)格組裝等方面，均實(shí)現(xiàn)了并行化處理?；谶@些方法，對(duì)旋翼多個(gè)狀態(tài)下的流場和氣動(dòng)力進(jìn)行了計(jì)算，這些計(jì)算工況對(duì)應(yīng)了槳葉表面基本為附著流動(dòng)，槳葉表面附近空間發(fā)生小范圍和大面積的氣流分離情況。討論了對(duì)于旋翼流場和氣動(dòng)力預(yù)測問題，不同工況下計(jì)算結(jié)果與傳統(tǒng)RANS方法的計(jì)算差異，比較的內(nèi)容既包括流場結(jié)構(gòu)，也包含宏觀的氣動(dòng)力差異。計(jì)算結(jié)果體現(xiàn)出RANS/LES混合方法在槳葉表面附近存在氣流分離情況下的特定優(yōu)勢，采用此方法在一定程度上可以提升旋翼流場和氣動(dòng)力的計(jì)算精度。

1 計(jì)算網(wǎng)格系統(tǒng)

為有效捕捉旋翼槳葉表面附近流動(dòng)和遠(yuǎn)場尾跡發(fā)展，在本文的嵌套網(wǎng)格系統(tǒng)中，貼體網(wǎng)格部分使用了純六面體網(wǎng)格，以提高網(wǎng)格質(zhì)量，背景網(wǎng)格使用了可自適應(yīng)的直角網(wǎng)格。

背景直角網(wǎng)格按照以下步驟生成：

步驟1根據(jù)預(yù)定的流動(dòng)區(qū)域外場大小，各方向均勻布置網(wǎng)格點(diǎn)，生成一個(gè)較為稀疏的初始網(wǎng)格。

步驟2根據(jù)槳葉運(yùn)動(dòng)掃掠過的空間區(qū)域范圍，進(jìn)行一次網(wǎng)格加密。

步驟3對(duì)預(yù)估的尾跡區(qū)域空間網(wǎng)格，進(jìn)行適度加密調(diào)整。

步驟4根據(jù)流場計(jì)算結(jié)果，使用渦量探針，加密流場中渦量集中區(qū)域的網(wǎng)格。

采用這樣的方法，貼體與背景網(wǎng)格插值單元的尺度比較接近，有利于維持?jǐn)?shù)值格式的精度；網(wǎng)格自適應(yīng)的第4步可以多次逐級(jí)進(jìn)行，用于有效捕捉遠(yuǎn)場尾跡發(fā)展，且使得計(jì)算網(wǎng)格的數(shù)量可控。

在直角網(wǎng)格每次自適應(yīng)加密后，還要進(jìn)行“洞”和“層差”的判斷，消除網(wǎng)格中的不光滑區(qū)域。實(shí)際旋翼流場計(jì)算時(shí)，可以根據(jù)槳葉的弦長預(yù)估可能的渦核的尺寸，從而確定初始和最終的網(wǎng)格加密后尺度。

在數(shù)據(jù)格式構(gòu)造方面，均以非結(jié)構(gòu)形式對(duì)貼體和背景網(wǎng)格進(jìn)行數(shù)據(jù)存儲(chǔ)。由于貼體和背景網(wǎng)格采用同一數(shù)據(jù)結(jié)構(gòu)[21]，在空間格式構(gòu)建等方面實(shí)現(xiàn)了統(tǒng)一處理，計(jì)算程序得到了極大的簡化。

2 流場數(shù)值計(jì)算方法

2.1 控制方程形式

采用格心格式的有限體積法，地面坐標(biāo)系下積分形式的雷諾平均Navier-Stokes方程可以寫為

(1)

式中：W為守恒變量；S為積分面面積；Ω為控制體體積；Fc和Fv分別為對(duì)流和黏性通量；p為壓強(qiáng)；ρ為密度；V為對(duì)流速度；u、v、w為3個(gè)方向速度；nx、ny、nz為網(wǎng)格交界面單位外法矢；H為單位質(zhì)量總焓；τ為應(yīng)力；T為溫度；K為傳熱系數(shù)。

不考慮與轉(zhuǎn)捩有關(guān)項(xiàng)，Spalart-Allmaras湍流模型[22]可以寫為

(2)

fw定義為

g=r+cw2(r6-r)

渦黏性υt定義為

模型常數(shù)定義為：cb1=0.135 5；σ=2/3；cb2=0.622；Von Karman常數(shù)κ=0.41；cw1=cb1/κ2+(1+cb2)/σ；cw2=0.3；cw3=2.0；cυ1=7.1。 采用RANS/LES混合方法時(shí)，長度尺度定義為

(3)

式中：d為空間網(wǎng)格到最近壁面的距離；Δ=max(Δx,Δy,Δz)為當(dāng)前單元和它所有鄰居單元中心間最大距離；CDES=0.65。

對(duì)于典型的RANS網(wǎng)格，在近壁面區(qū)域，RANS/LES混合模型中的RANS模式被激活。因?yàn)樵谠搮^(qū)域，壁面法向距離通常小于各方向網(wǎng)格間距的最大值。如果流向和展向的網(wǎng)格尺度在近壁面區(qū)域和法向尺度近似，LES模型則在邊界層內(nèi)被激活，導(dǎo)致邊界層過早分離。為克服這一問題，Spalart等又提出了DDES模型，長度尺度改寫為

(4)

fd=1-tanh(8rd)3

2.2 求解方法

非定常流場計(jì)算方面，主控和湍流方程以松耦合的方式，采用雙時(shí)間[23]結(jié)合LU-SGS[24]方法進(jìn)行時(shí)間推進(jìn)，物理時(shí)間步長選擇為旋翼旋轉(zhuǎn)周期的1/1 440。在最終的計(jì)算網(wǎng)格上，完成11個(gè)旋轉(zhuǎn)周期數(shù)值模擬后，再進(jìn)行3個(gè)周期的計(jì)算取結(jié)果進(jìn)行平均，獲得最后的計(jì)算結(jié)果。

空間格式方面，采用二階迎風(fēng)格式[25]計(jì)算無黏通量。由于網(wǎng)格采用無結(jié)構(gòu)化處理，物理量的梯度需要采用重構(gòu)方法實(shí)現(xiàn)(見文獻(xiàn)[26])，黏性通量采用中心格式進(jìn)行計(jì)算。

3 數(shù)據(jù)并行方法

3.1 網(wǎng)格分區(qū)方法

分布式并行計(jì)算時(shí)，首先要對(duì)網(wǎng)格進(jìn)行分區(qū)，該項(xiàng)工作在網(wǎng)格每次自適應(yīng)調(diào)整后進(jìn)行一次，使得各處理器上的負(fù)載盡量均衡。

旋翼流場計(jì)算采用多塊嵌套網(wǎng)格系統(tǒng)，這里按照處理器線程數(shù)對(duì)每塊網(wǎng)格進(jìn)行分區(qū)，網(wǎng)格按照非結(jié)構(gòu)的方式進(jìn)行存儲(chǔ)，因而可以使用metis[27]工具來進(jìn)行網(wǎng)格的切分。

設(shè)網(wǎng)格塊的標(biāo)號(hào)為B1,B2,…;網(wǎng)格B1分區(qū)后的小塊網(wǎng)格標(biāo)號(hào)分別為B1(1),B1(2)，…;網(wǎng)格B2分區(qū)后的小塊網(wǎng)格標(biāo)號(hào)分別為B2(1),B2(2)，…。在線程1上分布的網(wǎng)格塊為B1(1)，B2(1)，…，在線程N(yùn)上分布的網(wǎng)格塊為B1(N)，B2(N)，…。網(wǎng)格分區(qū)完成后，為構(gòu)建數(shù)值格式方便，每塊網(wǎng)格的分區(qū)邊界單元向外擴(kuò)展2層虛擬單元，形成處理器之間網(wǎng)格信息傳遞邊界。

3.2 網(wǎng)格組裝分布式并行方法

網(wǎng)格組裝工作分為洞切割、貢獻(xiàn)單元搜索和數(shù)據(jù)傳遞區(qū)域構(gòu)建等幾步進(jìn)行。為增加整個(gè)流場解算的并行段比重，提高計(jì)算效率，這里對(duì)這些步驟進(jìn)行并行化處理。

洞切割方面，首先在所有網(wǎng)格塊中依據(jù)距離物面的遠(yuǎn)近，根據(jù)已有網(wǎng)格面形成一個(gè)洞構(gòu)造表面，并把該構(gòu)造表面的網(wǎng)格信息廣播到所有的處理器上。這樣在各個(gè)處理器上，均可按照該構(gòu)造面進(jìn)行洞切割。

貢獻(xiàn)單元搜索方面，對(duì)每個(gè)處理器上洞切割后形成的插值邊界，均發(fā)送給其他處理器，進(jìn)行貢獻(xiàn)單元搜索，然后返回相應(yīng)結(jié)果。這樣就實(shí)現(xiàn)了貢獻(xiàn)單元搜索的并行化。

對(duì)每個(gè)處理器，在其他多個(gè)處理器上進(jìn)行貢獻(xiàn)單元搜索的結(jié)果回收后，在自身處理器上進(jìn)行結(jié)果匯總，并最終確定正確的數(shù)據(jù)交換關(guān)系。

3.3 流場計(jì)算相關(guān)分布式并行方法

在分布式并行情況下，采用隱式時(shí)間推進(jìn)時(shí)，由于在每個(gè)處理器上無法獲得全場的流動(dòng)信息更新情況，因而原始的LU-SGS方法需要改進(jìn)。這種改進(jìn)主要體現(xiàn)在各處理器分區(qū)邊界的流動(dòng)變量更新計(jì)算方面。實(shí)際計(jì)算時(shí)，可以先用退化的方式(雅克比迭代)計(jì)算出傳值邊界元的流動(dòng)變量更新量，然后各個(gè)處理器進(jìn)行數(shù)據(jù)交換，再進(jìn)行局部的LU-SGS迭代求解，具體的執(zhí)行方法見文獻(xiàn)[28]。

對(duì)于一些流場計(jì)算中需要的幾何信息，也可以進(jìn)行并行化處理。如壁面距離確定方面，可以通過壁面坐標(biāo)的廣播實(shí)現(xiàn)在每個(gè)處理器上的并行計(jì)算。

4 計(jì)算方法驗(yàn)證

4.1 翼型振蕩氣動(dòng)力計(jì)算

為驗(yàn)證本文方法在非定常氣動(dòng)力計(jì)算方面的能力，選取了NACA0015翼型強(qiáng)迫振蕩狀態(tài)進(jìn)行數(shù)值模擬。翼型弦長為0.304 8 m，來流馬赫數(shù)為0.29。翼型迎角α隨時(shí)間t的變化規(guī)律為

α=4.0+4.2sin(20πt)

(5)

圖1給出了計(jì)算得到的升力系數(shù)CL、阻力系數(shù)CD和力矩系數(shù)Cm隨迎角α變化的遲滯曲線與試驗(yàn)值[29]的對(duì)比。計(jì)算結(jié)果表明，采用RANS和DDES方法計(jì)算出的升力與試驗(yàn)結(jié)果符合得更好，阻力的結(jié)果相對(duì)差些，這可能與湍流模型的選取等多種因素有關(guān)；兩者計(jì)算的結(jié)果均能夠有效模擬出氣動(dòng)力隨迎角變化的趨勢，其相互間差異隨著翼型迎角的增大而增大。

圖1 NACA0015翼型振蕩狀態(tài)非定常氣動(dòng)特性Fig.1 Unsteady aerodynamic characteristics of oscillation of NACA0015 airfoil

4.2 槳尖附近區(qū)域激波誘導(dǎo)氣流分離

為驗(yàn)證本文方法在旋翼氣動(dòng)力計(jì)算方面的能力，選擇了有試驗(yàn)值可供對(duì)比的“Caradonna & Tung旋翼”。該旋翼槳葉數(shù)為2片，展長為1.143 m，展弦比為6，無扭轉(zhuǎn)尖削，槳葉截面段翼型為NACA0012。選取的計(jì)算狀態(tài)為槳尖馬赫數(shù)0.794，12°總距角。計(jì)算發(fā)現(xiàn)，在該狀態(tài)下，在靠近槳尖的槳葉上表面附近空間區(qū)域出現(xiàn)較強(qiáng)激波，在激波后方形成一定范圍的氣流分離。本文對(duì)比了采用非定常RANS和DDES方法的計(jì)算結(jié)果。圖2為計(jì)算得到的不同槳葉方位角ψ時(shí)槳葉上表面壓力的對(duì)比?？梢钥闯?，采用RANS方法時(shí)，計(jì)算的表面壓力幾乎不隨時(shí)間(槳葉方位角)變化，而采用DDES方法時(shí)，則預(yù)測到了槳葉上表面壓力分布時(shí)變特性(靠近槳尖區(qū)域)。這體現(xiàn)出DDES方法在分離流動(dòng)方面的預(yù)測能力。

圖2 槳葉上表面壓力對(duì)比(左：RANS，右：DDES)Fig.2 Comparison of blade upper surface pressure(left: RANS; right: DDES)

圖3 槳葉表面壓力系數(shù)計(jì)算值與試驗(yàn)值的對(duì)比Fig.3 Comparison of calculated and experimental surface pressure coefficients of blade

圖3為計(jì)算得到的靠近槳尖的部分截面時(shí)均表面壓力系數(shù)Cp計(jì)算值和試驗(yàn)值[30]的對(duì)比，r為槳葉截面段所在位置半徑，R為旋翼半徑,c為弦長?？梢园l(fā)現(xiàn)，DDES方法預(yù)測了激波后方的非定常分離流動(dòng)，激波后上表面的壓力值更小；相對(duì)RANS方法，DDES方法預(yù)測的激波位置也更靠近槳葉前緣。對(duì)于槳葉內(nèi)側(cè)截面，由于其流動(dòng)為附著流動(dòng)，兩種方法預(yù)測結(jié)果一致，這里不再給出。

5 算例分析

5.1 雙旋翼干擾流場

本算例中，計(jì)算模型為南京航空航天大學(xué)的一副試驗(yàn)共軸雙旋翼[31]，該旋翼的上下旋翼相距0.175 m，各由兩片槳葉組成，槳葉半徑為0.945 m，無扭轉(zhuǎn)，翼型為NACA0012，弦長為0.076 m，槳葉r=0.95R處向外有尖削，尖削比為1/3，試驗(yàn)狀態(tài)槳尖速度為118 m/s，上下旋翼槳葉槳距角分別為9°和10.27°，預(yù)錐角為5°。

計(jì)算網(wǎng)格如圖4所示，旋翼旋轉(zhuǎn)軸方向?yàn)閅軸，左圖是平行Y軸的網(wǎng)格截面圖，右圖為旋翼下方垂直Y軸的網(wǎng)格截面圖。網(wǎng)格尺度方面，旋翼下方尾跡區(qū)的網(wǎng)格尺度L設(shè)定為0.1C(C為槳葉根弦長)，按照渦量自適應(yīng)的網(wǎng)格尺度為0.05C，空間網(wǎng)格尺度隨著與旋翼距離增加逐漸放寬。

貼體網(wǎng)格單元為純六面體形式，每片槳葉貼體網(wǎng)格單元數(shù)約為430萬，背景直角網(wǎng)格單元數(shù)約為6 900萬。計(jì)算采用700個(gè)CPU并行，每次流場內(nèi)迭代時(shí)間約為2.8 s。在該計(jì)算狀態(tài)下，由于槳葉截面相對(duì)來流速度較低，迎角不高，槳葉表面觀察不到明顯的分離流動(dòng)。

圖4 共軸旋翼流場計(jì)算網(wǎng)格示意圖Fig.4 Schematic of grid for coaxial rotor flowfield calculation

圖5 時(shí)均軸向誘導(dǎo)速度Fig.5 Time-averaged axial induced velocity

圖5給出了計(jì)算得到的槳盤下方不同截面處(圖中y為槳盤下方截面的位置)軸向誘導(dǎo)速度時(shí)均值Vy與試驗(yàn)值的對(duì)比。在有試驗(yàn)結(jié)果的區(qū)域，計(jì)算與試驗(yàn)值較吻合(兩種計(jì)算方法得到的軸向誘導(dǎo)速度重合)。采用DDES方法與RANS方法的計(jì)算差別主要體現(xiàn)在靠近槳根和槳尖的區(qū)域，在靠近槳根區(qū)域，DDES方法預(yù)測出稍大一些的下洗速度。這可能是由于在槳葉靠近兩端的區(qū)域，氣流繞過槳葉，展向流動(dòng)效應(yīng)加強(qiáng)，流動(dòng)復(fù)雜度增加。

圖6為計(jì)算得到的下旋翼槳葉整體氣動(dòng)力(拉力系數(shù)CT和扭矩系數(shù)CQ)的比較。在槳葉方位角ψ=100°～130°(定義上下旋翼俯視重合時(shí)的槳葉方位角為0°)，下旋翼受到上旋翼尾跡影響較為嚴(yán)重，無論是槳葉拉力和扭矩都出現(xiàn)波動(dòng)。兩種方法計(jì)算得到的拉力差別體現(xiàn)在峰值大小上，而扭矩差別僅體現(xiàn)在出現(xiàn)氣動(dòng)干擾的方位上。總體來說，兩種方法計(jì)算出的宏觀氣動(dòng)力差別較小。這可能與槳葉表面空間附近無明顯分離有關(guān)。對(duì)于上旋翼槳葉，兩種計(jì)算方法得到的結(jié)果基本重合，這里不再給出。

圖6 下旋翼的拉力和扭矩系數(shù)Fig.6 Trust and torque coefficients of lower rotor

圖7 截面渦量等值圖Fig.7 Contour of vorticity in cutting plane

圖7為采用DDES方法計(jì)算得到流場空間特征截面上的渦量等值圖，上下旋翼的尾跡在經(jīng)過下旋翼的下方后，開始融合并一起向下發(fā)展。采用本文的網(wǎng)格系統(tǒng)和數(shù)值方法，計(jì)算得到了較為精細(xì)的流動(dòng)結(jié)構(gòu)。

5.2 下降狀態(tài)旋翼流場

下降狀態(tài)是直升機(jī)的一個(gè)重要工作狀態(tài)，該狀態(tài)下流場結(jié)構(gòu)更為復(fù)雜，容易發(fā)生氣動(dòng)干擾現(xiàn)象，引發(fā)氣動(dòng)和噪聲問題。

在該算例中，仍然選擇了“Caradonna & Tung旋翼”構(gòu)型作為研究對(duì)象。計(jì)算狀態(tài)為槳尖馬赫數(shù)0.439， 8°總距角。軸向誘導(dǎo)速度Vy分別選擇為20 m/s和40 m/s。

貼體網(wǎng)格單元為純六面體形式，每片槳葉貼體網(wǎng)格單元數(shù)約為300萬，背景直角網(wǎng)格單元數(shù)約為3 400萬。計(jì)算采用700個(gè)CPU并行，每次流場內(nèi)迭代時(shí)間約為1.2 s。

圖8為下降速度20 m/s時(shí)，采用DDES方法計(jì)算得到的不同槳葉方位角時(shí)流場特征截面上的渦量等值圖，可以發(fā)現(xiàn)計(jì)算的渦量具備明顯的時(shí)變特征。需要注意的是，由于來流方向?yàn)闃~下方，在槳葉的向下作用力和向上來流共同作用下，槳葉脫出的尾跡在經(jīng)過180°方位角時(shí)(到達(dá)另一片槳葉空間位置附近)并不能迅速遠(yuǎn)離槳盤平面。

圖9為下降速度20 m/s時(shí)，采用DDES方法和RANS方法計(jì)算得到的特征截面(x=0 m)上的流線圖(0°槳葉方位角)?？梢园l(fā)現(xiàn)，采用RANS方法計(jì)算的流場具備良好的對(duì)稱性，而采用DDES方法結(jié)果則存在明顯的左右不對(duì)稱。

圖8 不同槳葉方位角時(shí)的流場截面渦量等值圖(x=0 m)Fig.8 Contour of vorticity in flowfields cutting plane with different blade azimuthal angles (x=0 m)

圖9 截面流線(x=0 m)Fig.9 Streamlines in cutting plane (x=0 m)

由于DDES方法相比RANS方法在分離流動(dòng)與旋渦結(jié)構(gòu)捕捉能力方面更強(qiáng)些，因而其計(jì)算結(jié)果在流場中存在更多復(fù)雜渦結(jié)構(gòu)的情況下更接近于真實(shí)物理情況。

圖10給出了該狀態(tài)下，特征截面上的渦捕捉情況對(duì)比。可見DDES方法解析了更多的流場渦結(jié)構(gòu)，且其集中渦量的計(jì)算強(qiáng)度普遍高于RANS方法結(jié)果。

圖11和圖12分別為Vy=20，40 m/s，槳葉方位角為0°時(shí)，采用DDES方法和RANS方法計(jì)算得到的不同槳葉展向截面上的流線圖(已經(jīng)處理成相對(duì)槳葉截面速度)?？梢园l(fā)現(xiàn)在不同槳葉截面段上，Vy=20 m/s時(shí)，兩種方法計(jì)算流動(dòng)情況類似，即均無明顯分離流動(dòng)現(xiàn)象。在Vy=40 m/s時(shí)，兩種方法均計(jì)及槳葉表面大范圍的分離流動(dòng)，但采用DDES方法和RANS方法計(jì)算出不同的流動(dòng)狀態(tài)，尤其是靠近槳尖附近的截面，兩種方法計(jì)算得到的結(jié)果差異更大。

圖13給出了采用不同方法計(jì)算得到的槳葉升力展向變化情況。在Vy=20 m/s狀態(tài)下，RANS方法計(jì)算出時(shí)均的槳葉升力大于DDES方法的結(jié)果?？梢园l(fā)現(xiàn)，即使槳葉截面段上觀察不到明顯分離流動(dòng)(圖11)，RANS和DDES方法計(jì)算的結(jié)果差異也很明顯。分析其原因可能是，在該狀態(tài)下，旋翼整體尾跡和槳葉干擾作用較強(qiáng)，一方面旋翼仍提供正向升力，另一方面，旋翼尾跡向上方發(fā)展，不能迅速遠(yuǎn)離槳盤平面，由于兩種預(yù)測方法計(jì)算出的尾跡發(fā)展情況不同，因而尾跡與槳葉的作用情況也必然不同。還可以看出，在該狀態(tài)下，兩種方法預(yù)測的氣動(dòng)力都隨槳葉方位角變化明顯，這表明在這種情況下，即便是軸流入流狀態(tài)，旋翼的氣動(dòng)力也具備較強(qiáng)的時(shí)變特性。

圖10 截面渦量 (x=0 m)Fig.10 Vorticity in cutting plane (x=0 m)

圖11 Vy=20 m/s時(shí)槳葉截面流線(左：DDES，右：RANS)Fig.11 Streamlines in blade cross-sections for Vy=20 m/s (left: DDES; right: RANS)

圖12 Vy=40 m/s時(shí)槳葉截面流線(左：DDES，右：RANS)Fig.12 Streamlines in blade cross-sections for Vy=40 m/s (left: DDES; right: RANS)

圖13 槳葉展向升力比較Fig.13 Comparison of spanwise lift of blade

在Vy=40 m/s狀態(tài)下，RANS計(jì)算出時(shí)均的槳葉升力小于DDES方法的結(jié)果。這應(yīng)該與兩種方法預(yù)測出的槳葉表面附近的分離形式差異有關(guān)。從圖13還可以看出，采用RANS方法計(jì)算出的氣動(dòng)力相對(duì)穩(wěn)定，槳葉不同方位角時(shí)的氣動(dòng)力比較接近。采用DDES方法計(jì)算的結(jié)果則體現(xiàn)出較強(qiáng)的時(shí)間特性，這體現(xiàn)出DDES方法在物體表面附近存在強(qiáng)分離時(shí)的計(jì)算特性。

6 結(jié) 論

在可自適應(yīng)嵌套網(wǎng)格系統(tǒng)下，基于DDES方法，進(jìn)行了不同狀態(tài)下復(fù)雜旋翼的流場和氣動(dòng)力計(jì)算，并與相應(yīng)的RANS結(jié)果進(jìn)行比較。

1) 在槳葉表面無明顯分離，且不存在較強(qiáng)的旋翼/尾跡間相互作用情況下，DDES和RANS方法計(jì)算出的槳葉宏觀氣動(dòng)力方面差異較小。

2) 旋翼/尾跡間存在較強(qiáng)的相互作用情況下，即便槳葉截面流場觀察不出明顯分離流動(dòng)，兩種方法計(jì)算結(jié)果也存在明顯差異，這可能與尾跡預(yù)測的強(qiáng)度等因素有關(guān)。在該情況下，即便是在軸流入流狀態(tài)，無論是RANS方法還是DDES方法，計(jì)算出的旋翼氣動(dòng)力都體現(xiàn)出較強(qiáng)的時(shí)變特性。

3) 在槳葉表面附近空間存在大范圍分離流動(dòng)的情況下，采用DDES方法計(jì)算出更多的流場結(jié)構(gòu)，氣動(dòng)力的時(shí)變差異也更為顯著，體現(xiàn)出經(jīng)典的混合RANS/LES方法的計(jì)算優(yōu)勢。

4) 本文采用的網(wǎng)格系統(tǒng)和計(jì)算方法，能夠較好地預(yù)測流場中的流動(dòng)細(xì)節(jié)和旋翼氣動(dòng)力。

值得注意的是，對(duì)存在復(fù)雜分離流動(dòng)的旋翼流場預(yù)測來說，除去本文考慮的DDES方法研究外，提高計(jì)算格式的數(shù)值精度、網(wǎng)格數(shù)量和質(zhì)量，進(jìn)行湍流模型的選擇及轉(zhuǎn)捩研究等都是進(jìn)一步提升旋翼流場計(jì)算精度的關(guān)鍵因素，也是未來工作的發(fā)展方向。