王明哲，王建華，萬德成

（上海交通大學a.船海計算水動力學研究中心（CMHL）；b.船舶海洋與建筑工程學院，上海 200240）

0 引 言

在進行螺旋槳與其他結構物相互作用的數(shù)值模擬中，螺旋槳模型因為復雜的旋轉運動以及復雜的幾何結構，往往會需要較長的計算時間和較大的存儲誤差成本。體積力模型通過將螺旋槳與流體的相互作用力?；癁榉植荚诹鲌鲋械捏w積力，作為源項直接參與動量方程的求解，是代替真實螺旋槳模型的減小計算成本的有效手段，在船槳相互作用、槳舵相互作用以及船舶操縱性能模擬等方面得到了大量的應用。

作為一種簡化模型，體積力模型的復雜性與準確性是影響其應用性的關鍵因素。葉素動量理論模型作為一種復雜性與準確性適中的體積力模型，具有較大的應用潛力[1-3]。

使用葉素動量理論計算得到的螺旋槳性能與真實螺旋槳模型性能接近，但難以在不同工況下保證相同的準確性。Benini[4]研究了使用葉素動量理論模擬船用螺旋槳的適用性，其中葉元參數(shù)由基于勢流理論的二維葉元性能計算軟件XFOIL 計算，結果表明使用葉素動量理論計算螺旋槳性能的最小誤差為2%，但誤差對進速系數(shù)的變化較為敏感；Phillips[1]使用基于葉素動量理論的體積力模型進行了KVLCC 的PMM 試驗數(shù)值模擬，忽略了螺旋槳載荷分布的周向不均勻性，最后發(fā)現(xiàn)在小舵角工況下船體側向力與艏搖力矩預測誤差在（2～3）%以內，而當舵角增加時預報誤差增大。種種研究表明，可以通過修正手段提高葉素動量理論體積力模型的準確性。

葉元參數(shù)準確性與葉元水動力計算方法是影響葉素動量理論體積力模型的關鍵因素。對于葉元參數(shù)，已存在部分經驗公式方法用于修正非設計工況下的螺旋槳葉元參數(shù)，但無法保證此修正方法對于不同槳型及工況的適用性[5-6]。葉元水動力計算方法包括直接使用當?shù)厮俣纫约暗嬎阏T導因子。馮大奎等[7-10]基于KP505 真實螺旋槳模型的敞水試驗模擬結果建立了葉元參數(shù)與當?shù)厮俣鹊年P系，分析了不同葉數(shù)下葉元的性能變化，進行了KCS 模型尺度與實尺度在設計航速下的自航試驗模擬，模擬結果表明體積力模型與真實螺旋槳得到的自航因子結果吻合；Tokgoz[11]通過葉元參數(shù)的一般表達式與直接使用當?shù)厮俣确椒ㄟM行了AU螺旋槳敞水與斜流試驗模擬，結果表明敞水條件下得到的螺旋槳扭矩誤差較大，且進速系數(shù)越低準確性越低。

傳統(tǒng)的誘導因子迭代方法基于理想流體的動量理論推導得出，未考慮粘性求解器中的粘流環(huán)境[2]。Villa[12]擬合了單個進速工況下槳葉周向平均載荷的徑向分布函數(shù)，其中槳葉載荷由敞水實驗得到，并將此載荷分布作為體積力分布與RANS 求解器結合，結果表明該體積力模型與真實螺旋槳模型作用下的尾流速度分布吻合程度高，證明了按照真實螺旋槳載荷分布的體積力模型可以得到預期的誘導速度場。但此種方法需要已知工況，無法處理來流速度變化的非定常及非均勻流場。傅慧萍[13]使用HO模型模擬了KCS船自航，發(fā)現(xiàn)進速與槳盤面當?shù)厮俣冉咏€性關系。此研究結果表明在得到當?shù)厮俣鹊那闆r下，可以通過此對應關系根據槳盤面當?shù)厮俣确赐七M速，進而得到粘流條件下槳盤面各處誘導因子大小。

本文通過不同進速下螺旋槳敞水試驗的CFD 數(shù)值模擬，得到在粘性作用與葉間干擾作用下的螺旋槳各葉元參數(shù)，并依據此葉元參數(shù)得到符合真實螺旋槳載荷徑向分布的投影盤，進而得到此投影盤作用下各進速的敞水誘導速度場以及槳盤面處各位置的當?shù)厮俣?，最終通過當?shù)厮俣扰c進速計算誘導因子分布。此誘導因子分布將植入葉素動量理論體積力模型，避免了誘導因子的迭代求解，并通過KP505螺旋槳真實模型與體積力模型的敞水試驗驗證本文的體積力改進方法。

1 數(shù)值方法

1.1 控制方程

