一種耦合CFD修正的螺旋槳快速設(shè)計方法

郭佳豪，周洲，范中允

西北工業(yè)大學(xué) 航空學(xué)院，西安 710072

螺旋槳在飛行器中有著廣泛的應(yīng)用，為提高飛行器巡航效率，應(yīng)設(shè)計合適且高效的螺旋槳。關(guān)于螺旋槳的設(shè)計方法，國內(nèi)外學(xué)者已有大量研究。Larrabee[1]提出了最小誘導(dǎo)損失螺旋槳的設(shè)計方法。Angelo等[2]提出了根據(jù)給定工作狀態(tài)和槳葉剖面翼型，快速設(shè)計槳葉弦長和扭轉(zhuǎn)角分布的高效率螺旋槳設(shè)計方法。項松等[3-4]在設(shè)計時對槳葉截面弦長進(jìn)行限制，改進(jìn)了Angelo設(shè)計方法獲得的弦長可能過大的缺點，并對設(shè)計的螺旋槳進(jìn)行了風(fēng)洞試驗。翟若岱[5]應(yīng)用項松等改進(jìn)的方法，對臨近空間太陽能飛機(jī)的螺旋槳進(jìn)行了多點設(shè)計。

為進(jìn)一步提高螺旋槳的效率，眾多學(xué)者引入優(yōu)化設(shè)計的思想，對螺旋槳槳葉的弦長及扭轉(zhuǎn)角分布進(jìn)行優(yōu)化。許建華等[6]以求解雷諾平均Navier-Stokes方程為螺旋槳?dú)鈩有阅苡嬎惴椒ǎ云拭嬉硇蛪毫Ψ植紴槟繕?biāo)，對螺旋槳平面形狀進(jìn)行反設(shè)計。梁撐剛等[7]將Betz條件和遺傳算法相結(jié)合，對無人機(jī)螺旋槳槳葉弦長及扭轉(zhuǎn)角分布進(jìn)行優(yōu)化。焦俊等[8]對高空飛艇螺旋槳進(jìn)行了弦長和扭轉(zhuǎn)角優(yōu)化設(shè)計，并進(jìn)行車載試驗。劉芳等[9]采用Xfoil對翼型性能快速求解，首先進(jìn)行螺旋槳翼型優(yōu)化，之后對螺旋槳扭轉(zhuǎn)角分布和弦長分布進(jìn)行優(yōu)化，得到高空低雷諾數(shù)高效螺旋槳。而Kwon等[10]提出一種螺旋槳多層次優(yōu)化設(shè)計方法：首先對螺旋槳弦長和扭轉(zhuǎn)角進(jìn)行優(yōu)化，得到初始構(gòu)型；接著對槳葉不同站位的翼型進(jìn)行優(yōu)化設(shè)計，以提高效率；翼型優(yōu)化后，再優(yōu)化弦長和扭轉(zhuǎn)角，使得效率進(jìn)一步提高。

此外，Morgado等[11]基于JBLADE軟件，分別采用最大L/D及L3/2/D進(jìn)行螺旋槳設(shè)計，結(jié)果表明采用最大L3/2/D設(shè)計處的螺旋槳能夠產(chǎn)生更多的升力，且效率更高，其中:L為升力，D為阻力。Yonezawa等[12]針對火星飛行器，對低雷諾數(shù)螺旋槳進(jìn)行設(shè)計，并討論了螺旋槳半徑、槳葉個數(shù)及前進(jìn)比對設(shè)計結(jié)果的影響。

上述螺旋槳設(shè)計方法均基于葉素動量理論及渦流理論，不能真實且準(zhǔn)確地體現(xiàn)如徑向流動、槳葉槳根及槳尖效應(yīng)等。而對于太陽能無人機(jī)等飛行器應(yīng)用的螺旋槳，其主要面臨低雷諾數(shù)的問題，低雷諾數(shù)下槳葉葉素氣動力計算的準(zhǔn)確性及槳葉真實三維狀態(tài)與設(shè)計間的差異，均會對設(shè)計的準(zhǔn)確性產(chǎn)生較大的影響。一方面會導(dǎo)致螺旋槳真實拉力偏離設(shè)計拉力；另一方面，設(shè)計狀態(tài)下雖然能夠保證螺旋槳效率達(dá)到最高，但在真實情況下螺旋槳的高效率并不能夠保證。而若采用CFD求解對螺旋槳進(jìn)行優(yōu)化設(shè)計，由于單次螺旋槳數(shù)值模擬求解便需一定的計算時間，導(dǎo)致設(shè)計周期過長。

