劉 宇,易 賢,王 強,李維浩
(中國空氣動力研究與發展中心,綿陽 621000)
飛機結冰是飛行安全的主要威脅之一[1],為了評估結冰對飛行性能的影響,需要對飛機進行結冰風洞試驗。若結冰風洞受試驗段尺寸限制無法進行全尺寸模型試驗,或參考工況超過結冰風洞的模擬能力范圍,需要對試驗工況進行相似轉換。相似參數是描述結冰過程的重要無量綱數,相似轉換是指利用相似參數計算合適的試驗工況,使試驗可在風洞能力范圍內進行,并獲得與參考工況相同的冰形。但是在實際試驗過程中發現,現有相似轉換方法的結果很多時候超出風洞的實際能力范圍:在參考工況速度較高或液態水含量較小時,往往轉換結果的風速更高或液態水含量更小,導致轉換后的參數在風洞中難以實現。因此需要開展更具實用性的相似轉換方法研究。
國外的相似轉換研究始于20世紀50年代,對相似參數的研究已經相對完備,通過不同相似參數的組合形成了多種相似轉換方法。目前常用的轉換方法為Ruff方法[2](或稱AEDC方法),以及Ruff與Anderson提出的改進Ruff方法[3]。法國宇航局的Charpin等提出了ONERA方法[4]。還有LWC×Time方法[5]和Olsen方法[6],可以有效對液態水含量(LWC)進行轉換。單獨對溫度、速度和壓力進行的轉換方法[7-9]由于沒有試驗驗證很少被采用。國內對相似轉換的研究多在已有方法上拓展現應用場景,如針對過冷大水滴[10-11]轉部件[12-13]的相似轉換方法,對轉換方法的優化和改進較少。
國外現有的相似轉換方法對霜冰條件的轉換已經成熟,明冰條件的結冰機理還缺乏理解[14],轉換效果有限。對于轉換速度如何確定還存在一些討論,主要集中在韋伯數We特征長度的定義上,可能作為特征長度的有:翼型前緣半徑[15]、水滴直徑[16]、水膜厚度[17]等,對于如何合理確定轉換速度仍未得到合理解釋。其中Kind[17]和Feo[18]等從表面水膜流動和水滴表面張力等方面對Ruff方法進行了修正。在實際應用上國外以Ruff方法為主,近年關于ONERA方法的改進和驗證的文獻相對較少。
國外對結冰相似準則的研究隨應用需求轉入過冷大水滴[19]、后掠翼[20-21]、發動機冰晶結冰[22]和防冰試驗相似轉換[23]等研究領域,對傳統的相似轉換理論研究逐漸降溫。
為拓展相似轉換方法的適用范圍,提高其工程實用價值,發展適合我國大型結冰風洞的相似轉換方法,本文對相似參數及其變化規律進行了深入分析,根據相似參數和來流條件的關聯性,提出了一種混合相似轉換方法。該方法根據滿足的相似參數不同分為四種模式,每種模式的轉換參數選擇范圍不同。使用混合轉換方法對某機翼結冰試驗的工況進行了相似轉換,并通過數值模擬[24]計算了試驗工況、參考工況以及被普遍使用的Ruff方法轉換工況的冰形,并對試驗和計算結果進行了對比。
為保證水收集系數相似,需要滿足水滴運動軌跡相似。在流場相似的前提下滿足水滴慣性系數K0相等即可認為滿足水滴軌跡相似。其定義如下[25]:

其中,λ/λstokes是平均阻力比,K是未做修正的慣性系數[26]。
模型表面水收集總量受水收集系數、LWC和結冰時間影響,用積聚系數Ac表示:

其中,τ為結冰時間,ρi為冰密度,D為前緣直徑,V為來流速度。
根據Messinger模型[27]的能量守恒方程,可以推導出模型表面凍結系數:

其中,Λf為結冰相變潛熱,cp,ws為水的定壓比熱容。另外有描述熱平衡過程的三個參數:相對熱系數b、水滴能量傳遞勢φ和空氣能量傳遞勢θ。定義如下:
其中,ts為壁面溫度,tf為水膜溫度,pw為水蒸汽分壓力,pww為飽和水蒸氣壓力,Λf為水的蒸發潛熱,hc為對流換熱系數,hG為對流傳質系數。
翼型駐點位置凍結系數n0也可通過公式(3)計算。通常認為在流場和水收集系數相似的前提下,滿足n0相等則其余部分的冰生長過程也相似。
常見的相似轉換方法有ONERA方法、Ruff方法和Olsen方法。以上方法均有一定試驗數據進行支撐,但每種方法均有一定局限性。其中ONERA方法根據Modane風洞的硬件條件設計,沒有將壓力作為可調整變量[28]納入轉換方法中。Ruff方法建議使用WeL,s=WeL,r條件[19-20]確定Vs, 使轉換后速度增加,不利于開展試驗。WeL的定義為:

