基于流固耦合的螺旋槳水動(dòng)力性能數(shù)值仿真

黃 勝 白雪夫 孫祥杰 陳廣杰

（哈爾濱工程大學(xué) 船舶工程學(xué)院 哈爾濱150001）

引 言

復(fù)合材料螺旋槳作為一種新型的推進(jìn)器，在船舶節(jié)能減排方面的作用日益為人們所重視［1］。利用復(fù)合材料螺旋槳的可設(shè)計(jì)性，通過合理安排槳葉的纖維方向和鋪層順序，可以使螺旋槳槳葉在水動(dòng)力載荷作用下產(chǎn)生有利于水動(dòng)力性能的變形。由于復(fù)合材料螺旋槳的特性，造成傳統(tǒng)的模型水池試驗(yàn)成本高、周期長，傳統(tǒng)的螺旋槳理論設(shè)計(jì)與計(jì)算建立在勢(shì)流理論基礎(chǔ)之上，未能全面考慮粘性的影響且不考慮旋度，因而無法準(zhǔn)確預(yù)測(cè)槳葉邊界層、螺旋槳尾流場(chǎng)的結(jié)構(gòu)及槳葉梢渦的形成等真實(shí)情況下的流動(dòng)特性。基于RANS方程的粘性流場(chǎng)計(jì)算螺旋槳的流場(chǎng)特性的方法日趨完善，黃勝等［2］分析螺旋槳在不同工作狀態(tài)下的水動(dòng)力性能。關(guān)于螺旋槳流固耦合計(jì)算方法的研究，LIN H J等［3］采用升力面法和九節(jié)點(diǎn)退化殼單元耦合算法，實(shí)現(xiàn)了求解復(fù)合材料螺旋槳的水動(dòng)力性能的算法。Young Y L［4-5］研究面元法和軟件ABAQUS耦合的螺旋槳流固耦合計(jì)算方法，但這些方法均是基于勢(shì)流理論的螺旋槳水動(dòng)力計(jì)算。

本文主要根據(jù)螺旋槳理論，結(jié)合計(jì)算流體力學(xué)方法和結(jié)構(gòu)有限元方法［6］，基于ANSYS Workbench并調(diào)用CFD求解器Fluent和有限元Mechanical APDL，分別用于螺旋槳水動(dòng)力模型和有限元模型兩者的數(shù)值模擬計(jì)算，同時(shí)將兩個(gè)求解器通過System coupling模塊進(jìn)行數(shù)值傳遞，完成螺旋槳流固耦合數(shù)值方法的構(gòu)建。

1 數(shù)學(xué)模型

1.1 流體控制方程

考慮螺旋槳在粘性湍流中旋轉(zhuǎn)，其連續(xù)性方程可表示為：

動(dòng)量方程可表示為：

式中：p是靜態(tài)壓力，Pa；μ是湍流粘度；ρ是液體密度為Reynolds應(yīng)力，Pa。

1.2 湍流模型的選取

不同的湍流模型有不同的適用性。在運(yùn)用湍流模型進(jìn)行計(jì)算時(shí)，首先要綜合考慮不同湍流模型的模擬能力以及計(jì)算所需系統(tǒng)資源，在比較分析之后選擇合適的湍流模型進(jìn)行計(jì)算。常用的湍流模型有標(biāo)準(zhǔn)k-ε模型、可實(shí)現(xiàn)的k-ε模型、標(biāo)準(zhǔn)k-w模型和SSTk-w模型。本文通過比較分析幾種常用湍流模型，通過求解時(shí)間和求解精度方面的對(duì)比以及流固耦合本身對(duì)計(jì)算機(jī)性能的要求等方面進(jìn)行綜合考慮。螺旋槳水動(dòng)力性能計(jì)算模型選擇SSTk-w湍流模型，該模型是利用混合函數(shù)將k-ε和k-w方程相結(jié)合而構(gòu)建的湍流模型，在近壁區(qū)采用k-w方程，其他區(qū)域則采用k-ε方程以獲得湍流粘性作用，考慮了k-w方程近壁區(qū)模擬時(shí)的有效性及遠(yuǎn)場(chǎng)區(qū)無法準(zhǔn)確模擬的不足［7］，最終選擇SSTk-w模型求解螺旋槳的流固耦合特性［8］。

