曾志波，姚志崇，王瑋波，劉潤聞，韓用波，洪方文

（中國船舶科學研究中心，江蘇 無錫 214082）

1 引 言

先進復合材料于20世紀60年代問世。復合材料在航空、航天結(jié)構(gòu)上的應用帶來了突出的減重效果和綜合性能的顯著提高，使其成為軍用飛機四大結(jié)構(gòu)材料之一[1]。在船舶方面，自20世紀80年代以來，復合材料在水面艦船和水下潛艇方面的應用開始并迅速增長，在艦艇上層建筑、天線系統(tǒng)、螺旋槳、槳軸、舵、泵、管道及機座等方面得到了大范圍的應用[2]。

盡管復合材料螺旋槳在商用和軍用船艇上得到了相當應用，然而由于缺乏設計方法、系統(tǒng)的經(jīng)驗數(shù)據(jù)庫和可靠的計算工具，大多數(shù)螺旋槳仍然完全套用了金屬螺旋槳的設計與計算方法。最早的復合材料螺旋槳數(shù)值計算技術(shù)出現(xiàn)于1991年，Lin等[3-4]采用渦格法（VLM）計算水下無空泡螺旋槳的定常水動力性能，提取槳葉外載荷后，采用有限元（FEM）進行了槳葉應力和變形分析；1994年，Kane和Dow[5]也進行了類似的計算，結(jié)果表明:與傳統(tǒng)金屬槳相比，復合材料螺旋槳具有5～10倍的葉梢變形，但是以上計算結(jié)果和實船試航相差較大[2]。在早期的計算中，并未進行流固耦合迭代，因而沒有涉及流體、結(jié)構(gòu)的相互作用。1996年，Lin等[6]發(fā)展了考慮流固耦合的3D FEM/VLM方法，并用于螺旋槳水動力性能分析，隨后，Lin等進一步完善計算方法，2005年增加了應力評估和疲勞標準[7]。Young等[8]和Chen等[9]于2006年進行了復合材料螺旋槳設計與分析，在他們的研究中較全面地考慮流固耦合特性以及非均勻伴流的影響。Young等[10-11]采用面元法和有限元法開展了復合材料螺旋槳流固耦合特性數(shù)值計算分析,并與試驗結(jié)果做了比較，取得了較好的精度。近兩年來，Motley等[12-13]采用流固耦合方法針對復合材料節(jié)能性能做了分析研究。與國際水平相比，我國的復合材料螺旋槳研究屬起步階段，研究水平停留在小型螺旋槳的單方面的性能或者復合材料的性能上[14-15]。

本文基于面元法和有限元法，采取流固耦合中的弱耦合方法，開展了復合材料螺旋槳流固耦合分析方法研究。首先探討了復合材料槳葉有限元建模以及水動力外載荷和結(jié)構(gòu)變形的流固耦合交界面數(shù)據(jù)傳遞問題；然后在此基礎上集成為一項復合材料螺旋槳分析技術(shù)，測試了其收斂性；最后進行了比較計算驗證，為復合材料螺旋槳性能分析提供了必要的工具。

2 復合材料螺旋槳流固耦合數(shù)值分析

本文基于流固耦合中的弱耦合處理方法進行復合材料螺旋槳力學性能分析,其中弱耦合方法的基本思路是在每個時間步內(nèi)依次對流體力學方程和結(jié)構(gòu)力學方程求解，通過耦合交界面交換固體域和流體域的計算結(jié)果數(shù)據(jù)，反復交替計算從而實現(xiàn)耦合求解，該方法是研究氣固耦合（氣動彈性）問題的主流方法[16]。在弱耦合數(shù)值計算中，采用面元法進行流體域求解，采用有限元法進行固體域求解。對于復合材料螺旋槳的弱耦合分析流程包括:

1）每一時間步對槳葉水動力載荷的面元法計算；

2）水動力載荷向交界面中有限元計算網(wǎng)格的傳遞；

3）每一時間步對復合材料槳葉位移（變形）等信息的有限元計算；

