王昌林，李 磊，王全紅

(中國飛機(jī)強(qiáng)度研究所，陜西 西安 710065)

1 引 言

機(jī)翼是與飛機(jī)氣動(dòng)和結(jié)構(gòu)這兩部分關(guān)系最密切、最復(fù)雜的部分。在對(duì)機(jī)翼結(jié)構(gòu)進(jìn)行靜力、疲勞試驗(yàn)及結(jié)構(gòu)有限元分析的過程中，需要將機(jī)翼表面產(chǎn)生的氣動(dòng)載荷轉(zhuǎn)換到結(jié)構(gòu)上。氣動(dòng)載荷在設(shè)計(jì)階段通常采用CFD的方法獲得，在計(jì)算之前必須對(duì)計(jì)算區(qū)域進(jìn)行離散化處理，將空間上的計(jì)算區(qū)域劃分為多個(gè)子區(qū)域，并在每個(gè)子區(qū)域中確定節(jié)點(diǎn)，生成網(wǎng)格，即氣動(dòng)網(wǎng)格點(diǎn)壓力離散場(chǎng)的形式出現(xiàn)。此方法需要?jiǎng)澐忠?guī)模巨大的網(wǎng)格來進(jìn)行計(jì)算，所產(chǎn)生的氣動(dòng)載荷數(shù)據(jù)量也很巨大，而進(jìn)行結(jié)構(gòu)試驗(yàn)或有限元分析時(shí)所使用的加載節(jié)點(diǎn)或計(jì)算網(wǎng)格規(guī)模較小，此二者之間的數(shù)量相差幾個(gè)甚至幾十個(gè)量級(jí)，此時(shí)需要將氣動(dòng)載荷從氣動(dòng)網(wǎng)格節(jié)點(diǎn)轉(zhuǎn)換到試驗(yàn)加載節(jié)點(diǎn)或有限元分析所用的網(wǎng)格節(jié)點(diǎn)上。此外，氣動(dòng)力加載會(huì)引起結(jié)構(gòu)變形，而結(jié)構(gòu)變形會(huì)引起氣動(dòng)外形的變化，造成氣動(dòng)力分布的改變，需要在變形后的兩個(gè)網(wǎng)格之間產(chǎn)生位移和力的互相傳遞。因此，加載到結(jié)構(gòu)上的載荷能否精確地反映氣動(dòng)特性非常重要。

常用的二者之間載荷轉(zhuǎn)換與傳遞的方法有以下幾種：三點(diǎn)排[1]、四點(diǎn)排[2]、多點(diǎn)排[1]。這三種方法在原理上只支持氣動(dòng)載荷向有限元節(jié)點(diǎn)載荷的轉(zhuǎn)換與傳遞，并不支持載荷引起的結(jié)構(gòu)變形信息向氣動(dòng)外形的反饋傳遞。

本文使用MpCCI作為載荷轉(zhuǎn)換方法，分別與以上幾種方法進(jìn)行過程與結(jié)果的對(duì)比，通過算例說明這種方法在實(shí)際應(yīng)用中的意義，以期為進(jìn)行機(jī)翼結(jié)構(gòu)靜力、疲勞試驗(yàn)及有限元分析等過程中載荷轉(zhuǎn)換方法的選擇提供參考。

2 常用載荷轉(zhuǎn)換與分配方法

2.1 三點(diǎn)排方案

三點(diǎn)排(如圖1所示)表述了如何將氣動(dòng)節(jié)點(diǎn)A上的載荷PA按照靜力等效的原則分配轉(zhuǎn)換到鄰近的3個(gè)結(jié)構(gòu)節(jié)點(diǎn)上[3]。

圖1 三點(diǎn)排示意圖

假設(shè)選定點(diǎn)1、2、3作為結(jié)構(gòu)節(jié)點(diǎn)，那么氣動(dòng)點(diǎn)A分配到這3個(gè)有限元節(jié)點(diǎn)上的載荷計(jì)算公式如下：

(1)

式中，PA為氣動(dòng)節(jié)點(diǎn)A上的載荷；A、A1、A2、A3分別為△123、△A23、△A13、△A12的面積。

對(duì)所有的氣動(dòng)節(jié)點(diǎn)載荷都依據(jù)式(1)計(jì)算得到對(duì)應(yīng)的結(jié)構(gòu)節(jié)點(diǎn)載荷[4]，然后把同一結(jié)構(gòu)節(jié)點(diǎn)上分配到的載荷疊加起來，即為這一結(jié)構(gòu)節(jié)點(diǎn)上的總載荷。

2.2 四點(diǎn)排方案

