孫海濤，熊 鷹，時(shí)立攀

(1.國防科學(xué)技術(shù)大學(xué)指揮軍官基礎(chǔ)教育學(xué)院，湖南長沙 410072;2.海軍工程大學(xué)艦船工程系，湖北武漢 430033)

螺旋槳通常采用金屬材料制作，其設(shè)計(jì)和性能預(yù)報(bào)在剛性假設(shè)基礎(chǔ)上進(jìn)行，不考慮結(jié)構(gòu)變形的影響.但對于具有較大側(cè)斜角的螺旋槳或采用復(fù)合材料制作的螺旋槳，彈性的增加引起變形增大.變形后的螺旋槳將影響周圍流場及其承受的水動(dòng)力載荷，變化后的載荷又將影響槳的變形，這是流體與結(jié)構(gòu)相互作用的問題.此種情況下，槳葉變形對螺旋槳性能的影響不可忽略.

彈性變形槳的數(shù)值計(jì)算最早出現(xiàn)于1991年，文獻(xiàn)［1］采用渦格法計(jì)算了水下無空泡螺旋槳的定常水動(dòng)力性能，并利用有限元軟件 ABAQUS進(jìn)行應(yīng)力分析，未考慮流體結(jié)構(gòu)的相互作用.文獻(xiàn)［2-7］在后來的研究中，針對復(fù)合材料螺旋槳耦合變形的特點(diǎn)改進(jìn)上述方法，將渦格法計(jì)算程序與ABAQUS聯(lián)合迭代求解，并進(jìn)行了強(qiáng)度評估和模型試驗(yàn)研究，研究表明:材料纖維及鋪層設(shè)計(jì)對螺旋槳變形和性能有顯著影響，但僅靠改變纖維鋪層方式難以最大限度地發(fā)揮復(fù)合材料螺旋槳的優(yōu)勢，須配合調(diào)整螺旋槳的外形，進(jìn)行預(yù)變形設(shè)計(jì).文獻(xiàn)［8］以燃料消耗最小化為目標(biāo)進(jìn)行纖維的鋪層優(yōu)化設(shè)計(jì)，采用邊界元法與ANSYS軟件對設(shè)計(jì)結(jié)果進(jìn)行耦合計(jì)算.文獻(xiàn)［9-12］將邊界元法與有限元法結(jié)合起來，建立了流固耦合算法，并研究了槳葉變形特性及變形對水動(dòng)力性能的影響.流體部分采用基于速度勢的低階面元法，運(yùn)用Newton-Raphson方法對格林第3公式進(jìn)行離散求解，得出槳葉面元的速度勢進(jìn)而求得螺旋槳水動(dòng)力性能及槳葉表面壓力分布.結(jié)構(gòu)部分利用ABAQUS軟件的二次開發(fā)功能，編寫實(shí)體單元的生成程序和水動(dòng)力載荷的加載程序，以實(shí)現(xiàn)流體面元與結(jié)構(gòu)單元表面的對應(yīng)和載荷的正確加載.得到節(jié)點(diǎn)位移之后，將其導(dǎo)入流體計(jì)算程序以更新面元節(jié)點(diǎn)信息進(jìn)行新的求解.文獻(xiàn)［13-16］在復(fù)合材料螺旋槳的流固耦合算法方面進(jìn)行了大量的研究，建立了較為完善的流固耦合控制方程和求解算法，給出了復(fù)合材料螺旋槳預(yù)變形設(shè)計(jì)的原則.

上述數(shù)值算法在求解螺旋槳水動(dòng)力載荷時(shí)，均采用基于流體為無黏假定的勢流方法.隨著計(jì)算機(jī)技術(shù)的發(fā)展，基于雷諾時(shí)均Navier-Stokes方程的計(jì)算流體力學(xué)方法(CFD)逐漸成為一種預(yù)報(bào)螺旋槳性能的新興方法［17］，該方法綜合考慮了流體黏性和旋度等因素的影響，使獲得螺旋槳黏性流場的全場信息成為可能.筆者采用求解RANS方程的方法對螺旋槳3維流場進(jìn)行模擬，利用有限元方法計(jì)算螺旋槳結(jié)構(gòu)響應(yīng)，結(jié)合ANSYS軟件中的流體載荷和結(jié)構(gòu)響應(yīng)數(shù)據(jù)實(shí)時(shí)雙向傳輸技術(shù)，對螺旋槳進(jìn)行流固耦合特性分析.

