張建剛 孫仁俊 唐長紅

(中國航空工業(yè)集團公司第一飛機設計研究院，西安710089)

大型飛機氣動載荷向有限元節(jié)點等效分配的方法

張建剛1)孫仁俊 唐長紅

(中國航空工業(yè)集團公司第一飛機設計研究院，西安710089)

將氣動載荷分配到有限元節(jié)點上是工程實際中的一項重要而繁瑣的工作.對于二維的翼面氣動載荷，根據原始的氣動壓力點的壓力值，采用樣條曲面擬合的方法，擬合得到翼面壓力分布曲面，由該曲面得到有限元節(jié)點上的壓力值，再在有限元模型單元上積分得到有限元節(jié)點載荷供強度設計使用.大型飛機具有復雜的增升裝置，增升裝置的氣動載荷可能是三維的，對于三維的翼面載荷，直接在氣動網格上積分得到氣動載荷的小塊集中力，然后按照沿某一方向投影的方法，找到該集中力作用的單元，最后按照二次規(guī)劃的方法，將其分配到有限元節(jié)點上.

氣動載荷,有限元節(jié)點,薄板樣條

在飛機的結構強度設計中，氣動載荷是非常重要的設計輸入.在國內的航空院所里，氣動載荷由載荷專業(yè)計算，結構強度專業(yè)將該載荷施加在有限元模型上作為輸入進行下一步的設計.為了確保載荷輸入的準確性，一方面載荷專業(yè)要努力使其計算結果更準確，另一方面，強度專業(yè)也要確保將氣動載荷正確地施加在計算模型上.在工程中，將氣動載荷等效分配到有限元節(jié)點上，是一項非常繁瑣但又很重要的工作.而目前航空航天行業(yè)廣泛使用的強度計算軟件如MSC.Patran/Nastran,Abaqus等，對于不均勻面載荷加載還不能自動實現.如果人工手動地對節(jié)點進行加載，滿足不了工程實際的要求，因為實際的有限元模型具有大量的節(jié)點和工況.為了解決這一難題，出現了很多的算法.

國外相關研究常見于氣動彈性專業(yè)計算，氣彈計算每一迭代步需要將氣動載荷施加在結構節(jié)點上：首先在氣動網格上積分得到氣動載荷集中力，然后根據虛功原理，用樣條插值矩陣將空氣動力網格上的氣動力變換為結構網格點上的等效值[1-3].目前國內航空領域常用的算法是“多點排”方法[4]，該方法以應變能最小及靜力等效為約束條件，將每一個氣動載荷分配到若干個有限元節(jié)點上.王專利[5]、徐建新等[6]通過計算證明了該方法的正確性.高尚君等[7]提出了基于彈簧-懸臂梁模型最小變形能的氣動載荷分配方法，解決了氣動點與部分有限元點重合而導致的難以求解的問題.這些算法能有效解決工程實際中的一些問題，但也還存在著需要改進的地方.首先，載荷分配的起點是積分后氣動載荷的集中力.現代大型飛機翼面結構復雜，氣動分區(qū)和結構分區(qū)往往不一致，直接在氣動網格上積分成集中載荷有可能會跨區(qū)域傳遞載荷，導致傳力路線不真實；其次,這些方法分配載荷時人為地指定了一些節(jié)點去分配相應的氣動載荷，也可能導致傳力路線的不真實，實際上應該充分利用有限元模型的單元和節(jié)點信息，遵循就近分配的原則以保證載荷傳遞的正確；最后，分配的載荷都假設為二維載荷，這對翼面盒段等扁平結構是適用的，但對現代大型飛機的增升裝置及整流罩、雷達罩等三維氣動面則顯得較為牽強.

現代的大型飛機，比如大型運輸機，為了良好的氣動及起降性能，大多都布置有前緣縫翼、后緣多縫襟翼等復雜的增升裝置.這些裝置導致了大型飛機的氣動載荷分布極為復雜，并且這些翼面本身的結構也較為復雜，傳力路徑復雜.強度設計專業(yè)必須保證將這些氣動載荷以最小的代價，真實、準確地施加到結構有限元模型上，這是強度設計的首要條件，也是飛機精細化設計的必然要求.

