基于Excel開發(fā)的飛機(jī)有限元后處理強(qiáng)度計算平臺

樊建超

【摘 要】本文介紹了一種基于Excel軟件采用VBA語言開發(fā)的強(qiáng)度計算平臺，用于處理飛機(jī)結(jié)構(gòu)有限元內(nèi)力解。利用該計算平臺可實現(xiàn)模型導(dǎo)入、數(shù)據(jù)篩選、任意工程方法批量計算、工況循環(huán)找出最小裕度等功能，實現(xiàn)飛機(jī)結(jié)構(gòu)在多工況下受載后的強(qiáng)度分析，幫助工程師快速、方便、準(zhǔn)確的從龐大的數(shù)據(jù)信息中提煉出有價值的結(jié)果。

【關(guān)鍵詞】分析平臺；有限元；數(shù)據(jù)處理；VBA

The Stress Analysis System Used in the Postprocessing of Aircraft FEA Developped by Excel

FAN Jian-Chao

（ Shanghai Aircraft Design And Research Institute， COMAC， Shanghai 200232，China）

【Abstract】The stress analysis system based on Excel developed by the VBA language is introduced. Which is used to deal with the aircraft FEA result. It include model imported ，data selected ，formula formed， subcase cycled to find the dangerous results， it can be used to analyse the serious stress of structures at different subcases， which can help engineers to abstract the useful informatiom quickly and exactly from a mass of data

【Key words】Analysis system； FEA； Data analysis； VBA

0 引言

在飛機(jī)結(jié)構(gòu)設(shè)計中利用Nastran對全機(jī)有限元模型進(jìn)行求解，并采用各類工程方法處理有限元求解結(jié)果是一項重要的工作，飛機(jī)結(jié)構(gòu)的強(qiáng)度分析基本都圍繞著處理有限元結(jié)果而展開。

為提高效率，數(shù)十年來，工程師根據(jù)不同的結(jié)構(gòu)形式、不同的工程方法編制了一些程序，用于計算典型結(jié)構(gòu)，如普通框、普通壁板。但實際飛機(jī)結(jié)構(gòu)形式復(fù)雜多樣，工程方法五花八門，無法窮盡所有的工程方法都編寫一套程序。遇到具體問題，有時只能憑經(jīng)驗挑選某一種嚴(yán)重的工況進(jìn)行分析。

如何開發(fā)一套統(tǒng)一，標(biāo)準(zhǔn)、高效、且計算能力強(qiáng)大的強(qiáng)度分析平臺一直是各科研院所不斷追求的夢想。

目前國內(nèi)尚無此計算平臺問世，國外，空中客車公司委托SAMTECH通過多年實踐，開發(fā)了一款稱之為ISAMI的強(qiáng)度計算平臺，雖然ISAMI具有可視化、標(biāo)準(zhǔn)化等優(yōu)勢，但它本質(zhì)上是各種常用工程方法有限的集合，遇到新的工程方法還需專業(yè)人員編寫新方法的程序納入平臺，普通工程師更無法進(jìn)行二次開發(fā)，對于飛機(jī)上各類特殊計算無能為力。此外，ISAMI計算過程猶如一個黑匣子，對工程師理解問題也較為不利。

考慮到Excel可以通過VBA實現(xiàn)編程，如果能用Excel編程代替ISAMI的部分功能，同時利用Excel本身強(qiáng)大的計算功能供用戶進(jìn)行數(shù)據(jù)處理，那么就可以實現(xiàn)一個適用于任意工程方法，自由處理數(shù)據(jù)的真正的計算分析平臺。

1 強(qiáng)度計算平臺總體框架

飛機(jī)結(jié)構(gòu)強(qiáng)度計算工作分為有限元建模、全機(jī)有限元求解、有限元數(shù)據(jù)后處理三個過程。本文所介紹的是有限元數(shù)據(jù)后處理，即在已有全機(jī)有限元

模型及其內(nèi)力計算結(jié)果的情況下利用工程方法對具體結(jié)構(gòu)進(jìn)行強(qiáng)度分析的過程。基本流程如圖1所示。

圖1中“需要的單元”、“挑選單元內(nèi)力”、“多工況計算”、“求最小裕度”是計算的主線；“全機(jī)有限元模型”、“全機(jī)內(nèi)力解結(jié)果”、“結(jié)構(gòu)參數(shù)”是計算的輸入，“各類工程方法”是計算的準(zhǔn)則。

圖1 強(qiáng)度計算基本流程

Fig.1 The basic flowchart of stress analysis

本文介紹的方法利用VBA編程，按照計算主線實現(xiàn)各個功能，最后求得最小裕度，程序功能界面如圖2所示。