4）位移（變形）等信息向交界面中面元法計算網(wǎng)格的傳遞。

本文采用文獻[17]中的螺旋槳面元法和文獻[18-19]建立的螺旋槳結(jié)構(gòu)力學有限元分析方法建立復合材料螺旋槳弱耦合分析方法。

2.1 復合材料螺旋槳有限元模型與鋪層

在ANSYS結(jié)構(gòu)有限元分析中，選取8節(jié)點三維層狀結(jié)構(gòu)實體單元SOLID46對復合材料槳葉進行網(wǎng)格劃分，復合材料螺旋槳有限元模型如圖1所示。

螺旋槳采用正交各項異性碳纖維復合材料，材料鋪層結(jié)構(gòu)如圖2所示，左圖顯示的鋪層分布是[30/-30/90/-90/30]s對稱形式，右圖為根部鋪層結(jié)構(gòu)放大圖，共10層，各層之間是等厚度堆棧而成。

圖1 復合材料螺旋槳有限元模型Fig.1 The FE model of composite propellers

圖2 復合材料鋪層結(jié)構(gòu)（左:鋪層分布；右:根部鋪層結(jié)構(gòu)）Fig.2 The layer structure of composite materials(Left:Layer distribution;right:Layer structure in the root)

2.2 復合材料螺旋槳水動力外載荷傳遞

面元法計算外載荷分布結(jié)果通常采取壓力分布表示，如圖3所示，其中壓力系數(shù)CP定義如下:

式中:P、P0、ρ和 V 分別為槳葉表面壓力（Pa）、參考壓力（Pa）、水密度（kg/m3）和槳前方來流水速（m/s）。

采取插值的方式將面元法計算的水動力外載荷從面元法網(wǎng)格傳遞到有限元計算網(wǎng)格節(jié)點上。由于ANSYS復合材料結(jié)構(gòu)分析模塊中有限元外載荷由各節(jié)點的節(jié)點力構(gòu)成，如圖4所示，其中根部施加固定約束。因此，本文先在槳葉面元網(wǎng)格中由壓力分布結(jié)合各面元的面積和法向矢量計算各面元上的三維作用力，然后進行三維作用力的傳遞。

任一面元（圖5所示面元ABCD）法向方向的作用力計算公式如下:

上式中，cosα，cosβ，cosγ為方向余弦，采用螺旋槳面元法中的雙曲面元方向[17]。

圖3 面元法計算壓力分布Fig.3 Pressure distribution from PM

圖4 有限元節(jié)點力分布Fig.4 Force distribution for FEM

圖5 面元ABCDFig.5 A panel ABCD

將計算的各面元三維作用力（圖6所示）根據(jù)面元網(wǎng)格節(jié)點和有限元網(wǎng)格節(jié)點之間的位置關系（最小距離原則）傳遞到有限元網(wǎng)格節(jié)點上（圖4所示），其中葉背和葉面分別獨立進行。傳遞的三維節(jié)點力作為有限元計算邊界條件。

圖6 計算節(jié)點力分布（左:X方向力；中:Y方向力；右:Z方向力）Fig.6 The distribution of calculated forces(Left:X direction;middle:Y direction;right:Z direction)

2.3 復合材料螺旋槳幾何變形位移傳遞

在流固耦合分析中變形后的槳葉需重新進行外載荷計算，因此有限元中計算的復合材料槳葉各節(jié)點三向位移需傳遞給面元網(wǎng)格。本文采用如下插值方法[20-22]:

對于三維槳葉網(wǎng)格坐標（X，Y，Z），其中:X為軸向方向，向后為正，Y為槳葉參考線方向，向上為正，首先將槳葉投影到（X’,Y）平面，其中X’和X軸夾角選取0.7R半徑螺距角，如圖7所示，使得位移插值量的空間分布在（X’,Y）平面變得平坦，從而提高插值精度；然后在（X’,Y）平面中進行有限元網(wǎng)格點F向目標面元網(wǎng)格點P位移量的插值工作，步驟如下:

1） 在（X’，Y）平面中選出離目標P點最近的15個有限元網(wǎng)格源點Fi）（i=1，2…，15），如圖8所示；

