歐陽星， 余雄慶， 王 宇

(南京航空航天大學(xué)飛行器先進(jìn)設(shè)計(jì)技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室， 江蘇 南京 210016)

采用等效剛度有限元模型的復(fù)合材料機(jī)翼顫振分析

歐陽星， 余雄慶， 王 宇

(南京航空航天大學(xué)飛行器先進(jìn)設(shè)計(jì)技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室， 江蘇 南京 210016)

采用等效剛度方法，研究了一種適用于復(fù)合材料機(jī)翼初步設(shè)計(jì)階段的動(dòng)力學(xué)和顫振分析的簡化結(jié)構(gòu)有限元模型。該方法首先計(jì)算不同布局形式的加筋壁板的剛度矩陣，然后將其賦予與加筋壁板平面形狀相同的光板(等效板)上，使加筋壁板和等效板具相同的力學(xué)性能。以某客機(jī)概念方案的機(jī)翼為例，建立了反映實(shí)際結(jié)構(gòu)詳細(xì)有限元及其等效剛度有限元模型。通過編程實(shí)現(xiàn)了自動(dòng)化生成等效模型，極大地減少了建模和分析時(shí)間。計(jì)算結(jié)果和對比分析表明，兩種模型的固有頻率、振動(dòng)模態(tài)和顫振分析結(jié)果吻合得很好，從而驗(yàn)證了等效剛度方法在機(jī)翼結(jié)構(gòu)動(dòng)力學(xué)和顫振分析方面的準(zhǔn)確性。

機(jī)翼； 顫振； 復(fù)合材料； 有限元； 等效剛度

引 言

機(jī)翼主承力結(jié)構(gòu)采用復(fù)合材料可有效減輕重量，是現(xiàn)代飛機(jī)先進(jìn)性的重要標(biāo)志之一。機(jī)翼結(jié)構(gòu)必須保證在整個(gè)飛行包線不能發(fā)生顫振。這一設(shè)計(jì)要求對機(jī)翼結(jié)構(gòu)方案設(shè)計(jì)有重要影響，因此在機(jī)翼的初步方案設(shè)計(jì)中，就必須對機(jī)翼的顫振特性進(jìn)行分析。

結(jié)構(gòu)動(dòng)力學(xué)分析是顫振分析的基礎(chǔ)。適用于復(fù)合材料機(jī)翼顫振分析的結(jié)構(gòu)動(dòng)力學(xué)分析模型可概括為三類：1)等效梁模型[1-2]；2)等效板模型[3-5]；3)有限元模型[6-8]。等效梁模型是將機(jī)翼簡化為梁，主要適用于展弦比較大的機(jī)翼。等效板模型是將機(jī)翼等效至與機(jī)翼相同平面形狀的一組梯形板上，梯形板的力學(xué)性能逼近真實(shí)機(jī)翼，多用于小展弦比機(jī)翼的動(dòng)力學(xué)分析。有限元方法能詳細(xì)描述機(jī)翼結(jié)構(gòu)特征，具有通用和計(jì)算精度高等優(yōu)點(diǎn)。

目前飛機(jī)設(shè)計(jì)的一個(gè)重要趨勢是在飛機(jī)早期設(shè)計(jì)階段采用高精度分析方法[9]。這是因?yàn)椋?)高精度分析模型可用于新概念飛機(jī)方案(如翼身融合或聯(lián)接翼方案)；2)能提高分析結(jié)果的可信度，減少設(shè)計(jì)方案的風(fēng)險(xiǎn)。

有限元模型是結(jié)構(gòu)動(dòng)力學(xué)分析的高精度模型。但詳細(xì)的有限元模型復(fù)雜，其建模過程繁雜，所需計(jì)算時(shí)間長。由于在飛機(jī)初步設(shè)計(jì)階段，機(jī)翼外形和結(jié)構(gòu)方案還在不斷調(diào)整之中，需要一種快速的顫振評估方法，而有限元模型的上述缺點(diǎn)阻礙了其在初步設(shè)計(jì)階段的應(yīng)用。解決這個(gè)問題的一個(gè)措施是采用等效剛度方法對結(jié)構(gòu)有限元模型進(jìn)行簡化。

