大迎角下俯仰振蕩對細長體氣動特性的影響

葉 瑞，李 棟，夏 明

(西北工業大學翼型葉柵空氣動力學國防科技重點實驗室，陜西西安 710072)

0 引言

為了提高戰斗機和導彈等飛行器的機動性，大迎角下細長體非對稱渦分離流動的研究得到了很多學者的重視。在現代飛行器設計中，對細長體非對稱渦的控制與利用的一個難點就是其形成的機理尚不清楚。盡管研究固定迎角下細長體的非對稱流動很重要，但依然不能忽視固定迎角下細長體本身非定常運動對流場的影響。

對于非對稱渦的控制，我們必須考慮到細長體的運動和流場隨時間的變化。Ericsson［1］認為利用“動壁效應”，細長體可以采用某種運動消除掉其背風面分離渦的非對稱性。Gadel-Hak［2］等人在實驗中觀察了頭部長細比為3的細長體從0°到30°的俯仰運動，在中等頻率到高頻率作用下，雷諾數的影響可以被忽略掉，而非定常運動對流場影響是占主導作用。Stanek M J［3］等人通過數值方法，發現細長體從0°到20°的俯仰運動過程中，當俯仰運動速度降低時，細長體表面會形成較厚的剪切層。同時細長體表面和渦的相互作用將主渦從表面推開，對表面產生一個很高的吸力。Hoang N T［4－5］等人通過對比6∶1橢球體的非定常運動實驗發現，相對于定常運動，非定常運動可以改變背風面分離的位置并延遲分離發生迎角。但是前人的研究主要集中在細長體做振幅較大的俯仰運動，且很少有給出流場中渦系結構和其非定常變化的結果。此外，關于俯仰運動的實驗并沒有提供詳細非定常流動的數據(譬如速度和壓強)，因此很難用于工程研究。

本文運用DES(Detached Eddy Simulati-on)方法對細長體模型非定常俯仰運動進行數值模擬。結合了RANS(Reynolds-averaged Navier-Stokes)與LES(Large Eddy Simulation)的DES方法在模擬大分離流動時更具有優勢。Morton S［6］等分別運用DES與RANS方法對三角翼的分離渦破裂進行預測，結果表明DES具有明顯優勢;在大分離流動下，李棟［7］對幾種翼型進行了DES與RANS的對比計算，相同的網格DES計算結果與實驗更為吻合，并捕獲到了更細致的渦結構。文獻［8］在頭部未施加擾動塊的情況下，采用DES方法模擬出大迎角下細長旋成體背渦的非對稱性，并且計算結果與實驗符合較好。對于非定常俯仰運動，文獻［9］對頭部未施加擾動塊細長體計算發現，小振幅的振蕩可以抑制住背渦的非對稱。

為了驗證“動壁效應”在細長體俯仰運動時是主導渦系結構的主要因素，在振幅1°、無量綱頻率0.032俯仰振蕩作用下，本文分別通過對頭部無擾動和有擾動的兩個模型進行了對比計算，分析各個流場中渦系結構和氣動力非定常變化過程，并對比兩者之間的差別，對大迎角下細長體自身俯仰運動對其流場的影響進行探討。

1 數值方法

考慮到分離流動的非定常和三維特性，數值計算中使用3D網格和非定常的計算方法。計算中對NS方程和S-A湍流模型方程同時使用了偽時間步長技術，使用LU-SGS方法求解離散的S-A方程。

1.1 計算模型與網格

計算模型為尖拱-圓柱細長體，底部直徑為D，細長體總長度7.5D，頭部長細比等于3.5，后部圓柱段長細比等于4。頭部有擾動的模型是在頭部的尖端布置一微小擾動塊來代替頭部存在的微小幾何不對稱。擾動塊的長度 x/D=0.014，高度 h/D=0.0004，最大寬度0.004D，周向位置90°。圖1和圖2分別顯示擾動塊的表面示意圖和其截面附近的網格分布。

圖1 細長體與其頭部的擾動Fig.1 Slender body with geometrical disurbance

圖2 擾動附近的網格Fig.2 Grid around the disturbance

計算域徑向前外延5D，向后外延12D，法向外延10D，物面第一層法向網格間距為1×10－5D，網格數:160×100×120(流向×法向×周向)。細長體表面采用無滑移絕熱壁面，遠場邊界采用一維Riemann不變量。

1.2 DES 方法