2）對于源點 Fi）和位移量 Di（i=1，2…，15）進行二元三次多項式擬合，采用最小二乘法求得三次多項式系數(shù)ai（i=1，2…，10），得到擬合公式如下:

4）依次進行面元網(wǎng)格點的插值。

圖7 原坐標（X）與投影坐標（X’）Fig.7 Original coodinate(X)and projective coodiante(X’)

圖8 插值點（P）與源點（F）Fig.8 Interpolated point(P)and source points(F)

圖9為變形量插值前后的比較，其中左圖是有限元計算的總變形量分布，右圖是插值到面元法網(wǎng)格上的總變形量；從比較圖可以看出，變形量實現(xiàn)了很好的傳遞。

圖9 變形量插值前后比較（左圖:插值前，單位m；右圖:插值后，單位mm）Fig.9 The comparison of deformation(Left:Before interplolation,m;right:After interplolation,mm)

2.4 復合材料螺旋槳參數(shù)反解

由于槳葉形狀復雜，為了避免新槳葉面元網(wǎng)格產(chǎn)生畸變現(xiàn)象[22]，面元網(wǎng)格有時必須針對新的槳葉幾何進行重新劃分；另一方面，復合材料槳葉經(jīng)過流固耦合作用后，需從變形后的槳葉模型中提取螺距、拱度等特征參數(shù)進行分析。以上兩個方面需要開展槳葉模型的參數(shù)反解工作。

槳葉參數(shù)反解按以下兩個步驟進行:首先從槳葉三維模型中提取各半徑剖面葉背和葉面的三維坐標，可在UG建模軟件中實現(xiàn)；然后由各剖面的三維坐標計算槳葉參數(shù)，求解螺距、弦長、厚度、側(cè)斜和拱度等構(gòu)成新槳葉參數(shù)型值表。根據(jù)以上步驟反解的槳葉參數(shù)可達到較高的精度，如圖10所示，采用反解的參數(shù)建模與原模型的比較，結(jié)果十分吻合。

圖10 反解槳葉模型與原槳葉模型比較Fig.10 The comparison between reversely extracted blade model and original blade model

2.5 復合材料螺旋槳流固耦合分析集成

2.5.1 分析流程

復合材料螺旋槳流固耦合分析流程如下:

1）對給定的復合材料螺旋槳幾何型值，采用UG建模，并輸入到ANSYS有限元分析軟件，進行網(wǎng)格劃分；

2）對給定運行工況，采用面元法進行槳葉網(wǎng)格劃分，進而離散面元法方程，求解槳葉表面壓力分布，計算各面元三維作用力；

3）將面元三維作用力傳遞到（1）中有限元劃分的網(wǎng)格節(jié)點上，進行復合材料屬性及鋪層結(jié)構(gòu)設置，對結(jié)構(gòu)力學分析方程進行有限元求解，得到槳葉位移。

4）將有限元計算的各網(wǎng)格點位移插值到原始槳葉面元法中的面元網(wǎng)格節(jié)點上，得到變形后的槳葉幾何，提取新的槳葉幾何型值；

5）判斷新的槳葉幾何是否收斂，如果收斂，則結(jié)束；如果不收斂，則將新的槳葉幾何帶回（2），進行（2）～（5）的迭代計算，直至收斂。

上述分析流程見圖11。

2.5.2 收斂性

螺旋槳根據(jù)側(cè)斜方向不同可分為前側(cè)斜螺旋槳和后側(cè)斜螺旋槳，圖12給出了從槳后向槳前方向看的兩種側(cè)斜分布螺旋槳槳葉模型。目前，船舶螺旋槳大多采用后側(cè)斜螺旋槳，對于前側(cè)斜螺旋槳也有學者研究其對梢渦空泡起始性能的影響[23]。