等效剛度法是將加筋壁板等效到一塊相同形狀的光板(無加強(qiáng)筋的板)上，并使光板與加筋壁板具有相同的力學(xué)性能?；诘刃偠确ǖ挠邢拊Ｐ?以下簡稱等效剛度有限元模型)能模擬各種加強(qiáng)筋的力學(xué)特性，同時(shí)避免了對加強(qiáng)筋的有限元建模，大大簡化了有限元模型，易于實(shí)施參數(shù)化機(jī)翼結(jié)構(gòu)有限元建模。該方法已應(yīng)用于結(jié)構(gòu)靜強(qiáng)度和穩(wěn)定性分析。例如，Bradley[10]采用一種等效平板法對翼身融合布局飛機(jī)進(jìn)行有限元靜力分析。趙群[11]采用等效剛度法對復(fù)合材料加筋板總體失穩(wěn)性能進(jìn)行分析。Wenzel[12]采用等效剛度方法開展了以靜強(qiáng)度和屈曲作為約束的復(fù)合材料垂尾結(jié)構(gòu)重量計(jì)算方法的研究。但應(yīng)用等效剛度有限元模型來分析復(fù)材機(jī)翼的顫振特性，還鮮見有關(guān)資料的報(bào)導(dǎo)，動(dòng)力學(xué)特性和顫振特性的預(yù)測精度還不甚明了。

為給機(jī)翼初步設(shè)計(jì)階段提供一種高精度的、快速的結(jié)構(gòu)動(dòng)力學(xué)分析和顫振分析方法，本文應(yīng)用剛度等效方法簡化復(fù)合材料機(jī)翼結(jié)構(gòu)有限元模型，驗(yàn)證機(jī)翼等效剛度有限元模型的計(jì)算精度。

1 等效剛度方法

1.1 剛度等效的思路

現(xiàn)代飛機(jī)通常采用桁條和蒙皮組成的半硬殼式加筋結(jié)構(gòu)。桁條相互平行或者接近平行，其兩端有翼肋支持。桁條有兩個(gè)主要的作用：一是與蒙皮一起承受軸向拉伸/壓縮載荷；二是使加筋壁板的等效中面偏移蒙皮中面，提高了蒙皮局部抗彎剛度，提升了蒙皮局部抗失穩(wěn)性能。

等效剛度方法的目的就是簡化薄壁加筋結(jié)構(gòu)的有限元模型。它的主要思路是：將加強(qiáng)筋離散為一系列的板條元，然后計(jì)算加筋壁板的等效中面，根據(jù)復(fù)合材料經(jīng)典層合板理論，分別計(jì)算蒙皮的剛度系數(shù)和各板條元?jiǎng)偠认禂?shù)，疊加而成等效剛度矩陣；如果是夾芯結(jié)構(gòu)，則是面板剛度系數(shù)和夾芯結(jié)構(gòu)剛度系數(shù)疊加而成等效剛度矩陣，最后將等效剛度矩陣賦予等效板的屬性參數(shù)。圖1顯示了用等效剛度法對不同的加筋結(jié)構(gòu)方案的等效過程。等效剛度方法帶來的益處是：當(dāng)加筋結(jié)構(gòu)方案改變時(shí)，只需調(diào)整加強(qiáng)筋板條元的外形參數(shù)(寬和高)和加強(qiáng)筋布置參數(shù)(間距)，無需重新建立有限元模型，從而增強(qiáng)了有限元模型的參數(shù)化功能。

圖1 各種加筋/夾層壁板等效為光板Fig.1 Modeling approach of different stiffened panels

將剛度矩陣A中的系數(shù)定義薄膜材料(Membrane material)和剪切材料(Shear material)，將剛度矩陣B中的系數(shù)定義耦合材料(Coupling material)，將剛度矩陣D中的系數(shù)定義彎曲材料(Bending material)，然后將這些材料賦予有限元模型中等效板的屬性。對于不同布局的加強(qiáng)壁板，計(jì)算出不同的剛度矩陣，因此賦予光板不同的力學(xué)屬性。

1.2 等效剛度參數(shù)

(1)

式中Aij，Bij和Dij分別代表層板的面內(nèi)、耦合以及彎曲剛度系數(shù)。

圖2 層合板示意圖Fig.2 Coordinate locations of plies in a laminate

