葉禮裕 王 超 孫 帥 孫文林

（哈爾濱工程大學船舶工程學院 哈爾濱 150001）

0 引 言

隨著現代艦船逐步向大型化、高速化發(fā)展，以及大功率主機的應用，引起螺旋槳槳葉表面負荷增大；由于減振降噪方面有特殊的優(yōu)勢，大側斜螺旋槳在船舶領域得到廣泛的應用［1］，但特殊的幾何外形使其強度問題較為突出.另外，為了使螺旋槳實際應用中更加節(jié)能化和環(huán)?；?，設計槳應當符合輕量化的要求.因此，開展螺旋槳強度的研究不僅在于槳葉強度符合規(guī)范要求，更應當能夠準確地預報槳葉應力的分布.

目前，螺旋槳槳葉應力的預報通常采用有限元法和模型試驗［2］.對于有限元法，開展較多的是采用CFD計算和有限元分析軟件結合的方法預報槳葉應力分布［3－5］，也有學者基于面元法和有限元法耦合的方法預報槳葉應力分布［6］.在模型試驗方面，Boswell［7］對大側斜螺旋槳葉片進行了靜態(tài)應力測量試驗；趙波［8］開展了大側斜螺旋槳的靜態(tài)應力試驗和動態(tài)應力試驗，并與理論計算結果進行對比；楊向暉等［9］通過在大側斜槳模表面粘貼應變片的方法研究不同工況下的槳葉應變和應力分布.有限元法能夠較為準確預報螺旋槳槳葉應力分析，但是需要復雜的模型建立、網格劃分過程，不利于開展螺旋槳的設計.模型試驗方法成本較高、實驗難度大且需要花費大量時間.螺旋槳運轉時主要受到周向力、軸向力、離心力，以及偶然外力作用，從而產生彎曲、轉矩和拉伸變形，可將槳葉簡化為變截面扭曲的懸臂梁［10］.懸臂梁法也是一種比較方便且可行的槳葉應力預報方法.以往學者開展懸臂梁預報螺旋槳應力，因螺旋槳幾何外形復雜且橫截面為非對稱的，計算截面所受彎矩、離心力，以及面積、抗彎模量等的過程中運用了一些經驗公式，造成只能預報槳葉有限個點的應力或者應力預報結果，比較保守.

為準確預報槳葉各半徑切面所受的彎矩和離心力，將懸臂梁法理論中的由螺旋槳推力、旋轉阻力產生的彎矩、離心力及其產生的彎矩的積分計算公式進行離散處理.槳葉切面的形心和慣性矩可借助于工程力學中截面靜矩和慣性矩等計算公式來求解，以準確求解槳葉各半徑切面上各個點的抗彎剖面模數.本文采用FORTRAN語言編譯懸臂梁法槳葉應力預報程序，將其與螺旋槳的定常面元法性能預報程序對接.以無側斜槳和大側斜螺旋槳為例，通過本文計算結果和有限元計算結果對比分析來驗證本文所提方法的有效性.

1 懸臂梁離散方法

懸臂梁法預報螺旋槳強度是將螺旋槳葉梢視為自由端，而葉根固定于槳轂的一根變截面的懸臂梁，一般按正車最大航速工況預報螺旋槳靜載荷，將作用在螺旋槳上的外力載荷視為不變.以平切面作為校核面，該校核面為垂直于槳葉參考線的面且與槳葉相交所截的平面，對大側斜槳有更好的適用性［11］.

槳葉切面各點抗彎彈性模量計算過程如下：首先，對該切面進行分割，計算出其面積和靜矩，進而求出該剖面的形心坐標；然后，計算出該剖面的慣性矩；最后，由該剖面的形心坐標和慣性矩，可求處其各點的抗彎彈性模量.具體的切面形心、靜矩、慣性矩和抗彎彈性模量求解方法可參考文獻［12］.

離散法可準確地計算槳葉各半徑切面的彎矩.根據離散法原理，需對求解的螺旋槳幾何模型進行分割.為此，將槳葉沿展向分割成有限個區(qū)域，半徑r處取一葉元體，厚度為dr.為更好對接面元法預報程序，懸臂梁法的槳葉展向劃分方式與面元法相同，即采用余弦分割，節(jié)點為

1.1 推力和旋轉阻力引起的彎矩積分公式離散方法

對于常規(guī)螺旋槳，推力和旋轉阻力是產生彎矩的主要成分.設該葉元體所受的推力為dT和旋轉阻力dF，則槳葉推力和旋轉阻力對半徑處的彎矩分布為

為獲得上式的數值解，將該式離散處理，沿展向分割成Nr份，則

為計算推力和旋轉阻力產生的彎矩，必須計算出ΔTj和ΔFj大小，而ΔTj和ΔFj可通過面元法計算得到.面元法計算螺旋槳總的推力和轉矩的離散公式為［13－14］