1.3 數(shù)值離散方式

計(jì)算過程中對(duì)壓力速度耦合方程采用SIMPLEC算法進(jìn)行離散；動(dòng)量方程、湍流動(dòng)能、湍流耗散系數(shù)采用二階迎風(fēng)格式進(jìn)行離散；壓力離散格式采用PRESTO。

2 計(jì)算模型

2.1 槳的基本參數(shù)

本文流固耦合方法的構(gòu)建與驗(yàn)證均以DTMB P4381模型槳為研究對(duì)象，其水動(dòng)力模型主要采用多重參考模型（MRF模型），保持螺旋槳轉(zhuǎn)速n=600 r/min，選擇不同的進(jìn)速系數(shù)進(jìn)行水動(dòng)力計(jì)算。設(shè)置與傳統(tǒng)的CFD計(jì)算時(shí)稍有不同。表1為螺旋槳主要參數(shù)。

表1 螺旋槳主要參數(shù)

2.2 幾何模型的建立與網(wǎng)格劃分

圓柱體流場(chǎng)域使用與結(jié)構(gòu)部分模型統(tǒng)一的坐標(biāo)，圓柱體軸線與螺旋槳軸線一致。同時(shí)將整個(gè)流域分為大域與小域兩部分，有利于網(wǎng)格劃分與計(jì)算的設(shè)置。其中，大域圓柱體直徑為3D，進(jìn)口邊界距離螺旋槳2D，出口邊界距離螺旋槳5D；小域圓柱體直徑為1.1D。將入口邊界設(shè)置為速度入口，出口邊界設(shè)置為自由出流，圓柱面、槳葉面和槳轂面都設(shè)置為無滑移固壁條件。采用混合網(wǎng)格方法：大域采用全結(jié)構(gòu)化網(wǎng)格，以減少網(wǎng)格數(shù)目，加快計(jì)算速度；由于動(dòng)網(wǎng)格設(shè)置的需要，小域采用自適應(yīng)網(wǎng)格，并在葉片處對(duì)網(wǎng)格進(jìn)行加密，以提高計(jì)算精度。劃分完的網(wǎng)格如圖1所示。

圖1 流場(chǎng)網(wǎng)格劃分

本文使用的是雙向流固耦合的方法進(jìn)行計(jì)算，故壁面會(huì)產(chǎn)生位移。通過System coupling模塊將有限元求解器計(jì)算所得結(jié)構(gòu)位移從FSI-solid傳遞到FSI-fluid，在Fluent中完成邊界移動(dòng)。

為實(shí)現(xiàn)因邊界位移所產(chǎn)生的網(wǎng)格動(dòng)態(tài)變化，需要運(yùn)用動(dòng)網(wǎng)格技術(shù)。本文將Fluent提供的彈簧光順與網(wǎng)格重構(gòu)兩種動(dòng)網(wǎng)格方法相結(jié)合，成功完成了網(wǎng)格動(dòng)態(tài)變化的過程，保證了動(dòng)網(wǎng)格的質(zhì)量和計(jì)算精度。如圖2-圖3所示。

結(jié)構(gòu)部分網(wǎng)格劃分主要以四面體網(wǎng)格為主，在槳葉片處加密，以提高計(jì)算精度。劃分完成的網(wǎng)格如圖4所示，總網(wǎng)格數(shù)為120萬。

圖2 原邊界及網(wǎng)格

圖3 變化后的邊界及網(wǎng)格

圖4 結(jié)構(gòu)網(wǎng)格劃分

2.3 螺旋槳結(jié)構(gòu)有限元模型

本文不考慮螺旋槳材料的各向異性，選擇典型的金屬材料，其主要參數(shù)如表2所示。