復(fù)合材料蒙皮和加強(qiáng)筋可視為各向異性板，可按照復(fù)合材料經(jīng)典層合板理論計(jì)算蒙皮和各板條元的剛度系數(shù)，載荷/變形關(guān)系如下[13]

(2)

1.3 加筋壁板等效中面計(jì)算

計(jì)算蒙皮和各板條元的彎曲剛度需要參考中性面的位置。桁條的存在使加強(qiáng)區(qū)域的結(jié)構(gòu)中性面偏離蒙皮中面，而桁條之間的非加強(qiáng)區(qū)域的中性面仍然與蒙皮中面重合。因此，加筋壁板實(shí)際中性面并不是平面而是曲面，如圖3所示。

圖3 等效中面Fig.3 Equivalent neutral surface

當(dāng)加筋壁板承受彎矩時(shí)，蒙皮和桁條均以中性面產(chǎn)生形變，但由于該實(shí)際中性面為曲面，難以確定各板條元以其為參考面的本構(gòu)關(guān)系，所以采用一個(gè)等效中性面(圖3中xy平面)替代曲面的實(shí)際中性面，以便于計(jì)算。基于純彎曲假設(shè)，彎曲變形時(shí)等效中性面的曲率為кx，整個(gè)截面沿筋條方向的合力為零[14]，即

(3)

由于Δhxi=-Δzxc+Δzxi，因此式(3)可寫為

(4)

Δzxc表示等效中性面到蒙皮中面的距離，由式(4)可以確定其值為

(5)

1.4 剛度系數(shù)平移計(jì)算

由于加筋壁板的等效中面偏移蒙皮中面，蒙皮的彎曲剛度得到增強(qiáng)，需要進(jìn)行平行移面計(jì)算。蒙皮和加強(qiáng)筋相對于自身中面剛度系數(shù)與等效中面的剛度系數(shù)之間有如下關(guān)系[15]

(6)

圖4 中面平移Fig.4 Translation of neutral surface

1.5 蒙皮剛度矩陣計(jì)算

復(fù)合材料層合板為各向異性板，設(shè)計(jì)復(fù)合材料機(jī)翼蒙皮時(shí)，為避免固化時(shí)因耦合效應(yīng)引起翹曲變形，蒙皮采取對稱鋪層。對于正交各向異性對稱層合板，可以忽略拉彎耦合剛度，所以蒙皮的載荷/變形關(guān)系可簡化為

(7)

對于復(fù)合材料加筋板，等效中面偏移出蒙皮中面，這時(shí)蒙皮的彎曲剛度由兩部分組成，一是蒙皮相對于自身中面彎曲剛度，二是中面偏移增加的彎曲剛度，可由式(6)計(jì)算。然而采用對稱鋪層的蒙皮沒有拉彎耦合剛度，故平行移面時(shí)它不會(huì)對彎曲剛度產(chǎn)生貢獻(xiàn)。由于桁條軸向?yàn)閤向，提高蒙皮在zx平面內(nèi)的抗彎剛度，蒙皮和加強(qiáng)筋板條元只需在z向進(jìn)行平移計(jì)算。平移后的蒙皮剛度矩陣為

(8)

1.6 桁條剛度矩陣計(jì)算

對于復(fù)合材料桁條，腹板和上下凸緣可采用不同的鋪層，故需要分開處理。以工字型加強(qiáng)筋為例，將其離散為5個(gè)板條元，各板條元寬度為wi，厚度為hi，如圖5所示。圖中蒙皮厚度為tsk，桁條間距為bstr，蒙皮中面到組合截面等效中性面的距離為Δzxc。通過設(shè)置板條元的參數(shù)，很容易地模擬工程中常見不同截面形狀的加強(qiáng)筋。例如，將板條元5的寬度和高度參數(shù)設(shè)置為0時(shí)，加強(qiáng)筋則為J形截面。將板條元4，5的參數(shù)同時(shí)設(shè)為0時(shí)，則為T形截面，將板條元2，4的參數(shù)設(shè)置為0時(shí)，則是Z形截面。通過加筋壁板的參數(shù)化定義，可以快捷地描述不同加筋板的形狀。

圖5 開口加筋壁板的典型剖面Fig.5 A typical profile of stiffened panel with open stiffeners

1.6.1x向桁條剛度矩陣