使用時圖2中“模型導(dǎo)入功能”對應(yīng)圖1的“需要的單元”，圖2中“單元內(nèi)力挑選功能”對應(yīng)圖1的“挑選單元內(nèi)力”，這兩個步驟在用戶準(zhǔn)備好模型文件和結(jié)果文件后完全由程序自動完成。

“多工況計算”是計算的核心部分，即分析不同工況下結(jié)構(gòu)各部位應(yīng)力水平、裕度結(jié)果，利用前兩個步驟的數(shù)據(jù)，用戶按照一定規(guī)則在表格中編寫計算方法，借助圖2中“工況循環(huán)功能”，用戶輸入任意“工況號”代碼，即可得到任意工況下的計算結(jié)果。

強(qiáng)度計算通常關(guān)心結(jié)構(gòu)各種載荷工況下的“最小裕度”，利用圖2中“工況循環(huán)功能”，用戶在“工況號”中置空，可在最終的“總裕度”中得到最小裕度及其對應(yīng)的應(yīng)力等信息。

2 模型導(dǎo)入功能

由于全機(jī)有限元模型單元數(shù)量較多，約10萬個，且全機(jī)求解的工況數(shù)也較多，達(dá)上千種，由此得到了包含上億條單元內(nèi)力信息的文件。為提高計算效率，首先按需“切割”出分析部位結(jié)構(gòu)對應(yīng)的單元，如圖3、圖4所示。

切割后的有限元模型被保存成一個體量較小的模型文件，利用平臺中“模型導(dǎo)入”功能，程序自動把該模型文件中的單元編號及對應(yīng)的屬性信息讀入Excel建立的表“單元信息”Sheet中，如圖5所示。

導(dǎo)入模型后，平臺利用Excel自身的統(tǒng)計功能，自動計算各類單元數(shù)量，存放于圖5中B2、E2、H2......中。為便于后續(xù)計算，同時導(dǎo)入bar、rod單元兩側(cè)所對應(yīng)的殼quad單元，以框結(jié)構(gòu)為例，當(dāng)框結(jié)構(gòu)對應(yīng)的bar單元編號被導(dǎo)入后，程序分析并導(dǎo)入該框單元兩側(cè)相應(yīng)的蒙皮quad單元，quad單元編號存放于該bar單元后的表格C、D兩列內(nèi)。

3 單元內(nèi)力挑選

通過Nastran對全機(jī)有限元模型求解后通常輸出*.xdb或*.op2以及*.f06等文件，前兩者是二進(jìn)制文件，需用Patran、Hyperview等商用軟件查看后處理結(jié)果，*.f06為文本文件，存放著各單元、節(jié)點在各工況下的內(nèi)力、應(yīng)力、節(jié)點平衡力、位移等結(jié)果。

雖然內(nèi)力解文件*.f06為文本文件，但內(nèi)含信息量巨大，通常有數(shù)個G之多，即使查看1個單元在1個工況下的內(nèi)力已經(jīng)相當(dāng)困難，更不用說處理一批單元在不同工況下的內(nèi)力，因此，必須依靠程序把這些需要的單元內(nèi)力從*.f06文件中挑選出來使用。

需要的單元已經(jīng)在第2章中被導(dǎo)入，逐行讀取*.f06文件，對比單元編號，程序挑選出需要的單元內(nèi)力，根據(jù)單元類型不同存放于不同的表格內(nèi)，如rod、bar、quad、node等。

以最簡單的rod單元為例，如圖6所示，平臺把圖5中rod單元包含的內(nèi)力信息從內(nèi)力解文件中挑選出來，并按照單元編號、工況編號、單元內(nèi)力這樣固定的格式排列，此格式中每個單元的單元編號是唯一的，但每個單元可對應(yīng)多個工況編號及其相應(yīng)的內(nèi)力。此排列格式由程序自動完成，是后續(xù)多工況計算的基礎(chǔ)。

對于bar單元，被挑選輸出的內(nèi)力有：1平面A端彎矩、2平面A端彎矩、1平面B端彎矩、2平面B端彎矩、1平面剪切力、2平面剪切力、軸向力，如圖7所示。

圖7 bar單元導(dǎo)出的內(nèi)力信息

Fig.7 The force of bar element have been selected

對于quad單元，被挑選輸出的內(nèi)力有：X方向面線力，Y方向的線力，XY平面內(nèi)的剪切力，如圖8所示。

圖8 quad單元導(dǎo)出的內(nèi)力信息

Fig.8 The force of quad element have been selected

對于node節(jié)點，被挑選輸出的是節(jié)點的位移，包括：X、Y、Z三個方向的平動位移及Rx、Ry、Rz三個方向的轉(zhuǎn)動位移，如圖9所示。

圖9 node節(jié)點導(dǎo)出的位移信息