四點(diǎn)排示意圖如圖2所示，橫縱向?qū)嵕€分別表示長桁與肋的布置線，四點(diǎn)排就是將氣動(dòng)載荷節(jié)點(diǎn)C上的氣動(dòng)力轉(zhuǎn)換到結(jié)構(gòu)的結(jié)構(gòu)布線的節(jié)點(diǎn)上。

圖2 四點(diǎn)排示意圖

(2)

2.3 多點(diǎn)排方案

多點(diǎn)排方案(如圖3所示)的基本思路為：在分配過程中，距離氣動(dòng)節(jié)點(diǎn)近的結(jié)構(gòu)節(jié)點(diǎn)分配權(quán)重大一些，反之小一些。假設(shè)有一根無形的懸臂梁處于結(jié)構(gòu)節(jié)點(diǎn)和氣動(dòng)節(jié)點(diǎn)之間，EL為假想懸臂梁的抗彎剛度[5]。這根梁在氣動(dòng)節(jié)點(diǎn)的一端為固支，自由端上的結(jié)構(gòu)節(jié)點(diǎn)分配到載荷Pj時(shí)，其變形能力為：

(3)

圖3 多點(diǎn)排示意圖

那么，對(duì)于由所有結(jié)構(gòu)節(jié)點(diǎn)與某氣動(dòng)節(jié)點(diǎn)組成的所有假想懸臂梁系統(tǒng)而言，其變形能力為：

(4)

根據(jù)最小勢(shì)能原理計(jì)算分配到結(jié)構(gòu)節(jié)點(diǎn)上的載荷，在滿足連續(xù)條件和邊界條件的位移中，滿足平衡條件的位移其總勢(shì)能最小，那么欲使其系統(tǒng)變形能量最小，同時(shí)還應(yīng)滿足靜力等效條件：

(5)

(6)

(7)

式中，n為結(jié)構(gòu)節(jié)點(diǎn)數(shù)。

采用拉格朗日(Lagrange)乘子法建立拉格朗日函數(shù)：

(8)

為了使F(λλxλz)能夠取最小值，令：

(9)

得：

(10)

整理公式得：

由上式解出λ、λx、λz后，代入式(10)中即可得到結(jié)構(gòu)節(jié)點(diǎn)所分配到的載荷Pj。對(duì)整個(gè)氣動(dòng)節(jié)點(diǎn)網(wǎng)群上的每一個(gè)氣動(dòng)載荷，按上述方法分配，就得到了所期望的載荷轉(zhuǎn)換結(jié)果。

幾種常用載荷轉(zhuǎn)換方法的基本思路如上文所述，都遵循了靜力等效和能量原理。值得注意的是，在實(shí)際應(yīng)用過程中，往往需要解決的是成千上萬氣動(dòng)節(jié)點(diǎn)向成百上千結(jié)構(gòu)節(jié)點(diǎn)轉(zhuǎn)換載荷的問題，通過手工計(jì)算解決上述問題是不切實(shí)際的，此時(shí)就需要借助計(jì)算機(jī)以及有關(guān)軟件的幫助。此外，氣動(dòng)載荷的取得途徑是多樣的，常用的CFD軟件也種類繁多。上述3種載荷轉(zhuǎn)換的方法，針對(duì)不同來源的氣動(dòng)載荷和不同形式的結(jié)構(gòu)模型，處理起來也不夠友好、通用性也不強(qiáng)，在處理過程中需要編寫不同的接口程序來滿足轉(zhuǎn)換過程中上、下游數(shù)據(jù)流通的需要。

2.4 MpCCI方法

采用MpCCI進(jìn)行載荷轉(zhuǎn)換的基本過程是：直接調(diào)用CFD軟件對(duì)結(jié)構(gòu)的外流場(chǎng)進(jìn)行計(jì)算或調(diào)用已有的氣動(dòng)載荷，在耦合界面處提取氣動(dòng)載荷，進(jìn)行載荷轉(zhuǎn)換輸出結(jié)果或者直接將結(jié)果提交給用戶程序。由于氣動(dòng)節(jié)點(diǎn)組成的網(wǎng)格和結(jié)構(gòu)節(jié)點(diǎn)組成的網(wǎng)格的計(jì)算要求不同，在力的傳遞界面上網(wǎng)格不匹配，節(jié)點(diǎn)是不重合的，從而不能直接進(jìn)行數(shù)據(jù)交換，需要在傳遞界面上采用映射的方法。這里用Nf和Ns分別表示對(duì)應(yīng)氣動(dòng)網(wǎng)格和結(jié)構(gòu)網(wǎng)格上的節(jié)點(diǎn)，如果Nf和Ns不重合，可使用一定的算子映射來進(jìn)行計(jì)算，表示為：

