張龍, 賈普榮, 王波, 徐斌

(1.西北工業(yè)大學 力學與土木建筑學院, 陜西 西安 710072; 2.西北工業(yè)大學 航空學院, 陜西 西安 710072)

風扇葉片是航空發(fā)動機、船用推進器以及風力發(fā)電機的重要組成部分[1-3]。傳統(tǒng)風扇葉片是由金屬材料制成，但對于發(fā)動機而言，減輕結(jié)構重量能有效提高發(fā)動機的性能，因此為滿足日益增長的這一需求，采用復合材料替換傳統(tǒng)金屬材料是一種有效且可行的途徑[4-5]。相比于金屬材料，復合材料具有質(zhì)量輕、比模量高、比強度高、耐腐蝕、耐疲勞以及抗振動等優(yōu)點。同時復合材料具有各向異性這一特點，在受載時會產(chǎn)生拉伸、剪切、彎曲和扭轉(zhuǎn)之間互相的耦合效應，這使得復合材料結(jié)構設計具有更多靈活性，但也使設計更加復雜[6]。

由于傳統(tǒng)金屬葉片模量高，受載時葉片變形較小，因此在設計過程中一般忽略葉片氣動外形的變化對結(jié)構性能的影響[7]。然而復合材料葉片的變形較大，并且不同鋪層方案的變形差異也較大，過大的變形會導致葉片不再滿足氣動彈性外形的設計，進而導致葉片結(jié)構效率下降。因此在鋪層設計時必須考慮變形對葉片氣動彈性結(jié)構效率的影響[8]。周邢銀等人對復合材料懸臂梁不同區(qū)域上的彎扭耦合效應對結(jié)構變形的影響進行了研究[9]。彭峰和Xiao等人分別對船用復合材料葉片和發(fā)動機葉片進行了研究，結(jié)果表明鋪層角度、鋪層厚度、鋪層順序?qū)秃喜牧先~片結(jié)構的性能有顯著影響[10-11]。Abdul等人通過不同鋪層順序的試驗與仿真研究表明鋪層順序?qū)θ~片結(jié)構變形的影響很大[12]。同時相比于彎曲變形，葉片的扭轉(zhuǎn)變形會改變?nèi)~片的攻角，而攻角的改變會顯著降低葉片結(jié)構的效率[13]。為了提高葉片結(jié)構效率，必須減小葉片氣動外形的改變，尤其是減小扭轉(zhuǎn)變形。因此可以利用復合材料對稱層合板結(jié)構中的彎扭耦合效應，通過調(diào)整鋪層角度和鋪層順序，使結(jié)構中由彎矩通過耦合效應產(chǎn)生的扭轉(zhuǎn)變形，與扭矩產(chǎn)生的扭轉(zhuǎn)變形相互抵消，減小葉片的整體扭轉(zhuǎn)變形，從而提高葉片氣動彈性結(jié)構的效率[14-15]，達到鋪層優(yōu)化設計的目的。

本文基于經(jīng)典層合板理論，提出了2個分析參數(shù)：剛度權值和載荷系數(shù)。通過設計彎扭耦合試驗，研究了2個分析參數(shù)對彎扭耦合效應的影響，并結(jié)合有限元模擬對試驗數(shù)據(jù)與理論結(jié)果進行了對比驗證。進一步采用剛度權值和載荷系數(shù)對彎扭耦合效應進行了定量分析，將求解滿足最小曲率目標函數(shù)時的最優(yōu)鋪層順序轉(zhuǎn)化為求解最優(yōu)剛度權值，得到隨載荷系數(shù)變化的最優(yōu)剛度權值曲線。在給定任意彎扭載荷系數(shù)下，通過對最優(yōu)剛度權值曲線的逆向計算，能夠直接得到滿足結(jié)構設計條件的最優(yōu)鋪層順序，為復合材料葉片的鋪層優(yōu)化設計提供了參考和依據(jù)。

