直升機旋翼翼型的非定常氣動特性計算方法與驗證研究

孟微，胡和平，周云

中國直升機設計研究所 直升機旋翼動力學重點實驗室，江西 景德鎮 333001

直升機旋翼技術是直升機設計的關鍵，其振動載荷的準確分析對旋翼動力學設計與直升機減振有著至關重要的作用。旋翼振動載荷預估技術包括非定常氣動力建模、結構動力學建模和兩者之間的耦合分析研究等，其中翼型的非定常氣動特性建模是該技術的關鍵和難題之一。翼型非定常氣動特性的準確分析，對提高旋翼振動載荷分析水平、提升直升機旋翼設計能力有極為重要的作用。

翼型來流在準定常狀態下，隨著迎角的增加會從附著流狀態轉為氣流分離，在非定常狀態下，隨著迎角的增大則會產生動態失速。該現象是限制直升機飛行性能、引起直升機振動的因素之一，也是確定旋翼總升力、推力和使用限制的主要因素[1]。然而，準確預測動態失速現象對旋翼載荷和性能的影響是十分困難的，目前國內外對該現象的研究仍以基于試驗的數值分析為主[2，3]。雖然用計算流體力學（CFD）模擬動態失速問題已取得了一定的突破，但因其對計算機的硬件要求高、計算量大、計算周期長等原因仍不適用于目前的工程計算。對于工程中旋翼的設計和分析來說，國際上通常采用與試驗相結合的半經驗模型，該類方法的優點是計算效率高且準確性好[3，4]。

本文以Leishman-Beddoes動態失速模型為基礎，進行適用于旋翼振動載荷分析的二維翼型非定常氣動特性建模與驗證研究。將原動態失速模型中的10個經驗系數簡化為4個，并分別通過翼型靜態試驗和動態試驗數據獲得非定常模型參數表。最后通過試驗數據對本模型進行驗證。通過與試驗的相關性分析表明，本模型具有較高的計算精度，可最終應用于旋翼的振動載荷分析中。

1 翼型的非定常氣動模型

本文二維翼型非定常氣動模型以Leishman-Beddoes動態失速模型為基礎對其進行了簡化。本模型的最大特點在于適用于時域求解，對于整個非定常氣動問題的物理表示更加完善，方法簡單直觀，相比于其他半經驗模型涉及的經驗系數較少（僅為4個）。整個模型系統中使用的18個非定常氣動模型參數，均從翼型特性的靜態、動態試驗數據中獲得。

翼型來流在準定常狀態下，小迎角時為附著流狀態，迎角的進一步增大則轉為氣流分離（隨翼型厚度的增加分為薄翼型分離、前緣分離和后緣分離）狀態；在非定常狀態下，迎角在超過靜失速迎角后會產生動態失速。本模型根據翼型來流的物理過程，基于翼型的準定常非線性建模進行非定常氣動建模，包括翼型附著流、氣流分離和動態失速三部分，即三個子模型的綜合。對三個子模型附著流、氣流分離和動態失速及模型間的綜合與參數確定進行詳細介紹。

1.1 附著流計算

在本子模型中，采用經典的指數響應計算方法進行求解，其指數函數分為隨時間衰減的脈沖函數（非環量部分）和經過幾個弦長時間歷程后趨于穩態解的漸進函數（環量部分）兩類。目前的激勵為迎角和變距率兩項，如果假設翼型系統是線性的，那么根據疊加原理，通過Duhamel積分就可以得到本激勵下的氣動力響應。

本模型中，法向力系數Cn為：

其中，各項遞推公式為：

其中，CNα為法向力系數曲線斜率，c為弦長，a為聲速，V為翼型剖面速度。力矩系數Cm為：

遞推公式為：

遞推公式中變量為

模型簡化：本模型中將以上模型計算所需的Kan、Kqn、Kam、Kqm統一簡化為 Ki：

有 Kan=Ki，Kqn≈ Ki，法向力系數 Cn不變，力矩系數 Cm則轉化為：

遞推公式為：

因此，本模型指數函數計算中所需的10個經驗系數A1～A5，b1～b5簡化為以下4個，取值為：A1=0.3，A2=0.7，b1=0.14，b2=0.53。

阻力系數Cd為：

式中：Cd0為零升阻力系數，弦向力系數Cc為：

1.2 氣流分離

