二維射流中翼型氣動特性計算與分析

李旭, 周洲, 郭佳豪, 薛臣

(西北工業大學 航空學院, 陜西 西安 710072)

隨著固定翼類垂直起降飛機的發展,動力噴流中機翼的氣動特性越來越受到人們關注。在零來流或小來流狀態下,傳統的固定翼類飛機由于升力不足,無法正常起飛。然而,固定翼類垂直起降飛機憑借動力噴流對機翼的作用,可產生額外的升力,從而實現垂直或短距起降[1-2]。因此,清楚地認識噴流中機翼的氣動特性對此類飛機的安全起降來說非常重要。

動力噴流與升力面的相互作用是較為復雜的空氣動力學問題,國內外學者已經進行了不少研究。龔志斌等[3]采用CFD(computational fluid dynamic)技術對外吹式襟翼的氣動特性進行了計算,指出襟翼環量的增加與動力噴流對襟翼下表面的沖擊和上表面的吹氣加速有關。楊小川等[4]采用激勵盤方法來等效螺旋槳,對比了不同組合形式下分布式螺旋槳對機翼的影響,結果表明不同形式的滑流均能使機翼的升阻力增大。

CFD類方法可以獲取精細的流場細節,但需要求解整個流域,計算花費較大。因此,不少學者提出了一些簡化方法來研究動力噴流對機翼的作用。Cole等[5]采用高階自由尾跡方法對螺旋槳/機翼干擾進行了計算,結果表明該方法對輕載荷的槳有效。Veldhuis[6]采用滑流管模型來模擬滑流,研究了巡航階段螺旋槳對機翼氣動力的影響。但滑流管模型屬于固定尾跡的一種,忽略了滑流的收縮和偏轉,同樣只能適用于自由來流速度大、槳盤載荷小的情況。還有一些學者采用半經驗的噴流速度計算方法[7-8],但應用場景比較局限。

以上研究大多針對巡航階段動力/機翼間的干擾,Shollenberger[9-10]較早對小來流情況下,二維無黏射流中翼型的氣動特性進行了計算。其采用離散渦面來等效射流,基于射流內外總壓差不變假設,進行射流邊界位置的更新,對飛機起飛階段噴流/機翼干擾進行了初步研究。與Shollenberger的求解思路類似,Lewis[11]采用自由流線模型,利用半無限長渦面和有限長離散渦面來等效射流,對零自由來流情況下二維射流的偏轉進行了模擬。Bontempo等[12]利用自由流線模型,對均勻載荷槳盤的收縮效應進行了計算,進一步驗證了自由流線模型的有效性。

真實動力產生的噴流不僅有軸向速度,往往還存在旋轉速度,流動狀態比較復雜。直接進行實際動力的三維計算,時間花費較大,且不易對比分析各部分速度對機翼影響的差異。另外,三維機翼的氣動特性與二維翼型的氣動特性密切相關,工程中常基于二維翼型的氣動數據進行三維機翼特性的計算,如螺旋槳計算采用的片條理論。因此,為了解三維射流對機翼的影響機理,對二維射流中翼型的氣動特性要有較清楚的認識,目前這方面的研究還有待深入。

面元法可看作是一種半解析的方法,能排除黏性的影響,方便變化網格并進行不同影響因素的研究,有利于得到一些規律性的結論。因此,參考文獻[9]的求解思路,基于自由流線模型,利用渦面對射流進行建模,結合一階面元法,本文建立了二維射流中翼型氣動特性的計算方法,分析了射流中影響翼型氣動力的不同因素,以期為進一步認識動力增升原理提供幫助。

1 射流建模

射流內的流體一般比環境流動具有更大的速度,射流與周圍環境流動的邊界屬于自由剪切層,在此邊界上切向速度不連續,這與渦面的特性一致,因此可以用渦面來等效射流流動。

本文建立了二維射流中翼型氣動特性的計算方法,主要假設有:

1) 忽略黏性的影響,不考慮射流的卷吸[13];

2) 流動是二維且不可壓的;

3) 射流內的總壓保持不變。

1.1 無限長渦面