本文使用改進的葉素動量理論體積力模型與真實螺旋槳進行螺旋槳敞水試驗的數(shù)值模擬，使用的基礎數(shù)值求解器為開源平臺OpenFOAM 中的pimpleFoam 與pimpleDyMFoam，其中pimpleDyMFoam用于真實螺旋槳模型數(shù)值模擬，pimpleFoam 則用于實現(xiàn)體積力模型代碼在粘流求解器中的植入。對于粘性不可壓縮流體，控制方程為

式中，p為壓力，ρ為流體密度，μ為粘性系數(shù)，ui為速度分量為雷諾應力。本文使用湍流模型封閉含有雷諾應力項的方程。對于最適用于螺旋槳敞水試驗數(shù)值模型的湍流模型仍無定論[14]，文中所使用的為k-ωSST 湍流模型，這種湍流模型針對根據各區(qū)域到壁面距離的變化在對應位置計算時激活不同的湍流變量控制方程形式，是一種應用性較廣的方法[15]。(fε)i為體積力源項，由單獨的葉素動量理論體積力代碼計算得到[16]，在使用體積力模型進行數(shù)值模擬時激活此項。體積力模型與真實螺旋槳模型均使用非定常模擬，用以驗證體積力模型在處理非定常問題時的計算能力。速度場與壓力場的解耦使用PIMPLE 方法，體積力源項的計算采用顯式方法，即首先利用上一時間步的速度場計算體積力分布，再使用得到的體積力分布計算當前時間步的速度場。

1.2 葉素動量理論

每次完成控制方程的求解后，將進行一次葉素動量理論程序的計算以求取體積力分布。文中各葉元體作用于流場中稱為致動點的位置，致動點的位置不隨流場的迭代而發(fā)生改變，即螺旋槳轉速僅為計算葉元體水動力性能的參量。致動點處的葉元體水動力在每個時間步根據瞬時當?shù)厮俣扔嬎?，且各致動點互不影響，因此對于船槳配合的自航問題中螺旋槳載荷徑向與周向分布不均勻的問題，此體積力模型也可以進行有效處理，使用致動點的軸向當?shù)厮俣萔x與周向當?shù)厮俣萔θ依照下式計算葉元體的水動力，圖1展示了葉元體水動力計算各相關物理量的矢量關系：

圖1 葉素動量理論Fig.1 Blade element momentum theory

式中：a與b分別表示軸向誘導因子和切向誘導因子，表征螺旋槳引起的誘導速度與來流速度的比值；n為螺旋槳轉速；r為葉元體徑向位置,；Vx與Vθ分別表示當?shù)剌S向速度和當?shù)厍邢蛩俣?，即螺旋槳擾動產生的誘導速度與來流速度的和；Va與Vt分別表示大地坐標系下的軸向來流速度和切向來流速度，表征無窮遠處流體相對于葉元體的速度；β為進角；φ為葉元體螺距角；α為幾何攻角,；CL與CD分別表示葉元體的升力系數(shù)和阻力系數(shù)；cr為對應徑向位置處的葉元體弦長；dr為葉元體展向長度；dFL與dFD分別代表葉元體的升力和阻力。

1.3 葉元參數(shù)

葉元參數(shù)為葉元體升阻力系數(shù)與攻角的對應關系，在使用體積力模型前存儲到特定文件中，并由體積力計算程序實時讀取。本文通過對敞水試驗CFD 模擬結果中的葉元壁面進行壓力和粘性力積分，可以分別得到葉元體的升力和阻力系數(shù)，進而確定各葉元體升阻力系數(shù)與攻角的對應關系。在執(zhí)行體積力程序時，葉元性能由分段線性插值得到對應攻角下葉元體的升阻力系數(shù)。具體的實施流程如下：首先進行螺旋槳敞水試驗的CFD 數(shù)值模擬，并將螺旋槳壁面沿徑向等距分為多個葉元，各葉元體對推力與扭矩的貢獻dT與dQ通過表面壓力積分得到（圖2），進而計算葉元體的升阻力系數(shù)，并將所得力系數(shù)與攻角的對應關系保存為數(shù)據向量，由體積力計算程序通過分段線性插值方法讀取。在體積力程序的執(zhí)行過程中，首先通過誘導因子與當?shù)厮俣鹊玫饺~元來流，進而計算葉元攻角α，最終通過攻角數(shù)據向量α與升阻力系數(shù)向量CL、CD經分段線性插值得到式（8）與式（9）中的升阻力系數(shù)。分段線性插值格式如式（10）、式（11）所示，此種插值方式可盡量減小數(shù)據處理造成的誤差。

圖2 螺旋槳葉元體Fig.2 Blade element of propeller

1.4 真實螺旋槳誘導因子分布