假設(shè)桁條的軸向?yàn)閤向，垂直于桁條向右為y向，垂直蒙皮向下為z向，如圖5所示。蒙皮和桁條上的載荷按剛度分配，在單位寬度等效剛度蒙皮上，桁條上的軸向壓縮載荷與彎矩表示為各離散板條上載荷的疊加：

(9)

(10)

對于共固化的復(fù)合材料蒙皮和桁條，蒙皮接觸的板條元可以承受面內(nèi)剪切載荷，因此需計(jì)算該板條元剪切剛度

(11)

式中i=1,2，即只計(jì)算板條元1和2的剪切剛度。

同時(shí)桁條的扭轉(zhuǎn)剛度也需考慮，對于開口加強(qiáng)筋的扭轉(zhuǎn)剛度計(jì)算公式如下[15]

(12)

而對于閉口加強(qiáng)筋的扭轉(zhuǎn)剛度計(jì)算公式如下

(13)

式中i=2,3,4,6；Ω1為閉口加筋條壁厚中心線與蒙皮中面所圍成的面積；GJ1為閉室Ω1的抗扭剛度。

圖6 閉口加筋壁板的典型剖面Fig.6 A typical profile of stiffened panel with closed stiffeners

平行移面計(jì)算后，x向桁條的載荷/變形關(guān)系為

(14)

1.6.2y向桁條剛度矩陣

同理可得，平移移面計(jì)算后的y向桁條的載荷/變形關(guān)系為

(15)

而雙向加筋布局的桁條剛度矩陣為x向和y向桁條剛度矩陣的疊加。

綜上所述，統(tǒng)一的桁條剛度矩陣可寫為

(16)

1.7 組裝蒙皮桁條剛度矩陣

分別計(jì)算蒙皮和離散的加強(qiáng)筋板條元的剛度矩陣后，疊加即可組裝加筋壁板的等效剛度矩陣，簡化為

(17)

組裝后的蒙皮和桁條耦合剛度正負(fù)抵消，即等效剛度中的耦合剛度Beq=Bsk+Bst=0。

對于夾芯加強(qiáng)壁板，則是計(jì)算上、下面板和夾芯結(jié)構(gòu)對等效中面剛度矩陣，然后疊加而成等效剛度矩陣。

2 算例驗(yàn)證

2.1 機(jī)翼模型描述

本文以如圖7所示客機(jī)概念方案為研究背景，以其機(jī)翼作為例子，驗(yàn)證等效剛度有限元模型的結(jié)構(gòu)動(dòng)力學(xué)和顫振計(jì)算精度。

圖7 某客機(jī)概念方案Fig.7 A notional civil jet

該機(jī)翼展弦比為10.0，參考面積為95.75 m2，1/4弦線后掠角為24.5°，機(jī)翼根部相對厚度0.14。機(jī)翼的詳細(xì)有限元模型和等效剛度有限元模型見圖8。兩種模型的載荷和邊界條件相同。桁條(加強(qiáng)筋)截面為矩形，平行于后梁布置，與復(fù)合材料蒙皮共同固化，參見圖8(a)所示。詳細(xì)模型包含了53 697個(gè)單元，251 003個(gè)自由度。等效模型包含了6 205個(gè)單元，26 364個(gè)自由度，參見圖8(b)所示。可見，等效模型中單元數(shù)只有詳細(xì)模型的單元數(shù)量的1/9，模型規(guī)模大為降低。

復(fù)合材料機(jī)翼的蒙皮和桁條鋪層方向和比例如表1所示。沿展向劃分9個(gè)區(qū)域，每個(gè)區(qū)域具有相同的厚度。對于詳細(xì)有限元模型，由于其包含很多細(xì)節(jié)信息，有限元建模需要8 h左右，各設(shè)計(jì)區(qū)域復(fù)合材料定義和屬性的賦予也需要7 h，在Patran的Flightloads模塊中定義顫振計(jì)算工況需要5 h左右，最后提交Nastran進(jìn)行顫振計(jì)算需要20 min左右。而對于等效有限元模型，本文采用了Patran的PCL編程實(shí)現(xiàn)了有限元建模、等效材料計(jì)算、屬性定義和工況提交、顫振計(jì)算等所有步驟自動(dòng)化進(jìn)行，總運(yùn)行時(shí)間大概為4 min，顫振計(jì)算時(shí)間約為0.5 min。