首先對有限長渦面的誘導速度計算進行介紹,設渦面的渦強γ保持不變,規定γ順時針轉為正,逆時針為負，位置在y=0,即x軸上,如圖1所示。

圖1 有限長渦面

根據文獻[14],有限長等強度渦面對平面上一點的誘導速度的積分表達式如下所示：

當渦面長度無限大時,即x1→-∞,x2→+∞,積分可得

(3)

(4)

圖2 上下強度相反的無限長渦面

可以看出,相距為h,上下強度相反的無限長渦面誘導的速度場等價于寬度為h的射流。因此,可采用上下強度相反的無限長渦面來等效射流。當渦強不變的情況下,射流內部的速度不隨寬度h變化。

1.2 右半無限長上下渦面

沿x軸方向對稱分布上下2條半無限長渦面, 半無限長上渦面起點為(xR,yR),下渦面起點為(xR,-yR),渦面軸向與x軸平行,向右直到+∞。上渦面渦強為-γ,下渦面渦強為γ,進行積分計算,可得右半無限長上下渦面對給定點(x,y)的誘導速度如下所示：

(5)

(6)

當-yR

當y>yR時

當y<-yR時

同理可得左半無限長上下渦面誘導的速度場。可以驗證,無限長上下渦面誘導的速度場就等于左半和右半無限長上下渦面誘導速度場的疊加。

1.3 本文采用的射流模型

實際翼型一般距射流出口為有限值,引入有限長無厚度壁面來刻畫射流出口的壁面。另外,翼型存在于射流中,會對射流造成偏轉,影響偏轉部分渦面的渦強,因此引入有限長上下離散渦面來計算射流的偏轉。基于以上考慮,最終采用的射流模型如圖3所示。

圖3 采用的射流模型

射流模型包括左半無限長上下渦面+上下壁面+離散渦面+右半無限長上下渦面。基于此射流模型,本文將對射流中翼型的氣動特性進行計算。

2 計算模型介紹

2.1 射流邊界條件

圖4 射流偏轉示意

對于單獨射流,給定射流內部速度Vw,外部自由來流速度V∞,均沿x軸正向。在不受翼型影響時,射流不發生偏轉,其內部的壓力始終等于環境壓力,可求得射流內外總壓差

翼型存在的情況下,射流內部的壓力和速度將發生變化。基于射流內外總壓差不變假設,對于射流邊界上的任一點均有

由于自由渦面不能承力,在渦面上任一位置有pin=pout,則總壓差可寫為

ΔH=ργsheetVsheet

(7)

(7)式就是本文射流邊界滿足的計算條件,可用于離散渦面渦強的更新。

2.2 自由流線模型

自由流線模型認為射流邊界是一條流線,對于一個離散渦面,在起點坐標固定的情況下,根據當前渦面上各單元的速度,就可以對渦面節點的位置進行更新[11]。

設當前渦面上有N個單元,N+1個節點,每個單元的長度為Δsj(j=0,…N-1),每個單元中心的速度為(uj,vj),起始節點位置固定,則從第二個節點開始,各個節點的坐標可按(8)式進行更新

(8)

采用(8)式進行迭代,當計算過程中渦面節點的位置不再變化時,可認為計算收斂。

2.3 物面計算模型

參考文獻[14],無厚度壁面采用lumped渦模型進行求解,基于法向速度邊界條件,得到壁面上各單元的環量。有厚度翼型則采用等強度渦面元進行計算,采用切向速度邊界條件,求解得到翼型上各單元的渦強。線性方程組形式如下所示

AX=b

A為影響系數矩陣,只和物面本身有關,本文計算中存在3個物面,矩陣階數為3個物面單元數之和。X為方程解,包括壁面上各單元的環量和翼型上各單元的渦強。b為源項,指速度沿物面法向或切向的投影。對自由來流中的翼型而言,源項b只包含自由來流的影響。對本文計算的射流而言,源項b除了自由來流的速度外還需包含射流渦面對物面的誘導速度。

具體的系數矩陣組建過程可參考文獻[14],整個線性方程組采用Gauss消元法求解。

2.4 廣義Kutta-Joukowski定理

