戴 昱， 陳德龍， 信志強(qiáng)

（河海大學(xué)力學(xué)與材料學(xué)院，南京 211100）

隨著時(shí)代的進(jìn)步，微型飛行器的機(jī)翼、飛機(jī)的機(jī)翼以及風(fēng)機(jī)葉片都在飛速發(fā)展，這對(duì)能源的利用和信息的收集都有著巨大的幫助.因此對(duì)翼型進(jìn)行深入研究和改良創(chuàng)新是十分重要的.

如何提高翼型的氣動(dòng)性能是研究的核心問題.目前已有眾多研究采用主動(dòng)控制的方法使翼型發(fā)生變形來提升翼型的氣動(dòng)性能.張敬斌等［1］通過主動(dòng)控制技術(shù)研究翼型尾緣部分的變形，發(fā)現(xiàn)翼型尾緣部分發(fā)生適當(dāng)變形能夠有效地提高翼型的升力系數(shù)和升阻比.而Zhuang等［2］的研究發(fā)現(xiàn)，在一定范圍內(nèi)，隨著翼型尾緣部分主動(dòng)變形的范圍增大，翼型的氣動(dòng)性能也會(huì)隨之提高.對(duì)于俯仰振蕩翼型，張俊偉等［3］研究的結(jié)果表明，在翼型上仰過程中翼型主動(dòng)變形量的增大可以抑制渦的分離.Kan等［4］則對(duì)翼型后緣的主動(dòng)變形進(jìn)行了進(jìn)一步的研究，發(fā)現(xiàn)變形的角度如果發(fā)生周期性的變化能夠增大升力系數(shù)并且可以延遲失速.尾緣主動(dòng)變形的翼型和普通的襟翼相比存在一些差異，在這之前，普通襟翼的研究已有一些，例如賈亞雷等［5-6］的研究，而由于無縫光滑的特征，翼型尾緣的變形具有更優(yōu)越的氣動(dòng)性能［7］.通過這些研究可以看出，對(duì)翼型施加主動(dòng)控制使翼型變形從而提高氣動(dòng)性能的方法具有較為優(yōu)越的效果.雖然主動(dòng)變形控制的方法對(duì)氣動(dòng)性能的提升已經(jīng)呈現(xiàn)出較為優(yōu)越的效果，但它需要主動(dòng)控制系統(tǒng)，這往往導(dǎo)致其結(jié)構(gòu)較為復(fù)雜，被動(dòng)控制則相對(duì)簡(jiǎn)單.被動(dòng)控制的方式主要是對(duì)翼型采用柔性的結(jié)構(gòu)，使其可以產(chǎn)生自適應(yīng)變形從而提升翼型的氣動(dòng)性能.自適應(yīng)變形的潛能是巨大的［8］.Macphee和Beyene［9］采用被動(dòng)控制，利用流固耦合的方法對(duì)柔性翼型的氣動(dòng)性能展開了研究.研究表明，在垂直軸風(fēng)機(jī)中柔性翼型的功率系數(shù)要大于剛性翼型，即使是低彈性模量的柔性翼型，相比于剛性翼型依然可以提高翼型的氣動(dòng)特性［10］.柔性翼型在低于風(fēng)機(jī)的額定風(fēng)速下所表現(xiàn)出來的性能要優(yōu)于剛性翼型［11］.前人對(duì)柔性翼型的研究中還未考慮到柔性翼型的動(dòng)態(tài)失速特性，而動(dòng)態(tài)失速現(xiàn)象在風(fēng)力機(jī)和各種飛行器中又是普遍存在的［12-16］.因此對(duì)于柔性翼型在動(dòng)態(tài)失速現(xiàn)象中所表現(xiàn)出來的氣動(dòng)特性進(jìn)行研究是十分有必要的.