本文采用基于S-A湍流模型的DES方法，S-A方程的具體求解方法參見文獻［10］。

對于 DES 方法的求解，按照 Spalart［11］的說法，對典型的脫體分離為主導的流動，使用DES方法時，空間濾波尺度和非定常時間尺度在計算中是否適當，是可以根據一定的標準來判斷的，如果兩個尺度過大，或者格式精度過低導致出現低通濾波效應，LES的計算將退化相當于求解雷諾平均的NS方程(RANS)。因此對比DES和RANS的計算結果是一種驗證是否實現大渦模擬的方法。在細長體迎角固定的計算條件下Ma=0.2，α=40°，Re=3 ×106，文獻［8］給出了本文采用的 DES 與RANS計算的對比結果，見圖3。

圖3 DES與RANS計算細長體渦量云圖Fig.3 Vorticity contour using DES and RANS

從圖3中對比可以看出，在RANS方法計算的細長體流場中，細長體背風面兩側形成了一對對稱的分離渦。而在由DES方法計算的流場中，背風面形成了非對稱多渦系結構，且渦結構的非對稱性沿軸向而增強，這與前人的實驗結果［12］是一致的。表明本文DES方法是可行的。

1.3 數值方法驗證

為進一步驗證本文數值方法的可靠性。針對本文采用頭部無擾動的尖拱—圓柱細長體模型，計算了Ma=0.2，α =40°，Re=3 ×106定常狀態下的流場，與文獻［11］實驗結果進行了對比驗證。

圖4為DES計算的細長體不同截面處的周向壓力系數分布與實驗結果的對比，其中φ為周向角。在定常狀態下，可以從圖中看出采用DES方法模擬出了流場的非對稱。并且在細長體不同的截面上計算壓力系數分布和實驗吻合得比較好。這也進一步說明本文所采用的DES方法是可行的。

圖4 細長體不同截面的周向壓力系數分布Fig.4 Circumferential surface pressure coefficient on different cross section of the body

2 俯仰振蕩對細長體流場影響的數值分析

計算條件為 Ma=0.2，α =40°，Re=3 ×106。非定常計算中，取物理時間步長為0.1，每個物理時間步中子迭代步數為20步。

圖5為細長體俯仰振蕩的示意圖。圖6為細長體的迎角隨無量綱時間的變化過程。我們可以看出，俯仰非定常運動并不是在t=0時刻開始的。而是在t=40.0，α=40°時對細長體施加一個小振幅高頻率俯仰振動。此后迎角隨時間成正弦規律變化，實時迎角α=40°+Asin(ωt)，A=1°為振幅，ω =0.2 為圓頻率。

圖5 細長體俯仰振動示意圖Fig.5 Pitching motion of slender body

圖6 迎角隨無量綱時間的變化Fig.6 Angle of attack vs non-dimensional time

很多數值和實驗研究顯示［13－14］，在固定迎角狀態下，細長體背風面渦系結構的發展是其流場變化的主要特征。而渦形態的改變是影響大迎角流動的直接因素，因此研究背風面渦的強度和渦量輸入的變化，也是分析非定常運動下細長體渦結構的有效方法。

在細長體上取截面x/D=5來觀察細長體渦系結構的變化。從圖6中可以看出，整個計算過程包含了5個振蕩周期。而由于在第一個振蕩周期內，渦系從非對稱向對稱的轉變過程最為顯著。因此把第一個周期作為研究對象。為了便于研究，將俯仰振蕩的第一個周期分為四個時間段:上升①、下降①、下降②、上升②。圖7、圖8、圖9、圖10為兩個細長體模型截面x/D=5在四個時間段渦量圖隨時間的變化，左邊為頭部無擾動模型的截面，右邊為頭部有擾動模型的截面。

t=40，α=40°細長體振蕩開始，無擾動流場顯示為四渦系結構:細長體左側主渦VL1、右側主渦VR1、左側二次渦VL2、右側二次渦VR2。而由于頭部存在擾動，有擾動流場與無擾動流場的左右主渦強度正好相反，并且有擾動流場非對稱更明顯，已存在飄起主渦VR1’，新生主渦VR1開始形成。迎角上升至α=40.565°，雖然渦量圖沒有明顯的變化，但是從圖11中可以看到，無擾動和有擾動的截面法向力CN有顯著的增加。當α=40.985°時，截面法向力CN回復到原來的大小。