表2 螺旋槳金屬材料主要參數(shù)

2.4 約束與載荷

由于螺旋槳的槳轂和槳葉為一體化模型，故而在設(shè)置結(jié)構(gòu)部分約束時(shí)，只需在槳轂處設(shè)置固定約束（Fixed Support）。

在進(jìn)行螺旋槳的流固耦合分析時(shí)，以槳葉片處所受載荷為主。設(shè)置流固耦合面（Fluid Solid Interface），使用System coupling模塊，將葉片所受載荷通過FSI-fluid傳遞到FSI-solid，完成葉片動(dòng)態(tài)載荷的加載。

3 計(jì)算驗(yàn)證

運(yùn)用已經(jīng)建立的螺旋槳的流固耦合瞬態(tài)方法，對(duì)金屬材料螺旋槳進(jìn)行均勻來流情況下的水動(dòng)力性能計(jì)算，并與文獻(xiàn)中的實(shí)驗(yàn)值和傳統(tǒng)的CFD方法計(jì)算結(jié)果進(jìn)行比較，結(jié)果如圖5-圖7所示。

圖5 螺旋槳推力系數(shù)

圖6 螺旋槳轉(zhuǎn)矩系數(shù)

由上面3組圖可以看出，傳統(tǒng)的CFD方法及本文構(gòu)建的FSI穩(wěn)態(tài)方法計(jì)算所得結(jié)果與實(shí)驗(yàn)值較為吻合。在均勻來流情況下的水動(dòng)力性能曲線趨勢(shì)一致，且傳統(tǒng)CFD方法與FSI方法計(jì)算所得結(jié)果差異非常小。這是由于本章是以金屬材料螺旋槳進(jìn)行方法驗(yàn)證，流固耦合作用對(duì)螺旋槳槳葉的變形及敞水性能的影響極其微小，故無法直接通過敞水性能曲線準(zhǔn)確分析FSI方法與CFD方法的差異。

將各個(gè)進(jìn)速系數(shù)J下的CFD方法與FSI方法計(jì)算所得的KT、KQ及η0同實(shí)驗(yàn)值的誤差進(jìn)行比較分析，見下頁表3。

由表3可以看出，CFD方法及FSI方法計(jì)算所得的KT、KQ及η0與實(shí)驗(yàn)值的誤差均在±5%以內(nèi)，各項(xiàng)數(shù)據(jù)誤差整體趨勢(shì)隨著進(jìn)速系數(shù)J的增加而增大。FSI方法計(jì)算誤差較CFD方法誤差略小，且在低進(jìn)速情況下，F(xiàn)SI方法計(jì)算誤差小的優(yōu)勢(shì)相對(duì)更加明顯。隨著進(jìn)速的增加，兩種方法的誤差趨于一致。這是由于在低進(jìn)速情況下，螺旋槳處于重載狀態(tài)，水動(dòng)力載荷作用于螺旋槳葉面所產(chǎn)生的變形相對(duì)于高進(jìn)速輕載狀態(tài)下所產(chǎn)生的變形更大，而FSI方法能夠更好的捕捉槳葉的變形對(duì)金屬材料螺旋槳水動(dòng)力性能所帶來的微小變化，模擬更加接近實(shí)際情況。

表3 傳統(tǒng)CFD方法與FSI方法誤差分析

圖8顯示的是不同進(jìn)速系數(shù)下螺旋槳在軸向、徑向及周向的變形，選取的進(jìn)速系數(shù)J分別為0.3、0.7和1.0。由圖8可知，在不同進(jìn)速系數(shù)下，螺旋槳葉片均產(chǎn)生與來流方向相反的變形，變形沿著葉根向葉梢不斷增大，類似懸臂梁變形，可認(rèn)為螺旋槳產(chǎn)生了縱傾分布；同時(shí)在相同進(jìn)速系數(shù)下，任意半徑處葉元體變形量沿導(dǎo)邊向隨邊逐漸變小，可認(rèn)為螺旋槳各半徑處的螺距角變小。