本文的研究是在前人對(duì)柔性翼型研究的基礎(chǔ)上，進(jìn)一步考慮柔性翼型在動(dòng)態(tài)失速現(xiàn)象中是否還具有一定的價(jià)值.因此本文基于流固耦合的方法，研究柔性翼型在不同減縮頻率下的動(dòng)態(tài)失速特性.通過對(duì)比剛?cè)嵝砸硇偷臍鈩?dòng)特性、流場(chǎng)特點(diǎn)以及分析柔性翼型的結(jié)構(gòu)響應(yīng)，得到柔性翼型的動(dòng)態(tài)失速特性及其優(yōu)勢(shì).

1 問題描述

本文采用的柔性翼型的模型如圖1所示，以NACA0012為基礎(chǔ)翼型，弦長(zhǎng)c=1 m，并在距離翼型前端0.25 m處設(shè)置了一個(gè)0.02 m×0.06 m的矩形孔洞作為旋轉(zhuǎn)軸，整個(gè)翼型除矩形孔洞外，其余部分均為柔性材料.

圖1 柔性翼型模型Fig.1 Flexible airfoil model

來流介質(zhì)為空氣，密度取為1.225 kg/m3，來流速度為4.5 m/s，動(dòng)力黏度取1.837 5×10-5Pa·s，雷諾數(shù)為30 萬.參考前人對(duì)于柔性翼型材料參數(shù)的選擇，本文中的柔性材料密度為100 kg/m3，泊松比為0.37，這種材料可用聚氨酯泡沫實(shí)現(xiàn)，具有低密度和廣泛的彈性模量的特征［9］.在本次動(dòng)態(tài)失速研究中，通過翼型繞矩形旋轉(zhuǎn)軸做俯仰運(yùn)動(dòng)來模擬攻角的非定常變化，攻角的變化運(yùn)動(dòng)規(guī)律如式（1）所示：

其中：α0為初始攻角；αt為翼型俯仰過程中最大攻角；T為俯仰運(yùn)動(dòng)周期；n為俯仰運(yùn)動(dòng)周期的個(gè)數(shù).前3 s 翼型固定不動(dòng)以獲取穩(wěn)定流場(chǎng)，3 s后做角速度大小不變的俯仰運(yùn)動(dòng).

減縮頻率是翼型動(dòng)態(tài)失速的重要影響因素［17］，其表達(dá)式為k=，其中ω表示振蕩角速度，c表示弦長(zhǎng)，u為來流速度.本文針對(duì)剛?cè)醿煞N翼型分別設(shè)置了三種減縮頻率，具體工況設(shè)置如表1所示.

表1 工況設(shè)置Tab.1 Parameters of various working conditions

2 數(shù)值方法及驗(yàn)證

2.1 控制方程

2.1.1 流體控制方程 對(duì)于不可壓縮流體而言，流體密度為常數(shù)，質(zhì)量守恒方程可寫為

式中：i為下標(biāo)，表示坐標(biāo)軸的方向；u代表流體速度.

對(duì)于湍流情況，雷諾時(shí)均方程可表示為

式中：j與前者i同為下標(biāo)；uˉ為時(shí)均流速；u′為脈動(dòng)流速；t代表時(shí)間；f表示質(zhì)量力；ρ為流體密度；p代表壓強(qiáng)；μ為流體的動(dòng)力黏度.而雷諾時(shí)均方程相比于N-S方程多出了一項(xiàng)雷諾應(yīng)力，該項(xiàng)為對(duì)稱的二階張量，與速度梯度有關(guān).這導(dǎo)致了雷諾時(shí)均方程無法封閉，需要建立湍流模型.本文將采用SST K-omega湍流模型［9，18］，該模型結(jié)合了k-omega模型和k-epsilon模型，能夠精確地捕捉近壁面的流動(dòng)分離.

2.1.2 固體控制方程 柔性翼型采用的是齊次和各向同性結(jié)構(gòu)，可使用以下方程［19］：

其中：u為位移場(chǎng)；ρs為固體的密度；P為皮奧拉-基爾霍夫第一應(yīng)力張量；S表示皮奧拉-基爾霍夫第二應(yīng)力張量；F代表變形梯度；E代表格林應(yīng)變張量；I表示特征張量；μ和λ為拉梅參數(shù).