1 剛度權值和載荷系數(shù)的定義

1.1 剛度權值

為了分析對稱層合板結(jié)構中相同鋪層角度所在鋪層位置對彎扭耦合效應的影響，基于經(jīng)典層合板理論(CLT)，將面內(nèi)剛度矩陣A和彎曲剛度矩陣D中相同的鋪層角度進行合并,剛度矩陣表述如下:

(1)

式中，H為層合板總厚度,r為鋪層角度種類的數(shù)量;i為鋪層角度的編號,Qi為對應角度的偏軸剛度矩陣,αi和δi分別為該角度的面內(nèi)和彎曲剛度權值,由公式(2)給出:

(2)

式中，mi為第i個角度的鋪層數(shù)量,hi, t和hi, t-1分別為該角度下第t層的上下表面的Z坐標值。根據(jù)公式(2)的定義,αi與鋪層含量意義相同,且αi和δi滿足如下關系

(3)

1.2 載荷系數(shù)

為了分析不同載荷比例對耦合效應的影響,采用單位化對載荷進行變量代換,如公式(4)所示:

(4)

由此分別引入了體現(xiàn)鋪層因素的分析變量:剛度權值α和δ,以及體現(xiàn)載荷因素的分析變量:比例系數(shù)ω。

2 試驗與有限元模擬分析研究

2.1 彎扭耦合試驗方案

圖1為彎扭耦合試驗加載方式示意圖。圖2為復合材料試件尺寸示意圖。如圖所示,通過在復合材料懸臂板加載端5個不同位置施加載荷分別為LN2,LN1,P0,LP1和LP2。實現(xiàn)扭矩由負到正的5種彎扭載荷比例,試驗采用Instron1196電子萬能試驗機加載,試驗速度為2 mm/min,當載荷達到60 N時,結(jié)束試驗,保存數(shù)據(jù)。層合板選用M40J/CYCOM970高模量碳纖維復合材料,單層厚度0.125 mm,其單層屬性見表1。

圖1 試驗加載方式示意圖

圖2 復合材料試件尺寸示意圖(單位: mm)

E11/GPaE22=E33/GPaν12=v13ν23G12=G13/GPaG23/GPa179.1111.320.2540.455.262.76

表2為4種試件鋪層順序和對應鋪層角度的δ。通過改變試件鋪層順序,實現(xiàn)δ的變化。在設計試件偏軸鋪層角度時,為提高試驗數(shù)據(jù)的可分析性,即在改變ω和δ時,提高曲率和扭率的變化量,應選擇耦合效應較強的偏軸角度,同時考慮偏軸角度鋪設的難易,偏軸角度以整數(shù)為宜,最終根據(jù)本文所用復合材料屬性,選擇30°作為偏軸鋪層角度。

表2 試件鋪層順序和彎曲剛度權值δ

2.2 有限元分析模型

由于試驗只能通過計算應變得到測點曲率,無法測得內(nèi)力矩,因此不能直接對CLT理論進行計算驗證。通過引入有限元(FEA)仿真計算,可以同時得到測點曲率和內(nèi)力矩。通過對比試驗與FEA的曲率數(shù)據(jù),驗證FEA的計算結(jié)果的正確性。再將試驗曲率通過CLT理論計算得到的內(nèi)力矩與FEA的內(nèi)力矩進行對比分析,并以此來驗證采用CLT理論對彎扭耦合效應進行分析計算是否可行。

本文采用ABAQUS商業(yè)有限元件對試驗進行仿真模擬。模型尺寸和加載方式與圖2、圖1一致,材料屬性和鋪層順序與表1、表2一致,采用C3D20單元,單元數(shù)為720 000。模型邊界條件及加載方式如表3所示,均與試驗條件盡可能保持一致。

表3 有限元模型邊界條件及加載方式

2.3 試驗與有限元結(jié)果對比分析

圖3為4種試件試驗和FEA的曲率κ-加載點曲線。結(jié)合表2中的權值與圖中曲線可以看出:

