基于CFD/CSD耦合的全動(dòng)平尾氣動(dòng)彈性特性研究

羅文莉，陸琪

(上海飛機(jī)設(shè)計(jì)研究院，上海 201210)

0 引言

對于多數(shù)常規(guī)飛行器建模，采用適當(dāng)?shù)膭傮w假設(shè)是合理的。但是對于相對厚度較小的活動(dòng)面，如水平尾翼來說，很容易出現(xiàn)氣動(dòng)/結(jié)構(gòu)耦合的現(xiàn)象，因此剛體假設(shè)是不可行的。尤其是隨著飛行速度的不斷增大，氣動(dòng)彈性的影響已不容忽視。因此，對于全動(dòng)平尾的氣動(dòng)彈性研究是十分必要的。

關(guān)于氣動(dòng)彈性問題的研究從20世紀(jì)50年代開始興起，由于計(jì)算能力的限制，那時(shí)的研究主要集中于風(fēng)洞試驗(yàn)以及非定常氣動(dòng)理論。LAUTEN W T等[1]在氣動(dòng)彈性方面對X-15的全動(dòng)水平尾翼縮比模型進(jìn)行了風(fēng)洞試驗(yàn)研究?；诨钊碚撚?jì)算得到的顫振速度大約是試驗(yàn)得到速度的4倍。HEEG J等[2]對6種全動(dòng)尾翼模型進(jìn)行了風(fēng)洞試驗(yàn)，主要針對翼型、翼軸剛度以及翼面質(zhì)量分布進(jìn)行了研究，全部試件的試驗(yàn)顫振速度均低于使用2階活塞理論計(jì)算得到的值。

隨著計(jì)算能力的不斷提高，非定常氣動(dòng)力的計(jì)算方法逐漸由依靠理論分析轉(zhuǎn)向了與CFD計(jì)算相結(jié)合的方法。MCNAMARA J J等[3-4]基于活塞理論、牛頓理論、激波膨脹波等非定常氣動(dòng)力方法計(jì)算了一種雙楔形翼型的氣彈穩(wěn)定性，并與基于求解N-S方程的CFD方法進(jìn)行了對比，結(jié)果表明使用1階、2階活塞理論的誤差較大，并說明了黏性效應(yīng)對于二維楔形薄翼型的影響基本可以忽略。

此后，集合了CFD流場分析和結(jié)構(gòu)有限元計(jì)算(FEA)的CAE分析由于其在氣彈分析中的先進(jìn)性成為了一股研究趨勢。MCNAMARA J J等[5]對X-33進(jìn)行了氣動(dòng)彈性分析，使用MSC.NASTRAN進(jìn)行模態(tài)分析。GUPTA K K等[6]分別采用了基于活塞理論、CFD以及系統(tǒng)識(shí)別等方法計(jì)算非定常氣動(dòng)力，結(jié)構(gòu)采用有限元建模，得到的結(jié)果表明飛行器在工作環(huán)境中定常區(qū)域內(nèi)不會(huì)發(fā)生顫振現(xiàn)象。

采用常規(guī)的理論算法或是基于CFD的非定常氣動(dòng)力算法無法考慮翼型的變形。對于結(jié)構(gòu)模型的求解，國內(nèi)外大多理論方法均假設(shè)結(jié)構(gòu)在靜變形平衡位置附近作小幅振動(dòng)，采用線性結(jié)構(gòu)振動(dòng)方程，一些流固耦合降階方法也是基于線性考慮[7]。然而近年來，國外很多學(xué)者提出了HISSS/NASTRAN、CFL3D/GFEC、ZAERO/NASTRAN等CFD/CSD(Computational Fluid Dynamics/Computational Structural Dynamics) 耦合方法，開始將其應(yīng)用在常規(guī)亞/超飛行器氣彈分析中。CFD采用精確的流動(dòng)控制方程，CSD建立非線性結(jié)構(gòu)模型，基于CFD/CSD耦合分析具有高精度，可用于復(fù)雜問題的研究。國內(nèi)也有學(xué)者開展了相關(guān)研究，張華等[8]耦合FLUENT/NASTRAN，詳細(xì)研究了機(jī)翼結(jié)構(gòu)彈性對氣動(dòng)特性的影響，曾憲昂等[9]基于CFD/CSD耦合方法進(jìn)行了某機(jī)翼的顫振分析。

1 計(jì)算方法