1 數(shù)值方法

RANS方程是黏性流體運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)的普適性方程，與連續(xù)性方程一起構(gòu)成流體控制方程，將其作為求解螺旋槳水動(dòng)力性能的基本方程.RANS方程是通過將瞬時(shí)N-S方程中的速度、壓強(qiáng)、質(zhì)量力和密度等流體變量進(jìn)行時(shí)歷平均化后得到，其形式如下:

式中:ρ為流體密度;ui，uj，ul為速度分量;p0為靜壓;fi為單位質(zhì)量力.

槳葉在外載荷作用下將產(chǎn)生彎扭耦合變形，對其運(yùn)動(dòng)學(xué)方程進(jìn)行有限元離散可得有限元方程:

槳葉質(zhì)量可由密度和外形等確定，剛度可由彈性模量和泊松比確定.槳葉結(jié)構(gòu)阻尼較難計(jì)算，通常將其簡化為質(zhì)量矩陣和剛度矩陣的線性組合.

流場計(jì)算采用ANSYS中的CFX求解器，結(jié)構(gòu)響應(yīng)計(jì)算采用ANSYS中的結(jié)構(gòu)有限元求解器.為實(shí)現(xiàn)流場與結(jié)構(gòu)響應(yīng)計(jì)算數(shù)據(jù)的實(shí)時(shí)雙向傳遞，將有限元求解器的設(shè)置輸出為結(jié)構(gòu)求解程序，形成CFX可接受的INP文件，而后在CFX的分析類型中打開外部求解器功能，選擇ANSYS多場分析并讀入生成的結(jié)構(gòu)求解程序.這樣做的優(yōu)勢:可以充分發(fā)揮流體計(jì)算與結(jié)構(gòu)計(jì)算的已有經(jīng)驗(yàn)，并且所有流體與結(jié)構(gòu)的分析結(jié)果都將保存在流體模塊的結(jié)果文件中，便于查看.雙向流固耦合分析的設(shè)置相對較為復(fù)雜，設(shè)置不合理將導(dǎo)致計(jì)算的中斷，主要集中的問題:① 流體分析與結(jié)構(gòu)分析的時(shí)間步必須統(tǒng)一;②結(jié)構(gòu)變形導(dǎo)致的流場網(wǎng)格問題，需要采用網(wǎng)格移動(dòng)技術(shù)方可模擬變形對流場的影響.

2 網(wǎng)格劃分

若要求解較為準(zhǔn)確的流場信息，就需要在計(jì)算域中劃分出能正確反映螺旋槳邊界形狀且具有一定分辨率的計(jì)算網(wǎng)格.螺旋槳幾何形狀復(fù)雜，葉面呈螺旋狀扭曲結(jié)構(gòu)，槳葉剖面沿半徑方向先逐漸擴(kuò)大而后急劇收縮，對其流場劃分結(jié)構(gòu)化網(wǎng)格十分困難.并且結(jié)構(gòu)化網(wǎng)格以六面體網(wǎng)格為基礎(chǔ)，即使設(shè)法對包含有螺旋槳這樣復(fù)雜幾何形狀的計(jì)算域進(jìn)行了結(jié)構(gòu)化網(wǎng)格劃分，也會產(chǎn)生大量質(zhì)量較差的扭曲網(wǎng)格，影響計(jì)算的速度和精度，而非結(jié)構(gòu)網(wǎng)格就可以避免這一問題.非結(jié)構(gòu)化網(wǎng)格以四面體為基礎(chǔ)，對網(wǎng)格結(jié)點(diǎn)沒有嚴(yán)格的限制，比較適用于像螺旋槳流場這樣邊界幾何形狀復(fù)雜的計(jì)算域，但非結(jié)構(gòu)化網(wǎng)格的不足在于同等網(wǎng)格尺度條件下填充率較低，所需要的網(wǎng)格數(shù)量較大.