1) 對比4種試件在P0純彎曲加載時,隨δ0增大,|κx|和|κy|減小;隨|δ30-δ-30|增大,|κxy|增大;在δ30-δ-30>0的ABD型中κxy<0,反之在δ30-δ-30<0的C型中κxy>0。

2) 在δ30-δ-30>0的ABD型中,隨Mxy增大κx減小,反之在δ30-δ-30<0的C型中κx增大,隨|δ30-δ-30|增大,κx和κxy的斜率增大。

3) 對應5種不同彎扭載荷比例時,|κxy|最小值對應不同的鋪層順序。A型和B型|κxy|最小值對應在LP2加載點,其中A型|κxy|更小,C型和D型|κxy|最小值分別對應LN1和LP1加載點。

從以上分析可以發(fā)現(xiàn),①δ0控制抗彎曲性能;②|δ30-δ-30|的大小控制耦合效應的強弱,δ30-δ-30的正負控制耦合效應的正負。③在不同彎扭載荷比例下,最小|κxy|對應的鋪層也不相同。由此可將非連續(xù)的鋪層順序轉(zhuǎn)化為連續(xù)的剛度權值,對復合材料彎扭耦合效應進行定量分析和研究。

圖3 4種試件試驗和有限元的曲率κ-加載點曲線

圖4為4種試件試驗和FEA的力矩M-加載點曲線,其中試驗的合內(nèi)力矩M是由試驗的曲率數(shù)據(jù)通過CLT理論計算得到的。從圖中可以看出4種試件的受載情況基本一致。對4種試件的M求平均并計算得到平均ω,如表4所示。結(jié)合表4可以看出,在5種加載方式下:①Mx數(shù)值最大,且大小基本相同;②My數(shù)值很小,幾乎為0;③在P0加載點下,力矩Mxy幾乎為0,基本為純彎曲狀態(tài);加載點由LN2至LP2,Mxy由負變正。

圖4 4種試件試驗和有限元的力矩M-加載點曲線

Test/FEAωxωyωxyLN20.094/0.0920.987/0.9861.093/1.091 LN10.045/0.0460.995/0.9931.045/1.046 P00.002/0.0040.997/0.9961.000/1.000 LP10.050/0.0480.996/0.9940.950/0.952 LP20.096/0.0940.987/0.9860.905/0.906

為比較試驗和FEA 2組曲線的數(shù)據(jù)誤差,采用積分絕對值誤差準則(IAE)對數(shù)據(jù)進行分析,計算方法如公式(5)所示,式中DTest表示試驗數(shù)據(jù),DFEA表示有限元數(shù)據(jù),n為數(shù)據(jù)的數(shù)量。誤差計算結(jié)果如表5所示。

(5)

從表5中可以看出,試驗與FEA的κx和κy誤差較小,κxy誤差相對較大,這可能是由于κxy是通過應變轉(zhuǎn)軸公式計算得到,誤差被進一步放大。表5中試驗和FEA的Mx和Mxy誤差相對較小,My的誤差很大,但從圖4b)中可以看出My基本接近0,因此數(shù)據(jù)差異可能不大,但得到的誤差會較大。通過誤差分析可以看出,試驗與FEA的曲率數(shù)據(jù)基本一致,并且通過CLT理論計算得到的力矩與有限元仿真計算得到的力矩也相對一致。因此間接通過FEA模擬對CLT理論與試驗數(shù)據(jù)進行了驗證,可認為采用CLT理論對彎扭耦合效應的分析是可行的。

3 彎曲剛度權值的最優(yōu)路徑研究

3.1 問題描述

本文以工程中常用的由0°,90°,±45°組成的對稱層合板為例,以優(yōu)先滿足min(|κxy|),其次滿足min(|κx|)為目標函數(shù),進一步定量研究在不同ω下,δ的比例關系。由公式(1)可得層合板彎曲剛度矩陣如下