本文以工程實際需要為出發(fā)點，研究了如何準確、快速地將氣動載荷分配至有限元節(jié)點上，以供同行借鑒參考.

1 氣動載荷的提供形式及翼面結構的分類

載荷專業(yè)給出的氣動載荷一般形式是：根據試驗或理論計算給出一系列離散點上的壓力分布及翼面的總載、總矩.這些壓力點和有限元模型節(jié)點是不重合的.

根據翼面的結構及其受力特點，可以將其大致分為兩類：一類是翼面盒段等比較扁平的結構，這類結構主要承受垂直于弦平面的氣動載荷，另外兩個方向上的氣動載荷相比而言都是小量，并且在這兩個方向上翼面的剛度又很大，因此可以忽略掉這兩個方向的氣動載荷而不會引入較大的誤差，因此可認為這一類翼面載荷的分布是沿展向和弦向的二維分布.另一類是前緣縫翼、后緣襟翼等曲率較大的翼面，這類翼面有兩個以上方向的載荷具有相同的數量級，因此不能忽略任何一個分量，必須按照空間分布的三維載荷來對待.并且這類翼面上的氣動載荷變化梯度大，給出的氣動載荷的網格較二維的更為細密.

本文使用的載荷分配的基本思想是：對于二維載荷，根據已知氣動點上的壓力值，通過薄板樣條插值函數，擬合出壓力分布曲面，由該曲面得到每一個有限元節(jié)點上的壓力值，再在單元上對壓力進行積分，得到有限元節(jié)點上的載荷.插值方法采用航空航天工程中廣泛使用的薄板樣條曲面插值.對于三維載荷，直接在氣動網格上進行積分得到集中力，然后遵照就近分配的原則分配到有限元模型上.

2 薄板樣條插值函數

在分配載荷之前需要由氣動點的壓力值求出有限元節(jié)點上的壓力值，本文采用薄板樣條插值函數[8].設在N維歐氏空間內域D上有n個向量(每個向量代表一個節(jié)點坐標)及其對應的函數值，即已知

設向量X的函數W=W(X)(稱為真實函數)是單值連續(xù)的，則對D上的任一向量Xk(稱為插值點)來說，必有一個且只有一個確定的函數Wk與之對應.由于在插值問題中真實函數是不知道的，任一向量Xk所對應的函數值無法直接求出.為此，必須利用已知條件(節(jié)點信息)構成一個逼近函數，用它來代替真實函數，借以計算任一向量Xk所對應的函數之近似值.在這里，逼近函數的形式如下

式(2)中的待定系數由無窮遠處的邊界條件及n個已知節(jié)點函數值來確定

將上述方程組寫成矩陣形式如下

記為

式中S為待定系數組成的向量；f為已知函數值和0組成的向量；bji=r2jiln(r2ji+ε)；m=N+1+n.

式(6)是一個線性代數方程組，對其求解即可得到式(2)中的各待定系數，進而得到逼近函數.對于任一待求節(jié)點，將節(jié)點坐標代入式 (2)即可求得對應的函數值.

3 二維載荷的分配方法

二維載荷分配方法適用于機翼盒段等扁平翼面，載荷方向垂直于機翼的弦平面.氣動專業(yè)提供的載荷壓力分布是一系列離散點上的壓力，這里的壓力指的是當地靜壓減去來流靜壓后的剩余量，也就是壓力系數和速壓的乘積.在式(6)中，令N=2，根據氣動點的坐標構造出左邊的A矩陣，根據已知點的壓力值構造向量f，求解方程得到向量S，代入式(2)后得到有限元點上的壓力值.圖1給出了某機翼的已知壓力值和插值得到的壓力值的比較.圖中的x和y坐標是沿著機翼展向和弦向的坐標，z坐標是壓力.