為解決上述問題，本文在設(shè)計過程中耦合數(shù)值模擬求解，并以CFD計算結(jié)果為參考進(jìn)行槳葉氣動力反解，以此進(jìn)行修正與再設(shè)計，建立一種螺旋槳快速設(shè)計方法，在保證拉力要求下，提高螺旋槳效率。

1 傳統(tǒng)螺旋槳設(shè)計方法

傳統(tǒng)螺旋槳設(shè)計方法基于最小能量損失[13]。圖1展示了螺旋槳槳葉葉素受力分析，其中，α為葉素實際迎角，β為干涉角，φ0為幾何入流角，φ為實際入流角，W為合速度，W1為誘導(dǎo)合速度，Va1為軸向誘導(dǎo)速度，Vt1為環(huán)向誘導(dǎo)速度，V′為誘導(dǎo)螺距。這里假設(shè)合速度W與誘導(dǎo)合速度W1垂直。

圖1 螺旋槳葉素受力分析

給定設(shè)計高度H、來流速度V0、密度ρ、設(shè)計拉力T及設(shè)計轉(zhuǎn)速Ω，通過設(shè)計得到各葉素合適的翼型、弦長b及實際迎角α。設(shè)計過程如下：

1) 確定槳葉數(shù)NB，螺旋槳半徑R，將槳葉劃分為N個截面，并確定各截面翼型。

2) 根據(jù)設(shè)計拉力T及設(shè)計狀態(tài)求解誘導(dǎo)螺距V′。

最小能量損失下誘導(dǎo)螺距V′沿槳葉為定值[13]。根據(jù)渦流理論可知，槳盤徑向r處環(huán)向誘導(dǎo)速度與當(dāng)?shù)丨h(huán)量的關(guān)系為

(1)

寬度為dr的葉素產(chǎn)生的拉力dT為

dT=ρΓ(r)(Ωr-Vt1)dr

(2)

由式(1)、式(2)和圖1的幾何關(guān)系可得

(3)

則總拉力T為

(4)

其中：f為Prandtl修正因子，即

(5)

故可根據(jù)設(shè)計狀態(tài)及設(shè)計拉力，采用Newton迭代求解式(4)得到誘導(dǎo)螺距V′，從而確定各個截面的受力幾何關(guān)系，進(jìn)而得到考慮修正的螺旋槳最佳環(huán)量分布：

(6)

3) 求解r處截面的弦長b及扭轉(zhuǎn)角θ。

根據(jù)庫塔-茹科夫斯基定理，葉素升力與當(dāng)?shù)丨h(huán)量關(guān)系為

dL=ρΓ(r)Wdr

(7)

而升力又可表示為

(8)

聯(lián)立式(7)和式(8)可得

(9)

式中：CL為升力系數(shù)。

葉素工作在最大升阻比下，螺旋槳的效率最高。而根據(jù)當(dāng)?shù)叵议Lb，可得到對應(yīng)的雷諾數(shù)Re與馬赫數(shù)Ma，進(jìn)而求得葉素截面翼型最大升阻比對應(yīng)的迎角α及氣動力。根據(jù)式(9)迭代求解得到弦長b。對應(yīng)截面的扭轉(zhuǎn)角為

θ=α+φ

(10)

限制弦長范圍取為0.02R～0.3R，采用“黃金分割法”求取最接近最佳環(huán)量分布的弦長。弦長超過最大限制，則考慮增加槳葉數(shù)目。

4) 重復(fù)求解，得到所有截面幾何信息。

5) 對弦長b及扭轉(zhuǎn)角θ分布進(jìn)行光順處理。

采用四階貝塞爾樣條曲線對弦長b及扭轉(zhuǎn)角θ分布進(jìn)行擬合:

p=(1-t)4p0+4(1-t)3tp1+