LWC×Time方法僅適用于溫度較低的霜冰工況,Olsen方法的部分試驗結果顯示其轉換結果欠佳[29]。
通過上節中相似參數定義分析,可以定性的得出以下推論:1)K0主要受 δ和V影響,溫度和壓力變化對K0影響相對較小;2)當K0,s=K0,r條件滿足時,n0僅受到LWC和Tst影響,而三個熱平衡參數相等不是n0,s=n0,r的必要條件;3)相對熱系數b是LWC和V的函數;4)水滴能量傳遞勢φ是Tst和V的函數。
由以上推論可將相似參數分為三組:一是K0,通過相等條件K0,s=K0,r限制δs和Vs;二是n0,通過n0,s=n0,r限制LWCs和Tst,s;三是b、φ、θ,對來流條件的限制較弱,即使這三個參數均不相等或部分相等,相似轉換依然有效。以下將針對這三組參數進行分析。
慣性系數K0的經驗計算公式如下,該公式在同類型研究中得到廣泛應用:

其中,k和kl均為常數,取值方法可參見文獻[2]。實際試驗以水滴中值體積直徑(MVD)代替公式中的單粒徑水滴直徑δ。忽略溫度和壓力對空氣粘性和密度的影響,根據公式(8)繪制在滿足K0,s=K0,r條件下,K0,s關于V和MVD的變化曲面如圖1所示,參考工況已在圖中標注。從曲線變化趨勢可以看出,MVD隨V增大而減小,在K0,s取不同值時曲線變化趨勢一致。

圖1 K0,s 關于速度和MVD的曲面關系Fig. 1 Curvature plane of correlation between velocity and MVD in constantK0,s
根據3.1節的推論(2)可知,n0,s=n0,r是 較弱的約束條件,對任意狀態,理論上有無數LWC和Tst組 合滿足該條件。將公式(4 ~ 6)帶入到公式(3)并化簡,可以得到:

由于計算n0的公式展開后有若干物性參數,使用A1、A2、A3簡化展開公式產生的常數項,其中部分如空氣密度和靜壓與溫度直接相關的物性參數在簡化過程被作為常數處理。
根據公式(9)繪制在滿足n0,s=n0,r條件下,n0,s關于LWC和Tst變化的曲面如圖2所示,參考工況已在圖中標注,其中n0,s坐標軸取對數。從曲線的趨勢可以看出,Tst隨LWC增加而減少,n0越大可供選擇的LWC范圍越小。

圖2 n0,s 關于LWC和靜溫的曲面關系Fig. 2 Curvature plane of correlation between static temperature and LWC in constantn0,s
為了保證相似轉換前后的冰形盡可能相似,Ruff方法和ONERA方法均在n0相等的基礎上約束了φ或b相等,以期獲得與參考工況相似的熱平衡過程。為了進一步分析n0關于來流條件的變化規律,需要分別對b、φ、θ進行分析。
相對熱系數b是在n0,s=n0,r條件的基礎上對LWC的進一步約束,展開對流換熱系數,公式(5)可以寫成:

其中,ka為空氣導熱系數,Pra為普朗特數,通常取0.7。滿足了水滴慣性系數相等條件時式中β0可作為常數,即式中除LWC和V以外均為常數。
依據公式(10)繪制在滿足bs=br條件下,bs關于LWC和V的空間曲面如圖3所示,參考工況條件已在圖中標注。從曲線的規律可以看出,LWC隨V的增大而單調減小。

圖3 bs 關于LWC和速度的曲面關系Fig. 3 Curvature plane of correlation between velocity and LWC in constantbs
同理根據公式(5)做出φs關于Tst和V的空間曲面如圖4所示,參考工況已在圖中標注。從方程和投影曲線均可以看出,Tst和V是單調的線性關系,斜率與φs的值無關。