CFD/CSD耦合方法遵守基本守恒原則，在流體與固體耦合交界面處，滿足流體與固體應(yīng)力、位移、溫度等變量的相等或守恒。通過分離解法分別求解流體和固體控制方程，再通過流固耦合交界面進(jìn)行計(jì)算結(jié)果的傳遞。只要流固耦合面完全對應(yīng)，可以保證交界面上的參數(shù)從局部到全局精確傳遞。由于各自的物理屬性，在進(jìn)行數(shù)據(jù)傳遞時(shí)，并不是所有的變量都要進(jìn)行傳遞。流體將力傳遞給固體，固體將位移傳遞給流體。使用動(dòng)網(wǎng)格模擬流場形狀由于邊界運(yùn)動(dòng)而隨時(shí)間的改變，從而捕捉平尾的結(jié)構(gòu)變形。

本文采用ANSYS軟件，CFD計(jì)算基于有限體積法求解非定常雷諾平均N-S(RANS)方程。采用標(biāo)準(zhǔn)k-ε湍流模型結(jié)合壁面函數(shù)法，選擇密度基隱式求解，對于空間離散方法，擴(kuò)散項(xiàng)的無黏通量項(xiàng)選取AUSM通量差分分裂格式，而黏性通量項(xiàng)選取中心差分格式，對流項(xiàng)則使用2階迎風(fēng)格式。時(shí)間積分使用高斯賽德爾(gauss-seidel)方法。CSD計(jì)算基于有限元法，對固體結(jié)構(gòu)進(jìn)行離散，將連續(xù)結(jié)構(gòu)的無限維問題轉(zhuǎn)化為與節(jié)點(diǎn)位移相關(guān)聯(lián)的有限維問題。固體運(yùn)動(dòng)方程可由拉格朗日方程得到：

(1)

其中：M為質(zhì)量矩陣；C為阻尼矩陣；K為剛度矩陣；F為廣義力；x為節(jié)點(diǎn)位移向量。

2 建模

2.1 氣動(dòng)模型

本文計(jì)算的全動(dòng)平尾和坐標(biāo)系定義如圖1所示，z軸垂直于x軸和y軸，按右手定則給定正方向。翼根弦長6.309m，平尾寬度2.22m，翼型相對厚度0.04m，翼面積15.52m2。

圖1 平尾平面形狀

計(jì)算域邊界采用pressure far-field條件，高度為30 km，壓力、溫度等氣體參數(shù)由標(biāo)準(zhǔn)大氣表查得。對平尾周圍流場進(jìn)行網(wǎng)格劃分，如圖2所示。對稱面上采用對稱邊界條件，且在動(dòng)網(wǎng)格設(shè)置中，對稱面設(shè)置為可變形面，將變形約束在面內(nèi)。平尾上下表面、前后緣以及翼端面為設(shè)置流固耦合面，設(shè)置計(jì)算時(shí)間步長為0.001 s，最大迎角取40°。

圖2 平尾附近流場網(wǎng)格劃分

2.2 結(jié)構(gòu)模型

采用直軸式全動(dòng)平尾，即平尾與轉(zhuǎn)軸組成整體相對垂直于機(jī)身的轉(zhuǎn)軸旋轉(zhuǎn)。將轉(zhuǎn)軸位置定于翼根65%處，換算為平均氣動(dòng)弦長的43%。蒙皮厚度取2mm，翼軸長度取1m。翼軸選擇TC4鈦合金材料，翼面蒙皮為鎳合金GH4169，平尾內(nèi)部為全高度鎳合金蜂窩結(jié)構(gòu)。對平尾結(jié)構(gòu)進(jìn)行網(wǎng)格劃分，如圖3所示，網(wǎng)格數(shù)為28 000。其中翼根處翼軸末端面固定，翼根面變形約束在面內(nèi)，平尾上下表面、前后緣以及翼端面設(shè)置為流固耦合面。同樣設(shè)置時(shí)間步長為0.001 s。

圖3 平尾結(jié)構(gòu)網(wǎng)格

對平尾進(jìn)行模態(tài)分析，得到前4階模態(tài)如圖4所示，其中黑色線框表示未變形的平尾，云圖表示位移大小。各階模態(tài)、頻率見表1。

圖4 平尾模態(tài)分析

表1 模態(tài)頻率

3 氣彈特性分析

3.1 氣彈響應(yīng)分析