2.1.3 流固耦合 在耦合面上，流體模型和固體模型須滿足運(yùn)動(dòng)學(xué)條件和動(dòng)力學(xué)條件進(jìn)行耦合.運(yùn)動(dòng)學(xué)條件是指耦合面上的速度和位移必須連續(xù)，可用以下表達(dá)式進(jìn)行描述：

式中：Vf、uf、Vs、us分別為液體的速度、位移，固體的速度、位移.

動(dòng)力學(xué)條件指的是耦合面上的力必須處于平衡狀態(tài)，可用以下表達(dá)式描述：

式中：n為法線方向；σf、σs分別為液體、固體的應(yīng)力.

流體側(cè)和固體側(cè)通過在耦合面上進(jìn)行變量交換來完成耦合.在每一個(gè)時(shí)間步內(nèi)，都有流體部分提供壓力和黏性力給固體，固體部分提供位移增量和速度給流體［20］.一個(gè)時(shí)間步內(nèi)的耦合求解流程如圖2所示.

圖2 一個(gè)時(shí)間步內(nèi)耦合求解流程Fig.2 The process of fluid-structure interaction solution within one time step

2.2 數(shù)值方法驗(yàn)證

為了驗(yàn)證本文中流固耦合的流場(chǎng)及固體變形計(jì)算的準(zhǔn)確性，采用Krawczyk等［11］的研究，定量對(duì)比不同來流速度下柔性翼型變形的大小和升阻力系數(shù).計(jì)算的模型、參數(shù)等均與參考文獻(xiàn)保持一致.即采用C型計(jì)算區(qū)域，柔性翼型以弦長(zhǎng)1 m的NACA4412為基礎(chǔ)翼型.本文取文獻(xiàn)中彈性模量為2.5 MPa下的五種來流速度，分別為3.79、5.04、6.309、7.59、10.22 m/s，對(duì)比其變形大小和升阻力系數(shù).

計(jì)算的結(jié)果與文獻(xiàn)中的數(shù)據(jù)對(duì)比如圖3所示.圖中給出的變形大小及升阻力系數(shù)都是經(jīng)過歸一化處理的.以柔性翼型的尾緣部分在來流速度為7.59 m/s下的變形大小為0 m，其余工況下的變形大小為與其的差值，其中變形為正代表尾緣部分相對(duì)變形方向朝上，變形為負(fù)代表尾緣部分相對(duì)變形方向朝下；而升阻力系數(shù)則為與剛性翼型在來流速度為7.59 m/s工況下所得的升阻力系數(shù)的比值.可以看出本文的計(jì)算結(jié)果與文獻(xiàn)相比整體吻合較好，可以證明本文對(duì)柔性翼型的流固耦合計(jì)算是準(zhǔn)確的.

圖3 驗(yàn)證算例結(jié)果對(duì)比Fig.3 Comparison of calculation results of verification cases

2.3 計(jì)算域網(wǎng)格劃分及邊界條件

計(jì)算域的左半邊為半圓形，其半徑為3.5 m，右半邊為長(zhǎng)8.5 m、寬7 m 的矩形.翼型的尾部與半圓形的圓心重合，整個(gè)計(jì)算域的C形網(wǎng)格劃分如圖4所示.左端設(shè)定為速度入流邊界條件：u=4.5 m/s，v=0 m/s；上下設(shè)定為對(duì)稱邊界條件：?u/?y=v=0；右端出口設(shè)定為紐曼邊界條件：?u/?x=?v/?x=0；翼型則設(shè)定為流固耦合邊界條件.為了使邊界層區(qū)域內(nèi)的流場(chǎng)求解更加精確，壁面第一層網(wǎng)格y+設(shè)置為0.7，增長(zhǎng)率取1.15，同時(shí)，為了驗(yàn)證網(wǎng)格的無關(guān)性，對(duì)計(jì)算域劃分了三種不同數(shù)量的網(wǎng)格，網(wǎng)格數(shù)量分別為15 700、27 890和43 560并進(jìn)行了計(jì)算.結(jié)果表明，后兩者的升阻力系數(shù)最大僅相差0.5%，固體最大變形相差0.7%，而前兩者升阻力系數(shù)相差23%，固體最大變形相差18%，考慮到計(jì)算資源和時(shí)間成本，本文采用網(wǎng)格數(shù)量為27 890的網(wǎng)格進(jìn)行計(jì)算.對(duì)于時(shí)間步長(zhǎng)，取庫(kù)朗數(shù)為0.1時(shí)計(jì)算得出的時(shí)間步長(zhǎng)為1.5×10-5s，滿足計(jì)算所需的精度要求.