(6)

式中,Qi為i=0°,90°,±45°時單層板的偏軸剛度矩陣。根據(jù)試驗分析結(jié)果,取My=0對問題進行簡化,由公式(4)可得公式(7),獨立載荷因子剩下ωx,簡記為ω。

(7)

當鋪層角度確定后,κx和κxy分別是關于δ和ω的函數(shù)。根據(jù)對稱性分析可知,將ω最優(yōu)鋪層中的±θ°層的角度互換可得-ω的最優(yōu)鋪層。因此在公式(7)中只需討論ω∈[0, 1]的情況。

為滿足最小的|κx|和|κxy|條件,須滿足δ-45取值最小,根據(jù)層合板鋪層含量的約束,假設滿足α0≥0.6,其余鋪層角度含量α≥0.1時,根據(jù)公式(3)中條件,4個δ中只有3個相互獨立,此時權值可行域滿足如下關系:

(8)

根據(jù)公式(8)可知,δ的可行域如圖5所示,其中可行域邊界點的權值坐標如表6所示。

圖5 權值δ的可行域

可行域邊界點權值δ坐標δ0δ90F10.9730.025 F20.9730.001 F30.5100.001 F40.2160.295 F50.2160.488 F60.5100.488

3.2 權值最優(yōu)路徑的求解計算

本文采用Maple數(shù)學工程計算軟件對所需求解的方程進行計算。令H3/12=1，對公式(6)進行簡化，由此可得κx和κxy分別是關于δ和ω的函數(shù)κx(δ0,δ90,ω)和κxy(δ0,δ90,ω)。令κxy=0可得不同ω下δ關于κxy=0的等高曲線方程,如公式(9)所示

(9)

根據(jù)公式(9)作圖可得等高曲線δ90-δ0,如圖6所示。從圖中可以看出,ω存在一個臨界值ω1:當ω∈[0,ω1]時,存在等高曲線δ90-δ0使得κxy=0,當ω=ω1時,等高曲線δ90-δ0與可行域相交于點F4;將點F4的權值帶入公式(9)可得ω1=0.144;當ω∈[ω1, 1]時,κxy恒大于0。

圖6 不同ω下可行域上κxy=0的δ等高曲線δ90-δ0

1) 當ω∈[0,ω1]時:

求解公式(9)分別可得δ90和δ0的函數(shù)表達式:

δ90|κxy=0=f(δ0,ω),ω∈[0,ω1]

(10)

δ0|κxy=0=g(δ90,ω),ω∈[0,ω1]

(11)

將公式(10)帶入κx(δ0,δ90,ω)可得:

κx|κxy=0=κx(δ0,f(δ0,ω),ω),ω∈[0,ω1]

(12)

根據(jù)公式(12)可得在滿足κxy=0的條件下,κx-δ0關于ω的曲線族,如圖7所示。結(jié)合圖6與圖7可知,當ω∈[0,ω1]時,部分等高曲線關于δ0分為上下2個部分,部分等高曲線只有上半部分。且δ0取最大值時κx取得最小值。

圖7 不同ω下κxy=0時的κx-δ0曲線

為求得等高曲線上δ0的最大值,對圖6進行分析可以看出,受可行域的約束,隨ω增大,等高曲線上δ0的最大值被分為3個部分,分別在可行域邊界F1F2F3,可行域內(nèi)和可行域邊界F3F4上。令方程(11)關于δ90求導為0,如公式(13)所示:

(13)

將方程(13)帶入方程(11)并分別與F2F3和F3F4曲線方程聯(lián)立求解,可得對應臨界的ω分別為0.002和0.008。

分別對ω∈[0, 0.008],[0.008, 0.1]和[0.1,ω1]3個區(qū)間求解δ0的最大值,可得δ0(ω)函數(shù)方程。進一步帶入公式(11)可得δ90(ω)函數(shù)方程,再根據(jù)公式(3)可得δ45(ω)函數(shù)方程。通過以上計算可得ω∈[0,ω1]時,滿足κxy=0,且κx最小條件時,δ(ω)的最優(yōu)路徑。

