石力王，胡明輝,2，陳廣艷，陳倫國，石萬凱,2，曾 利

(1.重慶大學(xué) 機械與運載工程學(xué)院，重慶 400044；2.重慶大學(xué) 機械傳動國家重點實驗室， 重慶 400044；3.中國航發(fā)湖南動力機械研究所，湖南 株洲 412000；4.直升機傳動技術(shù)國防科技重點實驗室，湖南 株洲 412000)

0 引言

相比于常規(guī)直升機旋翼，共軸式剛性雙旋翼可以顯著提高直升機運行時的工作拉力，無需安裝力矩平衡尾槳，具有結(jié)構(gòu)緊湊、氣動效率高等特點，在高速直升機上應(yīng)用廣泛。但共軸式雙旋翼相較于單旋翼，氣動特性更為復(fù)雜，上下兩副旋翼的旋向相反，旋翼周圍流場為非定常流場[1]，上下兩副旋翼的間距小，上下旋翼間產(chǎn)生嚴重氣動干擾。氣動干擾不僅影響直升機的飛行平衡及其操縱穩(wěn)定，也使得共軸式雙旋翼氣動載荷計算難度驟增[2]。在正向設(shè)計直升機旋翼的時候，動力傳動系統(tǒng)的匹配是關(guān)鍵步驟之一，旋翼的拉力和力矩是選擇直升機發(fā)動機和減速器的重要指標之一。因此，建立精準的共軸式雙旋翼氣動載荷分析模型，分析旋翼設(shè)計參數(shù)對其氣動載荷的影響，探究平衡上下旋翼力矩方法，對于指導(dǎo)共軸式剛性雙旋翼直升機的正向設(shè)計具有重要意義。

國內(nèi)外眾多研究人員致力于開發(fā)精確有效的共軸式雙旋翼氣動載荷數(shù)值計算方法。自由尾跡、動量理論以及葉素理論等理論方法是開展共軸式雙旋翼氣動建模的常用方法[3-6]。此外，部分學(xué)者將自由尾跡法與升力面法以及二階升力線法相結(jié)合，進行雙旋翼懸停工況下的建模分析[7-8]。上述方法多屬于勢流法，雖計算量小，卻無法準確捕捉激波失速、槳尖渦流與氣流分離等現(xiàn)象，因此計算精度較低。為準確描述流體域內(nèi)非定?，F(xiàn)象，考慮渦量運輸與氣體粘性對流場的影響，提高共軸式雙旋翼氣動載荷計算精度，基于嵌套網(wǎng)格與Euler/Navier-Stokes方程的CFD方法得到研究人員的青睞。Lakshminarayan[9-10]建立了一套運用滑移與嵌套網(wǎng)格的混合網(wǎng)格方法，研究了共軸式雙旋翼氣動特性變化規(guī)律，但滑移網(wǎng)格會導(dǎo)致部分通量的損失。葉靚等[1]基于非結(jié)構(gòu)嵌套網(wǎng)格與網(wǎng)格自適應(yīng)技術(shù)，提出了一套共軸式雙旋翼氣動特性數(shù)值求解器。朱正等[11]則運用基于運動嵌套網(wǎng)格的CFD方法，對雙旋翼非定常流動干擾進行機理分析。

上述工作為共軸式雙旋翼的研究提供了重要參考，卻缺乏上下旋翼力矩平衡方法的研究，由于上下旋翼氣動干擾嚴重，在無平衡尾槳的情況下，上下旋翼間的力矩差嚴重影響共軸式剛性雙旋翼直升機的飛行平衡與操縱穩(wěn)定。鑒于此，本文基于先進的計算流體力學(xué)，建立了適用于共軸式剛性雙旋翼氣動性能求解的方法，利用雷諾時均Navier-Stokes(RANS)方程、雙時間方法和非結(jié)構(gòu)運動嵌套網(wǎng)格等技術(shù)，針對共軸式剛性雙旋翼懸停狀態(tài)進行了氣動載荷仿真分析，并通過改變上旋翼的總距角，研究雙旋翼總距差對上下旋翼氣動載荷的影響規(guī)律，獲得了雙旋翼的最佳總距角，在保證拉力達標的情況下，實現(xiàn)上下旋翼的力矩平衡，從而可以指導(dǎo)后期的關(guān)于直升機發(fā)動機的選型和動力傳動系統(tǒng)的正向設(shè)計。

1 運動嵌套網(wǎng)格與數(shù)值方法

1.1 運動嵌套網(wǎng)格生成