圖4 計(jì)算域網(wǎng)格Fig.4 Grid of the computational domain

3 計(jì)算結(jié)果與分析

3.1 不同彈性模量下的固定攻角柔性翼型氣動(dòng)特性對(duì)比

為了便于研究減縮頻率對(duì)翼型動(dòng)態(tài)失速特性的影響，首先對(duì)四組不同彈性模量下的固定攻角柔性翼型進(jìn)行模擬，分析不同彈性模量下的升阻比，選取使翼型具有最佳氣動(dòng)性能的彈性模量作為動(dòng)態(tài)失速時(shí)柔性翼型的彈性模量.

如圖5 所示為四種彈性模量下的升阻比.可以看出，在2.5°～17.5°攻角范圍內(nèi)，四種彈性模量的柔性翼型的升阻比均高于剛性翼型，這表明柔性翼型相對(duì)于剛性翼型普遍具有較優(yōu)的氣動(dòng)性能.在本次模擬的四種彈性模量中，彈性模量為3×107Pa 的柔性翼型在5°～17.5°攻角范圍內(nèi)所具有的氣動(dòng)特性都是最優(yōu)的，因此，對(duì)于俯仰振蕩中的柔性翼型，彈性模量均選用3×107Pa.

圖5 柔性翼型在不同彈性模量下的升阻比Fig.5 Lift-to-drag ratios of flexible airfoil under different elastic moduli

3.2 不同減縮頻率下剛?cè)嵋硇偷臍鈩?dòng)特性分析

柔性翼型在不同減縮頻率下的氣動(dòng)特性具有明顯的區(qū)別，且同一減縮頻率下，剛性翼型和柔性翼型也存在區(qū)別.

圖6 所示為剛性翼型和柔性翼型在不同減縮頻率下的升阻力系數(shù)遲滯回線圖.柔性翼型的升阻力系數(shù)遲滯回線除了局部具有波動(dòng)外，整體變化趨勢(shì)與剛性翼型保持一致，均有以下規(guī)律：當(dāng)k=0.08 時(shí)，升力系數(shù)在上仰過程中一直保持較低的狀態(tài)，在下俯過程中一直保持較高的狀態(tài)，阻力系數(shù)在整個(gè)過程中沒有出現(xiàn)大幅度的改變；當(dāng)k=0.024時(shí)，在攻角大于7°后升力系數(shù)一直大于k=0.08時(shí)的值，但在攻角變小的過程中卻突然驟降，直到18.5°時(shí)才回升，阻力系數(shù)在20°之后明顯提升；當(dāng)k=0.009 7時(shí)，上仰過程中升力系數(shù)在22°時(shí)發(fā)生了大幅度下降，在24°時(shí)出現(xiàn)回升，下俯過程中由25°至24°時(shí)有短暫的上升，隨后一直處于下降趨勢(shì)，阻力系數(shù)的變化趨勢(shì)則與升力系數(shù)相同.

圖6 剛性翼型和柔性翼型升阻力系數(shù)遲滯回線Fig.6 Hysteresis loops of the lift-drag coefficients for rigid airfoil and flexible airfoil