Msf:{fi|i∈Nf}→{fi|i∈Ns}

Mfs:{fi|i∈Ns}→{fi|i∈Nf}

上式中，fi和si分別為傳遞邊界上氣動(dòng)、結(jié)構(gòu)節(jié)點(diǎn)的矢量。由于氣動(dòng)域和結(jié)構(gòu)域的網(wǎng)格不同，即Msf≠M(fèi)fs-1。

在邊界上氣動(dòng)節(jié)點(diǎn)的位移如下：

uf=Msfus

力的平衡條件是將氣動(dòng)節(jié)點(diǎn)力先映射到結(jié)構(gòu)節(jié)點(diǎn)，然后通過積分得到結(jié)構(gòu)節(jié)點(diǎn)力，在此過程中采用固體邊界的插值函數(shù)。

式中，HST是記錄結(jié)構(gòu)單元沿著邊界的插值函數(shù)，這樣就可以計(jì)算出結(jié)構(gòu)模型的節(jié)點(diǎn)力，在有限元模型計(jì)算的過程中作為邊界條件使用。

3 算例對(duì)比

3.1 二維梯形機(jī)翼沿展向橢圓分布弦向三角形分布

氣動(dòng)載荷沿梯形機(jī)翼的展向呈橢圓分布、沿弦向呈三角形分布。弦向分布方程：在0到0.4以內(nèi)為0.8x+0.2，在0.4到2為-0.56x+1.22。其中，氣動(dòng)節(jié)點(diǎn)的數(shù)目為189，結(jié)構(gòu)節(jié)點(diǎn)的數(shù)目為78，氣動(dòng)分布圖和網(wǎng)格節(jié)點(diǎn)示意圖見圖4、圖5。

圖4 氣動(dòng)分布圖

圖5 網(wǎng)格節(jié)點(diǎn)示意圖

經(jīng)過三點(diǎn)排、四點(diǎn)排、多點(diǎn)排、MpCCI等方式進(jìn)行載荷轉(zhuǎn)換后，所得結(jié)果如圖6、圖7所示。

圖6 展向節(jié)點(diǎn)載荷分配示意圖

圖7 結(jié)構(gòu)力云圖

從圖中可以看出，三點(diǎn)排和四點(diǎn)排方案所取得的載荷有著很大幅度的跳躍，且最大分配載荷明顯比多點(diǎn)排和MpCCI所取得的載荷大。與三點(diǎn)排和四點(diǎn)排相比，多點(diǎn)排方案有了很大的改善，但與MpCCI取得的結(jié)果相比還是有差距的。MpCCI取得的結(jié)果呈橢圓分布，基本沒有出現(xiàn)跳躍，分布的規(guī)律更接近于氣動(dòng)載荷，曲線也最為光滑。表1給出了各方法所取載荷與原氣動(dòng)載荷的相關(guān)系數(shù)。

表1 各方法所取載荷與原載荷相關(guān)系數(shù)

3.2 二維平直機(jī)翼沿展向橢圓分布弦向三角形分布

平直機(jī)翼的氣動(dòng)分布與上文的梯形機(jī)翼的氣動(dòng)分布一致，氣動(dòng)節(jié)點(diǎn)數(shù)為189，結(jié)構(gòu)有限元節(jié)點(diǎn)數(shù)為78，氣動(dòng)分布圖和網(wǎng)格節(jié)點(diǎn)示意圖見圖8、圖9，各載荷轉(zhuǎn)換方式的展向節(jié)點(diǎn)載荷分配及結(jié)構(gòu)力云圖如圖10、圖11所示。

圖8 氣動(dòng)分布圖

圖9 網(wǎng)格節(jié)點(diǎn)示意圖

圖10 展向節(jié)點(diǎn)載荷分配示意圖

圖11 結(jié)構(gòu)力云圖

從圖中可以看出：無論是中線還是前緣，這幾種載荷轉(zhuǎn)換方式下，MpCCI取得的結(jié)果最接近真實(shí)氣動(dòng)載荷，跳躍小，過渡最為光滑。表2給出了各方法所取載荷與原氣動(dòng)載荷的相關(guān)系數(shù)。

表2 各方法所取載荷與原載荷相關(guān)系數(shù)

3.3 小 結(jié)

從以上兩個(gè)算例可以看出，MpCCI在氣動(dòng)載荷轉(zhuǎn)換為結(jié)構(gòu)有限元載荷方面，比其他轉(zhuǎn)換方法精度更高。

4 三維機(jī)翼載荷轉(zhuǎn)換實(shí)例

