李新凱，戴麗萍，康 順，2，梁思超

（1.華北電力大學電站設備狀態(tài)監(jiān)測與控制教育部重點實驗室，北京102206；2.西安現(xiàn)代控制技術研究所，西安710065）

隨著風力機單機功率的不斷提高，風力機葉片越來越長，為保證強度要求，葉片根部往往用大厚度鈍尾緣翼型［1-2］.葉片運行過程中根部會發(fā)生較大的流動分離，降低風輪對風能的捕獲功率.所以研究大厚度翼型的氣動特性，為改善葉片根部流動狀態(tài)提供理論依據(jù)非常有必要.進行大厚度翼型風洞實驗較為困難，目前高雷諾數(shù)風洞實驗數(shù)據(jù)較少［3-7］，而傳統(tǒng)的雷諾時均（RANS）方法對產生較大分離之后翼型氣動性能的預測效果較差，所以需要一種可以較準確地預測大厚度翼型氣動性能的數(shù)值計算方法.

Li等［8］利用大渦模擬（LES）方法對2.5D的NACA0018翼型在全攻角下進行了數(shù)值計算，結果表明在大部分攻角下LES方法均能較好地預測翼型的升力系數(shù)與阻力系數(shù)，效果明顯優(yōu)于非定常雷諾方法（URANS）方法.李棟等［9］采用分離渦模擬（DES）計算方法對3種小厚度的NACA 翼型進行了數(shù)值計算，通過與實驗結果比較，顯示該方法可以有效地預測翼型的失速特征.Kravchenko等［10］利用LES方法分析了不同網格密度及展向寬度的圓柱繞流，通過計算發(fā)現(xiàn)圓柱展向寬度及展向網格分辨率對計算結果有很大影響.高偉等［11］采用勢流方程與邊界層耦合方法對不同厚度翼型邊界層轉捩進行了研究，結果表明轉捩位置對翼型升阻力系數(shù)有一定影響.

目前，國內外學者針對風力機專用大厚度鈍尾緣翼型進行的數(shù)值模擬研究還很少，筆者以具有實驗數(shù)據(jù)的DU97-W-300和DU00-W2-401翼型為研究對象，對其進行了數(shù)值模擬研究.主要研究內容包括以下2個方面：（1）利用RANS方法對翼型進行計算，討論不同湍流模型（包括是否考慮轉捩）以及網格分布對RANS方法計算結果的影響；（2）分別采用RANS 和DES 方法對DU00-W2-401翼型在多個攻角下進行數(shù)值計算，對比RANS 與DES 方法的計算結果.

1 幾何模型及數(shù)值方法

1.1 幾何模型

圖1為DU97-W-300與DU00-W2-401翼型的幾何模型.其中，DU97-W-300 翼型弦長C＝0.6 mm，相對厚度為30%C，鈍尾緣，尾緣厚度為1.74%C；DU00-W2-401翼型弦長C＝0.6 mm，相對厚度為40%C，鈍尾緣，尾緣厚度為1%C.

圖1 幾何模型Fig.1 Geometry model

1.2 數(shù)值方法

用ICEM 軟件進行網格劃分，生成全域結構網格，第一次網格高度0.001mm，網格數(shù)目分布見表1，翼型網格分布見圖2.

圖2 翼型周圍網格分布Fig.2 Mesh distribution around the airfoil

邊界條件：所有計算均采用圓形計算域，計算域半徑R＝30C，采用速度進口、壓力出口邊界條件，葉片表面設置為無滑移固體壁面.

采用商用軟件Fluent進行數(shù)值計算，定常計算采用RANS 方法，非定常計算采用DES 方法.RANS方法中湍流模型包括SA、SST 以及考慮轉捩的SST（T-SST）湍流模型，DES方法中湍流模型為SA湍流模型，具體方案設置見表1.采用有限體積方法對控制方程進行離散，壓力-速度耦合基于Simple算法.RANS方法中控制方程的各項均采用二階迎風格式，DES方法中壓力項采用二階迎風格式，動量項采用中心差分格式.

表1 計算方案Tab.1 Calculation plan

2 計算結果及分析

2.1 湍流模型及網格對RANS方法計算結果的影響

圖3 給出了不同湍流模型下DU97-W-300 翼型的計算結果.圖中，Cl為升力系數(shù)，Cd為阻力系數(shù)α為攻角.從圖3 可以看出，在翼型失速之前，TSST 湍流模型更能準確地預測翼型的氣動性能，升力系數(shù)、阻力系數(shù)及升阻比更接近實驗值.但在翼型失速之后，全湍流模型（SA 和SST 模型）比T-SST湍流模型計算結果更接近實驗值.對比可以發(fā)現(xiàn)，無論是失速之前還是失速之后，相比全湍流模型，T-SST湍流模型計算得到的升力系數(shù)偏大，阻力系數(shù)偏小，升阻比偏大.這是由于考慮轉捩之后，葉片前緣部分流動變?yōu)閷恿?，從而使得葉片升力系數(shù)偏大，阻力系數(shù)偏小，這與文獻［12］和文獻［13］中得出的結論一致.全湍流模型中SA 湍流模型計算效果要優(yōu)于SST 湍流模型，計算值更加接近實驗值.