由圖1可以看出，插值后的壓力曲面和插值前的壓力曲面吻合程度很好，可以由有限元節(jié)點上的壓力積分得到節(jié)點的集中力.分配氣動載荷的單元一般是模型表面的單元，大多數是四邊形單元，還有少量的三角形單元.對于三角形單元，3個節(jié)點坐標

圖1 氣動點上的壓力值和插值后有限元點上的壓力值

為(xi,yi),i=1,2,3,假設壓力分布在單元上呈平面分布，即

根據已知節(jié)點的壓力值即可求出式(8)中的常數A,B和C.

假設三角形三個節(jié)點上分到的載荷是fi,i= 1,2,3，則有

式(9)中的積分區(qū)域S是三角形的表面，經過推導可以得到

式中，p0=p1+p2+p3,p1,p2和p3是三角形三個節(jié)點上的壓力值；S為三角形面積.

由以上方程組即可求解得到三角形節(jié)點的集中載荷.四邊形單元可以化為兩個三角形單元，再做如上的處理.

4 三維載荷的分配方法

三維載荷的分配方法理論上可以和二維的處理方法類似，即根據已知點的壓力值采用薄板樣條函數插值出有限元節(jié)點處的壓力值，然后在單元上積分即可得到集中載荷.但是在工程實際中，承受三維載荷的翼面往往壓力梯度較大，氣動網格密集，相比二維載荷,氣動點的數目大為增加.一塊翼面往往布置有上萬個氣動點，這時如果仍然進行插值運算，會導致計算量急劇增大.在個人計算機上求解一個工況的時間以小時計，實際的工況數又成百上千，這種做法在工程實際中不現實.如果在氣動點上挑選一部分點插值，無形中又會降低載荷的精度.

為了解決這一問題，權衡利弊后采用如下做法：在氣動網格上直接積分得到氣動單元上的小塊集中力，再將這些力等效到有限元節(jié)點上.有了小塊集中力之后，向有限元節(jié)點等效分配載荷之前需要先準確找到該力作用的單元.一個簡便可行的方是：沿著某一個方向(一般可選氣動單元平面法向量絕對值最大的分量方向)，將集中載荷的作用點向翼面結構投影，找到分配載荷的單元.假設沿著z向投影，某一個四邊形單元的4個節(jié)點坐標為(xi,yi),i=1,2,3,4，則集中力作用點在四邊形投影坐標為(x,y)，如圖2所示.

圖2 集中力作用點向四邊形投影示意圖

假設四邊形面積是S，載荷作用點的投影點和四邊形4個節(jié)點連成4個三角形，若則認為載荷點投影到了四邊形單元內部，則應該將該集中載荷分到這個四邊形的節(jié)點上.三角形投影方法類似.

將集中載荷分配到三角形單元可以按照靜力等效原則直接分配.對于四邊形單元，因為約束方程的個數少于自變量的個數，存在多種解.根據應變能最小原理，應變能近似正比于4個力的平方和，可以按照優(yōu)化理論中的二次規(guī)劃方法進行分配.比如需要分配z方向的氣動載荷Fz，載荷作用點坐標是(X,Y)，假設四邊形4個節(jié)點分配到的載荷和節(jié)點坐標分別為fzi和(xi,yi)，i=1,2,3,4，則問題等效為以下二次規(guī)劃問題[9]

其中，H是4階單位矩陣

式 (10)中的目標函數的意義是 4個力的平方和最小.首先定義Lagrange函數

令?xL(x,λ)=0及?λL(x,λ)=0.

由此得到方程組

求解式(11)即可得到各個節(jié)點上分配到的z方向的載荷.x和y方向的載荷分配方法完全相同.圖3～圖5是按照以上方法分配的幾個飛機部件的載荷.

圖3 某大型運輸機前緣縫翼載荷分配結果

圖4 某大型運輸機后緣襟翼載荷分配結果

圖5 某飛機噴管載荷分配結果

5 結論