以HG5-1螺旋槳為研究對象，模型直徑D為240 mm，將計(jì)算域分為圓柱形的固定域和旋轉(zhuǎn)域.外圓柱體區(qū)域直徑為9.0D，入口取在螺旋槳上游3.0D處，出口取在下游7.0D處.內(nèi)圓柱體區(qū)域直徑為1.5D，入口取在螺旋槳上游0.8D處，出口取在下游1.2D處.內(nèi)部旋轉(zhuǎn)域采用單旋轉(zhuǎn)系模型對螺旋槳的旋轉(zhuǎn)流場進(jìn)行求解.這種模型是建立一個(gè)固連于螺旋槳的旋轉(zhuǎn)坐標(biāo)系，并將螺旋槳周圍的流域獨(dú)立出來，從而使得這個(gè)流域的流場在旋轉(zhuǎn)坐標(biāo)系下是定常的，而這個(gè)流域之外的流域在絕對靜止坐標(biāo)系下也可以近似認(rèn)為是定常的.為降低總體計(jì)算網(wǎng)格的數(shù)量，采用結(jié)構(gòu)化-非結(jié)構(gòu)化多塊混合網(wǎng)格劃分方法，即對緊鄰于螺旋槳的形狀復(fù)雜流域進(jìn)行非結(jié)構(gòu)化網(wǎng)格劃分，如圖1，2所示.

圖1 螺旋槳表面的非結(jié)構(gòu)化網(wǎng)格

圖2 流場內(nèi)域的非結(jié)構(gòu)化網(wǎng)格

為彌補(bǔ)非結(jié)構(gòu)網(wǎng)格填充率較低的不足，對螺旋槳流場的外域，劃分高質(zhì)量的結(jié)構(gòu)化網(wǎng)格，流場計(jì)算域的總網(wǎng)格數(shù)為150.0萬個(gè)，如圖3所示.槳葉變形計(jì)算時(shí)，采用四面體單元對其進(jìn)行有限元剖分，總單元數(shù)為2.9萬個(gè)，如圖4所示.

圖3 流場外域的結(jié)構(gòu)化網(wǎng)格

圖4 螺旋槳有限元網(wǎng)格

3 結(jié)果與分析

螺旋槳轉(zhuǎn)速設(shè)為900 r·min-1，通過變化進(jìn)流速度大小來改變進(jìn)速系數(shù).螺旋槳材料為銅合金，其楊氏模量為110 GPa，泊松比為0.34，密度為8 300 kg·m-3.為進(jìn)行比較分析，增加螺旋槳的彈性，將材料的楊氏模量降為15 GPa，密度設(shè)為1 500 kg·m-3，其他參數(shù)保持不變.記彈性增加后的槳為彈性槳，銅合金槳為剛性槳.流固耦合面設(shè)為槳葉表面，數(shù)據(jù)即在此面上交換傳遞.

剛性槳水動(dòng)力性能的預(yù)報(bào)值與試驗(yàn)值的比較如圖5所示，KF，KM，η分別為彈性槳的推力系數(shù)、轉(zhuǎn)矩系數(shù)和效率，由于變形較小，考慮結(jié)構(gòu)變形影響與否的水動(dòng)力性能變化甚微，均與試驗(yàn)結(jié)果較為接近.因而，對常規(guī)金屬槳進(jìn)行性能預(yù)報(bào)時(shí)，可以不考慮流固耦合的影響.而彈性槳的情況則不同，彈性的增加使得變形增大，進(jìn)而對周圍流場產(chǎn)生影響.彈性槳與剛性槳預(yù)報(bào)值比較如圖6所示，考慮變形的影響后，彈性槳的推力和轉(zhuǎn)矩均比剛性槳小.在低進(jìn)速時(shí)，由于轉(zhuǎn)矩小更多，效率略高于剛性槳.高進(jìn)速時(shí)，推力減小更多，效率略低于剛性槳.進(jìn)速系數(shù)J=0.5時(shí)，彈性槳的最大變形為0.81 mm，而剛性槳的最大變形僅為0.11 mm;彈性槳的最大應(yīng)變?yōu)?.03×10-3，剛性槳的最大應(yīng)變?yōu)?.42×10-4;彈性槳的最大等效應(yīng)力為11.03 MPa，剛性槳的最大等效應(yīng)力為11.16 MPa.對于彈性槳而言，槳葉變形不僅影響水動(dòng)力性能，還將引起槳葉幾何參數(shù)的變化，進(jìn)而導(dǎo)致整個(gè)運(yùn)轉(zhuǎn)狀態(tài)及流場的改變.在這種情況下，考慮流固耦合的影響是很有必要的.