圖3 DU97-W-300翼型計算結果Fig.3 Computation result of DU97-W-300airfoil

圖4（a）為DU97-W-300翼型表面壓力系數(shù)Cp分布曲線.從圖4（a）可以看出，在9.3°、12.4°攻角下，T-SST 湍流模型的計算結果與實驗值更加吻合，SST 湍流模型計算效果最差.轉捩是翼型繞流中重要的流動現(xiàn)象，轉捩流動現(xiàn)象復雜，一直是CFD 模擬中的難點，翼型轉捩位置的預測從一定程度上反映了CFD 的計算精度.隨著對翼型繞流流動細節(jié)的研究，邊界層轉捩成為重要的研究對象.圖4（b）為CFD 預測翼型轉捩位置與實驗值的對比曲線.從圖4（b）可以看出，隨著攻角的增大，翼型吸力面轉捩位置向翼型前緣移動，壓力面轉捩位置向翼型尾緣移動.從計算轉捩位置來看，T-SST 湍流模型計算得到的翼型壓力面和吸力面轉捩位置與實驗值吻合較好，說明考慮轉捩的T-SST 湍流模型能較好地預測翼型轉捩位置.

圖4 DU97-W-300翼型表面壓力系數(shù)及轉捩位置分布Fig.4 Distribution of Cpand transition position for DU97-W-300airfoil

由上面的分析可以看出，在翼型線性段，T-SST湍流模型能更好地預測葉片升力系數(shù)、阻力系數(shù)及轉捩位置.筆者通過改變翼型一周網格數(shù)（方案A）及第一層網格高度ΔY（方案B）來研究網格對預測轉捩位置準確性的影響.具體方案及計算結果見表2，計算攻角選為α＝12.4°.

圖5為不同網格數(shù)及網格高度下DU97-W-300翼型表面壓力系數(shù)及轉捩位置分布圖.從圖5可以看出，方案A 中1、2、3、4計算所得Cp均與實驗值吻合較好，且均能正確預測轉捩位置，所以翼型一周網格數(shù)對預測轉捩位置的影響較小.方案B 中5、6能較準確地預測出轉捩位置，且Cp與實驗值吻合較好，但方案B 中7、8未能計算出翼型表面的轉捩位置，且Cp與實驗值相差較大，由此可見對于預測轉捩位置，第一層網格高度比較重要，當ΔY／C≥1.66×10-4（即y＋≥18）時已經不能預測出轉捩位置.

表2 網格對預測轉捩效果的影響Tab.2 Influence of mesh distribution on the predicted transition effects

圖5 不同網格數(shù)及網格高度下DU97-W-300翼型表面壓力系數(shù)及轉捩位置分布Fig.5 Distribution of surface pressure coefficient and transition position for DU97-W-300airfoil at different grid number and height

2.2 DES與RANS方法計算結果的對比

圖6給出了DU00-W2-401翼型DES與RANS方法計算結果的對比.其中計算攻角為5°～23°.RANS方法中湍流模型選擇SA 湍流模型，DES方法中近壁區(qū)域湍流模型也選擇SA 湍流模型.從圖6可以看出，相比于RANS方法的計算結果，DES方法能更好地預測大厚度翼型的氣動特性.相比實驗值，在翼型失速段，RANS 方法計算所得升力系數(shù)偏大，阻力系數(shù)偏小，升阻比偏大.與DES方法計算結果相比，在翼型失速之后，RANS 方法計算所得升力系數(shù)偏大，阻力系數(shù)偏小，升阻比偏大.DES方法計算得到的翼型升力系數(shù)、阻力系數(shù)及升阻比與實驗值吻合更好.從計算結果來看，在翼型失速之后，DES方法的模擬結果明顯優(yōu)于RANS方法.

圖6 DU00-W2-401翼型計算結果Fig.6 Computation result of DU00-W2-400airfoil

圖7 為DU00-W2-401 翼型表面壓力系數(shù)分布.從圖7可以看出，在15°、19°攻角時，RANS計算結果的吸力峰值均比DES方法的計算結果要高，所以計算所得氣動力較DES方法的計算結果偏大.