(1)根據以往的工程應用的實際情況，二維載荷按照以上做法求得的整個翼面的總載荷可能和氣動專業(yè)測力測矩得到的結果略有差異，對于大的載荷工況(往往是翼面本身的設計情況)誤差一般在1%以內；對于小載荷工況(這些工況可能是別的部件的設計情況，為了全機平衡，該翼面上的載荷可能是配平載荷，非翼面本身的設計情況)，誤差略大一些，最大不超過5%.為了保證全機求解時的平衡，可將上述載荷按線性化略加修正即可.

(2)對于三維載荷，按照該方法處理后，氣動載荷基本上垂直于當地翼面，如果有限元網格和氣動網格的疏密程度及外形、邊界完全一致，則處理后的氣動載荷完全垂直于翼面，這也證明了方法的正確性.

(3)該方法數學理論清晰，編程簡單.通過以往工程型號的驗證，證明可以有效地解決工程實際中的問題.

1 Prananta BB,Veul RPG,Boelens OJ.Manoeuvre-loads computations using CFD-based static aeroelastic simulations for multiple store conf i gurations.27th International Congress of the Aeronautical Sciences,Nice,France,2010

2 Johannes KSD,Mostafa MA,Yasser MM,et al.Static aeroelastic stif f ness optimization of a forward swept composite wing with CFD correct aero loads.International Forum on Aeroelasticity and Structural Dynamic,Saint Petersburg Russia,2015

3 Rodden WP,Johnson EH.MacNeal-Schwendler Corporation. MSC.Nastran Version 68,Aeroelastic Analysis’s Guide.MacNeal-Schwendler Corp,Santa Ana,CA,2004

4 解思適.飛機設計手冊第9冊載荷、強度和剛度.北京：航空工業(yè)出版社，2001

5 王專利.翼面結構有限元模型節(jié)點氣動載荷計算.洪都科技，2007,(1)：7-14

6 徐建新,王彪.方向舵靜強度試驗載荷等效計算研究.裝備制造技術，2013,(1)：14-16

7 高尚君，于哲峰，汪海.基于彈簧-懸臂梁模型最小變形能的氣動載荷分配方法.飛機設計，2013，(6)：1-4

8 趙永輝.氣動彈性力學與控制.北京：科學出版社，2007

9 蔣寶林.優(yōu)化理論與應用.北京：清華大學出版社,2006

(責任編輯：劉希國)

THE METHOD OF TRANSFERRING AERODYNAMIC LOAD TO FEM
NODES FOR LARGE AIRPLANE

ZHANG Jiangang1)SUN Renjun TANG Changhong
(The First Aircraft Institute of the Aviation Industry Corporation of China，Xi’an 710089，China)

Transferring the aerodynamic load to FEM(f i nite element method)nodes is a hard work but is very important.In order to transfer the two-dimensional aerodynamic load,a thin plate spline surface of the pressure is constructed according to the existing pressure on the aerodynamic nodes.Then the pressure on the FEM nodes can be computed.At last the concentrated force on the FEM Nodes can be obtained by integration of the pressure on the FEM elements.The large airplane often has complicated high-lift devices and their aerodynamic load may be three-dimensional.For the three-dimensional aerodynamic load,with the integration of the pressure on the aerodynamic webs,the concentrated force on the aerodynamic nodes can be obtained, and these forces can be transferred to the FEM nodes using the quadratic program method.

aerodynamic load,FEM(f i nite element method)nodes,thin plate spline

V211.41

A

10.6052/1000-0879-16-170

2016-05-13收到第1稿，2016-07-18收到修改稿.

1)張建剛，碩士，高級工程師，主要從事飛機強度及載荷設計.E-mail：155246438@qq.com

張建剛,孫仁俊,唐長紅.大型飛機氣動載荷向有限元節(jié)點等效分配的方法.力學與實踐,2017,39(1)：25-29

Zhang Jiangang,Sun Renjun,Tang Changhong.The method of transferring aerodynamic load to FEM nodes for large airplane.Mechanics in Engineering,2017,39(1)：25-29