經典的Kutta-Joukowski定理(以下簡稱K-J定理)只適合于自由來流情況下單個翼型升力的計算。本文所研究的流場則包含多個自由渦和附著渦,經典的K-J定理不再適用。文獻[15]將傳統的K-J定理進行了推廣,使之能適用于多個自由渦和附著渦的情況。參考文獻[15],本文將采用廣義K-J定理對射流中翼型的氣動力進行計算。

自由來流V∞沿x軸正向,定義y軸正向為升力方向,x軸正向為阻力方向。設翼型上有N個單元,則翼型的升阻力可按(9)式計算

(9)

式中：Vi,x指除翼型外,流場中其他自由渦(射流邊界)和附著渦(無厚度壁面)對翼型第i個面元x方向的誘導速度；Vi,y指這些渦對第i個面元y方向的誘導速度；Γi是翼型上第i個面元的環量。

3 方法驗證

3.1 求解流程

確定了射流和翼型的計算模型后,方法的整個求解步驟如下所示:

1) 確定射流內外總壓差,給定翼型的位置,確定離散渦面的長度和單元數,初始化渦面的位置和渦強。

若射流內部的速度為Vw,自由來流速度為V∞,則渦面上的渦強初始化為

γ0=Vw-V∞

(10)

2) 根據當前渦面的位置和強度求渦面對物面上各單元的誘導速度。

3) 求解線性方程組,更新物面上的渦強和環量。

4) 求物面對渦面,渦面對渦面的誘導速度,基于公式(8),更新渦面節點的位置。

5) 根據渦面上的速度,利用公式(7),更新渦面單元的渦強。

6) 重復步驟2)～5),直到翼型環量或渦面位置不再變化。

圖5 本文方法框架

在計算過程中,左半無限長上下渦面的位置和渦強保持不變,右半無限長上下渦面的渦強不變,起點位置隨離散渦面最后一個節點的坐標變化,有限長離散渦面的渦強和位置均會變化,整個計算過程采用松弛迭代進行。

本文計算包括的幾何參數有射流的寬度h,無厚度壁面的長度L,翼型的弦長c,離散渦面的長度Lw。

選用NACA0012翼型作為研究對象,翼型的旋轉中心定義在四分之一弦長處,在計算坐標系中對應為(xc,yc),翼型的攻角為α。計算坐標系原點定義在射流出口的中心,如圖6所示。

圖6 射流中翼型位置示意

為減小離散渦面上的網格量,可采用節點非均勻分布的方式,在靠近射流出口和翼型的位置上,布置較密的網格。在翼型后方,渦面上節點的間距可逐漸加大。

射流中翼型氣動力系數按(11)式計算

(11)

式中：cl表示升力系數；cd表示阻力系數；cp表示壓力系數。ρ∞和p∞分別指無窮遠處空氣的密度和靜壓；Vw指給定的射流速度。

本文中取ρ∞=1.225 kg/m3,p∞=101 325 Pa,射流速度Vw=30 m/s,射流的寬度h=0.16 m。

3.2 物面求解方法

選取NACA0012翼型,對自由來流中翼型進行計算,并與二維翼型計算開源軟件XFoil的無黏解進行對比, 檢驗本文程序對無黏情況下翼型氣動力計算的有效性。

對于面元法,翼型表面的壓力系數根據文獻[14],可表示為

(12)

式中,γ表示為翼型單元的渦強。

翼型上的離散面元N=256,攻角范圍從-4°～12°,計算得到升力系數曲線如圖7所示。

圖7 無黏計算升力系數對比

本文計算得到的升力系數與XFoil計算的基本符合,表明本文翼型氣動力計算方法是有效的。

3.3 單獨射流模型

在沒有翼型的情況下,對單獨射流模型進行計算,驗證本文建立的射流模型的有效性。

取無厚度壁面長度L=0.32 m,單元數96。離散渦面的長度Lw=4 m,單元數分別取為Nw=180,200,250和300,檢驗離散渦面網格量對射流速度的影響。

自由來流速度取V∞=1 m/s,射流速度保持不變。定義中點在x=0.32 m處,與y軸平行,長度為0.9h的線段,對比線段上的速度隨網格量的變化如圖8所示。

圖8 單元數對射流軸向速度的影響