圖8為DU00-W2-401翼型中間對稱面上渦量云圖.從圖8 可以看出，DES方法可以模擬出更細致的漩渦結構.從DES計算結果可以看出在分離區(qū)的上游有組織、有規(guī)律地脫落出分離渦，而在翼型的尾緣處，翼型上表面分離區(qū)內是低壓區(qū)，翼型下表面是高壓區(qū)，翼型下表面的流體會繞過葉片尾緣卷入分離區(qū)，這樣分離渦周而復始向下游脫落.

圖7 DU00-W2-401翼型表面壓力系數(shù)分布Fig.7 Distribution of surface pressure coefficient

圖8 渦量云圖Fig.8 Cloud map of vorticity

圖9 為帶速度云圖的Q等值面（Q＝1 000 s-2）.從圖9可以看出，隨著攻角的增大，翼型上表面的湍流分離渦尺度越來越大，分離區(qū)尺度越來越大.由于大厚度翼型極易發(fā)生流動分離，分離之后RANS方法不能模擬出細致的漩渦結構，導致其預測氣動性能不準確，而DES方法能較好地模擬出細致的湍流脈動漩渦，所以DES 方法較RANS 方法更適合大厚度翼型的氣動計算.

圖9 帶速度云圖的Q 等值面（Q＝1 000s-2）Fig.9 Isosurface of Q with velocity contour

3 結 論

（1）對于RANS方法的計算，在翼型最大升力系數(shù)攻角之前，考慮轉捩的T-SST 湍流模型的計算效果要優(yōu)于全湍流模型，在翼型最大升力系數(shù)攻角之后，全湍流模型的計算效果要優(yōu)于T-SST 湍流模型，全湍流模型中SA 湍流模型的計算效果優(yōu)于SST 湍流模型.網格對計算結果及轉捩預測的影響中，第一層網格高度更加敏感，而翼型一周網格數(shù)的影響不大.

（2）對比RANS與DES方法的計算結果發(fā)現(xiàn)，當葉片存在大尺度分離時，RANS方法不能模擬出更細致的分離渦結構，而DES方法則能在一定程度上模擬出細致規(guī)律的分離渦結構.對比氣動性能計算結果，DES 方法計算結果與實驗值更加吻合.DES方法比RANS方法更適合大厚度翼型的氣動性能計算.

［1］TPI Composites.Innovative design approaches for large wind turbine blades-final report［R］.Sandia，USA：TPI，2004.

［2］STANDISH K J，van DAM C P.Aerodynamic analysis of blunt trailing edge airfoils［J］.Journal of Solar Energy Engineering，2003，125（4）：479-487.

［3］BAKER J，MAYDA E，DAM C.Experimental and computational analysis of thick fatback wind turbine airfoils［R］.Reno，USA：AIAA，2006.

［4］LAW S P，GREGNREK G M.Wind turbine evaluation of a truncated NACA 64-621airfoil for wind turbine applications［R］.Ohio，USA：NASA，1987.

［5］TPI Composites.Flatback airfoil wind tunnel experiment［R］.Sandia，USA：TPI，2008.

［6］POST M L，JONES R B.Characterization of a flatback airfoil for use in wind power generation［R］.Reno，USA：AIAA，2008.

［7］BARONE M F，BERG D.Aerodynamic and aeroacoustic properties of a flatback airfoil：an update［R］.Oriando，USA：AIAA，2009.

［8］LI Chao，ZHU Songye，XU Youlin，etal.2.5Dlarge eddy simulation of vertical axis wind turbine in consideration of high angle of attack flow［J］.Renewable Energy，2013，51（3）：317-330.

［9］李棟，焦予秦，IGOR Men＇shov，等.Detached-Eddy Simulation方法模擬不同類型翼型的失速特性［J］.航空學報，2005，26（4）：406-410.

LI Dong，JIAO Yuqin，IGOR Men＇shov，etal.Detached-Eddy Simulation for airfoil stall［J］.Acta Aeronautica et Astronautica Sinica，2005，26（4）：406-410.

［10］KRAVCHENKO A G，MOIN P.Numerical studies of flow over a circular cylinder atReD＝3 900［J］.Physics of Fluids，2000，12（2）：403-417.

［11］高偉，李春，高月文，等.幾何參數(shù)對風力機翼型轉捩特性的影響［J］.動力工程學報，2013，33（6）：490-496.

GAO Wei，LI Chun，GAO Yuewen，etal.Influence of geometric parameters on transition characteristics of wind turbine airfoils［J］.Journal of Chinese Society of Power Engineering，2013，33（6）：490-496.

［12］范忠瑤.風力機定常與非定常氣動問題的數(shù)值模擬研究［D］.北京：華北電力大學，2011.

［13］MICHEL R，ARNAL D，COUSTOLS E.Stability calculations and transition criteria on two-or three-dimensional flows［J］.Laminar-Turbulent Transition，1984（7）：455-461.