在眾多網(wǎng)格生成技術(shù)中，1983年Steger等開創(chuàng)的嵌套網(wǎng)格技術(shù)，利用多塊子域網(wǎng)格對計算域進行全局覆蓋，在不同網(wǎng)格重疊部分采用信息共享的方法來實現(xiàn)信息交換，并不強制要求網(wǎng)格之間進行無縫拼接，極大減輕了子域網(wǎng)格生成的難度，而且能夠保證生成的子域的網(wǎng)格品質(zhì)，因此嵌套網(wǎng)格技術(shù)在實際工程的復(fù)雜外形的數(shù)值模擬和研究實驗當中得到了廣泛的應(yīng)用。

本文將計算懸停狀態(tài)下上下旋翼不同總距下共軸式剛性雙旋翼直升機的氣動載荷。上下旋翼由2片NACA 0012翼型槳葉構(gòu)成，旋翼半徑為8.53 m，弦長0.41 m，槳葉無扭轉(zhuǎn)和尖削，NACA 0012翼型截面形狀如圖1所示，雙旋翼算例模型的具體參數(shù)如表1所示。

圖1 NACA 0012翼型截面示意圖

表1 NACA 0012翼型參數(shù)

本文生成的運動嵌套網(wǎng)格系統(tǒng)由兩個層次的網(wǎng)格組成：一是圍繞旋翼表面生成的近似正交的貼體非結(jié)構(gòu)網(wǎng)格，其與旋翼固定連接并隨其同步周期性旋轉(zhuǎn)運動，用于近場尾渦的捕捉和模擬旋翼表面附近的流動；二是將旋翼表面貼體網(wǎng)格包圍的背景網(wǎng)格，它是靜止不動的，屬于笛卡爾網(wǎng)格，模擬旋翼遠場流動特征及流場結(jié)構(gòu)。為了能夠準確地捕捉槳尖渦，對背景網(wǎng)格進行3層加密處理，并將背景網(wǎng)格擴展至較大的計算域以減小遠場邊界的影響。非結(jié)構(gòu)運動嵌套網(wǎng)格系統(tǒng)劃分結(jié)果如圖2所示。

圖2 雙旋翼的非結(jié)構(gòu)運動嵌套網(wǎng)格系統(tǒng)

1.2 控制方程與數(shù)值求解方法

為了高效模擬非定常共軸式剛性雙旋翼氣動載荷[12]，本文利用有限體積法，基于STAR-CCM+商業(yè)軟件求解非定常不可壓Reynolds時均Navier-Stokes方程(URANS)。湍流模型采用一方程Spalart-Allmaras模型，S-A模型為航空旋轉(zhuǎn)機械領(lǐng)域廣泛運用來計算URANS粘性系數(shù)的，比零方程模型計算量更大但結(jié)果更加精確，比兩方程模型計算量更小但穩(wěn)定性更好，不需要網(wǎng)格處理得特別精細，更適合于非結(jié)構(gòu)網(wǎng)格[13]。

固定坐標系下的N-S方程如下：

(1)

式中：W為流場守恒變量，可以表示為W=[ρρuρvρwρE]T；F(W)為無黏通量；G(W)為黏性通量，分別表示為：

(2)

(3)

(4)

式中：T為流體溫度；k為熱傳導(dǎo)系數(shù)。本文采用有限體積法求解上述主控方程，時間離散上采用雙時間法，偽時間步長為多步迭代格式，湍流模型為Spalart-Allmaras模型；邊界條件中旋翼表面定義為無滑移壁面邊界條件，背景邊界則為壓力出口。

懸停狀態(tài)下上下旋翼轉(zhuǎn)速相等為300 r/min，轉(zhuǎn)向相反，定義上下旋翼貼體網(wǎng)格區(qū)域為旋轉(zhuǎn)體，轉(zhuǎn)速和方向與旋翼保持一致。在STAR-CCM+對求解模型進行選擇設(shè)置完成后，即可開始計算求解，具體設(shè)置如表2所示。

表2 求解模型設(shè)置

2 算例及結(jié)果分析

共軸式剛性雙旋翼直升機具有上下布置的兩副旋翼，兩副旋翼轉(zhuǎn)速相等轉(zhuǎn)向相反，可以抵消反力矩，充分釋放前行槳葉的升力潛能，大幅提高直升機的飛行速度。但是共軸式雙旋翼流場中存在較單旋翼流場更嚴重的槳-渦干擾和渦-渦干擾等氣動干擾現(xiàn)象，上旋翼產(chǎn)生的下洗流改變了下旋翼迎風(fēng)角，在上下旋翼總距相同的條件下，上下旋翼氣動特性存在差異，導(dǎo)致上下旋翼間力矩不平衡，影響直升機飛行穩(wěn)定性。因此，探尋上下旋翼力矩平衡方法是開展共軸式剛性雙旋翼直升機正向設(shè)計的必要過程。