圖5 水動(dòng)力性能預(yù)報(bào)值與試驗(yàn)值的比較

圖6 彈性槳與剛性槳預(yù)報(bào)值的比較

不同工況下，彈性槳與剛性槳的結(jié)構(gòu)響應(yīng)特性對比如圖7所示，Δsmax為最大位移，σmax為最大應(yīng)力，εmax為最大應(yīng)變.隨著進(jìn)速系數(shù)降低，螺旋槳承受的載荷增加，變形及應(yīng)力應(yīng)變隨之增加.在進(jìn)速系數(shù)相同時(shí)，彈性槳的最大變形及應(yīng)變明顯高于剛性槳，而最大等效應(yīng)力略低于剛性槳.可見，材料特性對變形和應(yīng)變的影響是顯著的.而應(yīng)力則是表示單位面積上的受力，在螺旋槳載荷未發(fā)生明顯變化的情況下，其變化范圍較小.前面的計(jì)算顯示，彈性槳承受的載荷較剛性槳為小，因而其最大等效應(yīng)力亦略小于剛性槳.

最大結(jié)構(gòu)響應(yīng)從宏觀上反映了螺旋槳運(yùn)轉(zhuǎn)時(shí)的情況，在強(qiáng)度評估與分析時(shí)，往往需要研究更為細(xì)致的應(yīng)力應(yīng)變分布.2種材料螺旋槳在J=0.5時(shí)的結(jié)構(gòu)響應(yīng)云圖如圖8所示，槳葉最大變形出現(xiàn)在葉梢，變形量從葉梢到葉根依次遞減，類似于懸臂梁受力后的變形特性.設(shè)R為螺旋槳半徑，在(0.6～0.7)R之間出現(xiàn)了明顯的應(yīng)力集中，最小值則位于葉梢和葉根隨邊部位.應(yīng)變的分布規(guī)律與等效應(yīng)力相近，因?yàn)槎叩南嗖顑H在于彈性矩陣.

圖9，10分別為槳葉壓力分布和尾流的計(jì)算結(jié)果，其分布規(guī)律相似，從圖中難以看出區(qū)別.具體數(shù)值上，彈性槳壓力的負(fù)壓峰值增加，正壓峰值降低，葉面和葉背的壓力差減小.從流線分布情況可知:尾渦半徑有一定的收縮，反映了螺旋槳對流體的抽吸作用.

圖11，12分別為螺旋槳在J=0.5時(shí)不同斷面處壓力分布和渦量分布的計(jì)算結(jié)果，圖11a，d，圖12a，d 為盤面上游 0.20D位置，圖 11b，e，圖 12 b，e為盤面位置，圖 11c，f，圖 12 c，f為盤面 下游0.15D位置.渦量幅值和壓力在盤面處最大，下游斷面次之.彈性槳的渦量幅值和壓力在盤面處的增加，上游與下游斷面有所降低，分布形式無明顯區(qū)別.

圖7 彈性槳與剛性槳位移、應(yīng)力及應(yīng)變比較

圖8 剛性槳與彈性槳變形、應(yīng)力及應(yīng)變云圖

圖9 剛性槳與彈性槳壓力分布云圖

圖10 J=0.5時(shí)螺旋槳尾流特性的計(jì)算結(jié)果

圖11 不同斷面處的壓力分布

圖12 不同斷面處的渦量分布

4 結(jié)論

利用ANSYS軟件中流體載荷與結(jié)構(gòu)響應(yīng)數(shù)據(jù)的雙向傳輸技術(shù)，將流體與結(jié)構(gòu)計(jì)算耦合在一起，給出了求解螺旋槳流固耦合問題的3維數(shù)值模擬方法.設(shè)定2種材料特性，對螺旋槳的水動(dòng)力性能、槳葉變形及應(yīng)力應(yīng)變特性等進(jìn)行了數(shù)值模擬.結(jié)果表明:彈性槳的變形與應(yīng)變較剛性槳約高1個(gè)數(shù)量級.槳葉的變形又對水動(dòng)力性能及周圍流場產(chǎn)生了影響，彈性槳的推力和轉(zhuǎn)矩均比剛性槳小.當(dāng)轉(zhuǎn)矩減小更多時(shí)，槳的效率有所提高.這一特性可用于復(fù)合材料螺旋槳的優(yōu)化設(shè)計(jì)，可以實(shí)現(xiàn)螺旋槳整體性能的提升.