圖4 φs 關于靜溫和速度的曲面關系Fig. 4 Curvature plane of correlation of static temperature and velocity in constantφs
由以上分析可以知道,在相似轉換方法僅滿足ns=nr的弱約束條件時,LWC和Tst存在非線性的函數關系;在引入bs=br或φs=φr條件后,LWC和Tst實際上是兩個完全不相關的參數。因此在速度確定的情況下,可計算出同時滿足bs=br和φs=φr的LWCs和Tst,s。此時θs=θr成了滿足凍結系數相似的充要條件。
在沒有其他約束條件作為前提時,θ受到多個來流條件的影響,無法通過 θ約束單一來流條件。但在其他參數已經確定的情況,可求解Pst使 θ滿足相似條件,分析 θ關于Pst的變化規律具有一定的實用價值。
根據公式(6),當其余來流條件確定時,θ關于Pst的變化規律如圖5所示。需要額外注意的是,根據Ruff對 θ的分析,當Vs大于Vr時,滿足θs=θr的靜壓Pst,s小于Pst,r,反之亦然。這個規律意味著當參考工況的靜壓指定為當地壓力時,更小Vs會導致相似轉換后的試驗壓力高于當地壓力。

圖5 不同條件下θ 關于靜壓的變化曲線Fig. 5 Curve of correlation between static pressure and θ in different condition
經過對水滴慣性系數、凍結系數以及三個熱平衡參數的分析,已經確定了來流的主要參數,僅剩結冰時間。根據Acs=Acr計算相應的結冰時間即可滿足相似。
相似參數可分為三組:約束 MVDs和Vs的水滴慣性系數K0;約束LWCs和Tst,s的駐點凍結系數n0;最后一組為熱平衡參數b、φ和θ,對試驗參數沒有直接約束,即使轉換前后沒有滿足相等關系也可獲得相似的冰形。分析了三組相似參數關于來流參數的變化規律,并繪制了當滿足相似參數相等時來流參數之間的關系曲線。
一種混合相似轉換方法的流程如圖6所示,具體實施流程為:
1)指定速度Vs,根據K0條件計算MVDs,Vs可通過WeL條件計算獲得,若計算所得MVDs超出風洞能力包線,可參考圖1中MVD-V的曲線,重新調整Vs;
2)根據φ和b的相似條件,分別計算Tst,s和LWCS,當計算結果超出風洞能力包線,可參考圖3和圖4中兩個參數分別關于V的曲線關系,重新調整Vs并重復步驟1;

圖6 混合轉換方法流程圖Fig. 6 Flow chart for hybrid scaling method
3)a. 如果風洞具備連續調節壓力的能力,可根據公式(7)計算Pst,s,使之滿足 θ條件。由于Pst,s改變導致K0,s出現變化,進一步導致bs不相等,可再次調整MVDs和LWCs使之重新滿足。此時為混和相似轉換方法模式1;
b. 如果風洞無法連續調節壓力,則在b和φ條件中選擇一個,并通過n0條件計算LWCs或Tst,s。此時為混合相似轉換方法模式2或模式3;
c. 當模式2或模式3所計算的參數超出風洞能力范圍時改用Olsen方法,指定LWCs或Tst,s并通過n0條件計算另一個。此時為模式4;
d. 根據Ac條件計算結冰時間τs,并再次檢查各參數是否超出風洞能力包線,轉換結束。
根據滿足的相似參數不同,混合相似轉換方法又區分為四種模式。混合相似轉換方法與常見方法的對比見表1。