面元法是基于勢流理論的分析方法，計算中沒有考慮流體粘性的影響.為使結果更符合實際情況，計算水動力性能時，采用表面摩擦阻力的方法近似計入粘性的影響，用以下公式求解推力和扭矩的粘性修正量.

只要對面元法計算螺旋槳推力和扭矩的離散公式稍加處理，即可用于求解每個葉元體的推力和旋轉阻力，即

1.2 離心力及其產生的彎矩積分公式離散方法

螺旋槳高速運轉時，槳葉各徑向剖面將受到離心力的作用.離心力引起的彎矩相對于推力和旋轉阻力要小，但有側斜和縱傾的螺旋槳由于離心力作用線已不是徑向剖面幾何中心的連線，此時產生的離心力彎矩較大.

假設螺旋槳的轉速為n（r／s），材料的重量密度為γ（10kN／m3），計算半徑處的切面面積為S（m2），重力加速度為g（9.8m／s2），則取一小微元dr，它對半徑rP處的離心力為

與推力和旋轉阻力產生彎矩數值計算方法類似，將上述進行離散處理，得到

式中：Sj為切面的面積.

對于有縱傾的螺旋槳，將因槳葉縱傾產生離心力彎矩MC.假設在不同半徑處槳葉縱傾分布為h（r），下式縱傾引起的槳葉rP處切面彎矩MC及其離散結果

對于有側斜的螺旋槳，將因槳葉側斜引起離心彎矩MS.通常，螺旋槳側斜的方向是與其旋轉方向相反的.因此，側斜引起的離心力彎矩MS與縱傾引起的彎矩MF的方向相反.假設各半徑側斜分布為p（r），則側斜引起的槳葉rp處切面彎矩MC及其離散結果為

1.3 切面上各點應力求解

針對半徑rp處承受的彎矩情況進行分析，見圖1.將推力、旋轉阻力和離心力引起的彎矩沿ξ－ξ軸和η－η軸方向分解，可以得到合成彎矩.

圖1 槳葉剖面承受彎矩情況

校核切面的面積及切面上任何一點的抗彎剖面模數的求解方法已在以上作簡要介紹.假設P（ξ，η）為切面上的一點，以Wξ（η）和Wη（ξ）分別表示P點對于ξ－ξ軸和η－η軸抗彎截面模量.那么，合成彎矩Mξ和Mη在P點處的應力可表示為：

若P（ξ，η）點所在校核切面的面積為S，則離心力C在P點產生的應力為

綜述所述，合成彎矩和離心力所產生的總應力為

式中：σξ（η），ση（η）正負號取值需要根據具體情況而定，即ξ－ξ軸和η－η軸把校核剖面分割為4個區(qū)域，表1給出了每個區(qū)域受力狀態(tài)相同.若為拉應力，則取正值.反之，取負值.

表1 槳葉剖面各區(qū)域的應力性質

2 算例分析

根據以上介紹的懸臂梁離散方法，筆者編譯了懸臂梁法預報程序，并與定常面元法水動力性能預報程序進行對接，以便更好地把螺旋槳水動力性能預報傳入到懸臂梁法程序計算中.程序運行計算后，可將槳葉的應力預報結果傳入到Tecplot軟件中進行后處理，得到槳葉應力分布云圖.為了證明本文方法的可行性和有效性，將本文方法與有限元軟件的槳葉應力的預報結果進行比較.

這里需要強調的是，一般情況下，任意一點的主單元體的應力狀態(tài)是三向的應力.懸臂梁法是以最大拉壓應力σ1為校核條件的，將其他2項應力σ2和σ3忽略不計，而有限元法計算的以von Mises應力（等效應力）σε為校核條件的，見文獻［15］，σε計算方法見式（20）.根據螺旋槳的受力情況可知，槳葉上任意一點的最大拉壓應力σ1遠大于其他2項應力σ2和σ3.因此，槳葉任意一點的最大拉壓應力σ1的絕對值與von Mises應力σs的大小相差不大.

本文選了2種類型的螺旋槳進行算例分析，包括無縱傾無側斜的P4119、P4381槳，以及某大側斜調距槳.

2.1 P4119槳和P4381槳的模型算例分析

對于無縱傾無側斜槳的算例，本文選用22nd ITTC Propulsion Committee提供的P4119槳和P4381槳為算例，其幾何參數見表2.

表2 P4119槳和P4381的主要參數

為了驗證本文方法預報螺旋槳強度的準確性，將本文的應力預報結果與文獻［16］進行比較.預報P4119槳的應力分布的工況為：進速系數0.833，轉速600r／min.在該工況下.本文面元法程序和文獻［16］預報的P4119槳水動力性能見表3，可見2種的誤差較小.

