李玉蓮 袁勝弢

(哈爾濱飛機工業(yè)集團有限責(zé)任公司飛機設(shè)計研究所，哈爾濱150066)

應(yīng)力分析是直升機結(jié)構(gòu)強度設(shè)計最重要的工作內(nèi)容之一，隨著商業(yè)分析軟件的快速發(fā)展，有限元已經(jīng)成為復(fù)雜結(jié)構(gòu)應(yīng)力分析中最重要的工具[1]。結(jié)構(gòu)強度有限元分析過程包括有限元前處理、求解計算和后處理。對于直升機艙門結(jié)構(gòu)強度有限元分析前處理工作，氣動載荷是分析設(shè)計工作中非常重要的載荷設(shè)計輸入[2]。在國內(nèi)航空界大部分科研院所，氣動載荷分布由載荷專業(yè)計算，結(jié)構(gòu)強度分析專業(yè)在有限元的前處理工作中將該載荷模擬并施加在有限元模型節(jié)點(單元)上，以進行有限元的求解計算工作。在工程設(shè)計中，將不均勻分布氣動載荷施加到有限元節(jié)點上是一項非常繁瑣但又很重要的工作。截至目前，航空航天行業(yè)廣泛應(yīng)用的有限元分析軟件如Patran/Nastran，Abaqus等，盡管提供了強大的有限元分析計算功能，但對于不均勻分布面載荷還不能自動實現(xiàn)加載。人工手動對有限元節(jié)點或者單元進行逐個加載，滿足不了實際工程設(shè)計的需求，因為直升機結(jié)構(gòu)強度設(shè)計中的載荷工況較多，而且有限元細節(jié)模型的節(jié)點數(shù)目非常龐大，逐個單元或者節(jié)點的加載方式嚴(yán)重影響載荷施加的效率以及質(zhì)量。如何提高有限元分析前處理工作中加載的效率與質(zhì)量成為國內(nèi)外專家學(xué)者重視的問題。

為了解決有限元分析前處理中施加不均勻分布氣動載荷這一難題，國內(nèi)外出現(xiàn)了很多算法。目前航空領(lǐng)域針對翼面結(jié)構(gòu)的氣動載荷分布，常用的算法是“多點排”方法[3]，該方法以靜力等效以及應(yīng)變能最小為約束條件，將每一種氣動載荷分配到結(jié)構(gòu)有限元分析模型的節(jié)點上。林小夏等[4]提出基于特征函數(shù)分布對不均勻載荷進行加載。王專利[5]、徐建新等[6]使用了“多點排”的計算方法。高尚君等[7]提出基于彈簧懸臂梁模型最小變形能的氣動載荷分配方法。張建剛等[8]針對翼面結(jié)構(gòu)提出薄板樣條插值函數(shù)來計算節(jié)點壓力值，這些算法能解決實際中的問題，但需要在結(jié)構(gòu)有限元模型與氣動模型之間建立一種載荷數(shù)據(jù)傳遞關(guān)系，這種傳遞關(guān)系通過樣條矩陣來實現(xiàn)。這些方法都是基于數(shù)值插值分析，即每一次計算均要調(diào)用原始數(shù)據(jù)，增加了算法的時間以及復(fù)雜度，而針對直升機結(jié)構(gòu)氣動載荷如何實現(xiàn)快速加載、等效分配的方法比較少。

本文以直升機艙門強度設(shè)計的工程實際需求為出發(fā)點，根據(jù)已知氣動模型節(jié)點的壓力值，通過數(shù)學(xué)分析最小二乘法，擬合出壓力分布曲面，由該曲面擬合函數(shù)計算得到結(jié)構(gòu)有限元模型每一節(jié)點的壓力值，并對數(shù)值回歸性進行分析。通過該擬合壓力分布曲面函數(shù)的方法，可以準(zhǔn)確、快速地將直升機艙門的不均勻分布氣動載荷分配施加到有限元節(jié)點上，以供直升機結(jié)構(gòu)強度工程師參考。

1 直升機艙門氣動載荷