6(1-t)2t2p2+4(1-t)t3p3+t4p4

(11)

式中：t取值范圍為0～1。

沿槳葉徑向進(jìn)行積分，可得設(shè)計下的螺旋槳總拉力T、扭矩M及效率η分別為

(12)

(13)

(14)

(15)

(16)

式中：CD為阻力系數(shù)。

2 數(shù)值模擬方法驗證

對于螺旋槳?dú)鈩恿Φ那蠼猓啾扔趧恿咳~素理論的快速求解方法，CFD方法基于真實模型計算更加接近真實物理狀態(tài)，已發(fā)展出如多重參考系法、滑移網(wǎng)格及嵌套網(wǎng)格模擬等多種方法。其中多重參考系(Multiple Reference Frames，MRF)法作為一種準(zhǔn)定常方法，計算效率高，在螺旋槳數(shù)值模擬中應(yīng)用廣泛[14-16]。

這里以某型太陽能無人機(jī)螺旋槳為例，對其分別進(jìn)行地面試驗與數(shù)值模擬。該型螺旋槳槳葉半徑R=0.203 2 m，槳葉數(shù)NB=2。地面試驗在無來流狀態(tài)下進(jìn)行，通過改變轉(zhuǎn)速測量螺旋槳拉力及扭矩的變化。

數(shù)值模擬基于MRF法，采用SA湍流模型，通過混合網(wǎng)格技術(shù)進(jìn)行網(wǎng)格生成。其中旋轉(zhuǎn)域采用四面體非結(jié)構(gòu)網(wǎng)格，網(wǎng)格量為260萬，外流域生成結(jié)構(gòu)網(wǎng)格，網(wǎng)格量為50萬。局部網(wǎng)格如圖2所示。

拉力和扭矩測量結(jié)果(EXP)與數(shù)值模擬結(jié)果(CFD)對比如圖3和圖4所示。可以看出，在較小拉力下，數(shù)值模擬結(jié)果與試驗測試結(jié)果符合較好，驗證了數(shù)值模擬方法的可行性，且相比于葉素動量理論，拉力的計算更加準(zhǔn)確。

圖2 局部網(wǎng)格示意圖

圖3 拉力對比

圖4 扭矩對比

3 基于傳統(tǒng)設(shè)計方法的螺旋槳設(shè)計

3.1 設(shè)計狀態(tài)

對某型太陽能無人機(jī)的螺旋槳進(jìn)行設(shè)計，設(shè)計高度H=3 km，來流速度V0=13 m/s，拉力T=10 N，槳盤半徑R=0.27 m。取螺旋槳轉(zhuǎn)速Ω=2 700 r/min，槳葉數(shù)NB=2，槳葉采用Clark Y翼型。

3.2 設(shè)計結(jié)果

設(shè)計得到的螺旋槳弦長、扭轉(zhuǎn)角分布及三維幾何圖形如圖5所示。主要對0.2R以上的槳葉進(jìn)行設(shè)計，而槳轂至0.2R處做過渡處理。

采用上述螺旋槳方法對設(shè)計得到的螺旋槳進(jìn)行計算，設(shè)計結(jié)果與CFD計算結(jié)果如表1所示，可以看出實際的螺旋槳拉力及效率均小于設(shè)計結(jié)果，不能滿足設(shè)計要求。

圖5 螺旋槳設(shè)計結(jié)果

表1 螺旋槳性能對比

4 耦合CFD求解的螺旋槳快速設(shè)計

設(shè)計結(jié)果與真實情況間的差異最終體現(xiàn)在槳葉各葉素截面的氣動力之上，若能夠得到較為準(zhǔn)確的葉素截面氣動力，那么一方面設(shè)計的拉力較為準(zhǔn)確，另一方面設(shè)計更接近最小能量損失，可提高螺旋槳效率。因此，本文提出耦合CFD求解的螺旋槳快速設(shè)計方法，該方法的主要思想便是：采用CFD數(shù)值模擬對設(shè)計進(jìn)行修正，根據(jù)CFD的計算結(jié)果反解槳葉截面氣動力，之后根據(jù)新的氣動力重新設(shè)計新的螺旋槳。

4.1 給定葉素氣動力的螺旋槳設(shè)計