誘導因子體現(xiàn)了來流速度與當?shù)厮俣鹊年P系，如圖3 所示。傳統(tǒng)葉素動量理論中誘導因子的迭代求解基于理想流體理論[2，4]，不適用于粘性求解器中的體積力模型，因此需要一種方法來研究粘性流場螺旋槳載荷作用下各葉元體來流與當?shù)厮俣鹊年P系。為得到這種關系，本文利用1.3節(jié)得到的葉元參數(shù)得到了多組來流速度下與真實螺旋槳模型載荷分布相符的體積力分布，并將此體積力分布作用于敞水流場中得到誘導速度場，如圖4所示。在得到誘導速度場的過程中，葉元體的水動力大小不會隨流場的變化而迭代。在此種情況下，最終得到的誘導速度場與真實螺旋槳的誘導速度場吻合程度較高。各葉元體的當?shù)厮俣扔捎嬎惴€(wěn)定后的速度場在致動點位置處插值得到，并由同一徑向位置處各致動點當?shù)厮俣鹊钠骄Y果代表此徑向位置處的葉元體當?shù)厮俣取Ｗ罱K由預設來流速度V與各徑向位置處當?shù)厮俣萔x計算各徑向位置處的誘導因子大小。由于敞水試驗中橫向流動較弱，因此本研究中尚未考慮切向誘導因子。

圖3 理想流體理論[2]Fig.3 Ideal fluid theory[2]

圖4 基于真實螺旋槳載荷分布的誘導速度場Fig.4 Induced velocity field based on real propeller load distribution

經由上述過程，可以建立單個工況下誘導因子與徑向位置的關系，也可建立同一徑向位置處誘導因子與當?shù)厮俣鹊年P系。將各徑向位置處誘導因子數(shù)據記為向量a，當?shù)厮俣扔洖橄蛄縑x并存儲到特定文件中。在體積力模型計算過程中，通過插值得到致動點位置處的當?shù)厮俣群?，通過在此數(shù)據向量中進行分段線性插值可得到各徑向位置處的誘導因子，如式（12）所示。由于建立了分段線性插值關系，在流場迭代求解中插值數(shù)據也將經過反復迭代，體積力模型能否通過迭代收斂到螺旋槳CFD數(shù)值模擬的載荷分布是本文提出的體積力模型是否具有應用性的關鍵。

2 算例設置

2.1 模型與工況

本文進行了KP505 螺旋槳真實模型及體積力模型的敞水試驗數(shù)值模擬。圖5 為KP505 螺旋槳幾何模型。真實螺旋槳模型敞水試驗轉速為10 r/s，包含0.4、0.5、0.6、0.7、0.8 五個進速系數(shù)。擬合誘導因子的體積力模型工況參數(shù)與之相同。KP505螺旋槳主要尺度及擬合誘導因子工況參數(shù)如表1所示，體積力模型直徑及轂徑與之相同。

表1 KP505螺旋槳主要參數(shù)及敞水試驗工況參數(shù)Tab.1 Main parameters of KP505 propeller and open-water test

圖5 KP505螺旋槳Fig.5 Geometry of KP505 propeller

2.2 計算域及網格劃分

圖6為真實螺旋槳模型與體積力模型敞水試驗模擬計算域大小。真實螺旋槳模型及體積力模型的計算域大小相同，均為圓柱形計算域，直徑為5D，槳盤面至入口為3D、至出口為5D。圖7 為螺旋槳敞水試驗網格劃分，真實螺旋槳模型算例網格量為251 萬，時間步長取0.0005 s。體積力模型算例網格量為46 萬，時間步長取0.001 s。真實螺旋槳模型采用滑移網格方法，滑移面附近及螺旋槳壁面均進行逐級加密。圖7（b）中圓柱區(qū)域為體積力作用區(qū)域，徑向網格數(shù)量為38 個，周向網格數(shù)量為104個。網格無關性驗證與時間步長無關性驗證已由之前相關工作完成[17]。

圖6 KP505螺旋槳敞水試驗模擬計算域Fig.6 Computational domain of KP505 propeller open-water test simulation

圖7 KP505螺旋槳敞水試驗模擬網格Fig.7 Grids of KP505 propeller open-water test simulation

3 結果分析