圖8 有限元模型對比Fig.8 Comparison of two finite element models

表1 復(fù)合材料蒙皮和桁條鋪層順序和比例

2.2 固有頻率和振動(dòng)模態(tài)計(jì)算

應(yīng)用MSC.Nastran軟件對詳細(xì)有限元模型和等效剛度有限元模型的固有頻率和振動(dòng)模態(tài)進(jìn)行計(jì)算分析，其結(jié)果見表2所示。由表2可知，等效剛度模型的各階振動(dòng)模態(tài)與詳細(xì)模型的完全一致。兩種模型的各階頻率很接近，誤差均小于5%。

2.3 顫振分析

采用亞聲速偶極子法進(jìn)行顫振計(jì)算，求解方法使用p-k法[16]。在工程分析中，通常取機(jī)翼的前5階模態(tài)來做顫振分析。為了驗(yàn)證等效剛度方法的準(zhǔn)確性，分別從表2中取詳細(xì)模型和等效模型的前10階模態(tài)進(jìn)行顫振計(jì)算，計(jì)算獲得的速度-阻尼(V-g)圖和速度-頻率(V-ω)圖，分別見圖9和10所示。從圖9和10可以看到，兩種模型的速度-阻尼圖和速度-頻率圖的走勢一致。

表2 動(dòng)力學(xué)分析結(jié)果比較

圖9 詳細(xì)模型顫振分析結(jié)果Fig.9 Flutter analysis results of the detailed FE model

圖10 等效模型顫振分析結(jié)果Fig.10 Flutter analysis results of the equivalent FE model

根據(jù)詳細(xì)模型和等效模型計(jì)算獲得的顫振速度如表3所示。比較表3中的顫振計(jì)算結(jié)果，顫振速度誤差很小，僅為0.849%，而顫振頻率幾乎完全吻合。

表3 顫振計(jì)算結(jié)果比較

3 結(jié) 論

為簡化復(fù)合材料機(jī)翼有限元模型，本文研究了一種將加筋壁板結(jié)構(gòu)等效為具有相同力學(xué)屬性光板的等效剛度方法。通過對比詳細(xì)有限元模型和等效剛度有限元模型的計(jì)算結(jié)果，得出以下結(jié)論：

1)等效剛度方法大大減小了有限元網(wǎng)格規(guī)模和復(fù)雜程度，單元數(shù)量減至為詳細(xì)模型的1/9，且有利于機(jī)翼結(jié)構(gòu)有限元的參數(shù)化建模。

2)兩種模型計(jì)算出復(fù)合材料機(jī)翼的固有頻率和振動(dòng)模態(tài)一致，最大誤差小于5%。

3)兩種模型計(jì)算出的顫振速度和頻率非常一致，顫振速度的誤差不到1%。

4)等效剛度有限元模型具有簡單快速和準(zhǔn)確的特點(diǎn)，可用于復(fù)合材料機(jī)翼初步設(shè)計(jì)階段對顫振特性的評估。

[1] Zhao Yonghui, Hu Haiyan. Structural modeling and aeroelastic analysis of high-aspect ratio composite wings[J]. Chinese Journal of Aeronautics, 2005,18(1):25—30.

[2] Wright J R, Cooper J E. Introduction to Aircraft Aeroelasticity and Loads[M]. Chichester, England: John Wiley & Sons, Ltd, 2007:427—436.

[3] Giles G L. Equivalent plate analysis of aircraft wing box structures with general planform geometry[R]. NASA-TM-87697, 1986.

[4] Livne E, Navarro I. Nonlinear equivalent plate modeling of wing-box structures[J]. Journal of Aircraft, 1999,36(5)：851—865.

[5] 晏飛.復(fù)合材料翼面結(jié)構(gòu)綜合優(yōu)化設(shè)計(jì)技術(shù)研究[D].西安:西北工業(yè)大學(xué),2001.

Yan Fei. Study on the optimization technologies of composite wing structure[D]. Xi′an: Northwestern Polytechnical University, 2001．

[6] Cesnik C E S, Hodges D H, Patil M J. Aeroelastic analysis of composite wings[R]. AIAA 96-1469-CP. 1996:1 113—1 123.