圖9為不同進(jìn)速系數(shù)下的螺旋槳徑向變形分布。由圖9可知，在不同進(jìn)速系數(shù)下，螺旋槳葉片在徑向產(chǎn)生變形，徑向變形在導(dǎo)邊靠近葉梢處為最大值，在葉面上朝向隨邊逐漸減小，而邊緣處達(dá)到最小值。

圖8 不同進(jìn)速系數(shù)（J） 下的 螺旋槳槳葉軸向變形分布

圖9 不同進(jìn)速系數(shù)（J） 下的螺旋槳槳葉徑向變形分布

圖10 不同進(jìn)速系數(shù)（J） 下的螺旋槳槳葉周向變形分布

圖10為不同進(jìn)速系數(shù)下的螺旋槳槳葉周向變形分布。由圖10可知，在不同進(jìn)速系數(shù)下，螺旋槳葉片均產(chǎn)生與螺旋槳轉(zhuǎn)向相反的周向變形，變形沿著葉根向葉梢不斷增大，可認(rèn)為螺旋槳產(chǎn)生了側(cè)斜。

總體上看，相同轉(zhuǎn)速下，螺旋槳的軸向變形、徑向變形及周向變形在低進(jìn)速下均大于高進(jìn)速，軸向變形與周向變形的最大值比徑向變形大1～2個(gè)數(shù)量級(jí)。與之相對(duì)的螺旋槳的縱傾變化、側(cè)斜變化以及螺距角的變化均隨著進(jìn)速系數(shù)的增加而減小，這也與前面所作的分析一致。

4 結(jié) 論

本文以DTMB P4381模型槳為研究對(duì)象，對(duì)所構(gòu)建的螺旋槳流固耦合方法進(jìn)行了分析驗(yàn)證。由驗(yàn)證結(jié)果可知，當(dāng)計(jì)算同一金屬材料螺旋槳時(shí)，本文所構(gòu)建的FSI方法在低進(jìn)速下計(jì)算所得結(jié)果要較傳統(tǒng)CFD方法更為準(zhǔn)確，而在高進(jìn)速時(shí)，兩種算法結(jié)果相對(duì)差別較小，這也與螺旋槳低進(jìn)速重載狀態(tài)下，葉片變形較大的的計(jì)算結(jié)果相吻合，從而從多角度驗(yàn)證了本文構(gòu)建的螺旋槳流固耦合方法的可靠性。

［1］ Mouritz, Gellert E, Burchill P, et al.Review of advanced composite structures for naval ships and submarines ［J］ .Composite Structures, 2001, 53：21-41.

［2］ 黃勝，王超，王詩洋.不同湍流模型在螺旋槳水動(dòng)力性能計(jì)算中的應(yīng)用與比較 ［J］.哈爾濱工程大學(xué)學(xué)報(bào)，2009（5）：484-485.

［3］ Lin H J，Lin J J.Nonlinear hydro-elastic behavior of propellers using a finite-element method and lifting surface theory ［J］.Journal of Marine Science and Technology,1996（1）：114-124.

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

［5］ Young Y L.Time-dependent hydroelastic analysis of cavitating propulsors［J］.Journal of Fluids and Structures,2007, 23 （2）： 269-295.

［6］ 孔為平，王金寶，馮學(xué)梅.螺旋槳敞水性能數(shù)值預(yù)報(bào)方法研究 ［J］.船舶，2012（6）： 13-21.

［7］ 李巍，王國強(qiáng)，汪蕾.螺旋槳粘流水動(dòng)力特性數(shù)值模擬［J］.上海交通大學(xué)學(xué)報(bào)，2007（7）：1200-1203.

［8］ 張帥，朱錫，侯海量.船舶螺旋槳流固耦合穩(wěn)態(tài)求解算法 ［J］.哈爾濱工程大學(xué)學(xué)報(bào)，2012（5）： 615-621.