2) 當∈[ω1, 1]時:

κxy恒大于0,根據(jù)圖6可以看出,當ω=ω1時,κxy最小值在點F4。過點F4的κxy等高曲線方程如公式(14)所示

κxy(δ0,δ90,ω)=κxy|δ=F4,ω∈[ω1,1]

(14)

根據(jù)公式(14)可得不同ω下過邊界點F4的κxy等高曲線族,如圖8所示。從圖中可以看出,存在一個臨界值ω2,當ω∈[ω1,ω2]時,κxy最小值取在點F4處,當ω∈[ω2, 1]時,κxy最小值由點F4沿可行域邊界向F3移動。

圖8 不同ω下過邊界點F4的κxy等高曲線族

將公式(8)中可行域邊界F3F4帶入κxy,并令κxy關于δ0求導為0,如公式(15)所示

(15)

將點F4的δ0帶入公式(15)求解可得臨界值ω2=0.843。對公式(15)求解可得ω∈[ω2, 1]范圍內(nèi)κxy最小值的δ0(ω)函數(shù)方程,帶入公式(8)中可行域邊界F3F4可得δ90(ω)函數(shù)方程,進一步可得δ45(ω)函數(shù)方程。

3.3 權值最優(yōu)路徑計算結(jié)果與分析

通過以上分析可得,在可行域約束下的δ最優(yōu)路徑,如圖9所示。

圖9 不同比例系數(shù)ω下可行域內(nèi)權值δ的最優(yōu)路徑

圖中虛線為可行域邊界帶有形狀符號的曲線為對應不同ω時κxy=0的等高曲線,實曲線H1F2H2H3F4H4為整體δ最優(yōu)路徑δ90-δ0曲線。在該最優(yōu)路徑上,能夠優(yōu)先保證|κxy|最小,其次保證|κx|最小。由于受可行域邊界的約束,δ最優(yōu)路徑被分為6個部分。當ω∈[0.008, 0.1]時,最優(yōu)路徑在可行域內(nèi)的曲線H2H3上,當ω不在該范圍時,最優(yōu)路徑在可行域的邊界上。δ最優(yōu)路徑的拐點坐標與對應的ω如表7所示。

通過δ最優(yōu)路徑,可得δ(ω)曲線,如圖10所示。結(jié)合表7與圖10可以看出:

表7 權值δ最優(yōu)路徑拐點坐標與對應比例系數(shù)ω

圖10 可行域上最優(yōu)權值路徑的δ-ω曲線

1) 當ω∈[0, 0.002]時,隨ω增大,δ0不變,δ90減小,δ45增大。δ最優(yōu)路徑為可行域邊界H1F1。

2) 當ω∈[0.002, 0.008]時,隨ω增大,δ0減小,δ90不變,δ45增大。δ最優(yōu)路徑為可行域邊界F1H2。

3) 當ω∈[0.008,ω1]時,隨ω增大δ0減小,δ90和δ45增大,且δ45增幅大于δ90。以ω=0.1為分界,δ最優(yōu)路徑分為2個部分:可行域內(nèi)的H2H3和可行域邊界上的H3F4。

4) 當ω∈[ω1,ω2]時,隨ω增大,δ0,δ90和δ45保持不變,δ最優(yōu)路徑停留在可行域邊界點F4上。

5) 當ω∈[ω2, 1]時,隨ω增大,δ0增大,δ90減

小,δ45保持不變,δ最優(yōu)路徑由點F4沿可行域邊界向H4移動。

通過最優(yōu)路徑δ(ω),可以計算得到在最優(yōu)路徑上κ-ω曲線,如圖11所示。從圖中可以看出,在最優(yōu)路徑δ(ω)上:

1) 當ω∈[0,ω1]時,通過改變δ可以保證,隨ω增大,κxy=0,κx取得最小值。

2) 當ω∈[ω1, 1]時,隨ω增大,κxy開始增大,κx開始減小。

圖11 可行域上最優(yōu)權值路徑的κ-ω曲線

3.4 權值最優(yōu)路徑的算例驗證與討論

本文以16層對稱層合板為例,對最優(yōu)路徑δ(ω)進行驗證。以試驗中P0,LP1和LP2 3種情況的載荷因子為基礎,分別選取ω=0,0.05,0.1以及ω1進行分析計算。根據(jù)鋪層含量的約束關系可知4種角度的鋪層數(shù)分別為0°鋪層有10層,90°和±45°鋪層各2層。

由于實際情況中δ的取值是在三維空間中的離散點(δ0,δ90,δ45)上,不一定剛好取在理論最優(yōu)δ點上,因此實際最優(yōu)δ點是離散的三維空間點中包圍理論最優(yōu)δ點最小六面體8個頂點中的1個。通過計算實際與理論δ最優(yōu)點的距離可以計算得到六面體8個頂點的權值δ,從8種情況中選取滿足最優(yōu)條件的最優(yōu)實際δ,即為所求最優(yōu)實際最優(yōu)解。

表8 對應比例系數(shù)ω下權值δ,曲率κ的理論與實際最優(yōu)解及鋪層順序

表8為采用Maple編程計算得到對應ω的δ和κ的理論、實際最優(yōu)解,及對應的實際鋪層順序。從表中可以看出,最優(yōu)解的理論與實際的δ和κ非常接近。①當ω=0純彎曲時,0°層在最外層,±45°層在最內(nèi)層;②隨ω增大,45°層最先移動到外層,90°和-45°層相繼向外層移動。當ω=0.05時,δ45的實際值為0.330,已經(jīng)達到實際的最大值,但仍小于δ45的理論值0.341,隨ω增大,二者差異增大,κxy也開始增大。這是由于實際鋪層的非連續(xù)性以及受到鋪層含量的約束導致δ的實際最優(yōu)值不能完全與理論最優(yōu)值重合。

4 結(jié) 論

本文基于經(jīng)典層合板理論提出了2個分析參數(shù)：①將非量化的鋪層順序轉(zhuǎn)化為可量化分析的剛度權值；②將合內(nèi)力矩轉(zhuǎn)化為歸一化的載荷系數(shù)。并通過彎扭耦合試驗與有限元模擬討論了2個分析參數(shù)對彎扭耦合效應的影響。

以彎曲剛度權值的可行域為約束條件,結(jié)構曲率最小為目標函數(shù),對含有0°,90°,±45°的層合板進行了優(yōu)化分析,得到了彎曲剛度權值的最優(yōu)路徑。受到可行域的限制,彎曲剛度權值最優(yōu)路徑被分成6個部分。①在前4個部分中,即ω∈[0,ω1]時,存在優(yōu)化權值比例使得曲率κxy=0,|κx|取得最小值。②當ω∈[ω1, 1]時,受到可行域的限制,在彎曲剛度權值最優(yōu)路徑上,隨ω增大,曲率κxy開始增大,κx開始減小。

根據(jù)本文提出的剛度權值最優(yōu)路徑,在給定任意的彎扭載荷條件下,可以得到對應的理論最優(yōu)剛度權值,經(jīng)過逆向計算可以得到實際最優(yōu)鋪層順序。通過該方法能夠快速有效地找到滿足曲率設計目標的最優(yōu)鋪層順序,能夠?qū)秃喜牧先~片彎扭耦合效應的優(yōu)化設計提供一定的參考和依據(jù)。在本文的分析研究中,將葉片模型簡化假設為層合板,忽略曲面葉形結(jié)構的影響,因此在后續(xù)的工作中需要進一步研究在曲面葉形結(jié)構上彎扭耦合效應的優(yōu)化分析。