表3 P4119槳的推力和轉矩系數對比

圖2為本文方法預報DTRC4119槳的葉背和葉面的應力分布.文獻［16］是采用有限元計算軟件計算DTRC4119槳槳葉應力分布，預報結果見圖3.由圖2～3可知，本文計算的槳葉的應力分布趨勢與文獻［16］計算結果類似.本文方法預報槳葉最大拉應力為，而文獻［16］預報最大von Mises應力為，兩者的最大應力預報結果偏差較小.另外，注意到文獻［16］提到其最大應力發(fā)生在槳葉葉根弦向中部處，而本文預報結果最大應力也是發(fā)生在葉根弦向中部處.

圖2 本文方法預報的P4119槳槳葉應力分布

圖3 文獻［16］預報P4119槳槳葉應力分布

P4381槳的槳葉應力計算工況為進速2.711 45m／s、轉速600r／min.P4381槳的有限元軟件應力預報步驟如下：首先，在ICEM軟件中建立P4381槳的幾何模型并對其劃分網格，導入到Fluent軟件中進行水動力計算；將P4381槳的幾何模型導入到ANSYS結構力學模塊中進行結構網格劃分，將P4381槳的水動力載荷傳入到調距槳該槳的有限元模型中，計算出P4381槳葉面和葉背的應力分布.

在計算工況下，本文的定常面元法程序和CFD方法預報得到的該槳推力系數和轉矩系數見表4，2種方法得到的該槳推力系數和轉矩系數誤差較小.

表4 P4381槳的推力和扭矩系數對比

運行懸臂梁法預報程序，得到了P4381槳的葉背和葉面的應力分布情況，見圖4.有限元法預報P4381槳的葉背和葉面的應力分布見圖5.從葉背應力預報結果分析，本文方法和有限元法預報的P4381槳葉背的分布趨勢基本是一致的；本文方法預報葉背最大壓應力為3.21MPa，而有限元法預報葉背最大von Mises應力為3.61MPa，2種方法計算得到葉背最大應力值相近且均是發(fā)生于葉根弦向中部處.從葉面應力預報結果分析，本文方法和有限元法預報的P4381槳葉面的分布趨勢有一定的差別，但是葉根處的分布趨勢比較相近；本文方法預報葉面最大拉應力為3.16 MPa，而有限元法預報葉面最大von Mises應力為3.61MPa，2種方法計算得到葉面最大應力值相近且均是發(fā)生于葉根弦向中部處.

圖4 本文方法預報P4381槳槳葉應力分布

圖5 有限元法預報P4381槳槳葉應力分布

綜上所述，本文方法可用于預報常規(guī)槳的應力分布趨勢，能夠準確預報常規(guī)槳的最大應力值及其發(fā)生的位置.但是，在預報常規(guī)槳葉面的外半徑處的應力分布有所缺陷，主要原因是葉面外半徑處葉切面比較平坦，導致同一半徑處的葉切面處的應力值趨向相同.

2.2 大側斜槳算例分析

本文大側斜螺旋的計算模型選的是某大側斜調距槳，主要參數見表5.

表5 大側斜槳的主要參數

該槳的計算工況：進流速度為11.8m／s，轉速為3.4r／s.有限元法計算該槳的步驟與P4119槳和P4381槳的模型算例一樣.

從敞水性能預報的準確性分析，表6給出了計算工況下本文面元法和CFD方法計算得到的調距槳的敞水性能，兩者的誤差很小，即不會由于水動力性能計算結果的不同而引起懸臂梁法和有限元法預報該槳的槳葉應力分布有較大的差別.

表6 大側斜槳的敞水性能對比

本文方法預報該槳的槳葉應力分布云圖見圖6.由圖6可知，本文方法預報的槳葉最大壓應力為83.6MPa，發(fā)生于葉背中部的弦向中部處，槳葉壓應力也主要集中于葉背中部弦向中部處；最大拉應力為92.7MPa，發(fā)生于葉面中部的隨邊處，槳葉的拉應力主要集中與槳葉中部的導邊和隨邊處.有限元軟件預報該調距槳的槳葉應力分布云圖，見圖7.由圖7可知，有限元法預報的葉背最大von Mises應力為83.2MPa，發(fā)生于葉背中部的弦向中部靠近隨邊處，槳葉的拉應力也集中于此處；葉面最大von Mises應力為91.6 MPa，發(fā)生于葉面中部的弦向隨邊處.2種方法預報大側斜槳葉背和葉面的最大應力大小差別小，且發(fā)生位置比較相近，驗證了懸臂梁法在計算大側斜槳槳葉應力分布的可行性和準確性.但是，2種方法預報的槳葉應力分布趨勢有所差別，即懸臂梁法預報拉應力集中于葉面中部的導邊和隨邊處，而有限元法預報拉應力只集中于葉面中部的隨邊處，主要原因是：由于懸臂梁法的彎矩計算原理的局限性，導致計算得到同一徑向處不同弦向位置的彎矩是相同的，而有限元法計算得到同一徑向處不同弦向位置的彎矩是不同的.