由圖8可知,對渦面進行離散,會引入離散誤差,因此得到的射流速度與理論值會存在差異。但隨著網格數的增加,誤差逐漸減小。當Nw=200時,軸向速度與理論值的誤差不超過0.2%,以上結果表明本文建立的射流模型是有效的。

3.4 射流中翼型氣動特性計算

由于缺少實驗數據,本文將利用CFD計算結果作為參考,對不同偏角下射流中翼型的氣動力進行計算,檢驗本文方法的有效性。

選用對稱的NACA0012翼型進行計算,翼型的弦長c=0.2 m,旋轉中心取為xc=0.32 m,yc=0 m。離散渦面單元數Nw=300,其余參數與3.3節單獨射流計算相同。 翼型攻角從-4°～20°,角度間隔為4°。自由來流速度V∞=1 m/s,射流速度Vw=30 m/s,均沿x軸正向。

根據下離散渦面最后一個節點的位置變化來判斷收斂,當前后兩次迭代節點位置誤差絕對值小于0.000 1時,認為計算收斂。給出4°和12°偏角時,渦面位置誤差隨迭代次數的變化如圖9所示。

圖9 位置誤差收斂曲線

可以看出,經過多次迭代,渦面的位置最終不再變化,表明本文方法是收斂的。

計算收斂后,給出翼型4°和12°攻角時射流邊界的位置如圖10所示。可以看出,隨著翼型攻角增大,射流方向的偏轉也在變大。

圖10 不同攻角下的射流偏轉對比

CFD計算域遠場邊界分別采用壓力入口和壓力出口邊界。壓力入口給定自由來流的總壓,壓力出口給定自由來流靜壓。射流出口采用壓力入口邊界,根據射流速度計算出口處的總壓,上下壁面則采用無滑移邊界條件。

圖11 CFD計算域

計算域網格數為4萬,通過旋轉翼型,生成不同攻角下的翼型網格。12°攻角時翼型網格如圖12所示。

圖12 CFD計算網格

CFD計算自由來流速度V∞=1.0 m/s,射流速度Vw=30 m/s。翼型采用無滑移邊界條件,保持xcr=0.32 m,ycr=0 m。湍流模型分別選用S-A和SSTk-w模型,計算收斂后對不同方法得到的氣動力系數進行對比。由于無黏計算無法得到摩擦阻力,CFD的阻力將只取壓差阻力。

本文方法以下簡稱為fsm(free streamline model),將采用廣義K-J定理對翼型的升力和阻力進行計算。 fsm方法與CFD計算的氣動力對比如圖13所示。

圖13 不同方法翼型氣動力對比

可以看出,在20°攻角范圍內, fsm與CFD計算的氣動力符合較好,升力系數基本呈線性變化。

以上結果表明,本文的fsm方法可以反映出射流中翼型的氣動特性,可用于射流中翼型氣動力的計算。

4 射流中翼型氣動特性影響因素研究

在保持射流總壓不變的情況下,利用fsm方法,對翼型在射流中不同上下位置,射流寬度與翼型弦長比值等參數進行研究,分析這些參數對射流中翼型氣動特性的影響。

4.1 不同上下位置下翼型氣動力對比

保持旋轉中心xcr=0.32 m不變,改變ycr的取值,研究翼型在射流中的上下位置對其氣動力的影響。令ycr=0,±0.03,±0.05 m,取攻角α=0°,4°,8°,12°對應的升力系數和阻力系數進行對比。

圖14 不同上下位置翼型氣動力對比

總體上看,翼型位置越靠上,升力和阻力系數就越小。值得一提的是,雖然NACA0012是對稱翼型,但在射流中,攻角為0°時翼型的升力也不一定為0。

取攻角α=12°時,翼型在不同上下位置時射流上邊界的偏轉情況對比如圖15所示。

圖15 翼型不同上下位置射流偏轉對比

可以看出,相同攻角下,翼型處于射流中不同位置引起的射流偏角是不同的。翼型越靠下,射流偏轉越明顯,因此氣動力越大。

4.2 不同寬弦比下翼型氣動力對比