調(diào)整上下旋翼的總距角是減小上下旋翼力矩差的關(guān)鍵手段，共軸式剛性雙旋翼系統(tǒng)的本質(zhì)是利用轉(zhuǎn)向相反的上下旋翼來平衡反力矩，但是因為上下旋翼間氣動干擾的存在，上下旋翼受到的氣動力不同，使得反力矩?zé)o法平衡，這時只能通過上下旋翼的總距差動來平衡上下旋翼間由于氣動干擾帶來的力矩差。

下面算例中以不同的上下旋翼總距角來分析共軸式雙旋翼的氣動載荷，探究總距角和共軸式剛性雙旋翼上下力矩與拉力的關(guān)系，來獲得符合要求的總距差，各個算例旋翼的總距角情況如表3所示。

表3 算例情況

在對所有算例進行仿真分析前，本文基于算例1工況，共計開展了6套不同網(wǎng)格數(shù)量級的數(shù)值計算，單片槳葉網(wǎng)格數(shù)量范圍為60 712～110 856，且每套的加密區(qū)范圍和尺寸相同。發(fā)現(xiàn)單片槳葉網(wǎng)格量增加到92 057后，拉力系數(shù)已經(jīng)收斂，不再受到網(wǎng)格數(shù)量的影響。具體的上旋翼拉力系數(shù)隨網(wǎng)格數(shù)量的變化情況如圖3所示。因此，在確保網(wǎng)格無關(guān)的前提下提升數(shù)值仿真效率，同時考慮工作站的計算能力和計算時間的問題，所以最終選擇了單片槳葉為10萬左右的網(wǎng)格數(shù)量進行計算。

圖3 網(wǎng)格無關(guān)性驗證

下面依次對每個算例進行數(shù)值仿真計算，待迭代計算收斂，得到穩(wěn)定旋翼拉力和力矩周期解，取旋翼旋轉(zhuǎn)一周的數(shù)據(jù)進行比較分析。

圖4—圖6為不同總距差下，共軸式剛性雙旋翼在懸停狀態(tài)下上旋翼和下旋翼的拉力變化曲線對比圖。

圖4 上旋翼拉力變化曲線

圖5 下旋翼拉力變化曲線

圖6 上下旋翼拉力和變化曲線

由圖可見，在上下旋翼槳葉形狀大小完全相同的情況下，隨著上下旋翼總距差的改變，上下旋翼拉力的變化趨勢一致，呈周期性變化，其主波動周期大約為90°，拉力幅值和波動幅度有所改變，呈“載荷效應(yīng)”“厚度效應(yīng)”與“渦誘導(dǎo)效應(yīng)”[14]。

保持下旋翼總距角不變而增大上旋翼總距角，總距差增大，而翼型剖面存在一定的厚度和迎角，因此上旋翼下表面與下旋翼上表面間距隨上旋翼總距角的增大而縮小，雙旋翼間的氣體流速增大，從而使上旋翼的下表面與下旋翼的上表面壓強減小，導(dǎo)致上旋翼拉力減小、下旋翼拉力增大，如圖4與圖5所示，為雙旋翼的“厚度效應(yīng)”。此外，算例4相較于其他算例，上旋翼氣動環(huán)境更為良好，該算例下的上旋翼拉力有所提升。下旋翼浸沒于上旋翼產(chǎn)生的下洗流中，氣動環(huán)境較上旋翼更為惡劣，下洗流引起下旋翼的有效工作迎角減小，因此下旋翼拉力通常小于上旋翼拉力。上下旋翼拉力波動幅度差異，主要由旋翼間“載荷效應(yīng)”和“厚度效應(yīng)”共同影響造成[15]。上下旋翼在相對旋轉(zhuǎn)時，每個周期內(nèi)的相互接近和遠離過程中，旋翼槳葉前緣的附著渦對另一片旋翼槳葉產(chǎn)生的上洗作用會轉(zhuǎn)變?yōu)樾順~后緣的附著渦對另一片旋翼槳葉產(chǎn)生的下洗作用，從而導(dǎo)致上下旋翼的相對迎角和拉力表現(xiàn)為先增大后減小的周期性波動，這種跟載荷相關(guān)的干擾現(xiàn)象則被稱為“載荷效應(yīng)”。在“厚度效應(yīng)”的作用下，雙旋翼拉力波動可以互相抵消削弱一部分，因此上下旋翼拉力和的波動幅度會有所減小。圖6中，上下旋翼的拉力和也呈周期性變化，此模型的直升機重量為5 T，除了算例1和算例2，其他算例的拉力和均能滿足重量的要求，同時還能提供一部分余量。

