曾志波,姚志崇,王瑋波,劉潤聞,韓用波,洪方文
(中國船舶科學(xué)研究中心,江蘇 無錫 214082)
先進(jìn)復(fù)合材料于20世紀(jì)60年代問世。復(fù)合材料在航空、航天結(jié)構(gòu)上的應(yīng)用帶來了突出的減重效果和綜合性能的顯著提高,使其成為軍用飛機(jī)四大結(jié)構(gòu)材料之一[1]。在船舶方面,自20世紀(jì)80年代以來,復(fù)合材料在水面艦船和水下潛艇方面的應(yīng)用開始并迅速增長,在艦艇上層建筑、天線系統(tǒng)、螺旋槳、槳軸、舵、泵、管道及機(jī)座等方面得到了大范圍的應(yīng)用[2]。
盡管復(fù)合材料螺旋槳在商用和軍用船艇上得到了相當(dāng)應(yīng)用,然而由于缺乏設(shè)計(jì)方法、系統(tǒng)的經(jīng)驗(yàn)數(shù)據(jù)庫和可靠的計(jì)算工具,大多數(shù)螺旋槳仍然完全套用了金屬螺旋槳的設(shè)計(jì)與計(jì)算方法。最早的復(fù)合材料螺旋槳數(shù)值計(jì)算技術(shù)出現(xiàn)于1991年,Lin等[3-4]采用渦格法(VLM)計(jì)算水下無空泡螺旋槳的定常水動(dòng)力性能,提取槳葉外載荷后,采用有限元(FEM)進(jìn)行了槳葉應(yīng)力和變形分析;1994年,Kane和Dow[5]也進(jìn)行了類似的計(jì)算,結(jié)果表明:與傳統(tǒng)金屬槳相比,復(fù)合材料螺旋槳具有5~10倍的葉梢變形,但是以上計(jì)算結(jié)果和實(shí)船試航相差較大[2]。在早期的計(jì)算中,并未進(jìn)行流固耦合迭代,因而沒有涉及流體、結(jié)構(gòu)的相互作用。1996年,Lin等[6]發(fā)展了考慮流固耦合的3D FEM/VLM方法,并用于螺旋槳水動(dòng)力性能分析,隨后,Lin等進(jìn)一步完善計(jì)算方法,2005年增加了應(yīng)力評(píng)估和疲勞標(biāo)準(zhǔn)[7]。Young等[8]和Chen等[9]于2006年進(jìn)行了復(fù)合材料螺旋槳設(shè)計(jì)與分析,在他們的研究中較全面地考慮流固耦合特性以及非均勻伴流的影響。Young等[10-11]采用面元法和有限元法開展了復(fù)合材料螺旋槳流固耦合特性數(shù)值計(jì)算分析,并與試驗(yàn)結(jié)果做了比較,取得了較好的精度。近兩年來,Motley等[12-13]采用流固耦合方法針對(duì)復(fù)合材料節(jié)能性能做了分析研究。與國際水平相比,我國的復(fù)合材料螺旋槳研究屬起步階段,研究水平停留在小型螺旋槳的單方面的性能或者復(fù)合材料的性能上[14-15]。
本文基于面元法和有限元法,采取流固耦合中的弱耦合方法,開展了復(fù)合材料螺旋槳流固耦合分析方法研究。首先探討了復(fù)合材料槳葉有限元建模以及水動(dòng)力外載荷和結(jié)構(gòu)變形的流固耦合交界面數(shù)據(jù)傳遞問題;然后在此基礎(chǔ)上集成為一項(xiàng)復(fù)合材料螺旋槳分析技術(shù),測試了其收斂性;最后進(jìn)行了比較計(jì)算驗(yàn)證,為復(fù)合材料螺旋槳性能分析提供了必要的工具。
本文基于流固耦合中的弱耦合處理方法進(jìn)行復(fù)合材料螺旋槳力學(xué)性能分析,其中弱耦合方法的基本思路是在每個(gè)時(shí)間步內(nèi)依次對(duì)流體力學(xué)方程和結(jié)構(gòu)力學(xué)方程求解,通過耦合交界面交換固體域和流體域的計(jì)算結(jié)果數(shù)據(jù),反復(fù)交替計(jì)算從而實(shí)現(xiàn)耦合求解,該方法是研究氣固耦合(氣動(dòng)彈性)問題的主流方法[16]。在弱耦合數(shù)值計(jì)算中,采用面元法進(jìn)行流體域求解,采用有限元法進(jìn)行固體域求解。對(duì)于復(fù)合材料螺旋槳的弱耦合分析流程包括:
1)每一時(shí)間步對(duì)槳葉水動(dòng)力載荷的面元法計(jì)算;
2)水動(dòng)力載荷向交界面中有限元計(jì)算網(wǎng)格的傳遞;
3)每一時(shí)間步對(duì)復(fù)合材料槳葉位移(變形)等信息的有限元計(jì)算;
4)位移(變形)等信息向交界面中面元法計(jì)算網(wǎng)格的傳遞。
本文采用文獻(xiàn)[17]中的螺旋槳面元法和文獻(xiàn)[18-19]建立的螺旋槳結(jié)構(gòu)力學(xué)有限元分析方法建立復(fù)合材料螺旋槳弱耦合分析方法。
在ANSYS結(jié)構(gòu)有限元分析中,選取8節(jié)點(diǎn)三維層狀結(jié)構(gòu)實(shí)體單元SOLID46對(duì)復(fù)合材料槳葉進(jìn)行網(wǎng)格劃分,復(fù)合材料螺旋槳有限元模型如圖1所示。
螺旋槳采用正交各項(xiàng)異性碳纖維復(fù)合材料,材料鋪層結(jié)構(gòu)如圖2所示,左圖顯示的鋪層分布是[30/-30/90/-90/30]s對(duì)稱形式,右圖為根部鋪層結(jié)構(gòu)放大圖,共10層,各層之間是等厚度堆棧而成。
圖1 復(fù)合材料螺旋槳有限元模型Fig.1 The FE model of composite propellers
圖2 復(fù)合材料鋪層結(jié)構(gòu)(左:鋪層分布;右:根部鋪層結(jié)構(gòu))Fig.2 The layer structure of composite materials(Left:Layer distribution;right:Layer structure in the root)
面元法計(jì)算外載荷分布結(jié)果通常采取壓力分布表示,如圖3所示,其中壓力系數(shù)CP定義如下:
式中:P、P0、ρ和 V 分別為槳葉表面壓力(Pa)、參考?jí)毫Γ≒a)、水密度(kg/m3)和槳前方來流水速(m/s)。
采取插值的方式將面元法計(jì)算的水動(dòng)力外載荷從面元法網(wǎng)格傳遞到有限元計(jì)算網(wǎng)格節(jié)點(diǎn)上。由于ANSYS復(fù)合材料結(jié)構(gòu)分析模塊中有限元外載荷由各節(jié)點(diǎn)的節(jié)點(diǎn)力構(gòu)成,如圖4所示,其中根部施加固定約束。因此,本文先在槳葉面元網(wǎng)格中由壓力分布結(jié)合各面元的面積和法向矢量計(jì)算各面元上的三維作用力,然后進(jìn)行三維作用力的傳遞。
任一面元(圖5所示面元ABCD)法向方向的作用力計(jì)算公式如下:
上式中,cosα,cosβ,cosγ為方向余弦,采用螺旋槳面元法中的雙曲面元方向[17]。
圖3 面元法計(jì)算壓力分布Fig.3 Pressure distribution from PM
圖4 有限元節(jié)點(diǎn)力分布Fig.4 Force distribution for FEM
圖5 面元ABCDFig.5 A panel ABCD
將計(jì)算的各面元三維作用力(圖6所示)根據(jù)面元網(wǎng)格節(jié)點(diǎn)和有限元網(wǎng)格節(jié)點(diǎn)之間的位置關(guān)系(最小距離原則)傳遞到有限元網(wǎng)格節(jié)點(diǎn)上(圖4所示),其中葉背和葉面分別獨(dú)立進(jìn)行。傳遞的三維節(jié)點(diǎn)力作為有限元計(jì)算邊界條件。
圖6 計(jì)算節(jié)點(diǎn)力分布(左:X方向力;中:Y方向力;右:Z方向力)Fig.6 The distribution of calculated forces(Left:X direction;middle:Y direction;right:Z direction)
在流固耦合分析中變形后的槳葉需重新進(jìn)行外載荷計(jì)算,因此有限元中計(jì)算的復(fù)合材料槳葉各節(jié)點(diǎn)三向位移需傳遞給面元網(wǎng)格。本文采用如下插值方法[20-22]:
對(duì)于三維槳葉網(wǎng)格坐標(biāo)(X,Y,Z),其中:X為軸向方向,向后為正,Y為槳葉參考線方向,向上為正,首先將槳葉投影到(X’,Y)平面,其中X’和X軸夾角選取0.7R半徑螺距角,如圖7所示,使得位移插值量的空間分布在(X’,Y)平面變得平坦,從而提高插值精度;然后在(X’,Y)平面中進(jìn)行有限元網(wǎng)格點(diǎn)F向目標(biāo)面元網(wǎng)格點(diǎn)P位移量的插值工作,步驟如下:
1) 在(X’,Y)平面中選出離目標(biāo)P點(diǎn)最近的15個(gè)有限元網(wǎng)格源點(diǎn)Fi)(i=1,2…,15),如圖8所示;
2)對(duì)于源點(diǎn) Fi)和位移量 Di(i=1,2…,15)進(jìn)行二元三次多項(xiàng)式擬合,采用最小二乘法求得三次多項(xiàng)式系數(shù)ai(i=1,2…,10),得到擬合公式如下:
4)依次進(jìn)行面元網(wǎng)格點(diǎn)的插值。
圖7 原坐標(biāo)(X)與投影坐標(biāo)(X’)Fig.7 Original coodinate(X)and projective coodiante(X’)
圖8 插值點(diǎn)(P)與源點(diǎn)(F)Fig.8 Interpolated point(P)and source points(F)
圖9為變形量插值前后的比較,其中左圖是有限元計(jì)算的總變形量分布,右圖是插值到面元法網(wǎng)格上的總變形量;從比較圖可以看出,變形量實(shí)現(xiàn)了很好的傳遞。
圖9 變形量插值前后比較(左圖:插值前,單位m;右圖:插值后,單位mm)Fig.9 The comparison of deformation(Left:Before interplolation,m;right:After interplolation,mm)
由于槳葉形狀復(fù)雜,為了避免新槳葉面元網(wǎng)格產(chǎn)生畸變現(xiàn)象[22],面元網(wǎng)格有時(shí)必須針對(duì)新的槳葉幾何進(jìn)行重新劃分;另一方面,復(fù)合材料槳葉經(jīng)過流固耦合作用后,需從變形后的槳葉模型中提取螺距、拱度等特征參數(shù)進(jìn)行分析。以上兩個(gè)方面需要開展槳葉模型的參數(shù)反解工作。
槳葉參數(shù)反解按以下兩個(gè)步驟進(jìn)行:首先從槳葉三維模型中提取各半徑剖面葉背和葉面的三維坐標(biāo),可在UG建模軟件中實(shí)現(xiàn);然后由各剖面的三維坐標(biāo)計(jì)算槳葉參數(shù),求解螺距、弦長、厚度、側(cè)斜和拱度等構(gòu)成新槳葉參數(shù)型值表。根據(jù)以上步驟反解的槳葉參數(shù)可達(dá)到較高的精度,如圖10所示,采用反解的參數(shù)建模與原模型的比較,結(jié)果十分吻合。
圖10 反解槳葉模型與原槳葉模型比較Fig.10 The comparison between reversely extracted blade model and original blade model
2.5.1 分析流程
復(fù)合材料螺旋槳流固耦合分析流程如下:
1)對(duì)給定的復(fù)合材料螺旋槳幾何型值,采用UG建模,并輸入到ANSYS有限元分析軟件,進(jìn)行網(wǎng)格劃分;
2)對(duì)給定運(yùn)行工況,采用面元法進(jìn)行槳葉網(wǎng)格劃分,進(jìn)而離散面元法方程,求解槳葉表面壓力分布,計(jì)算各面元三維作用力;
3)將面元三維作用力傳遞到(1)中有限元?jiǎng)澐值木W(wǎng)格節(jié)點(diǎn)上,進(jìn)行復(fù)合材料屬性及鋪層結(jié)構(gòu)設(shè)置,對(duì)結(jié)構(gòu)力學(xué)分析方程進(jìn)行有限元求解,得到槳葉位移。
4)將有限元計(jì)算的各網(wǎng)格點(diǎn)位移插值到原始槳葉面元法中的面元網(wǎng)格節(jié)點(diǎn)上,得到變形后的槳葉幾何,提取新的槳葉幾何型值;
5)判斷新的槳葉幾何是否收斂,如果收斂,則結(jié)束;如果不收斂,則將新的槳葉幾何帶回(2),進(jìn)行(2)~(5)的迭代計(jì)算,直至收斂。
上述分析流程見圖11。
2.5.2 收斂性
螺旋槳根據(jù)側(cè)斜方向不同可分為前側(cè)斜螺旋槳和后側(cè)斜螺旋槳,圖12給出了從槳后向槳前方向看的兩種側(cè)斜分布螺旋槳槳葉模型。目前,船舶螺旋槳大多采用后側(cè)斜螺旋槳,對(duì)于前側(cè)斜螺旋槳也有學(xué)者研究其對(duì)梢渦空泡起始性能的影響[23]。
復(fù)合材料螺旋槳采取不同的側(cè)斜分布,呈現(xiàn)不同的流固耦合收斂過程,如圖13所示。圖13左圖表示后側(cè)斜復(fù)合材料螺旋槳的收斂過程,圖中螺旋槳推力系數(shù)和扭矩系數(shù)以及0.7R螺距比呈現(xiàn)交替收斂形式,最終收斂槳葉螺距變小,相應(yīng)的載荷變小。而對(duì)于前側(cè)斜復(fù)合材料螺旋槳的收斂過程,呈現(xiàn)單調(diào)收斂過程,如圖13右圖所示,最終收斂槳葉螺距變大,相應(yīng)的載荷變大。
圖11 復(fù)合材料螺旋槳流固耦合分析流程Fig.11 The flow chart of FSI analysis for composite propellers
圖12 不同側(cè)斜螺旋槳(左:后側(cè)斜;右:前側(cè)斜)Fig.12 Different skewed propellers(Left:Backward skew;right:Forward skew)
圖13 復(fù)合材料螺旋槳流固耦合收斂過程(左:后側(cè)斜;右:前側(cè)斜)Fig.13 The convergence process of FSI of composite propellers
采取本文建立的流固耦合分析技術(shù),針對(duì)文獻(xiàn)[12]中5474槳進(jìn)行了比較計(jì)算驗(yàn)證。在比較計(jì)算中采用相同的Hexcel IM7-8552碳纖維復(fù)合材料、鋪層方式,材料屬性如表1所示。
表1 Hexcel IM7-8552碳纖維復(fù)合材料屬性Tab.1 Material properties of Hexcel IM7-8552
圖14給出了復(fù)合材料螺旋槳變形前后面元法計(jì)算的推力系數(shù)和扭矩系數(shù),其中流固耦合迭代計(jì)算點(diǎn):進(jìn)速系數(shù)為J=0.66,轉(zhuǎn)速n=780 rpm。在圖14中也給出了文獻(xiàn)[12]的計(jì)算結(jié)果,本文計(jì)算的推力系數(shù)和扭矩系數(shù)與文獻(xiàn)[12]結(jié)果的差別在4%以內(nèi)。圖15為變形前后的螺距角分布比較,變形后從根部到梢部螺距角變形量逐漸增大,變形后梢部螺距角計(jì)算結(jié)果與圖中文獻(xiàn)[12]計(jì)算結(jié)果吻合良好。
圖14 5474復(fù)合材料螺旋槳水動(dòng)力比較Fig.14 The comparison of hydrodynamic performance of composite propeller 5474
圖15 5474復(fù)合材料螺旋槳螺距角分布比較Fig.15 The comparison of pitch angle distribution of composite propeller 5474
本文研究了復(fù)合材料螺旋槳流固耦合(弱耦合)分析方法。重點(diǎn)探討了復(fù)合材料槳葉有限元建模以及流固耦合中載荷及幾何變形量的傳遞問題,采用的方法對(duì)復(fù)合材料螺旋槳具有較好的適用性。集成的復(fù)合材料螺旋槳流固耦合分析技術(shù)對(duì)不同側(cè)斜螺旋槳進(jìn)行的迭代計(jì)算均可獲得收斂解,針對(duì)5474復(fù)合材料螺旋槳的計(jì)算結(jié)果與文獻(xiàn)的吻合較好。后續(xù)的工作是開展復(fù)合材料螺旋槳性能研究,剖析復(fù)合材料螺旋槳流固耦合變形規(guī)律。
[1]楊乃賓,章怡寧.復(fù)合材料飛機(jī)結(jié)構(gòu)設(shè)計(jì)[M].北京:航空工業(yè)出版社,2004.
[2]Mouritz A P,Gellert E.Review of advanced composite structures for naval ships and submarines[J].Composites Structures,2001,53:21-41.
[3]Lin G.Comparative stress-deflection analyses of a thick-shell composite propeller blade[R].David Taylor Research Center,DTRC/SHD-1373-01,1991a.
[4]Lin G.Three-dimensional stress analyses of a fiber-reinforced composite thruster blade[C]//Symposium on Propellers/Shafting Society of Naval Architects and Marine Engineers.Virginia Beach,VA,USA,1991b.
[5]Kane C,Dow R.Marine propulsors design in fibre reinforced plastics[J].J Defence Sci,1994,4:301-308.
[6]Lin H,Lin J.Nonlinear hydroelastic behavior of propellers using a finite element method and lifting surface theory[J].Journal of Marine Science and Technology,1996,1(2):114-124.
[7]Lin H,Lin J.Strength evaluations of a composite marine propeller blade[J].Journal of Reinforced Plastics and Composites,2005,17:1791-1807.
[8]Young Y L,Michael T J,Seaver M,Trickey S T.Numerical and experimental investigations of composite marine propellers[C]//26th Symposium on Naval Hydrodynamics.Rome,Italy,2006.
[9]BenjaminY-H Chen,Stephen K Neely,et al.Design,fabrication and testing of pitching-adapting(flexible)composite propellers[C].The SNAME,Propellers/Shafting Symposium,2006.
[10]Young Y L.Time-dependent hydroelastic analysis of cavitating propulsors[J].Journal of Fluids and Structures,2007,23(2):269-295.
[11]Young Y L.Fluid-structure interaction analysis of flexible composite marine propellers[J].Journal of Fluids and Structures,2008,24:799-818.
[12]Motley M R,Liu Z,Young Y L.Utilizing fluid structure interactions to improve energy efficiency of composite marine propellers in spatially varying wake[J].Composite Structures,2009,90(3):304-313.
[13]Motley M R,Young Y L.Performance-based design of adaptive composite marine propellers[C]//28th Symposium on Naval Hydrodynamics.Pasadena,California,USA,2010.
[14]楊傳勇.小型船舶復(fù)合材料船用螺旋槳的設(shè)計(jì)與分析[D].大連:大連海事大學(xué),2005.
[15]洪 毅.復(fù)合材料船用螺旋槳結(jié)構(gòu)設(shè)計(jì)研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2006.
[16]葉正寅,張偉偉,史愛明等.流固耦合力學(xué)基礎(chǔ)及其應(yīng)用[M].哈爾濱:哈爾濱工業(yè)大學(xué)出版社,2010.
[17]董世湯,唐登海,周偉新.CSSRC的螺旋槳定常面元法[J].船舶力學(xué),2005,9(5):48-62.Dong Shitang,Tang Denghai,Zhou Weixin.Panel method of CSSRC for propeller in steady flows[J].Journal of Ship Mechanics,2005,9(5):46-60.
[18]姚志崇,洪方文,丁恩寶.基于work bench的螺旋槳流固耦強(qiáng)度校核數(shù)值模擬[C].2007ANSYS-CFD年會(huì)文集,2007.
[19]姚志崇,曾志波.復(fù)合材料船用螺旋槳流固耦合性能分析[C].2010全國固體力學(xué)大會(huì),華中科技大學(xué),2010.
[20]李立州,王婧超,呂震宙,岳珠峰.學(xué)科間載荷參數(shù)空間傳遞方法[J].航空動(dòng)力學(xué)報(bào),2007(07):32-36.
[21]虞跨海,岳珠峰.渦輪冷卻葉片參數(shù)化建模及多學(xué)科設(shè)計(jì)優(yōu)化[J].航空動(dòng)力學(xué)報(bào),2007(08):144-149.
[22]岳珠峰等.航空發(fā)動(dòng)機(jī)渦輪葉片多學(xué)科設(shè)計(jì)優(yōu)化[M].北京:科學(xué)出版社,2007.
[23]Kuiper G,Van Terwisga T J C,Zondervan G J,Jessup S D,Krikke E M.Cavitation inception tests on a systematic series of two-bladed propellers[C]//26th Symposium on Naval Hydrodynamics.Rome,Italy,2006.