圖7 第一個周期上升①x/D=5截面渦量圖Fig.7 Up① of the first cycle x/D=5 vorticity contour

當迎角開始下降，從圖8的渦量圖可以看出兩側的主渦VL1、VR1強度都在變大，這說明兩側主渦都得到了渦量注入。圖11顯示截面側向力系數Cz都在減小，相對于無擾動模型，有擾動模型減小的更快。這也顯示相對于兩側較強的主渦，弱的一側主渦得到了更多的渦量注入，使得兩側渦強度差變小，這也阻止了渦的繼續飄起，使得其并沒有像固定迎角一樣繼續脫落。同時兩側二次渦強度也慢慢變強。

當迎角繼續下降，此時細長體模型的角加速度很大，流場變化很劇烈，兩側的二次渦開始遠離物面，并且慢慢開始合并成一個渦V2，新合并的渦位置處于背風面中軸位置。兩側飄起主渦VL1’、VL2’由于與各自渦VL1、VL2的聯系中斷，得不到渦量的注入，其渦強是一直降低。對于無擾動模型來說，圖11顯示截面側向力系數Cz都在減小，兩側的渦強趨于平衡。但對于有擾動模型，側向力系數Cz的大小是先減小，而后又有增大。這也反映兩側渦強度差的變化趨勢。

迎角又開始上升時，兩圖中飄起的渦VR1’強度都在減小，合并渦V2仍然處于背風面中軸位置。同時兩側主渦VL1、VR1隨著兩側剪切層渦量的注入，強度慢慢增大。不同于無擾動形成的一對平衡渦，有擾動左側的主渦的強度比右側要大，流場趨于向不對稱轉變。

圖8 第一個周期下降①x/D=5截面渦量圖Fig.8 Down① of the first cycle x/D=5 vorticity contour

在整個的第一個振蕩周期內，我們可以從細長體兩側渦的變化過程可以看出，非定常俯仰運動促使兩側原先不對稱的渦系向對稱轉變。但不同的是，頭部無擾動的細長體流場中的非對稱渦系向對稱轉變的趨勢在這一過程中始終是不變的。而對于頭部存在擾動的細長體，頭部存在這種人為的擾動增強了整個流場的非對稱性，而這種人為的擾動引起的流場中非對稱降低了俯仰非定常運動促使流場向對稱轉變的這一作用效果，使得第一個振蕩周期結束后，細長體兩側主渦分布仍然與開始振蕩時相似。第一個振蕩周期開始(t=40，α=40°)的有擾動的渦量圖顯示左側主渦強度占優，第一個振蕩周期結束(t=71，α=39.917°)有擾動的渦量圖仍然顯示左側主渦強度占優。這一現象支持了上面的說法。

流場渦的作用效果可以通過截面力系數的分布來體現。

圖9 第一個周期下降②x/D=5截面渦量圖Fig.9 Down② of the first cycle x/D=5 vorticity contour

法向力系數CN定義為

側向力系數Cz定義為

其中，FN為模型受到的總法向力，Fz為其總側向力。動壓為周向角，D 為底部直徑。

從圖11(a)中可以看出，盡管頭部有擾動模型的流場與無擾動的有很大的不同，但是在俯仰振蕩過程中，兩者截面法向力系數CN的大小和變化沒有很大的差別。這說明頭部是否存在微小擾動對于截面法向力系數CN影響不大。但對于截面法向力系數Cz，無擾動模型的Cz基本或總體趨勢是減小的，而有擾動模型的Cz變化則有一個明顯的波動過程。這也顯示了在第一個周期內，振蕩對無擾動細長體流場的非對稱抑制作用更明顯。

圖10 第一個周期上升②x/D=5截面渦量圖Fig.10 Up② of the first cycle x/D=5 vorticity contour

圖11 第一個周期x/D=5截面力系數變化Fig.11 x/D=5 sectional side force coefficient in the first cycle

圖12為俯仰振蕩中無擾動和有擾動細長體法向力系數CN隨時間的變化，可以看出它們的法向力系數CN隨著俯仰振蕩的進行都成周期性變化，并且CN大小與變化的趨勢相差不大，這進一步說明頭部是否存在微小擾動對于法向力是影響不大。但同時應該注意是在振蕩開始時(無量綱時間t=50)，兩者法向力系數CN與其側向力系數Cz在大小上相差很小。這說明此時側向力對細長體的影響是相當大的。另外經過4個周期振蕩后，圖13顯示在很短時間內，法向力系數CN出現的最大值(圖12中A點)與最小值(B點)相差幾十倍，而圖13顯示此時的側向力基本是保持在一個接近零的小量。這說明細長體在持續的振蕩下，法向力的變化是其主要的受力特征，并且其變化快、幅度大。