復合材料螺旋槳采取不同的側(cè)斜分布，呈現(xiàn)不同的流固耦合收斂過程，如圖13所示。圖13左圖表示后側(cè)斜復合材料螺旋槳的收斂過程，圖中螺旋槳推力系數(shù)和扭矩系數(shù)以及0.7R螺距比呈現(xiàn)交替收斂形式，最終收斂槳葉螺距變小，相應的載荷變小。而對于前側(cè)斜復合材料螺旋槳的收斂過程，呈現(xiàn)單調(diào)收斂過程，如圖13右圖所示，最終收斂槳葉螺距變大，相應的載荷變大。

圖11 復合材料螺旋槳流固耦合分析流程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 復合材料螺旋槳流固耦合收斂過程（左:后側(cè)斜；右:前側(cè)斜）Fig.13 The convergence process of FSI of composite propellers

3 復合材料螺旋槳流固耦合分析方法比較驗證

采取本文建立的流固耦合分析技術(shù)，針對文獻[12]中5474槳進行了比較計算驗證。在比較計算中采用相同的Hexcel IM7-8552碳纖維復合材料、鋪層方式，材料屬性如表1所示。

表1 Hexcel IM7-8552碳纖維復合材料屬性Tab.1 Material properties of Hexcel IM7-8552

圖14給出了復合材料螺旋槳變形前后面元法計算的推力系數(shù)和扭矩系數(shù)，其中流固耦合迭代計算點:進速系數(shù)為J=0.66，轉(zhuǎn)速n=780 rpm。在圖14中也給出了文獻[12]的計算結(jié)果，本文計算的推力系數(shù)和扭矩系數(shù)與文獻[12]結(jié)果的差別在4%以內(nèi)。圖15為變形前后的螺距角分布比較，變形后從根部到梢部螺距角變形量逐漸增大，變形后梢部螺距角計算結(jié)果與圖中文獻[12]計算結(jié)果吻合良好。

圖14 5474復合材料螺旋槳水動力比較Fig.14 The comparison of hydrodynamic performance of composite propeller 5474

圖15 5474復合材料螺旋槳螺距角分布比較Fig.15 The comparison of pitch angle distribution of composite propeller 5474

4 結(jié) 論

本文研究了復合材料螺旋槳流固耦合（弱耦合）分析方法。重點探討了復合材料槳葉有限元建模以及流固耦合中載荷及幾何變形量的傳遞問題，采用的方法對復合材料螺旋槳具有較好的適用性。集成的復合材料螺旋槳流固耦合分析技術(shù)對不同側(cè)斜螺旋槳進行的迭代計算均可獲得收斂解，針對5474復合材料螺旋槳的計算結(jié)果與文獻的吻合較好。后續(xù)的工作是開展復合材料螺旋槳性能研究，剖析復合材料螺旋槳流固耦合變形規(guī)律。

[1]楊乃賓,章怡寧.復合材料飛機結(jié)構(gòu)設計[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]楊傳勇.小型船舶復合材料船用螺旋槳的設計與分析[D].大連:大連海事大學,2005.

[15]洪 毅.復合材料船用螺旋槳結(jié)構(gòu)設計研究[D].哈爾濱:哈爾濱工業(yè)大學,2006.

[16]葉正寅,張偉偉,史愛明等.流固耦合力學基礎及其應用[M].哈爾濱:哈爾濱工業(yè)大學出版社,2010.

[17]董世湯,唐登海,周偉新.CSSRC的螺旋槳定常面元法[J].船舶力學,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的螺旋槳流固耦強度校核數(shù)值模擬[C].2007ANSYS-CFD年會文集,2007.

[19]姚志崇,曾志波.復合材料船用螺旋槳流固耦合性能分析[C].2010全國固體力學大會,華中科技大學,2010.

[20]李立州,王婧超,呂震宙,岳珠峰.學科間載荷參數(shù)空間傳遞方法[J].航空動力學報,2007(07):32-36.

[21]虞跨海,岳珠峰.渦輪冷卻葉片參數(shù)化建模及多學科設計優(yōu)化[J].航空動力學報,2007(08):144-149.

[22]岳珠峰等.航空發(fā)動機渦輪葉片多學科設計優(yōu)化[M].北京:科學出版社,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.