現(xiàn)代飛行器結(jié)構(gòu)建模和氣動模型普遍采用有限單元法[9]，兩種模型劃分是為了各自的分析任務(wù)而獨立進行的。由于兩種模型的分析目的不同，建模方式也存在很大的不同，以某型號直升機艙門為例，艙門的氣動模型與結(jié)構(gòu)分析有限元模型節(jié)點存在很大的差異。艙門氣動模型的節(jié)點數(shù)為1892個，氣動壓力分布計算結(jié)果為1892組離散數(shù)據(jù)；而有限元模型為了保留結(jié)構(gòu)邊界，區(qū)分玻璃區(qū)域與結(jié)構(gòu)區(qū)域，模型劃分方法以及節(jié)點數(shù)目與氣動模型均存在不同，兩種模型的有限元節(jié)點差異對比如圖1所示。

圖1 氣動節(jié)點與結(jié)構(gòu)有限元節(jié)點差異示意圖

在艙門的結(jié)構(gòu)強度分析有限元模型中，要實現(xiàn)施加不均勻分布的氣動載荷，需要在結(jié)構(gòu)有限元模型與氣動模型之間建立一種載荷數(shù)據(jù)傳遞關(guān)系，工程上一般采用數(shù)值插值法或數(shù)據(jù)擬合法。數(shù)據(jù)擬合法根據(jù)氣動載荷分布數(shù)據(jù)擬合出相對簡單的數(shù)學(xué)模型，不需要每一個有限元節(jié)點的計算均調(diào)用原始氣動分布數(shù)據(jù)，降低了算法時間和復(fù)雜度，擬合函數(shù)更逼近真實函數(shù)，本文中這種傳遞關(guān)系通過數(shù)學(xué)分析，采用最小二乘法擬合來實現(xiàn)[10]。

2 最小二乘法基本理論

最小二乘法是一種基于已知離散點而擬合出近似函數(shù)的方法，使用離散點數(shù)據(jù)擬合得到曲面。最小二乘法具有擬合精度高、通用性強的特點。已知離散點函數(shù)(稱為真實氣動載荷)是單值連續(xù)的，則對空間區(qū)域D上的任何一個向量坐標(biāo)點來說，理論上必有且只有一個確定的函數(shù)與之對應(yīng)。由于在計算過程中真實函數(shù)是不知道的，任一向量所對應(yīng)的函數(shù)值無法直接求出。為此，必須利用已知條件(有限元的節(jié)點以及氣動載荷信息)構(gòu)成一個逼近函數(shù)，用它來代替真實函數(shù)，借以計算任一向量所對應(yīng)的函數(shù)的近似值。

假設(shè)空間區(qū)域D內(nèi)，已知數(shù)值的n個數(shù)據(jù)點集p i，其中i=1,2,···,n，其坐標(biāo)可以表示為(x i,y i,z i)，在n個節(jié)點z i=(x i,y i)有已知值

已知數(shù)值用矢量表示為

設(shè)擬合函數(shù)φ(z)=[b1(z)b2(z)···b m(z)]。已知數(shù)值用矢量表示為

其中b m(z)為m個線形無關(guān)的基函數(shù)，λ(z)為m維系數(shù)矢量，是基函數(shù)的待定系數(shù)。最小二乘法即按照偏差平方和最小原則選取擬合曲線，在計算區(qū)域內(nèi)，由式(3),利用加權(quán)最小二乘法構(gòu)成二次形式

其中

基函數(shù)b m(z)是已知的，為求得待定系數(shù)，式(4)中對λ(z)一次求導(dǎo)等于0，求得

從而求得擬合近似函數(shù)為

3 艙門氣動分布擬合

直升機艙門的氣動載荷為不均勻的分布載荷，載荷專業(yè)提供的艙門壓力為氣動模型離散點的分布，在某型號直升機設(shè)計中，艙門某種狀態(tài)的氣動壓力分布等高線分布圖如圖2所示。從圖中可以看出，對于直升機艙門結(jié)構(gòu)，沿著艙門縱向坐標(biāo)X以及高度方向坐標(biāo)Z的不同，氣動壓力分布呈現(xiàn)為非均布載荷。

圖2 某型機艙門氣動分布

艙門氣動載荷擬合函數(shù)具有最小二乘法的性質(zhì)，得到直升機在空中各種不同飛行狀態(tài)時壓力分布的逼近曲面函數(shù)。在計算過程中，待定系數(shù)的快速求解采用數(shù)學(xué)分析軟件Mathematica中的程序進行[11]。

