基于附面層的摩擦阻力數值算法研究

彭 聰

(成都大學生物產業學院，四川 成都 610106)

基于附面層的摩擦阻力數值算法研究

彭 聰

(成都大學生物產業學院，四川 成都 610106)

利用附面層的速度型求解壁面的空氣剪切應力，進而求得摩擦阻力，由此計算了一維層流平板邊界層和二維層流NACA0012翼型的摩擦阻力.平板邊界層計算結果同布拉修斯理論解相比較，吻合性良好.翼型計算結果同實驗數據比較發現，小攻角氣流不產生分離的情況下，摩擦阻力值與實驗數據接近，隨著迎角增大，分離區的擴展，壓差阻力的比重增加，計算誤差明顯增加.

摩擦阻力;附面層;數值計算;空氣粘性

0 引言

從20世紀70年代末開始，在先進理論方法和成功工程實踐的基礎上，由于計算流體動力學(Computational Fluid Dynamics，CFD)計算方法的不斷突破和計算機技術的快速發展，CFD進入了蓬勃發展的新時期，開始越來越多地應用到航空航天飛行器的研究和設計中［1－4］.目前，阻力計算一直是CFD中的一大難題和研究熱點，根據其產生機理，阻力主要分為壓差阻力和摩擦阻力2種.壓差阻力是飛行器各部分表面受到的壓力在速度方向上的合力，又分成波阻、渦阻、型阻等;而摩擦阻力是飛行器各部分表面受到的剪切力在速度方向的合力，它是由空氣的粘性而引起的.目前的CFD方法對于流場中物體表面壓力的分布能得到很好的計算結果，在壓差阻力的計算方面精度較高，但對于摩擦阻力的計算卻存在一定難度，這是因為摩擦阻力的計算對CFD有較高的要求，需要計算出飛行器合理的流場結構，尤其需要正確捕捉壁面邊界層內的速度型以準確計算壁面處的剪切力［5－6］.平板繞流問題是附面層流動最簡單也是最重要的情況，而對于曲面繞流，只要不發生顯著分離，摩擦阻力與平板情況相差不多.本研究采用層流平板邊界層和NACA0012翼型邊界層來驗證摩擦阻力計算的有效性及精度.

1 摩擦阻力計算方法

1.1 控制方程

本研究利用商用Fluent軟件求解N-S方程，不考慮體積力和外部熱源，其形式［6］為，

式中，Q為守恒變矢量;f，g，h分別為x，y，z坐標方向的通量.

1.2 壁面摩擦阻力

摩擦阻力由切向應力引起，并跟附面層的速度型相關.切向應力為，

式中，μ為流體粘性系數，摩擦應力τ與附面層間的切向速度u的變化有關.

摩擦阻力為壁面摩擦應力xF的積分，其形式為，

平板摩擦阻力系數CF根據定義為，

1.3 壁面應力的精度

壁面應力τw的計算采用二階精度的中心差分格式，其形式為:

2 算例及分析

2.1 層流平板附面層及摩擦阻力計算

層流平板附面層的計算模型如圖1所示.網格數為106×50(流向 ×法向)，遠場為4倍平板長.為準確地模擬粘性效應，附面層內布置了至少4個網格點.自由來流條件為:來流速度Ma∞=0.3，雷諾數 ReL=1.0 ×105，攻角 α =0°，研究網格質量對流場及摩擦阻力的影響.表1為平板摩擦阻力計算結果，并與布拉修斯解進行了比較，圖2為平板各占位的速度型.

圖1 層流平板附面層計算網格

表1 摩擦阻力計算結果

圖2 層流平板各站位的速度型

從圖2可以看出，平板附面層內各站位的速度型與布拉休斯解曲線擬合得很好，由此計算得到的阻力系數值很接近.但由于阻力系數值本身為一個小量，網格的粗細和計算格式的粘性分辨率對計算精度都產生影響.從計算效率來說，對于本算例，附面層內布置3～4個網格點就能很好地捕捉速度型.此對于大型網格計算具有指導意義，但對于復雜外形，網格數量更大的情形，還需要進一步驗證.

2.2 NACA0012翼型附面層及摩擦阻力計算

NACA0012翼型附面層的計算模型如圖3所示.網格數為497×100(C型網格)，遠場為10倍弦長.邊界層內同樣布置了4個網格點.自由來流條件為:來流速度 V∞=21.1 m/s，雷諾數 ReL=1.44×106，攻角分別為 0°、2°、4°、8°時，著重研究翼型繞流情況對摩擦阻力計算的影響.

圖3 NACA0012翼型計算網格

圖4、圖5、表2為模型在不同攻角下的計算與實驗結果.

圖4 不同攻角下的流線圖