References)

［1］Lin G F.Comparative stress-deflection analyses of a thick-shell composite propeller blade［R］.Bethesda，MD:David Taylor Research Center DTRC/SHD-1373-01，1991.

［2］Lin H J，Lin J 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.

［3］Lin H J，Lai W M，Kuo Y M.Effects of stacking sequence on nonlinear hydroelastic behavior of composite propeller blade ［J］.Journal of Mechanics，2010，26(3):293-298.

［4］Lin C C，Lee Y J.Regression of the response surface of laminated composite structures［J］.Composite Structures，2003，62:91-105.

［5］Lin C C，Lee Y J.Stacking sequence optimization of laminated composite structures using genetic algorithm with local improvement ［J］.Composite Structures，2004，63:339-345.

［6］Lee Y J，Lin C C.Optimized design of composite propeller［J］.Mechanics of Advanced Materials and Structures，2004，11(1):17-30.

［7］Lin H J，Lin J J，Chuang T J.Strength evaluation of a composite propeller blade ［J］.Journal of Reinforced Plastics and Composites，2005，24(17):1791-1807.

［8］Blasques José PedroAlbergariaAmaral， Berggreen Christian，Andersen Poul.Hydro-elastic analysis and optimization of a composite marine propeller［J］.Marine Structures，2010，23:22-38.

［9］Sun Haitao，Xiong Ying.Fluid-structure interaction analysis of flexible marine propellers［J］.Applied Mechanics and Materials，2012，226:479-482.

［10］孫海濤，熊 鷹.槳葉變形對復(fù)合材料槳水動(dòng)力性能影響［J］.華中科技大學(xué)學(xué)報(bào):自然科學(xué)版，2013，41(6):81-85.

Sun Haitao，Xiong Ying.Effect of blade deformation of a composite propeller on its hydrodynamic performance［J］.Journal of Huazhong University of Science and Technology:Natural Science Edition，2013，41(6):81-85.(in Chinese)

［11］孫海濤，熊 鷹，黃 政.大側(cè)斜螺旋槳水動(dòng)力及變形特性研究［C］∥第十一屆全國水動(dòng)力學(xué)學(xué)術(shù)會議暨第二十四屆全國水動(dòng)力學(xué)研討會.北京:中國海洋出版社，2012:961-967.

［12］孫海濤，熊 鷹.考慮變形的螺旋槳水動(dòng)力及變形特性研究［J］.哈爾濱工程大學(xué)學(xué)報(bào)，2013，34(9):1108-1112.

Sun Haitao，Xiong Ying.Study on hydrodynamic and deformation performance of propellers considering the blade deformation［J］.Journal of Harbin Engineering U-niversity，2013，34(9):1108-1112.(in Chinese)

［13］Young Y L.Time-dependent hydroelastic analysis of cavitating propulsors［J］.Journal of Fluids and Structures，2007，23:269-295.

［14］Young Y L.Fluid-structure interaction analysis of flexible composite marine propellers［J］.Journal of Fluids and Structures，2008，24:799-818.

［15］Liu Z K，Young Y L.Utilization of bend-twist coupling for performance enhancement of composite marine propellers［J］.Journal of Fluids and Structures，2009，25:1102-1116.

［16］Motley M R，Liu Z K，Young Y L.Utilizing fluid-structure interactions to improve energy efficiency of composite marine propellers in spatially varying wake［J］.Composite Structures，2009，90:304-313.

［17］劉志華，熊 鷹，葉金銘.基于多塊混合網(wǎng)格的RANS方法預(yù)報(bào)螺旋槳敞水性能的研究［J］.水動(dòng)力學(xué)研究與進(jìn)展，2007，22(4):450-456.

Liu Zhihua，Xiong Ying，Ye Jinming.Study on the prediction of propeller open-water performance using RANS formula and multi-block hybrid meshes［J］.Chinese Journal of Hydrodynamic，2007，22(4):450-456.(in Chinese)