分別計(jì)算平尾迎角為0°～40°間隔20°的情況，依次分析平尾的氣動(dòng)參數(shù)以及結(jié)構(gòu)參數(shù)響應(yīng)，其中氣動(dòng)參數(shù)包括升阻力系數(shù)，結(jié)構(gòu)參數(shù)包括應(yīng)力和平尾前緣點(diǎn)位移，初始值為不考慮結(jié)構(gòu)彈性時(shí)的值。0°迎角時(shí)的氣彈響應(yīng)如圖5所示。

圖5 0°迎角時(shí)的氣彈響應(yīng)

可以看出升阻力系數(shù)均很快收斂，并無明顯波動(dòng)。由于采用對稱翼型，0°迎角下升力系數(shù)為0，阻力系數(shù)也僅在10-4量級。考慮結(jié)構(gòu)彈性之后升力系數(shù)幾乎不變，阻力系數(shù)增加30%，但幅值仍較小。結(jié)構(gòu)響應(yīng)結(jié)果同樣可以看出，應(yīng)力和位移兩者均很快收斂穩(wěn)定。其中最大應(yīng)力值為0.42MPa，最大位移值為0.03mm，從量級上看，幾乎可以忽略。

圖6所示分別為20°和40°迎角下的氣動(dòng)響應(yīng)結(jié)果?？梢钥闯?，升阻力系數(shù)均出現(xiàn)波動(dòng)，無相位差且呈收斂趨勢。20°迎角下升力系數(shù)初始值為0.282，收斂值約為0.272，減小約3.6%。阻力系數(shù)初始值為0.106，收斂值約為0.101，減小約4.7%。升阻力波動(dòng)半衰時(shí)長約0.726 s。40°迎角下升力系數(shù)初始值為0.666，收斂值約為0.632，減小約5.1%。阻力系數(shù)初始值為0.563，收斂值約為0.526，減小約6.6%。升阻力波動(dòng)半衰時(shí)長約0.478 s。

圖6 大迎角下的氣動(dòng)響應(yīng)

圖7所示分別為20°和40°迎角下的結(jié)構(gòu)響應(yīng)結(jié)果。可以看出20°迎角下平尾前緣點(diǎn)位移最大值約120mm，且隨時(shí)間逐漸衰減，收斂值約為60mm。在初始?xì)鈩?dòng)載荷的作用下，結(jié)構(gòu)的應(yīng)力響應(yīng)曲線迅速達(dá)到最大值，約115MPa。此后隨著時(shí)間的增長，響應(yīng)逐漸衰減，最后收斂到約63MPa。40°迎角下位移最大值約330mm，且隨時(shí)間逐漸衰減，收斂值約為97mm。在初始?xì)鈩?dòng)載荷的作用下，結(jié)構(gòu)的應(yīng)力響應(yīng)曲線迅速達(dá)到最大值約320MPa，此后隨著時(shí)間的增長，響應(yīng)逐漸衰減，最后收斂到約190MPa。

可見除0°迎角迅速收斂至平衡位置以外，在大迎角下平尾的升阻力系數(shù)以及前緣點(diǎn)位移、翼軸最大應(yīng)力點(diǎn)等效應(yīng)力曲線均出現(xiàn)波動(dòng)，隨時(shí)間變化逐漸衰減至平衡位置。迎角越大，初始振幅越大，升阻力減小的比例越大，但衰減得越快。

圖7 大迎角下的結(jié)構(gòu)響應(yīng)

3.2 流動(dòng)特性分析

圖8所示為0°迎角下考慮結(jié)構(gòu)彈性后翼根附近對稱面以及平尾上、下表面的壓力分布云圖。可在0°迎角下上下表面的壓力基本呈對稱分布，因此升力系數(shù)幾乎為0。

圖8 壓力分布云圖(0°迎角)

圖9所示為20°迎角時(shí)，初始流場與平衡位置所對應(yīng)的平尾下表面壓力分布云圖以及平尾變形示意圖。可以看出，結(jié)構(gòu)變形導(dǎo)致下表面的壓力分布發(fā)生變化。相比初始位置，平衡位置對應(yīng)的平尾下表面后緣膨脹區(qū)低壓區(qū)域擴(kuò)大，平尾前緣高壓區(qū)靠近前緣，使得翼尖處的高壓區(qū)域明顯縮小。