槳葉各截面采用相同的翼型進(jìn)行設(shè)計時，傳統(tǒng)螺旋槳設(shè)計方法得到的沿徑向各葉素截面最大升阻比對應(yīng)的氣動力如圖6所示。除在槳尖和槳根處存在一定的變化外，沿徑向最大升阻比對應(yīng)的迎角與升力系數(shù)變化不大，而槳尖和槳根處本身的環(huán)量偏小，對螺旋槳?dú)鈩恿Φ呢暙I(xiàn)較低。因此，設(shè)計時可視各葉素截面最大升阻比對應(yīng)的氣動力相同，采用固定的迎角α、升力系數(shù)CL及阻力系數(shù)CD進(jìn)行設(shè)計。

用固定的迎角、升力系數(shù)及阻力系數(shù)進(jìn)行設(shè)計時，主要步驟同傳統(tǒng)設(shè)計方法。首先根據(jù)設(shè)計拉力T及設(shè)計狀態(tài)求解誘導(dǎo)螺距V′。而在求解r處截面的弦長b及扭轉(zhuǎn)角θ時，由于升力系數(shù)CL確定，那么根據(jù)式(9)可直接求解弦長，不必進(jìn)行迭代求解。剩余設(shè)計步驟與傳統(tǒng)設(shè)計方法一致。

圖6 迎角及升力系數(shù)徑向分布

4.2 根據(jù)CFD結(jié)果反解槳葉素氣動力

假設(shè)各葉素截面最大升阻比對應(yīng)的氣動力相同，對于給定設(shè)計狀態(tài)及槳葉幾何的螺旋槳，可根據(jù)CFD計算得到的拉力及扭矩結(jié)果，反解最大升阻比對應(yīng)的氣動力。步驟如下：

1) 將槳葉劃分為N個截面，根據(jù)設(shè)計拉力T及設(shè)計狀態(tài)求解誘導(dǎo)螺距V′，具體步驟同傳統(tǒng)設(shè)計方法。

2) 求解最大升阻比對應(yīng)的迎角α。

根據(jù)誘導(dǎo)螺距V′，計算得到第i個截面的實際入流角φi，進(jìn)而根據(jù)式(10)求得迎角αi，故取最大升阻比對應(yīng)的迎角α為

(17)

3) 已知螺旋槳總拉力T、扭矩M，根據(jù)式(12)、式(13)，求得滿足要求的升力系數(shù)CL及阻力系數(shù)CD。

4.3 耦合CFD修正的螺旋槳設(shè)計方法

基于以上方法及分析，提出一種耦合CFD修正的螺旋槳快速設(shè)計方法，主要步驟如下：

1) 根據(jù)給定的設(shè)計狀態(tài)，采用傳統(tǒng)螺旋槳設(shè)計方法，計算得到初始螺旋槳幾何形狀(弦長及扭轉(zhuǎn)角分布)。

2) 由于設(shè)計采用的葉素氣動力與真實情況存在差異，為保證螺旋槳效率最高，需找到真實情況葉素截面最大升阻比對應(yīng)的迎角。因此，在初始幾何外形的基礎(chǔ)上，將槳葉的扭轉(zhuǎn)角整體調(diào)整“±θ0”，生成n個螺旋槳外形。接著，采用CFD分別計算其拉力及扭矩，并反解其對應(yīng)的槳葉葉素氣動力。選取升阻比最大的氣動力數(shù)據(jù)。

3) 根據(jù)選定的數(shù)據(jù)，通過給定截面氣動力的螺旋槳設(shè)計方法，設(shè)計新的螺旋槳。

設(shè)計流程如圖7所示。

圖7 螺旋槳快速設(shè)計流程

4.4 設(shè)計結(jié)果

在圖5所示傳統(tǒng)方法求得的螺旋槳基礎(chǔ)上，通過整體偏轉(zhuǎn)扭轉(zhuǎn)角，得到5個新的構(gòu)型，對其分別進(jìn)行CFD求解，并反解氣動力，結(jié)果如表2所示。

選取升阻比最大的數(shù)據(jù)，迎角α=5.819°，升力系數(shù)CL=0.864，阻力系數(shù)CD=0.041 2，采用給定葉素氣動力的設(shè)計方法重新設(shè)計。設(shè)計結(jié)果如圖8所示。