Fig.9 The displacement of node have been selected

由于Excel本身行數(shù)的限制，最大為1048576行，所以對于任意類型的單元，當(dāng)單元數(shù)×工況數(shù)量過多時，需適當(dāng)?shù)木珳p單元的數(shù)量，或減少計算的工況數(shù)，可根據(jù)實際情況靈活調(diào)整。

4 強(qiáng)度分析計算

飛機(jī)結(jié)構(gòu)的強(qiáng)度分析過程是利用有限元的內(nèi)力解結(jié)果，根據(jù)結(jié)構(gòu)形式，計算出剖面上的應(yīng)力，并把此數(shù)值與結(jié)構(gòu)的許用值進(jìn)行對比的過程，所計算的結(jié)果以裕度的形式得到，如公式（1）所示。

M.S.=■-1（1）

式中：[σ]為結(jié)構(gòu)的許用應(yīng)力；

σ為結(jié)構(gòu)的工作應(yīng)力。

若M.S.大于0，說明結(jié)構(gòu)是安全的，若M.S.小于0，說明結(jié)構(gòu)不滿足強(qiáng)度要求。

在計算工作應(yīng)力σ時，根據(jù)計算部位的單元編號，在“靜強(qiáng)度計算過程”Sheet中通過vlookup（對象單元，對應(yīng)查找區(qū)域，內(nèi)力結(jié)果所在列）函數(shù)，從rod、bar、quad等單元信息中引用當(dāng)前工況對應(yīng)的當(dāng)前內(nèi)力。如計算某段框結(jié)構(gòu)應(yīng)力時根據(jù)公式：

σ=■×y+■（2）

式中： 為結(jié)構(gòu)的工作應(yīng)力；

M為單元的1平面彎矩；

F為單元的軸向力；

I為框的慣性矩；

y為內(nèi)外緣距型心的高度；

A為框的面積。

其中M、F是變量，工況不同則數(shù)值不同；

I、y，A是常數(shù)，與結(jié)構(gòu)本身的特性相關(guān)。

變量M，F(xiàn)的值通過vlookup（）函數(shù)從“bar”sheet中引用。索引的對象為“bar”sheet中的第1列單元編號；對應(yīng)查找的區(qū)域為從A列開始到當(dāng)前內(nèi)力中的最后一列結(jié)束，對于bar單元；內(nèi)力結(jié)果所在列根據(jù)需要選擇不同的列，如需要1平面A端彎矩，則輸入12，如需要梁軸力，則輸入18。完成輸入后由vlookup（）函數(shù)自動得到對象單元所在工況的M和F，結(jié)合結(jié)構(gòu)常數(shù)即可算出應(yīng)力，進(jìn)而得到裕度M.S.，如圖10所示。在計算過程中可充分利用Excel本身的各種運(yùn)算功能，例如可通過下拉來進(jìn)行批量處理各個框段。

至此，得到了當(dāng)前工況（180006）下的結(jié)果，若需查看其它工況，則執(zhí)行圖2“工況循環(huán)功能”，輸入工況號，單擊開始計算即可。由于當(dāng)前工況的內(nèi)力以及靜強(qiáng)度計算中的數(shù)據(jù)都是索引而得，程序通過改變當(dāng)前的工況編號，所有結(jié)果隨之而變，瞬間即可得到任意工況下的結(jié)果。

通常強(qiáng)度分析不僅要得到任意工況下的結(jié)果，還需比較所有工況，找出最嚴(yán)重的結(jié)果。此時程序?qū)λ泄r進(jìn)行一輪循環(huán)計算，比較各工況的裕度，若當(dāng)前工況裕度小于之前計算的裕度則把當(dāng)前工況裕度結(jié)果保留，同時保留對應(yīng)的內(nèi)力和參數(shù)。

操作時在“總裕度”sheet中通過索引“靜強(qiáng)度計算過程”得到當(dāng)前工況的結(jié)果，如圖11中黃色區(qū)域，當(dāng)發(fā)現(xiàn)黃色區(qū)域的裕度值小于藍(lán)色區(qū)域的最小裕度值時，則把黃色區(qū)域的最小值及其對應(yīng)的內(nèi)力和參數(shù)粘貼至藍(lán)色區(qū)域。經(jīng)過一輪工況的循環(huán)，得到最終藍(lán)色區(qū)域為所有工況下計算部位各自的最小裕度。

圖11 當(dāng)前工況與裕度最小工況的比較

Fig.11 The compare of current subcase with dangerous subcase