再通過比較同一減縮頻率下，剛性翼型和柔性翼型的升阻力系數(shù)遲滯回線可以發(fā)現(xiàn)，在k較小時(shí)，二者無明顯區(qū)別；而k較大時(shí)，即k=0.08時(shí)，柔性翼型的升力系數(shù)整體略高于剛性翼型，阻力系數(shù)整體小于剛性翼型，且阻力系數(shù)的最大值明顯小于剛性翼型.由此可見，在高減縮頻率下，柔性翼型的氣動(dòng)特性要顯著優(yōu)于剛性翼型，具有高升阻比和低阻力的特點(diǎn).

3.3 柔性翼型的結(jié)構(gòu)響應(yīng)分析

在柔性翼型的俯仰過程中，柔性翼型會(huì)受到流體的作用力而發(fā)生自適應(yīng)變形，柔性翼型的變形具有一定的規(guī)律.

圖7 展示了柔性翼型在k=0.08 時(shí)整個(gè)俯仰過程中尾緣變形量的變化，其中上仰至25°攻角時(shí)的變形及應(yīng)力在圖8 中進(jìn)行展示.通過圖8 中展示的應(yīng)力可以看出高應(yīng)力主要集中在矩形固定孔洞附近及離孔洞最近的上下表面處，從孔洞至尾緣應(yīng)力逐漸減小.觀察圖8中的變形分布規(guī)律可以明顯看到，在整個(gè)過程中柔性翼型的變形主要集中在尾緣部分，其余部位的變形量接近零.變形量從翼型的中間至尾緣逐漸增大，在尾部達(dá)到最大值.從定性的角度來看，變形的分布規(guī)律與Macphee和Beyene［9］的研究結(jié)果不謀而合，但從定量的角度來分析，動(dòng)態(tài)失速中的柔性翼型的變形在數(shù)值及變化規(guī)律上又有其獨(dú)有的特征.通過圖7可以看出在整個(gè)俯仰過程中，變形的大小隨著攻角的增大（減小）快速地上下波動(dòng)，這表明柔性翼型的變形十分靈活，對(duì)氣動(dòng)荷載較為敏感.在上仰過程中，隨著攻角的增大，變形的數(shù)值整體趨勢(shì)是略有降低的，而在下俯過程中卻相反，隨著攻角的減小，其變形的整體趨勢(shì)也隨之減小.上仰階段和下俯階段的變形差距存在一定的規(guī)律，在小攻角范圍內(nèi)（0°～12.5°），上仰階段和下俯階段變形大小的差距并不明顯，而在大攻角范圍內(nèi)（12.5°～25°），下俯階段的變形在數(shù)值上大約是上仰階段的1.5倍至2倍，差距十分明顯.

圖7 柔性翼型俯仰過程中尾部變形量的變化Fig.7 Deformation variations of flexible airfoil trailing edges during pitching

圖8 柔性翼型上仰至25°攻角的變形及應(yīng)力Fig.8 Deformations and stresses of flexible airfoil pitching up to 25°angle of attack

3.4 不同減縮頻率下剛?cè)嵋硇偷牧鲌?chǎng)分析

本小節(jié)將通過分析流場(chǎng)來對(duì)柔性翼型動(dòng)態(tài)失速的機(jī)理以及剛?cè)嵋硇椭g氣動(dòng)特性存在差異的原因進(jìn)行研究.

圖9為柔性翼型在k=0.08時(shí)上仰過程的渦量和流線.圖10為柔性翼型在k=0.009 7時(shí)上仰過程的渦量和對(duì)應(yīng)的流線.通過對(duì)比圖9和圖10中的流場(chǎng)可以發(fā)現(xiàn)，上仰至15°攻角時(shí)，二者尾流均出現(xiàn)了波動(dòng)，尾緣上表面出現(xiàn)回流現(xiàn)象，且k=0.08的流場(chǎng)尾流部分已經(jīng)出現(xiàn)渦街；上仰至25°時(shí)，二者出現(xiàn)了明顯區(qū)別，k=0.08的流場(chǎng)尾流部分的渦街增強(qiáng)，翼型尾緣上表面出現(xiàn)流動(dòng)分離現(xiàn)象，而k=0.009 7的流場(chǎng)已經(jīng)出現(xiàn)渦旋從翼型表面大規(guī)模脫落的情況，翼型已經(jīng)處于失速狀態(tài)，這導(dǎo)致升力系數(shù)在上仰過程中突然下降.低減縮頻率下的翼型在上仰過程出現(xiàn)失速現(xiàn)象，表明了高減縮頻率下的翼型相比低減縮頻率具有失速延遲的特性.