由平尾變形示意圖可以看出，平衡時(shí)平尾發(fā)生扭轉(zhuǎn)。從正視圖可以看出50%展長以內(nèi)的部分前緣向上偏轉(zhuǎn)，相對迎角增大，從而使得激波偏折角增大，波后壓力增大，因此造成高壓區(qū)前移。而50%展長以外的部分前緣向下偏轉(zhuǎn)，相對迎角減小，相應(yīng)的激波后壓力減小。從后視圖可以看出平尾后緣均向上偏轉(zhuǎn)，相對迎角增加，從而增大偏折角，但此時(shí)前方激波經(jīng)過最大厚度處產(chǎn)生膨脹波，偏折角增大反而造成壓力減小。綜上，翼根前緣處壓力的增大無法彌補(bǔ)翼尖前緣以及整個(gè)平尾后緣壓力的降低，因此整體壓力減小，升力系數(shù)降低。

圖9 結(jié)構(gòu)變形后的流場(20°迎角)

圖10給出了一個(gè)典型的振蕩周期內(nèi)平尾結(jié)構(gòu)變形以及表面壓力分布云圖的變化過程。取0.03s～0.08s，間隔0.01s的6個(gè)時(shí)刻。可以明顯看出平尾前緣扭轉(zhuǎn)和翼尖彎曲變形。變形最大時(shí)升力系數(shù)最小，對應(yīng)地，變形最小時(shí)升力系數(shù)最大。

圖10 平尾結(jié)構(gòu)變形和表面壓力分布云圖

3.3 結(jié)構(gòu)特性分析

圖11所示為0°迎角時(shí)平尾的等效應(yīng)力分布，可以看出，翼軸處應(yīng)力分布較為集中，且應(yīng)力最大出現(xiàn)在翼軸的x向與平尾相連接處。根據(jù)前文壓力分布特點(diǎn)可以得出0°迎角下平尾結(jié)構(gòu)受到的力主要是阻力，集中在x軸方向，而在z向由于上、下表面對稱升力為0，因此應(yīng)力較小。

圖11 結(jié)構(gòu)應(yīng)力(0°迎角)

圖12為20°和40°迎角下變形最大時(shí)對應(yīng)的等效應(yīng)力分布云圖，可以看出，靠近翼軸處應(yīng)力較為集中，而翼軸上的應(yīng)力較小。迎角越大，最大應(yīng)力越大。平尾內(nèi)部最大應(yīng)力在40°時(shí)達(dá)到1.9 GPa，已經(jīng)達(dá)到所用鎳合金材料的

屈服強(qiáng)度極限，理論上說材料已經(jīng)破壞，由于采用線彈性材料，計(jì)算結(jié)果可能較實(shí)際結(jié)果偏大，但為了保險(xiǎn)起見，實(shí)際結(jié)構(gòu)設(shè)計(jì)時(shí)應(yīng)在翼軸與內(nèi)部蜂窩接觸部位附近加固，防止蜂窩結(jié)構(gòu)破壞。

圖12 大迎角下的等效應(yīng)力分布云圖

4 結(jié)語

本文采用CFD/CSD流固耦合法對一種典型全動(dòng)平尾進(jìn)行了氣動(dòng)彈性數(shù)值模擬。給出了不同迎角下平尾的氣動(dòng)彈性響應(yīng)，并分別針對流場和結(jié)構(gòu)響應(yīng)做出了分析，結(jié)果表明：

1) 除0°迎角迅速收斂至平衡位置以外，其余各迎角時(shí)的氣動(dòng)力和結(jié)構(gòu)響應(yīng)曲線均出現(xiàn)波動(dòng)，隨時(shí)間變化逐漸衰減至平衡位置。迎角越大，初始振幅越大，氣動(dòng)力減小的比例越大，結(jié)構(gòu)應(yīng)力越大，但隨時(shí)間衰減得越快。

2) 結(jié)構(gòu)變形導(dǎo)致下表面壓力分布發(fā)生變化。靠近翼根部分前緣上偏，壓力增大?？拷砑獠糠智熬壪缕?，壓力減小。后緣均上偏，壓力減小。整體壓力減小，升力系數(shù)降低。迎角越大，現(xiàn)象越明顯。

3) 靠近翼軸處應(yīng)力較為集中，而翼軸上的應(yīng)力較小。平尾存在彎曲/扭轉(zhuǎn)耦合現(xiàn)象，隨時(shí)間逐漸收斂至平衡位

置，但相對初始位置的變形隨迎角增大而增大。