圖6 本文方法預報的大側斜槳槳葉的應力分布

圖7 有限元法預報大側斜槳槳葉的應力分布

4 結 論

由于螺旋槳的幾何外形復雜，懸臂梁法理論螺旋槳推力、旋轉阻力產生的彎矩、離心力及其產生的彎矩的積分公式很難用解析的方法求出.為準確計算出槳葉切面的彎矩和離心力，本文建立了懸臂梁離散方法用于求解槳葉切面所受彎矩和離心力的.基于FORTRAN語言編譯懸臂梁法槳葉應力預報程序，并與螺旋槳定常面元法性能預報程序對接.以無縱傾無側斜槳和大側斜螺旋槳為例，預報槳葉應力分布，將計算結果與有限元法計算結果相比較，分析得到如下結論.

1）本文方法能較為準確的預報無縱傾無側斜槳的葉面和葉背最大應力值，且對于該類型的槳最大應力發(fā)生于槳葉葉根的弦向中部.除了預報常規(guī)槳葉面的外半徑處的應力分布有所缺陷，本文方法預報常規(guī)槳的槳葉的應力分布趨勢與有限元法預報趨勢是相同.

2）本文方法可較為準確預報大側斜槳的葉面和葉背最大應力值.對于本算例的平衡式側斜槳，最大拉應力一般發(fā)生于槳葉中部的隨邊處，最大壓應力發(fā)生于槳葉中部的弦向中部處.懸臂梁法預報槳葉應力分布趨勢有缺陷，即懸臂梁法預報螺旋槳的拉應力在葉面中部的導邊和隨邊處有集中現象，而實際上對于該槳拉應力只集中于葉面中部的隨邊處.

總而言之，應用本文方法可快速得找到在計算工況下的螺旋槳槳葉最大應力大小及其發(fā)生位置，便于對螺旋槳的強度校核.雖然本文方法預報在槳葉應力分布趨勢與有限元法還存在某些缺陷，但是總體上還是可行的.在后續(xù)的工作中，作者將開展將該方法的程序嵌入到螺旋槳的理論設計和優(yōu)化設計過程中，以保證設計槳有足夠的強度.

［1］CUMMING R A，MORGAN W B，BOSWELL R J.Highly skewed propellers［J］.SNAME Transactions，1984（1）：90：98.

［2］S?NTVEDT T.Propeller blade stresses，application of finite element methods［J］.Computers ＆ Structures，1974，4（1）：193－204.

［3］文學棟，王永生，李堅波.船用螺旋槳槳葉應力數值計算［J］.船海工程，2010，39（1）：27－30.

［4］柳章存.船用螺旋槳的數值模擬及流固分析［D］.鎮(zhèn)江：江蘇科技大學，2012.

［5］任 弘，李范春，杜 玲.流固耦合作用對螺旋槳強度影響的數值計算［J］.武漢理工大學學報：交通科學與工程版，2015，29（1）：144－147.

［6］劉竹青，陳奕宏，姚志崇.基于面元法及有限元法耦合的螺旋槳強度計算［J］.中國造船，2012，53（1）：25－30.

［7］BOSWELL R J.Static stress measurements on a highly－skewed propeller blade［R］.Naval Ship Research and Development Center Washington DC，1969.

［8］趙 波.大側斜螺旋槳強度研究［D］.武漢：華中科技大學，2003.

［9］楊向暉，程爾升.大側斜螺旋槳槳葉應力測試研究［J］.船海工程，2004（1）：6－9.

［10］盛振邦，劉應中.船舶原理：下冊［M］.上海：上海交通大學出版社，2004.

［11］高同兵，程爾升.大側斜螺旋槳設計中的考慮［J］.船海工程，2004（2）：1－3.

［12］楊在林.工程力學［M］.哈爾濱：哈爾濱工程大學出版社，2010.

［13］蘇玉民，黃 勝.船舶螺旋槳理論［M］.哈爾濱：哈爾濱工程大學出版社，2003.

［14］楊在林，宋天舒，楊 勇.材料力學：下冊［M］.哈爾濱：哈爾濱工程大學出版社，2011.

［15］KOYAMA K.Comparative calculations of propellers by surface panel method－workshop organized by 20th ITTC propulsion committee［J］.Papers of Ship Research Institute，1993，16（S）：57－66.

［16］顧鋮璋，鄭百林.船用螺旋槳敞水性能與槳葉應力的數值分析［J］.力學季刊，2011，32（3）：440－443.