圖9 柔性翼型上仰過程渦量、流線、變形及應(yīng)力，k=0.08Fig.9 Vorticity contours，streamlines，deformations and stresses of the flexible airfoil during pitching up，k=0.08

圖10 柔性翼型上仰過程渦量及流線，k=0.009 7Fig.10 Vorticity contours and streamlines of the flexible airfoil during pitching up，k=0.009 7

圖11所示為剛性翼型在k=0.08時(shí)上仰至25°攻角的渦量，與圖9（b）中相同條件下的柔性翼型相比，二者的區(qū)別十分明顯，剛性翼型此時(shí)已經(jīng)在上表面產(chǎn)生了集中強(qiáng)渦旋，而柔性翼型上表面僅僅是尾緣部分出現(xiàn)了流動(dòng)分離，并未出現(xiàn)渦旋.圖8展示了此時(shí)柔性翼型因?yàn)槭艿搅黧w的作用力所產(chǎn)生的自適應(yīng)變形，該變形會(huì)反作用于流場(chǎng)而改變流場(chǎng)的特性，此時(shí)柔性翼型的變形抑制了翼型表面渦旋的生成，對(duì)流場(chǎng)起到了一定的調(diào)控作用，使柔性翼型表面僅出現(xiàn)了流動(dòng)分離.通過對(duì)流場(chǎng)差異性的分析可以得出，剛?cè)嵋硇偷臍鈩?dòng)特性差異是由變形對(duì)流場(chǎng)的調(diào)控導(dǎo)致的.再通過對(duì)應(yīng)圖6（a）和（b）中上仰至25°時(shí)的升阻力系數(shù)發(fā)現(xiàn)柔性翼型的升力系數(shù)僅有微小提升，而阻力系數(shù)出現(xiàn)了大幅度降低.

圖11 剛性翼型上仰至25°攻角渦量及流線，k=0.08Fig.11 Vorticity contours and streamlines of the rigid airfoil pitching up to 25°of attack，k=0.08

4 結(jié)論

本文采用流固耦合數(shù)值模擬方法，研究了柔性翼型在不同減縮頻率下其升阻力系數(shù)和流場(chǎng)的變化規(guī)律，通過結(jié)構(gòu)響應(yīng)分析揭示了剛?cè)嵋硇蜌鈩?dòng)特性具有明顯差異的物理機(jī)理.研究發(fā)現(xiàn)即使在動(dòng)態(tài)失速中柔性翼型依然具有較優(yōu)的氣動(dòng)性能，這對(duì)未來翼型的發(fā)展具有一定的指導(dǎo)意義.具體的研究結(jié)論如下：

1）柔性翼型相比剛性翼型在氣動(dòng)特性上所展現(xiàn)的優(yōu)勢(shì)是大幅度減阻和小幅度增升.2）柔性翼型的被動(dòng)變形主要產(chǎn)生在尾緣部分，變形在動(dòng)態(tài)失速下的分布規(guī)律與前人研究的靜態(tài)中的變形分布規(guī)律相同，但區(qū)別在于動(dòng)態(tài)失速中變形的大小會(huì)隨攻角的變化上下波動(dòng).3）柔性翼型良好的氣動(dòng)特性取決于柔性翼型變形對(duì)流場(chǎng)調(diào)控的結(jié)果.相比于剛性翼型，該變形能夠抑制表面渦旋的產(chǎn)生.