圖13則顯示經過一個振蕩周期后，無擾動和有擾動細長體側向力系數Cz的大小都迅速降低，隨著振蕩的繼續進行，側向力幾乎保持在接近零的一個量級。表明對大迎角下的細長旋成體施加俯仰振蕩，會有效抑制其流場的非對稱。

圖12 法向力系數CN隨無量綱時間變化Fig.12 Normal force coefficient CNvs non-dimensional time

圖13 側向力系數Cz隨無量綱時間變化Fig.13 Side force coefficient Czvs non-dimensional time

3 結論

本文數值模擬了大迎角俯仰振蕩對頭部無擾動和頭部有擾動的細長體氣動特性的影響，計算研究表明:

(1)對頭部無擾動和頭部有擾動的細長體施加俯仰振蕩，在開始的第一個周期，細長體背風面渦系都由原來非對稱形態向對稱形態轉變。相對于有擾動的細長體模型，在第一個周期內非定常俯仰運動對頭部無擾動的細長體作用更明顯。計算結果表明，非定常運動產生的“動壁效應”此時主導著渦系的發展。

(2)在俯仰振蕩過程中，法向力系數隨著俯仰振蕩的進行都成周期性變化，頭部是否存在微小擾動對法向力影響作用很小。經過一段時間的持續振蕩后，法向力的大幅度的迅速變化是細長體主要的受力特征。

(3)在持續的俯仰振蕩下，頭部無擾動和頭部有擾動的細長體側向力始終保持在一個小的量級，表明流場非對稱被有效抑制。

［1］ERICSSON L E.The fickle effect of nose microasymmetry on the high-alpha aerodynamics［R］.AIAA-90-0067，1990.

［2］GAD-EL-HAK M，HO C.Unsteady flow around an ogive cylinder［J］.Journal of Aircraft，1986，23(6):520-528.

［3］STANEK M J，VISBAL M R.Investigation of vortex development in a pitching slender body of revolution［J］.Journal of Aircraft，1993，30(5):711-718.

［4］HOANG N T，WETZEL T G，SIMPSON R L.Unsteady measurement over a 6∶1 prolate spheroid undergoing a pitch-up maneuver［R］.AIAA-94-0197，1994.

［5］HOANG N T，WETZEL T G，SIMPSON R L.Surface pressure measurements over a 6∶1 prolate spheroid under-going time-development maneuver［R］.AIAA-94-1908，1994.

［6］MORTON S.DES and RANS simulations of delta wing vortical flows［R］.AIAA-2002-0587，2002.

［7］LI D.Numerical simulation of thin airfoil stall by using a modified DES approach［J］.Numerical Methods in FLUIDS，2007，54(3):325-332.

［8］XIA M，SONG B F，LI D，YE R.Comparison of predictive capabilities of DES and RANS simulations of the separated flow around a circular cylinder［R］.AIAA-2010-328.

［9］XIA M，LI D，YE R，SONG B F.Numerical study of pitching oscillation effect on aerodynamic characteristics over a slender body using DES［A］.The 8th Asian Computational Fluid Dynamics Conference［C］.Hong Kong，2010.

［10］SPALART P R，ALLMARAS S R.A one-equation turbulence model for aerodynamic flows［R］.AIAA-92-0439，1992.

［11］LAMONT P J.The complex asymmetric flow over a 3.5D ogive nose and cylindrical afterbody at high angles of attack［R］.AIAA-82-0053，1982.

［12］SPALART P R，JOU W H，STRELETS M，ALLMARAS S R.Comments on the feasibility of LES for wings，and on a hybrid RANS/LES approach［A］.First AFOSR International Conference on DNS/LES［C］.Ruston，LA，1997.

［13］CHEN X R，DENG X Y，WANG Y K，LIU P Q，GU Z F.Influence of nose perturbations on behaviors of asymmetric vortices over slender body［J］.ACTA Aeronouticaet Astronautica Sinica，2002，18(6):581-593.

［14］周乃春，葉正寅，胡漢東，馬明生，楊其德.細長尖頭旋成體大迎角非對稱渦系結構［J］.空氣動力學學報，2004，22(4):486-493.