表1 相似轉換方法滿足相似參數表Table 1 Scaling parameters fitted by different scaling method
模式1同時滿足了n0相關的三個熱平衡參數,約束最多,理論上可以獲得最接近參考冰形的轉換結果;模式2可視為改良的ONERA方法,模式3與修正的Ruff方法類似;模式4直接使用了Olsen方法,在模式4下參數選擇范圍最大,任意參考工況均能找到試驗能力范圍內的轉換結果。
根據需要轉換的參考工況和試驗設備的能力在四種模式中進行選擇:試驗設備具有高度模擬的能力可以考慮使用模式1;使用模式2或模式3的參數選擇范圍相對自由,同時也有大量的文獻和試驗結果支撐其轉換效果,可信度最高;模式4最具泛用性,考慮到有公開文獻驗證Olsen方法的轉換結果并不總是令人滿意,特定情況下模式4下轉換效果可能較差。
分別使用混合相似轉換方法的四種模式,對文獻[30]中的試驗工況進行了相似轉換,使用數值模擬方法計算了轉換工況的冰形,并與文獻中的試驗結果進行了對比。參考狀態使用弦長0.533 4 m的NACA0012翼型,轉換狀態的試驗和計算均使用0.266 7 m的縮比模型。算例如表2所示。
圖7展示了兩組不同速度的計算結果與試驗的對比結果,試驗與計算結果基本吻合,主要差異在于:轉換工況的駐點冰厚較試驗更低;下冰角角度和厚度基本一致,但上冰角吻合程度稍差;結冰極限相對試驗更靠近前緣??紤]到數值模擬采用單粒徑水滴和單步法計算結冰,與試驗結果的誤差在可接受范圍內。
從驗證結果來看,混合相似轉換方法的四種模式均能有效獲得與參考工況相近的轉換結果。速度不同時轉換結果的冰形存在一定差異,但區別并不明顯。
試驗在中國空氣動力研究與發展中心(CARDC)的3 m×2 m結冰風洞進行。該風洞試驗能力覆蓋液態水含量0.2~2.0 g/m3,水滴中值體積直徑15~50 μm的云霧參數范圍,主試驗段最大風速210 m/s。
使用混合轉換方法,對某型飛機的結冰試驗工況進行了相似轉換,具體工況見表3至表5。其中工況1是因模型存在縮比而進行相似轉換;工況2和工況3是因原定試驗參數超出結冰風洞試驗能力而進行相似轉換。試驗工況均采用混合轉換方法所得,在表中以#號標出。通過熱刀法截取的翼型中截面冰形。對試驗工況、原始工況和其他相似轉換工況進行了數值模擬,并對比了試驗結果和計算結果。

表2 混合相似轉換方法驗證算例Table 2 Validation cases by hybrid scaling method
工況1中對比了Ruff方法和混合相似準換方法,如圖8所示。從計算所得冰形結果分析,Ruff方法和混合轉換方法均能獲得與參考冰形類似的結果,且均與試驗結果對比較好。但Ruff方法的轉換結果風速高,液態水含量小,超過了試驗設備的能力范圍。

圖7 混合相似轉換方法計算驗證結果圖Fig. 7 Ice shape results comparison with validation cases by hybrid scaling method

表3 工況1參數表Table 3 Reference and scaled conditions for test 1

表4 工況2參數表Table 4 The reference and scaled conditions for test 2

表5 工況3參數表Table 5 Reference and scaled conditions for test 3
工況2和工況3來流條件相同,但所用翼型不同,對比結果分別如圖9、圖10所示。轉換前后的冰形基本相同,與試驗結果對比,冰形特征基本吻合。其中工況3計算冰形的冰角角度和下冰角厚度與試驗結果有表明顯的差異,其來源可能是沒有對冰增長過程進行多步法求解。

圖8 工況1試驗與計算結果冰形對比Fig. 8 Ice shape results comparison with validation cases for test 1

圖9 工況2試驗與計算結果冰形對比Fig. 9 Ice shape results comparison with validation cases for test 2

圖10 工況3試驗與計算結果冰形對比Fig. 10 Ice shape results comparison with validation cases for test 3
本文分析了相似參數隨試驗參數的變化規律,提出了一種混合相似轉換方法,通過試驗和數值模擬進行了對比驗證,得到以下結論:
1)使用數值模擬方法驗證了混合相似轉換方法的四種模式,對比文獻試驗結果,四種模式均與試驗吻合良好,驗證了本文方法的有效性。又以混合相似轉換方法對某機翼的結冰試驗進行了相似轉換,對比了試驗結果和計算結果,初步驗證了本方法的實用性;
2)混合相似轉換方法相比傳統的相似轉換方法,提供了更寬的適用范圍,并根據參考狀態的需求確定不同的初始轉換條件,選擇對應的轉換模式,能夠更好的解決實際工程問題;
3)混合相似轉換方法的四種模式中,模式1創新的提出了一種使三個熱平衡參數同時滿足相等的轉換方法,相比經典的Ruff方法和ONERA方法只能滿足三個熱平衡參數中的一個,其轉換結果在理論上更接近參考狀態的結冰熱力學過程;
4)后續工作一方面會繼續在理論上拓展相似轉換方法,通過考慮水膜流動對結冰過程的影響,對來流風速的選擇缺乏有力的理論約束;另一方面將積極開展結冰風洞標模試驗,進一步驗證混合相似轉換方法的有效性。