將射流寬度與翼型弦長的比值定義為寬弦比表示為RA=h/c。基于fsm方法,保持翼型的弦長和位置不變,通過改變射流寬度,研究寬弦比對射流中翼型氣動特性的影響。

仍然采用NACA0012翼型,翼型的旋轉中心固定為xcr=0.32 m,ycr=0。射流寬度分別取為h=0.16,0.32,3.2 m,對應寬弦比RA=0.8,1.6和16,得到翼型的升力系數對比如圖16所示。

圖16 不同寬弦比下翼型升力系數對比

其中,自由來流指面元法計算的,在自由來流中,NACA0012翼型對應的無黏升力系數曲線。可以看出,翼型在有限寬射流中的氣動特性與無界的自由來流中差異較大。隨著寬弦比的增大,翼型的升力線斜率不斷增大,二者間的差異逐漸變小。

取攻角α=12°,給出不同寬弦比下,射流上邊界的偏轉情況對比如圖17所示。

圖17 不同寬弦比下射流偏轉對比

由圖17可知,隨著寬弦比的增大,射流的偏轉程度逐漸變小,這主要是因為射流邊界與翼型間的距離變大,翼型對射流邊界的影響在逐漸變小。

繼續增大射流寬度,給出α=12°時,射流中翼型阻力系數隨寬弦比的變化如圖18所示。由圖可知,當寬弦比較小時,射流中翼型的阻力隨射流寬度的增大會先變大。當寬弦比RA≥16后,阻力則隨寬弦比增大而不斷減小。可以預見,當寬弦比趨于無窮時,翼型周圍的流動將接近于自由來流,翼型的阻力也將趨于零。

圖18 阻力系數隨寬弦比的變化

4.3 無黏射流中翼型氣動特性分析

為對射流中翼型氣動特性的變化有更深的認識,本文將利用動量理論對射流中翼型的受力特性進行分析。

L=Tsinθ,D=T(1-cosθ)

(13)

可以看出,在零來流情況下,當翼型對射流出口動量沒有影響時,參數F≡T,且與射流偏角θ無關。也就是說,翼型處于不同攻角下,計算出的參數F均是相同的。

(14)

圖19 不同攻角下的參數f對比

下面對這一等式進行驗證。將第3.4節fsm方法計算出的不同攻角下翼型的氣動力系數代入(14)式，計算得到參數f，如圖19所示。在第3.4節中,cT=1.6,可以看出fsm方法計算的f與cT值的相對誤差很小,不超過1%,表明動量理論分析的結果是有效的,射流中翼型氣動力的大小與射流的偏轉密切相關,本文fsm方法可以計算出這一特征。

再來看寬弦比對升力系數的影響,當射流寬度增加時,雖然射流的偏角θ在變小,但由于h增加,射流出口的動量是變大的,二者的乘積最終決定升力的變化。

對翼型阻力而言,從動量理論得到的阻力公式來看,只要射流偏角非零,翼型就會受到阻力。從廣義K-J定理阻力公式來看,射流偏轉后,翼型會受到向下的誘導速度(可看作下洗),因而產生阻力,這與三維有限展長機翼誘導阻力的產生類似。當寬弦比足夠大時,此時翼型受到的下洗速度接近于零,因此翼型阻力也趨于零。

5 結 論

基于自由流線模型,本文發展了一種無黏射流中翼型氣動力的計算方法。在射流內總壓不變的情況下,研究了不同因素對翼型氣動特性的影響,并根據動量理論對射流中翼型的受力特性進行了分析,得出的主要結論有:

1) 翼型的攻角、翼型在射流中上下位置對其氣動力的影響,均可通過射流偏角反映出來,本文的fsm方法能夠計算出這一特征。

2) 在二維有限寬射流中,翼型產生升力的同時也會帶來阻力。

3) 當射流寬度給定時,射流偏角越大,射流中翼型的升力和阻力就越大。

4) 二維有限寬射流中翼型的升力線斜率小于自由來流中翼型的升力線斜率。

5) 攻角不變時,翼型在射流中位置越靠下,射流偏角越大。

6) 當射流寬弦比無限大時,翼型在射流中的氣動特性就接近于自由來流中的情形。無黏情況下,二維翼型的阻力為零。