盡管當翼型迎角隨時間變化時，將臨界前緣壓力和壓力梯度作為判斷產生靜態失速的準則不再適用，但是仍然可以通過與前緣壓力直接相關的前緣氣流分離的臨界法向力系數CN1來判斷翼型是否產生氣流分離，該系數可由靜態翼型試驗得到。在非定常情況下，該模型對法向力系數做一階滯后補償處理，引入表明前緣分離的時間常數Tp，該參數是馬赫數的函數，隨翼型變化不大。分離點的確定則采用Kirchhoff-Helmholtz理論推導的后緣分離模型。該模型將翼型的法向力和迎角與后緣分離點聯系起來，與翼型的靜態升力失速特性試驗結果相結合，從而確定等效分離點f（相對于弦長無量綱化）的變化，引入判斷后緣分離的時間常數Tf。

本模型中，前緣分離修正的法向力系數為：

非定常情況下出現前緣氣流分離的條件是Cnp> CN1，而后緣分離的法向力系數Cnf可近似為：

其中：

系數S1和S2定義了翼型的靜態失速特性，有效迎角α1=Cn′/CNα，α1為分離點 f=0.7 時的迎角，f隨迎角及 S1、S2和α1隨馬赫數的變化關系都可以通過靜態升力數據得到。

力矩系數Cnf為：

式中：K0、K1、K2都可以通過靜態升力數據Cm/Cn曲線擬合得到。m對于NACA0012翼型取值為2，對于其他翼型，可以根據翼型數據選擇0.5或1。

弦向力系數Ccf為：

式中：η為根據翼型系數引入的恢復因子，一般η=0.95，對于無黏流η=1。

1.3 動態失速

在非定常狀態下，翼型來流隨迎角的增大會產生動態失速現象。動態失速最大的特點為氣流在翼型前緣產生集中渦并向后緣遷移至完全脫離后流入尾流中。雖然集中渦停留在翼型上時升力有一定的增加，但該集中渦極不穩定，很快脫離翼型使壓心迅速后移，產生非常大的低頭力矩，并增加了槳葉的扭轉載荷，為其影響分析帶來困擾。

本子模型采用公式模擬動態失速的物理過程，即渦在前緣的積累、形成集中渦、向后緣遷移和脫離翼型的整個過程。引入了渦衰減的時間常數Tv和渦沿弦長傳播的時間常數Tvl，當生成渦的時間參數τv=0時分離點產生，τv=Tvl時渦到達后緣。Tv和Tvl在較大的馬赫數范圍內變化不大，且一般認為翼型對兩者的影響不大。

在動態失速模型中，法向力系數為：

其中，動態失速集中渦產生的升力Cv為：

力矩系數Cmv為：

1.4 模型的綜合

對以上三個子模型中各項系數進行綜合，本方法的最終結果為：

（1）總法向力系數：

（2）總力矩系數：

（3）總弦向力系數：

式中：DF為弦向力修正系數。

（4）總阻力系數：

（5）總升力系數：

通過以上分析，可以看出三個子模型間互相耦合。附著流模型的輸出結果作為氣流分離模型的輸入，而氣流分離模型的輸出結果作為動態失速（集中渦流出）模型的輸入。附著流模型可以獨立求解，但氣流分離與動態失速模型之間因氣流是否附著、集中渦是否脫離而產生的相互影響則要通過對時間常數的修正加以體現。為了擴大對于其他翼型的適用性，一般這種時間常數的修正僅限于Tf和Tv兩個量值。為確保氣流分離模型與動態失速模型中使用的時間常數準確，在每一次迭代步中時間常數都要重新修正。本模型的流程圖如圖1所示。

1.5 模型參數的確定

本模型中使用的4個經驗系數在前面已經介紹。以下18個非定常氣動模型參數根據翼型試驗數據獲得（見表1）。其中，前1～14項通過靜態試驗數據（翼型C81表）獲得，后4個（15～18）表示動態失速特征的時間常數項通過動態試驗數據獲得，這4個時間常數項隨馬赫數變化較大，但隨翼型的變化較小，本模型中在NACA0012翼型的時間常數項的基礎上，根據本文中選取翼型的動態試驗結果做微小修正后確定。

圖1 非定常氣動模型流程圖Fig.1 Unsteady aerodynamic model fl ow

表1 非定常模型參數Table 1 The unsteady model parameters

2 計算驗證

本文中的試驗數據采用2011年俄羅斯中央空氣流體動力學研究院（TsAGI）在TsAGI SVS-2風洞中進行的翼型動態失速特性試驗數據，該翼型弦長為0.18m。在模型的驗證中，分別進行Ma=0.3和Ma=0.6兩個不同馬赫數狀態的計算與試驗結果對比分析，驗證本文的翼型非定常氣動模型具有較好的計算精度和可靠性。