圖8為改進的體積力模型與真實螺旋槳模型敞水曲線模擬結果，其中real、bf與exp分別代表真實螺旋槳模型和體積力模型的計算結果以及螺旋槳敞水性能的物理實驗結果。可以發(fā)現(xiàn)三者有較好的吻合程度，表明真實螺旋槳的CFD具有準確性。同時表明在應用了符合真實螺旋槳模型的葉元性能與誘導因子分布后，在不使用翼型理論經驗公式的情況下，體積力模型的準確度得到了有效提高。圖9 展示了各進速條件下體積力模型相對于真實螺旋槳敞水曲線計算結果的誤差，除J=0.4 時推力系數(shù)誤差為1.05%外，各進速條件下推力系數(shù)誤差與扭矩系數(shù)誤差均在1%以內，其中扭矩系數(shù)誤差均在0.8%以內，且誤差水平未隨進速系數(shù)出現(xiàn)較大變化，表明改進后的體積力模型誤差水平不再與進速系數(shù)敏感[4，7，11]，其準確性可適用于不同工況。兩種模型模擬條件下槳載荷在其中三種進速條件下的徑向分布如圖10 所示。其中axial 與tangential 分別代表軸向載荷與切向載荷，并進行了無因次化處理。兩種模型數(shù)據的吻合程度較高，表明改進后的葉素動量理論模型可以通過不同徑向位置處的當?shù)厮俣鹊玫椒险鎸嵚菪龢d荷分布的葉元水動力穩(wěn)定值。

圖8 敞水曲線模擬結果Fig.8 Open-water-curve simulation results

圖9 體積力模型敞水曲線誤差分布Fig.9 Error distribution of open-water-curve of body-force model

圖10 螺旋槳載荷分布模擬結果Fig.10 Simulation results of propeller load distribution

圖11 與圖12 分別為J=0.5 與J=0.7 條件下兩種模型的尾流場軸向速度對比，其中圖12 為槳盤面后0.5D、1D與1.5D處的位置數(shù)據，此分布范圍囊括了一般槳舵干擾工況下舵模型的位置[17]。真實螺旋槳模型的尾流場數(shù)據為一個旋轉周期內的時均數(shù)據。誘導速度峰值分布在0.6R～0.8R范圍內，與圖10中最大槳載荷位置相同?？梢园l(fā)現(xiàn)誘導速度分布峰值附近的較大范圍內兩種模型的速度分布在兩種工況下的吻合程度均較高，表明體積力模型可以模擬真實螺旋槳尾流不同截面中的動量輸運。在輪轂處兩組數(shù)據存在一定差距的原因為在進行體積力模型模擬時未考慮槳軸及輪轂壁面的建模，因此忽略了真實螺旋槳模型槳軸與輪轂的阻塞效應，導致體積力模型誘導速度場輪轂附近速度偏高。此差異隨流場的向后發(fā)展逐漸變小，且對于螺旋槳尾流動量輸運不起主要影響作用。同時從圖10中可以發(fā)現(xiàn)，槳轂附近的葉元體載荷為整個槳葉中的最低水平，對螺旋槳性能不起主要作用，因此槳轂存在對螺旋槳性能的影響可以忽略。兩種模型尾流分布存在細微差異的另一原因是真實螺旋槳壁面的擾動作用產生了大量復雜渦結構，而體積力模型由于忽略了葉片旋轉效應，渦結構更加平滑有序。圖13為兩種模型計算得到的渦結構（Q=50等值面）對比。體積力模型無法還原較為精細的流場信息，但可以通過替代真實螺旋槳模型有效減少螺旋槳尾流動量輸運模擬的計算成本。

圖11 尾流場模擬結果Fig.11 Simulation results of wake flow field

圖13 體積力模型與真實螺旋槳模型的渦結構Fig.13 Vortex structure of body-force model and real propeller model

4 結 語

本文通過擬合誘導因子的方式，將誘導因子與當?shù)厮俣鹊膶獢?shù)據存儲到葉素動量理論的執(zhí)行過程中，使用螺旋槳敞水CFD數(shù)值模擬計算葉元性能，并將改進后的體積力模型應用于五組進速下的KP505螺旋槳敞水試驗數(shù)值模擬?；诟倪M葉素動量理論體積力模型模擬的KP505螺旋槳敞水試驗結果表明，體積力模型可以從初始流場迭代到螺旋槳CFD 數(shù)值模擬的工作狀態(tài)，其螺旋槳性能、載荷分布以及尾流分布均與螺旋槳直接數(shù)值模擬的結果相近，且不同進速下均有較好的預報精度，敞水誤差在五種進速下均未超過1.2%。

以上研究表明，基于符合真實螺旋槳誘導因子分布改進后的葉素動量理論模型，適用于不同進速下的螺旋槳敞水性能預報，克服了以往葉素動量理論體積力模型的準確性隨模擬工況變化敏感度高的問題。在不均勻來流問題如船槳干擾問題中，槳盤面各處的來流速度不同，特別是在船舶操縱條件下螺旋槳的入流條件更加復雜，因此在較大進速范圍內保證葉素動量理論體積力模型的計算準確性對于解決此類問題具有較大的意義。本文提出的改進體積力模型預報得到的螺旋槳的誘導速度場與真實螺旋槳吻合程度較好，因此可以較為準確地代替真實螺旋槳進行槳與其他結構物的相互作用模擬，在降低計算時間與存儲空間需求基礎上還能較好地提升計算精度，為未來船槳舵干擾下的水動力性能高效評估提供了必要的技術支撐。