表2 CFD結(jié)果及反解得到的氣動力

圖8 新設(shè)計方法所得弦長與扭轉(zhuǎn)角分布

采用CFD方法計算其氣動力，并與傳統(tǒng)設(shè)計結(jié)果進(jìn)行對比，如表3所示。在設(shè)計點處，本文提出的設(shè)計方法，一方面能夠較好地滿足拉力要求，另一方面相比于傳統(tǒng)設(shè)計方法，設(shè)計得到的螺旋槳效率提高2.75%。在非設(shè)計點，大前進(jìn)比(λ)下設(shè)計得到的螺旋槳拉力及效率略小于傳統(tǒng)方法設(shè)計結(jié)果，而前進(jìn)比小于0.6時，螺旋槳的拉力及效率均優(yōu)于傳統(tǒng)方法設(shè)計結(jié)果，如圖9所示。驗證了設(shè)計方法的可行性。

表3 傳統(tǒng)和新設(shè)計方法螺旋槳設(shè)計點性能對比

圖9 傳統(tǒng)和新設(shè)計方法螺旋槳性能對比

5 考慮翼型優(yōu)化的螺旋槳設(shè)計方法

為進(jìn)一步提高螺旋槳效率，在前文提出方法的基礎(chǔ)上，對螺旋槳翼型進(jìn)行優(yōu)化。

5.1 優(yōu)化方法

采用“代理優(yōu)化”[17]方法對螺旋槳翼型進(jìn)行優(yōu)化。文中采用的代理優(yōu)化方法基于“Kriging代理模型”[18]，子優(yōu)化方法采用“遺傳算法+乘子法”。遺傳算法采用擬二進(jìn)制交叉、多項式變異、最佳保留策略及聯(lián)賽選擇；乘子法采用Rockfellar推廣的乘子法;并在優(yōu)化過程中采用“最小代理模型預(yù)測準(zhǔn)則[17](Minimum of Prediction, MSP)+改善期望準(zhǔn)則[19](Expected Improvement, EI)”兩種加點法則更新代理模型。

文獻(xiàn)[20]研究表明：Xfoil在低雷諾數(shù)下求解精度雖不及采用高級湍流模型的CFD方法，但其本身具有一定的精度。且對于翼型的優(yōu)化，更關(guān)注其相比于原翼型的改進(jìn)趨勢，因此對翼型采用Xfoil快速求解，這也能夠節(jié)省大量設(shè)計時間。而翼型采用“CST參數(shù)化方法”[21]。

5.2 優(yōu)化目標(biāo)及約束

對于上文采用耦合CFD修正設(shè)計得到的螺旋槳，選取0.2R、0.5R及0.8R這3個槳葉截面，優(yōu)化迎角α=5.819°狀態(tài)下的翼型氣動力，對應(yīng)的雷諾數(shù)分別取為50 000、70 000、90 000。

優(yōu)化目標(biāo)：平均阻力最小，即MinC′D，且

C′D=(CD1+CD2+CD3)/3

(18)

式中：CD1對應(yīng)第1個狀態(tài)的阻力系數(shù)。這里對阻力進(jìn)行優(yōu)化的原因是，低雷諾數(shù)下阻力的變化更加劇烈，優(yōu)化的空間更大。

約束條件：平均升力不減小，C′L≥C′L0，下標(biāo)0代表基礎(chǔ)翼型；翼型面積不減小，A≥A0；翼型最大厚度不減小，t≥t0。

5.3 翼型優(yōu)化結(jié)果

優(yōu)化過程收斂曲線如圖10所示。優(yōu)化前后翼型對比如表4和圖11所示。優(yōu)化后的翼型最大厚度和面積略有增加，后緣彎度增大。從圖11中可以看出，優(yōu)化后翼型前緣上表面的低壓降低，后緣下表面的高壓增大，因此阻力的降低隨之帶來翼型升力的提高，這也是采用“升力不減小”而非“升力不變”約束進(jìn)行翼型優(yōu)化的原因，圖中:Cp為壓力系數(shù)。

圖10 優(yōu)化收斂曲線