圖7是算例3、算例4和算例5的上旋翼槳葉上下表面在方位角為360°時的壓強分布云圖，從圖中可知，上旋翼上表面的壓力隨著總距差的增大而減小，由于算例4中上旋翼氣動條件更為良好，其上表面的負壓大于算例3和算例4，而下表面壓力場與其他算例無明顯區(qū)別，因此算例4上旋翼拉力更大。從圖7中還能看出上面分析的旋翼“厚度效應(yīng)”帶來的影響，隨著總距差的增大，上下旋翼間距減小，上旋翼的低壓區(qū)減小，使得旋翼上下表面壓差變大。

圖7 上旋翼槳葉上下表面壓力分布云圖

上下旋翼轉(zhuǎn)向相反，取上下旋翼力矩的絕對值分析旋翼力矩變化規(guī)律，如圖8—圖10所示。上下旋翼的力矩隨總距差的變化規(guī)律與旋翼拉力變化情況相近，在此則不加贅述。

圖8 上旋翼力矩變化曲線

圖9 下旋翼力矩變化曲線

圖10為不同總距差下，上下旋翼間力矩差(上旋翼力矩絕對值-下旋翼力矩絕對值)的變化情況，上下旋翼的力矩差值隨總距差的增大而先減小再增大，若總距差大于-0.5°這一臨界值，則力矩差為負。上述現(xiàn)象主要是“渦誘導(dǎo)效應(yīng)”引起，由于雙旋翼尾跡相互誘導(dǎo)，下旋翼受上旋翼產(chǎn)生的下洗流干擾，其有效工作迎角變小，若要達到上下旋翼力矩平衡，應(yīng)增大下旋翼的總距角或者減小上旋翼的總距角來彌補力矩差，上下旋翼總距差應(yīng)由翼型參數(shù)確定。

圖10 上下旋翼力矩差變化曲線

圖11給出了懸停狀態(tài)下算例3、算例4和算例5的上下旋翼槳葉在方位角為360°時的誘導(dǎo)速度場(Vj)圖，從中可以看到明顯的“渦誘導(dǎo)效應(yīng)”，其下旋翼槳盤的誘導(dǎo)速度明顯大于上旋翼槳盤的誘導(dǎo)速度，由于上、下旋翼的相互誘導(dǎo)作用，上旋翼受到的上洗流也明顯小于下旋翼受到的下洗流，在旋轉(zhuǎn)速度一定的情況下其下旋翼的實際工作迎角變小，誘導(dǎo)功率更大，導(dǎo)致功率損失更大，所以必須增大下旋翼總距[16]。

圖11 雙旋翼槳葉誘導(dǎo)速度場

一般的，在懸停狀態(tài)下共軸式剛性雙旋翼的上下旋翼總距角相差在1°以內(nèi)，且上旋翼總距角總小于下旋翼總距角[16-18]。根據(jù)之前的分析，在上下旋翼相同的功率下，上旋翼的下洗流對下旋翼的干擾比下旋翼的上洗流對上旋翼的干擾更加強烈，正常工作時，下旋翼的誘導(dǎo)速度等于自身產(chǎn)生的誘導(dǎo)速度與上旋翼產(chǎn)生的下洗速度之和，下旋翼的誘導(dǎo)速度增大導(dǎo)致下旋翼槳葉剖面的有效工作迎角減小，若上下旋翼采用相同的總距角，則下旋翼產(chǎn)生的拉力和力矩更小，所以，必須使下旋翼的總距角大于上旋翼的總距角，才能平衡直升機的反力矩[19]。所以圖12和圖13是在5個算例的基礎(chǔ)之上，插值擬合出的總距差在-1°～0°，上下旋翼平均拉力和力矩相關(guān)參數(shù)的變化情況。

上下旋翼平均拉力與力矩表達式：

(5)

(6)

式中：Fu(i)與Fd(i)分別為i時刻上下旋翼拉力；n為一個周期內(nèi)的步長數(shù)；Tu(i)與Td(i)分別為i時刻上下旋翼力矩。