[7] 李少華，楊智春，谷迎松，等.復(fù)合材料連翼的氣動(dòng)彈性剪裁研究[J].西北工業(yè)大學(xué)學(xué)報(bào),2008,26(3):292—296.

Li Shaohua, Yang Zhichun,Gu Yingsong, et al. A different aeroelastic tailoring of a composite joined-wing[J]. Journal of Northwestern Polytechnical University, 2008,26(3):292—296.

[8] 陳全龍,韓景龍，員海瑋.考慮發(fā)動(dòng)機(jī)推力影響的機(jī)翼顫振分析[J].振動(dòng)工程學(xué)報(bào)，2012,25(2):110—116.

Chen Quanlong, Han Jinglong, Yuan Haiwei. Flutter analysis of wings subjected to engine thrusts[J]. Journal of Vibration Engineering, 2012,25(2):110—116.

[9] Weck O De, Agte J, Sobieszczanski-Sobieski J, et al. State-of-the art and future trends in multidisciplinary design optimization[R]. AIAA 2007-1905, 2007.

[10]Bradley K R. A sizing methodology for the conceptual design of blended-wing-body transports[D]. George Washington: George Washington University, 2004.

[11]趙群，金海波,丁運(yùn)亮.加筋板總體失穩(wěn)分析的等效層合板模型[J].復(fù)合材料力學(xué)報(bào),2009,26(3):195—201.

Zhao Qun, Jin Haibo, Ding Yunliang. Equivalent laminates model for stiffened panel global buckling analysis[J]. Acta Materiae Compositae Sinica, 2009,26(3):195—201.

[12]Wenzel J, Sinapius M, Gabbert U. Primary structure mass estimation in early phases of aircraft development using the finite element method[J]. CEAS Aeronautical Journal, 2012,3(1):35—44.

[13]Jones R M. Mechanics of Composite Materials[M]. 2nd ed. Philadelphia: Taylor & Francis Inc, 1999:55—118.

[14]Timshenko S, Woinhowsky-Krieger S. Theory of Plates and Shells[M]. 2nd ed. New York: McGraw-Hill Book Company, 1959:36—45.

[15]李順林.復(fù)合材料工作手冊[M].北京:航空工業(yè)出版社,1988:334—335.

Li Shunlin. An Handbook of Composite Materials[M]. Beijing :Aviation Industry Press, 1988:334—335.

[16]Rodden W P, Harder R L, Bellinger E D. Aeroelastic addition to NASTRAN[R]． NASA CR-3094. Washington: NASA, 1979.

Flutter prediction of composite wing using finite element method with equivalent stiffness

OUYANGXing,YUXiong-qing,WANGYu

(Key Laboratory of Fundamental Science for National Defense-Advanced Design Technology of Flight Vehicle,Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China)

A method using the concept of equivalent stiffness was developed for the modal and flutter analysis of composite wing in preliminary design. The stiffness matrices calculated from different stringer-stiffened panels were applied to FEM shell properties of the same planform panels (i.e., unstiffened panels) that have the same mechanical properties as those of the original stringer-stiffened panels. The advantage of this method is that the complicated definition of composite material and modeling of stiffeners are avoided in the equivalent finite element model (EFEM) and the complicacy of finite element model is reduced, while the feature of structural layout of wing was still retained. A wing of the notional civil jet was used as an example to verify the method. Both the detailed finite element model (DFEM) and the EFEM of wings are created and analyzed for flutter prediction. Especially, the EFEM is automatically created by running scripts, and lots of time is saved from the process of modeling and analysis. The results indicate that the natural frequencies, mode shapes and flutter speed from the two models are in a very good agreement. The accuracy of the EFEM for the modal and flutter analysis of wings is verified. The EFEM provides an attractive method for wing flutter prediction in aircraft preliminary design phase due to its rapid modeling and satisfactory accuracy.

wing; flutter; composite; finite element; equivalent stiffness

2014-03-15；

2014-12-09

V211.41; V215.3+4

A

1004-4523(2015)03-0404-07

10.16385/j.cnki.issn.1004-4523.2015.03.009

歐陽星(1986—)，男，博士研究生。電話：(025)84892102；E-mail：ouyangxing@nuaa.edu.cn

余雄慶(1965—)，男，教授，博士生導(dǎo)師。電話：(025)84892102；E-mail：yxq@nuaa.edu.cn