以上，是以典型簡化的框結(jié)構(gòu)為例進(jìn)行應(yīng)力及最小裕度計算，實際飛機(jī)結(jié)構(gòu)類型多樣，各類計算方法也更加復(fù)雜多變，但本質(zhì)上都是一系列使用單元內(nèi)力及結(jié)構(gòu)參數(shù)的公式f（M，F(xiàn)，F(xiàn)x，F(xiàn)y，F(xiàn)xy......I，A，t......），其中單元內(nèi)力是根據(jù)工況不同而不同的變量，結(jié)構(gòu)參數(shù)是與結(jié)構(gòu)特性相關(guān)的常數(shù)。

因此，用戶可根據(jù)自身計算需要，編輯任意計算方法，并控制需要的輸出結(jié)果，對于計算方法相同的結(jié)構(gòu)又可以使用Excel本身的功能快速的進(jìn)行批量處理。強(qiáng)度計算平臺不同于普通程序，只能計算一種或幾種類型，它提供的是一種計算模式，由用戶在Excel表格內(nèi)自由編輯定義計算方法，因此是一個真正的計算平臺。

此外，在飛機(jī)設(shè)計分析時還采用一種稱為DFR的疲勞強(qiáng)度分析方法，疲勞分析需同時考慮幾十種載荷情況的組合，形成所謂的疲勞載荷譜，本程序在循環(huán)工況時記錄各個計算點的工作應(yīng)力，從而實現(xiàn)對應(yīng)力譜的獲得，最終把應(yīng)力譜用于疲勞計算表格得到疲勞裕度。疲勞計算在模型導(dǎo)入、單元內(nèi)力挑選上與靜強(qiáng)度計算方法完全相同，在計算應(yīng)力時與靜強(qiáng)度也相同，但是靜強(qiáng)度直接利用當(dāng)前工況的應(yīng)力值計算出當(dāng)前工況的裕度，而疲勞計算把所有的工況的應(yīng)力值組合在一起同時使用，得到一個裕度。

5 結(jié)論

強(qiáng)度計算平臺與普通的小程序以及空客公司的ISAMI平臺相比，他的最大特點是程序內(nèi)并不包含具體的計算方法，而是把計算方法放在程序外由用戶編寫。同時他又是基于Excel軟件來實現(xiàn)的，因此可以充分的利用Excel的功能，歸納起來，強(qiáng)度計算寄生平臺相較與普通的計算程序具有以下優(yōu)點：

（1）通用性。強(qiáng)度計算寄生平臺并不計算特定具體結(jié)構(gòu)，而是為計算分析提供一種分析環(huán)境，用戶根據(jù)需要可計算任意結(jié)構(gòu)。

（2）易用性。平臺不需要用戶懂得某種語言進(jìn)行編程，只需用戶在工作表中編輯計算過程，平臺即可幫助用戶實現(xiàn)靜強(qiáng)度和疲勞強(qiáng)度分析。

（3）靈活性。平臺“寄生”在Excel中，用戶可充分利用Excel的各種函數(shù)和計算功能，靈活方便的進(jìn)行計算，輸出清晰、規(guī)則的計算結(jié)果。

（4）集成性。平臺把原始數(shù)據(jù)、計算過程、計算結(jié)果集成在一個Excel文件內(nèi)，數(shù)據(jù)簡潔，便于保持、追溯，檢查，校對。

（5）高效性。平臺計算速度快，數(shù)據(jù)更新方便，只需數(shù)次點擊即可替換新一輪內(nèi)力解結(jié)果。

（6）透明。平臺可集成大量工程算法，也可由用戶自己定義計算方法，用戶可查看計算過稱，了解計算方法，避免了使用傳統(tǒng)程序，計算過程是“黑匣子”的不透明性。

（7）綠色。強(qiáng)度計算寄生平臺不需要特殊安裝，只需計算機(jī)安裝Excel軟件即可使用，走到哪、用到哪，對于出差工作，對外合作交流都較為便利。

6 展望

本強(qiáng)度計算平臺是基于飛機(jī)結(jié)構(gòu)的強(qiáng)度分析而誕生的，由于受專業(yè)限制，其應(yīng)用目前也僅限于飛機(jī)結(jié)構(gòu)強(qiáng)度，但是對于其他同樣采用有限元進(jìn)行求解的領(lǐng)域如汽車，船舶，建筑，地質(zhì)，大氣等也可采用此平臺的思路，通過增加單元類型，進(jìn)行分析。

【參考文獻(xiàn)】

[1]飛機(jī)設(shè)計手冊 9 載荷、強(qiáng)度和剛度，北京，航空工業(yè)出版社，2001.

[2]中國商用飛機(jī)有限責(zé)任公司上海飛機(jī)設(shè)計研究院.SADRI單元內(nèi)力挑選及強(qiáng)度計算軟件V1.0.： 2012SR011885.

[責(zé)任編輯：朱麗娜]