表4 翼型優(yōu)化前后氣動力對比

圖11 優(yōu)化前后翼型壓力分布對比(Re=50 000)

5.4 螺旋槳優(yōu)化設(shè)計流程

考慮翼型優(yōu)化的螺旋槳快速設(shè)計流程如圖12所示，優(yōu)化前后翼型外形對比如圖13所示。得到優(yōu)化翼型后，在耦合CFD修正的設(shè)計結(jié)果上替代基礎(chǔ)翼型，生成新的螺旋槳并進(jìn)行計算。因為對翼型優(yōu)化時只限制升力的下限，故會存在升力系數(shù)變化的情況，此時得到的拉力則會偏離設(shè)計拉力。并且新的翼型應(yīng)用在槳上后，其葉素截面最大升阻比對應(yīng)的迎角也無法保證與基礎(chǔ)翼型一致。因此，需再次通過改變槳葉扭轉(zhuǎn)角生成不同螺旋槳，根據(jù)CFD計算結(jié)果進(jìn)行反解氣動力及再設(shè)計，從而得到滿足設(shè)計要求的高效螺旋槳。

圖12 考慮翼型優(yōu)化的螺旋槳快速設(shè)計流程

圖13 優(yōu)化前后翼型外形對比

5.5 螺旋槳性能對比

反解得到葉素截面最大升阻比對應(yīng)的最佳氣動力特性包括：α=4.565°，CL=0.835，CD=0.036 4。再采用該氣動力重新設(shè)計，得到如圖14所示結(jié)果。

表5及圖15展示了3種設(shè)計方法得到的螺旋槳性能對比。考慮翼型優(yōu)化的設(shè)計方法，能夠在滿足拉力要求的同時，進(jìn)一步提高螺旋槳效率。在設(shè)計點處，相比于傳統(tǒng)方法，效率提高3.95%。

圖14 優(yōu)化翼型前后設(shè)計結(jié)果對比

表5 不同設(shè)計方法螺旋槳設(shè)計點性能對比

圖15 不同設(shè)計方法螺旋槳性能對比

6 結(jié) 論

本文根據(jù)CFD計算對螺旋槳的設(shè)計進(jìn)行修正，提出一種能夠較好滿足設(shè)計拉力要求的高效螺旋槳設(shè)計方法，得到的結(jié)論如下：

1) 傳統(tǒng)螺旋槳設(shè)計方法存在一定缺陷。首先，設(shè)計得到的螺旋槳其真實的拉力與目標(biāo)拉力存在一定偏差；其次，不能保證真實情況下槳葉葉素截面的升阻比最大，因此不能保證設(shè)計螺旋槳的高效率。

2) 將設(shè)計與實際情況葉素氣動力的差異，以及真實情況下如展向流動、槳葉槳根及槳尖處的三維效應(yīng)帶來的影響，歸結(jié)于螺旋槳槳葉氣動力的改變。提出了根據(jù)給定氣動力進(jìn)行螺旋槳設(shè)計的方法。假設(shè)槳葉葉素截面最大升阻比對應(yīng)的氣動力沿徑向不變，通過整體改變螺旋槳扭轉(zhuǎn)角，生成多個螺旋槳構(gòu)型，再通過CFD計算及反解氣動力，可得到升阻比最大的氣動力數(shù)據(jù)。根據(jù)該數(shù)據(jù)進(jìn)行螺旋槳的再設(shè)計，可得到高效率的螺旋槳。

3) 經(jīng)過實例的對比分析，文中所建耦合CFD求解的螺旋槳設(shè)計方法，一方面能夠較好地滿足設(shè)計拉力要求，另一方面相比于傳統(tǒng)設(shè)計方法，設(shè)計得到的螺旋槳效率提高了2.75%。而采用代理優(yōu)化的方法對螺旋槳翼型進(jìn)行優(yōu)化后，相比于傳統(tǒng)設(shè)計方法螺旋槳的效率可進(jìn)一步提高3.95%。此外，該方法只需進(jìn)行少量的CFD計算即可，計算量遠(yuǎn)小于直接采用數(shù)值模擬優(yōu)化螺旋槳弦長及扭轉(zhuǎn)角分布的方法。