圖12中，隨著上旋翼總距角的增大，上旋翼產(chǎn)生的下洗氣流方向也因此改變，下旋翼氣動條件改善，有效迎風(fēng)角增大，因此下旋翼拉力逐步提升。上下旋翼的總距差增大，在總距差為 -0.25°時，由于上旋翼周圍的氣動條件更優(yōu)，表現(xiàn)為算例4上旋翼產(chǎn)生的拉力高于其他算例。

圖12 上下旋翼拉力相關(guān)參數(shù)變化曲線

本文以上下旋翼力矩差系數(shù)來描述上下旋翼力矩在整個周期內(nèi)的差異程度，“力矩差系數(shù)”越大，說明上下旋翼力矩波動越明顯，雙旋翼系統(tǒng)越不穩(wěn)定。力矩差系數(shù)計算公式如下所示：

(7)

從圖13中可以發(fā)現(xiàn)，隨著總距差的變化，上下旋翼的力矩變化和拉力的變化情況類似，上下旋翼的力矩差系數(shù)存在最小的極值點。上下旋翼周期內(nèi)平均力矩有3個交點，分別是A(-0.6,23 140)、B(-0.4,24 925)和C(-0.1,27 370)點，代表周期內(nèi)的上下旋翼力矩平衡。3點中，B點對應(yīng)的力矩差系數(shù)P(420.67)最小，選取力矩平衡點B的上下旋翼總距角作為力矩平衡的總距角。

圖13 上下旋翼力矩相關(guān)參數(shù)變化曲線

以B點算例所代表的上下旋翼的總距角對上面的力矩平衡分析結(jié)論進行仿真驗證，由圖14和圖15可以看出，在上下旋翼總距角為5.6°和6.0°的情況下，基本與上述分析吻合，實現(xiàn)上下旋翼的力矩平衡。懸停狀態(tài)下，在一個穩(wěn)定的運行周期內(nèi)，上旋翼產(chǎn)生了更大的拉力，而上下旋翼產(chǎn)生的力矩相當，計算得到上旋翼的平均力矩為27 168 N·m，下旋翼的平均力矩為27 044 N·m，平均力矩差為124 N·m，上下旋翼的平均力矩差僅為上下旋翼平均力矩的0.46%，力矩差系數(shù)為470。

圖14 B點算例上下旋翼拉力變化曲線

圖15 B點算例上下旋翼力矩變化曲線

圖16給出了上下旋翼總距角分別為5.6°和6.0°時計算得到的氣動性能結(jié)果，圖中可見，上下旋翼的氣動載荷的計算值均較為正常，可以進一步將其作為直升機傳動系統(tǒng)正向設(shè)計的輸入。其中CT為共軸雙旋翼拉力系數(shù)、CQ為共軸雙旋翼扭矩系數(shù)，定義如下：

圖16 配平狀態(tài)下雙旋翼的氣動性能

(8)

(9)

式中：Fup和Flow分別為上、下旋翼拉力；Tup和Tlow分別為上、下旋翼力矩；ρ為海平面空氣密度；Ω為旋翼轉(zhuǎn)速；R為旋翼半徑。

3 結(jié)論

1) 在共軸式剛性雙旋翼懸停狀態(tài)下，保持下旋翼總距角不變，隨著共軸式雙旋翼的總距差的增大，由于旋翼的“厚度效應(yīng)”，上旋翼的拉力與力矩會先減后增再減小，下旋翼的拉力和力矩會一直增大。

2) 由于旋翼“載荷效應(yīng)”和“厚度效應(yīng)”的共同影響，雙旋翼在相對旋轉(zhuǎn)過程中，上下旋翼的拉力均會出現(xiàn)上下波動，且波動方向相反的現(xiàn)象，而且上旋翼的總距角越大，受“載荷效應(yīng)”的影響越大。

3) 由于“載荷效應(yīng)”與“渦誘導(dǎo)效應(yīng)”的存在，共軸式剛性雙旋翼拉力和力矩呈周期性變化，而且在變化過程中，下旋翼處在上旋翼的下洗流中，使下旋翼的實際工作迎角變小，下旋翼的氣動性能降低，導(dǎo)致下旋翼的拉力比上旋翼小。

4) 根據(jù)本文的仿真分析與力矩平衡計算，在懸停狀態(tài)下，上下旋翼力矩平衡時，下旋翼的總距角要大于上旋翼，本文中共軸式剛性雙旋翼模型懸停下的最佳總距角為上旋翼5.6°，下旋翼6.0°，此時不僅能夠?qū)崿F(xiàn)上下旋翼的力矩平衡，而且其力矩波動系數(shù)最小。