M6機(jī)翼為跨音速試驗(yàn)驗(yàn)證機(jī)翼，展長1.1963m，前緣后掠角30°，后緣后掠角15.8°，機(jī)翼面積0.7537451m2，蒙皮厚度0.0045m，梁腹板厚度0.011m，翼肋厚度0.006m，材料為鋁合金，彈性模量E=70GPa，泊松比μ=0.33。

外流場(chǎng)網(wǎng)格劃分使用的是ICEM軟件，總共劃分了73萬個(gè)網(wǎng)格。計(jì)算氣動(dòng)力的軟件為FLUENT，模型單位尺寸為毫米；三維基于密度隱式求解器，K-ωSST湍流模型，湍流強(qiáng)度1.86%，采用二階迎風(fēng)格式，流場(chǎng)進(jìn)口使用壓力遠(yuǎn)場(chǎng)邊界條件。氣動(dòng)網(wǎng)格如圖12所示。

圖12 氣動(dòng)網(wǎng)格

結(jié)構(gòu)模型網(wǎng)格劃分使用的是SimXpert，總共劃分網(wǎng)格數(shù)量為1075，其中所有節(jié)點(diǎn)均為梁肋布置線的交點(diǎn)，結(jié)構(gòu)有限元網(wǎng)格如圖13所示。本文對(duì)機(jī)翼結(jié)構(gòu)模型作了適當(dāng)簡化：

(1)只考慮了組成承力翼盒的主要部件：蒙皮、梁腹板、翼肋，長桁和梁緣條厚度“打扁”后計(jì)入蒙皮的厚度，以減少有限元模型的復(fù)雜程度和計(jì)算量，文獻(xiàn)[6]和文獻(xiàn)[7]均采用這種簡化方法；

(2)在實(shí)際的機(jī)翼結(jié)構(gòu)中，上下蒙皮厚度、前后梁腹板厚度是沿展向變化的，本文認(rèn)為蒙皮厚度、梁腹板厚度沿展向是不變的，同樣是為了減少模型的復(fù)雜程度；

(3)在實(shí)際情況中，機(jī)翼和機(jī)身的連接情況是復(fù)雜的，本文假設(shè)機(jī)翼根部為固支。

圖14為在MpCCI中載荷轉(zhuǎn)換過程。

圖14 轉(zhuǎn)換過程

載荷轉(zhuǎn)換結(jié)果分析：由于是將載荷轉(zhuǎn)換后提交Nastran進(jìn)行計(jì)算，所以根據(jù)Nastran不同的設(shè)置可以得到不同的輸出結(jié)果。圖15是輸出的機(jī)翼上下表面應(yīng)力云圖，可以看出，應(yīng)力云圖的分布與氣動(dòng)力的分布一致。此外，計(jì)算所得翼尖的位移為0.0065m，這基本與試驗(yàn)值相符合。

在以往的研究中，通過Fluent、Star-CD、CFX等CFD軟件實(shí)現(xiàn)了與結(jié)構(gòu)有限元軟件之間的單向流固耦合，這樣沒有真正地實(shí)現(xiàn)數(shù)據(jù)在迭代求解過程中的動(dòng)態(tài)傳輸，降低了求解精度，與實(shí)際情況有著較大的差異。借助MpCCI雙向流固耦合接口，將有限元軟件與CFD程序軟件進(jìn)行連接，實(shí)現(xiàn)了多相耦合，多維度耦合分析，以及動(dòng)態(tài)雙向流固耦合分析。同時(shí)，MpCCI能夠?qū)⒔^大多數(shù)的CFD軟件和有限元軟件以及熱力學(xué)軟件耦合起來，還可以使用簡單的代碼將Matlab等工具結(jié)合起來。

5 結(jié) 論

本文使用MpCCI作為氣動(dòng)節(jié)點(diǎn)載荷向結(jié)構(gòu)有限元節(jié)點(diǎn)

圖15 機(jī)翼上下表面應(yīng)力云圖

載荷轉(zhuǎn)換的平臺(tái)，通過算例，對(duì)比分析了其與傳統(tǒng)轉(zhuǎn)換方法之間的差異，突顯出了這種方法在使用上的優(yōu)越性。研究結(jié)果表明，MpCCI在轉(zhuǎn)換精度上優(yōu)于傳統(tǒng)方法。

同時(shí)，本文通過MpCCI對(duì)FLUENT和NASTRAN的調(diào)用，模擬并取得機(jī)翼結(jié)構(gòu)在氣動(dòng)加載條件下的響應(yīng)，顯示了MpCCI強(qiáng)大的通用性，除了用于載荷轉(zhuǎn)換外，在多物理場(chǎng)耦合方面也有強(qiáng)大的作用。此方法有一定的計(jì)算效率、精度及工程應(yīng)用價(jià)值。