從圖4中的流線可以看出，0°和2°攻角的流場情況較為相似，氣流均貼著翼型呈明顯的層流狀態;4°攻角時，翼型對來流產生擾動，上翼面1/2弦長位置出現局部氣泡，但氣流未出現明顯分離，至后緣氣流再附于翼型表面，基本保持層流狀態;8°攻角時，擾動進一步增強，上翼面氣泡出現位置向前緣移動，且出現明顯的分離趨勢，至后緣與下翼面繞流形成分離渦結構.

圖5 NACA0012翼型計算結果同實驗數據的比較

表2 NACA0012翼型在不同攻角下阻力系數計算值與實驗值的比較

此外，從圖4中可看出，在對稱流動時，駐點位于翼的前緣，翼型繞流為附體狀態，阻力的主要貢獻為空氣粘性影響的翼型表面很薄的附面層，表現為摩擦阻力，表面壓強分布對阻力影響很小.層流狀態下，摩擦阻力也很小，即使在很小的攻角情況下，只要繞流不產生分離，摩擦阻力的增加不明顯.隨著迎角的增大，駐點向下翼面移動，并逐漸后移，上翼面最大速度點向前緣移動.對于圓頭較厚翼型，在中等迎角下，上翼面靠近后緣的附面層在逐漸增大的逆壓梯度作用下發生局部分離.隨著攻角增加，分離區向前擴展，當迎角達到某個臨界值后上翼面的附體流動被徹底破壞，翼型表面產生很大的壓強差，導致升力降低，阻力顯著增大.由此可見，層流到完全破壞的過程，壓差阻力的貢獻越來越大.

同時，從圖5與表2可以看出，在小攻角的層流狀態下，根據摩擦阻力計算的阻力系數與實驗值接近，誤差約為10%.攻角增大后，隨著分離氣流的出現，壓差阻力的增加，阻力系數計算誤差也顯著增加，4°攻角時為15%，8°攻角時達40%.這主要有2方面原因:一是分離氣流的速度方向不再沿著翼型流向，法向速度分量的壓力貢獻產生的合力增加;二是分離氣流本身的不穩定性，計算模型在模擬分離流動的精度不夠而造成.采用式(2)的阻力計算方法建立在精確模擬流場特性的基礎上，對此，需尋求更好的計算模型和計算精度來進行模擬求解.

3 結論

本研究基于對一維平板邊界層和二維NACA0012翼型采用數值方法得到的摩擦阻力計算結果，可以看出:采用低階精度格式計算低速條件下的層流邊界層能夠得到精度很高的摩擦阻力值，與理論解和實驗數據相當吻合;采用壁面摩擦阻力計算方法在層流狀態下能夠得到準確的結果，但隨著流動分離的出現，壓差阻力的比重增加，采用本方法得到的結果誤差增大，適用性降低.基于上述2點，本研究下一步的工作思路是研究湍流狀態下的摩擦阻力計算方法，包括選擇合適的湍流模型，采用更高精度的計算格式，等等.

［1］Tinoco E N.An assessment of CFD prediction of drag and other longitudinal characteristics［C］//39th Aerospace Sciences Meeting and Exhibit.Reno，NV:AIAA，2011:1002.

［2］Peavey C.Drag prediction of military aircraft using CFD［C］//38th Aerospace Sciences Meeting and Exhibit.Reno，NV:AIAA，2000:383.

［3］葉建，林國華，鄒正平，等.低雷諾數下二維翼型繞流的流場特性分析［J］.航空動力學報，2003，18(2):38 －45.

［4］侯志勇，石磊，聶萬勝.對稱翼型低速繞流流場特性的數值分析研究［J］.科學技術與工程，2009，9(18):5610－5613.

［5］周丹，禹建軍，閆超.層流平板摩擦阻力的數值計算［J］.北京航空航天大學學報，2007，52(6):663－667.

［6］閆超.計算流體力學方法及應用［M］.北京:航空航天大學出版社，2007.

［7］Rice M S.Hand book of airfoil sections for light aircraft［M］.Milwaukee，WI:Aviation Publications，1976.

Study on Numerical Algorithm for Skin Friction Based on Boundary Layers

PENG Cong
(School of Biological Industry，Chengdu University，Chengdu 610106，China)

In order to obtain skin friction of the boundary layers，we adopt the air velocity model to solve the viscid stress.The skin frictions of two cases，1-D laminar flat-plate boundary layer and 2-D laminar NACA0012 airfoil，are calculated by this method.The result of the flat-plate boundary layer matches Blasius solution quite well.Comparing the result of the airfoil calculation with the experimental data，we find that the calculation precision is acceptable at a relatively small angle of attack without airflow separation.As the angle of attack increases and the separation region expands and the proportion of pressure resistance increases，the calculation errors of this method will increase significantly.

skin friction;boundary layer;numerical calculation;air viscosity

V211.3

A

1004－5422(2014)01－0029－04

2013－09－05.

彭 聰(1980—)，女，從事計算機軟件編程與算法研究.