2.1 馬赫數Ma=0.3

迎角α隨時間t的變化分為三段。第一段迎角范圍0°～6.9°，對應圖中試驗1和計算1；第二段迎角范圍6.5°～14.6°，對應圖中試驗2和計算2；第三段迎角范圍11.7°～22°，對應圖中試驗3和計算3。本狀態馬赫數Ma=0.3，折合頻率k為0.029。如圖2～圖4所示。

圖2 升力系數隨迎角變化（Ma=0.3）Fig.2 The lift coeff i cient varies with the attack angle（Ma=0.3）

圖3 阻力系數隨迎角變化（Ma=0.3）Fig.3 The drag coeff i cient varies with the attack angle（Ma=0.3）

圖4 力矩系數隨迎角變化（Ma=0.3）Fig.4The torque coefficient varies with the attack angle（Ma=0.3）

從圖中可以看出，在該狀態下計算結果與試驗數據符合地較為理想。升力系數的動態失速現象得以較好地捕捉；阻力系數反映出動態失速時阻力“8”字環的變化；力矩系數與試驗結果符合的非常理想。總體來看，計算結果較為理想，且可以從試驗數據判斷出在該狀態下動態失速迎角在12°附近，與計算結果基本一致。

2.2 馬赫數Ma=0.6

迎角α隨時間t的變化分為三段。第一段迎角范圍0°～6.6°，對應圖中試驗1和計算1；第二段迎角范圍7.4°～15.7°，對應圖中試驗2和計算2；第三段迎角范圍13.9°～23.8°，對應圖中試驗3和計算3。本狀態馬赫數Ma=0.6，折合頻率k為0.0146。如圖5～圖7所示。

圖5 升力系數隨迎角變化（Ma=0.6）Fig.5 The lift coeff i cient varies with the attack angle（Ma=0.6）

圖6 阻力系數隨迎角變化（Ma=0.6）Fig.6 The drag coeff i cient varies with the attack angle（Ma=0.6）

在Ma=0.6的狀態下，整體來看計算結果與試驗數據吻合地較好。動態失速產生的遲滯環小于計算結果，是由于折合頻率較小造成的。升力系數最大值計算結果略大于試驗結果，初步判斷為渦的耗散時間常數Tv過小造成的；力矩系數偏差略大，可能是由于計算結果的動態失速迎角小于試驗結果。根據試驗和計算結果可以判斷本計算中的翼型在該馬赫數下的動態失速迎角為10°。

圖7 力矩系數隨迎角變化（Ma=0.6）Fig.7The torque coefficient varies with the attack angle（Ma=0.6）

3 結論

本文對建立可適用于直升機旋翼的振動載荷分析的翼型非定常氣動模型進行研究，通過模型對模擬物理過程的闡述、模型參數的確定及模型計算與試驗結果的對比分析，最終得到的如下結論：

（1）本模型在保證一定的計算精度的前提下，將原有Leishman-Beddoes模型指數函數中的10個經驗系數簡化為4個，且較于其他模型物理過程表述清晰，方法簡單直觀。

（2）本模型可以捕捉升力、阻力、力矩系數因動態失速產生的遲滯環；且計算結果準確地反映了該翼型在Ma=0.3、Ma=0.6下的動態失速迎角分別為12°、10°。

（3）本模型在不同馬赫數下的計算精度都較為良好，與試驗結果在變化趨勢上一致性高；模型的可靠性也較高，可在旋翼的振動載荷預估分析中進行應用。

[1] Leishman J G. Principles of helicopter aerodynamics[M]. Second Edition. Cambridge： Cambridge University Press， 2006.

[2] Brain C. Modeling dynamic stall of SC-1095 airfoil athigh Mach Numbers[D]. Georgia Institute of Technology in Partial Fulfillment of the Requirements for the Degree Master of Science in the School of Aerospace Engineering， 2010.

[3] Ericsson L E， Reding J P. The difference between the effects of pitch and plunge on dynamic airfoll stall[C]// The 9thEuropean Rotorcraft Forum， 1983.

[4] Ashish B. Cobtributions to the mathematical modeling of rotor flow-fields using a pseudo-implicit free-wake analysis[D]. The University of Maryland in Partial of the Requirements for the Degree of Doctor of Philosophy， 1995.