分析計算中，有限元節(jié)點的坐標(biāo)值參考坐標(biāo)系為機身全局坐標(biāo)系。艙門機身軸向向后為X坐標(biāo)正向，機身高度向上為Z坐標(biāo)正向。最小二乘法擬合函數(shù)可以是一次、二次函數(shù)等所有初等函數(shù)，在艙門的壓力分布中，分別對平面分布和曲面分布進行多次擬合以及回歸性分析，以計算校核逼近函數(shù)的精度。其中線性基擬合時，式(3)中m取3，基函數(shù)表達式為式(9)；二次基擬合時，式(3)中m取6，基函數(shù)表達式為式(10)；立方基擬合時，式(3)中m取10，基函數(shù)表達式為式(11)。

當(dāng)m=3時，艙門的壓力分布載荷呈平面分布，基函數(shù)取二次以上時，艙門的壓力分布為曲面分布。根據(jù)已知氣動點的壓力值求出逼近函數(shù)，進而求得有限元節(jié)點壓力值。有限元軟件Patran中有Field(域)加載功能，通過最小二乘法擬合，可快速實現(xiàn)不均勻分布載荷的施加處理過程。

4 艙門氣動載荷擬合結(jié)果分析

4.1 數(shù)學(xué)分析與擬合

某型號直升機艙門的非均布氣動分布載荷，當(dāng)采用不同的基函數(shù)進行最小二乘法擬合時，確定系數(shù)如表1所示。圖3～圖6是不同基函數(shù)所計算得出的擬合函數(shù)與艙門原氣動點壓力值對比分析，從圖4與圖3的對比結(jié)果可以看出，隨著擬合函數(shù)系數(shù)的增高，擬合后的壓力曲面和擬合前的節(jié)點壓力分布趨勢一致，吻合程度逐漸提高，結(jié)合表1數(shù)據(jù)分析，從線性基到二次基的確定系數(shù)提高量最大，二次基后隨著擬合函數(shù)次數(shù)提高，確定系數(shù)的數(shù)值區(qū)域增加得較為平緩。

表1 擬合誤差分析表

圖3 艙門氣動載荷數(shù)據(jù)與一次擬合曲面

圖4 艙門氣動載荷數(shù)據(jù)與二次擬合曲面

圖5 艙門氣動載荷數(shù)據(jù)與三次擬合曲面

圖6 艙門氣動載荷數(shù)據(jù)與四次擬合曲面

根據(jù)以上分析，二次基的擬合結(jié)果較為優(yōu)異，當(dāng)采用二次基函數(shù)擬合時，對氣動節(jié)點進行擬合精度分析，計算結(jié)果見表2，擬合結(jié)果與原始數(shù)據(jù)的單點誤差最大值小于1%，最小二乘法的曲面擬合結(jié)果具有較高精度，能夠滿足工程上的設(shè)計要求。

表2 擬合結(jié)果與氣動節(jié)點壓強對比分析表

4.2 擬合結(jié)果分析

從圖6中可以看出，插值擬合的壓力曲面和插值前的氣動點壓力分布吻合程度很好，通過這種數(shù)據(jù)傳遞方式，可以快速實現(xiàn)由氣動模型的散點壓力分布求出有限元節(jié)點的壓力分布，進而通過逼近函數(shù)實現(xiàn)對有限元模型的快速加載。

5 結(jié)論

本文根據(jù)工程實際需求，基于氣動離散點壓力分布，通過最小二乘法擬合函數(shù)，可以實現(xiàn)直升機艙門壓力載荷擬合加載，提高有限元分析中加載的效率，擬合加載的精度較高。

Mathematica程序中的數(shù)學(xué)理論清晰，編程擬合插值簡單，可以有效解決工程中的實際問題。根據(jù)不同的工況，擬合逼近函數(shù)隨著氣動離散點壓力分布的改變而改變，尤其針對直升機的飛行狀態(tài)多變、壓力分布多變的情況，插值擬合函數(shù)具有較好的適應(yīng)性，對多工況可以實現(xiàn)快速方便地加載，對于主要承受不均勻分布的氣動載荷的結(jié)構(gòu)，提高了強度有限元分析的效率。