导管桨是在工程类船舶采用较多的一种特种推进方式，导管在保护螺旋桨的同时还能产生推力，提高推进系统效率，改善螺旋桨空化、噪声性能，尤其是在重载状态下性能更优。研究表明，导管推力的产生与导管的剖面形状及导管内表面附着涡相关，为进一步了解导管回收尾涡能量的作用机理，深入对导管桨内流场的流动特性(特别是涡结构特性)分析十分有必要^[1]。

现阶段导管桨内流场特性主要采用数值计算方法，早期开展的研究大多基于势流理论进行数值计算^[2-3]，苏玉民等^[4]、解学参等^[5]基于面元法对导管桨的水动力性能和内流场进行了计算，势流方法计算效率高但很难准确模拟复杂流动；随着计算机性能的不断发展，基于粘流的多种数值计算方法被应用到导管桨不同方向的数值模拟中^[6-9]，胡健等^[10]、时立攀等^[11]基于雷诺平均方法(Reynolds-averaged navier-stokes, RANS)对导管桨的周围流场及尾涡特性进行了预报，孙瑜等^[12]对导管结构改进前后的导管桨水动力性能进行了数值模拟；Bhattacharyya等^[13-15]基于RANS方法对导管桨的尺度效应问题进行了研究，并探讨了不同导管类型对导管桨水动力性能的影响；Gaggero等^[16]对导管桨的空化梢涡进行了数值预报，其网格划分及加密拓扑对拟捕捉流动细节的数值模拟研究方向具有较高参考价值。

本文基于分离涡模拟(detached eddy simulation, DES)方法进行导管桨敞水试验数值模拟，通过与模型试验结果相比验证数值计算精度与可信度，叙述数值模型的网格类型、湍流模型选择及计算方法等，重点分析导管桨内流场典型物理量分布特征，研究导管桨内流场涡系组成，分析导管内壁面附着涡具体形态与分布特性。

1 几何模型与计算区域

数值计算几何模型为某动力定位系统导管桨模型，螺旋桨为定螺距4叶桨，桨模直径D=0.240 m(半径R=0.120 m)，毂径比为0.190，盘面比为0.417，导管剖面为荷兰MARIN设计的No.19A导管(图 1(a))，叶梢间隙为2 mm。参考坐标系为笛卡尔坐标系O-xyz，坐标系O原点位于螺旋桨的几何中心，顺流方向为x轴正向，y轴、z轴正向按右手定则确立。

数值计算域为圆柱形，圆柱直径为6D，数值计算域向上游延伸2D，向下游延伸8D距离以避免产生限制效应。螺旋桨的旋转运动通过滑移网格技术实现，计算域被分为静态域和旋转域2个子域。计算域边界条件设定如下：进口及圆柱形外边界为速度入口，出口为压力出口，螺旋桨及导管模型表面为非滑移壁面，静态域和旋转域的交界面设定为interface以实现2个子域内数据交换和迭代。

2 数值模型建立 2.1 湍流模型与离散方法

数值模拟为三维黏性、不可压缩、常密度的非稳态数值计算，应用滑移网格技术实现导管桨的敞水性能试验数值模拟。分离涡模型(DES)是一种混合模型^[17]，其应用RANS方法处理近壁区湍流，应用大涡模拟方法(large eddy simulation, LES)求解距离壁面较远的大涡结构，因此DES方法能用相对较少的计算资源提高对小尺度涡结构的计模拟精度。所有的DES模拟中均应用Spalart & Allmaras一方程湍流模型用于流场中湍流结构求解。

基于Spalart & Allmaras一方程模型的DES方法简要介绍如下，湍流粘性系数通过设置辅助变量求解，代数方程为：

$ \nu_{t}=\tilde{\nu} f_{\nu 1} $

(1)

$ f_{\nu 1}=\chi^{3} /\left(\chi^{3}+C_{\nu 1}^{3}\right) $

(2)

$ \chi=\tilde{\nu} / \nu $

(3)

式中C_ν1³为常数。辅助变量根据涡粘性输运方程求解：

$ \begin{array}{c}{\frac{\partial \tilde{\nu}}{\partial t}+u_{j} \frac{\partial \tilde{\nu}}{\partial x_{j}}=\frac{1}{\sigma}\{\nabla \cdot[(\nu+\tilde{\nu}) \nabla \tilde{\nu}]+} \\ {C_{b 2}\left|\nabla_{\nu}\right|^{2} \}+C_{b 1}\left[1-f_{t 2}\right] \tilde{S} \tilde{\nu}-} \\ {\left[C_{w 1} f_{w}-\frac{C_{b 1}}{\kappa^{2}} f_{l 2}\right]\left(\frac{\tilde{\nu}}{d}\right)^{2}+f_{t 1} \Delta U^{2}}\end{array} $

(4)

式中：等号右侧波形括号内第1项为扩散项，第2项诱导层流区内诱导湍流向前传播，第3项和第4项代表湍流的生成项和耗散项，最后一项用于在特定位置触发湍流。诸多C项和σ项均为常数，冯卡曼常数κ=0.41，f是壁面距离d的函数，$ \tilde S$是d和$S = \sqrt {{\mathit{\Omega} _{ij}}{\mathit{\Omega} _{ij}}} $的函数，其中Ω_ij为旋转张量。

定义一个新的壁面距离$ {\tilde d}$来区分RANS求解区和LES求解区^[18]此来避免LES方法在近壁面区域难以求解的问题。$ {\tilde d}$的定义为：

$\tilde{d}=\min \left(d, C_{\mathrm{DES}} \mathit{\Delta}\right) $

(5)

式中：C_DES是常数；Δ=max(Δ_x, Δ_y, Δ_z)。

DES模拟采用SIMPLE算法完成压力-速度耦合方程的求解。其中，对流项通过二阶迎风格式进行离散，扩散项通过中心差分格式进行离散，时间项二阶隐式格式离散，均采用Spalart-Allmaras湍流模型封闭N-S方程组。

2.2 网格划分

计算域网格类型为混合网格，其中静态域为切割体网格，旋转域为多面体网格。数值计算域整体网格尺寸需逐渐过渡以避免发生较大的数值耗散^[19]叶梢和导管之间的计算属于近壁面相对运动问题，因此需要进行更为精细的离散处理以提升网格质量。叶梢区域内设定包含叶梢的环状圆柱进行区域体加密处理，同时，叶梢与导管间隙区域设定间隙加密区域；叶片表面通过设定导边和随边特征线加密以控制体网格生成；螺旋桨和导管物面设置8层棱柱层网格完成边界层模拟，为满足DES算法要求，壁面y⁺值控制在(-∞，1)^[20]。最终确认计算域网格总数约为951万(其中静态域736万，旋转域215万), 导管桨部件面网格划分结果见图 1(b)。

3 计算结果与讨论 3.1 敞水性能计算结果验证

导管桨的水动力系数定义如下：

$ \left\{\begin{array}{l}{J=V_{x} / n D} \\ {K_{t p}=T_{p} / \rho n^{2} D^{4}} \\ {K_{t d}=T_{d} / \rho n^{2} D^{4}} \\ {K_{q}=Q_{0} / \rho n^{2} D^{5}} \\ {\eta_{0}=\frac{J}{2 \pi} \frac{K_{t p}+K_{t d}}{K_{q}}}\end{array}\right. $

(6)

式中：J为进速系数；V_x为进流速度；n为转速；D为螺旋桨直径；T_p为螺旋桨推力；T_d为导管推力；Q₀为螺旋桨扭矩；ρ为水密度；K_tp、K_td为螺旋桨、导管的推力系数；K_q为导管桨扭矩系数；η₀为导管桨效率。

数值计算工况选定J为0.3, 0.4, 0.5这3个进速系数，转速设定为n=11 r/s，相应的雷诺数范围为3.5×10⁵~3.6×10⁵，时间步长设定为1°。

导管桨模型前期已经基于船模拖曳水池进行过敞水试验，数值计算结果与敞水试验结果对比如表 1所示。从对比结果来看，导管推力系数K_td的最大相对偏差达到4.42%，螺旋桨周期性旋转导致导管推力产生明显脉动；螺旋桨的推力系数、扭矩系数和试验结果吻合良好，偏差均在3%以内。上述结果说明数值计算模型的数值计算精度较高，后续关于内流场特性及涡特性的分析具备一定的可信度。

3.2 内流场特性分析

选取计算域内典型的剖面进行内流场特性分析，其中，横剖面(y-z平面)分别选定桨前盘面(x=-0.5R)和桨后盘面(x=0.3R)，剖面坐标系均基于螺旋桨半径进行无因次化。后文详细地对上述盘面内瞬态的轴向速度场分布、涡量分布、压力分布等进行分析，选定分析工况为进速系数J=0.5。

导管的存在增强了螺旋桨对桨前区域的诱导作用，桨前盘面内水流的速度场明显发生改变。从径向变化规律来看，导管桨桨前盘面(x=-0.5R)的轴向诱导速度(V_x-in=(V_x-V₀)/V₀)沿径向呈现先增加再减小的趋势(图 2)，这与导管壁面和桨毂对水流流动存在阻滞作用有关。一周范围内盘面不同半径处(0.4R, 0.6R, 0.8R)轴向诱导速度分布如图 3所示，从周向变化规律来看，V_x-in沿周向呈现出明显的周期性，轴向诱导速度峰、谷值数与桨叶数一致，图 3中可见，0.8R位置轴向诱导作用波动最为强烈，整体来看，越靠近盘面外侧轴向诱导速度沿周向的波动幅度越大。

桨前水流经由桨叶加速后，导管桨内流场产生剧烈变化且涡特征明显。从桨后盘面(x=0.3R)的瞬态轴向速度分布(图 4(a))来看，桨叶对尾流的加速作用集中体现在0.5~0.8R，无因次轴向速度V_x/V₀分布在1.4~1.8，由于粘性作用，导管壁面和桨毂附近区域流体以一定的速度梯度流动。桨后盘面内涡量分布(图 4(b))能体现出导管内流场涡系的组成成分：导管内壁面的剪切层流产生剪切层涡，剪切层涡在横剖面内体现为贴近导管内壁面的圆环形状，螺旋桨产生的涡系由外到内依次为叶梢梢涡(图 4(b)(Ⅰ))、叶片脱落涡(图 4(b)(Ⅱ))及叶根涡(图 4(b)(Ⅲ))，漩涡涡核分布于涡强最大值附近；湍流粘度反映涡扩散运动，从分析结果来看湍流粘度值在梢涡分布区、叶梢及叶根区较大，湍流粘度分布特性与涡分布特性保持一致。

导管表面的瞬态压力分布如图 5所示，导管内表面叶梢前后附近存在压力峰值与谷值点。在导管内壁面x=0截面处设定探点P₁(图 5)，监测在导管桨正常工作过程中P₁处压力脉动变化。从时域范围内分析结果来看(图 6(a))，螺旋桨旋转一周的过程中，探点P₁处经历4个压力峰值与谷值，导管内表面压力脉动幅度大，流场变化剧烈，探点处的压力值与桨叶的位置直接相关；从频域来看(图 6(b))，压力脉动在叶频(BPF)处峰值最高，整个频域范围内多倍叶频处连续出现峰值，说明导管内壁面压力脉动特征与叶频相关。

3.3 内流场涡特性分析

导管桨内流场涡特性及其分布规律直接影响导管桨的水动力性能，选定水平剖面(z=0)分析内流场涡量沿流向的分布规律。如图 7所示，导管与螺旋桨之间的强相互作用导致叶梢间隙附近及导管内壁面区域(Ⅰ框内)分布大量泄出涡，螺旋桨叶根涡(Ⅱ框内)和毂涡(Ⅲ框内)的分布相对清晰，导管剪切层涡与内流场泄出涡在导管尾缘产生相互复杂干扰，尾涡能量相互传递并向下游传输，导管桨近尾流场中涡量出现连续椭圆形分布(Ⅳ框内)。

通过体渲染方式对导管桨内流场三维空间涡量分布进行分析，图 8中可以清晰分辨出螺旋桨叶片随边处产生的片状脱落涡，脱落涡连接叶梢涡与叶根涡，涡量大小随着向尾流中扩散迅速降低，不同桨叶叶根涡呈螺旋状发展并逐步汇集于桨毂处；叶梢涡分布于导管内壁面区域，梢涡的形态明显产生变形，而且还与导管内壁面剪切层涡发生融合与干扰，因此梢涡涡能的分布发生明显变化，涡量更大面积地分布于导管内壁面，引起导管表面产生更强的环流进而产生导管推力。

通过局部流动特性分析可以对变形梢涡结构内部特征进行深入辨析，图 9是z=0平面内梢涡截面内流线结果。黑色虚线为涡量/s^-1=100 s^-1等值线，其外部轮廓与图 8中一致，在投影视角上呈现“W”形，流线分布可以看到梢涡管内部诱导速度存在2个漩涡中心，这与独立螺旋桨梢涡涡量的分布特征存在明显差异。进一步读取y/R=1位置处涡管范围(0.2 < x/R < 0.5)内的涡量大小(图 10)，发现2个漩涡中心涡量沿顺流方向逐渐下降。

4 结论

1) 导管内壁区域流动复杂，导管内壁面压力脉动特性与螺旋桨叶频相关；

2) 导管桨尾涡系包含导管剪切层涡、叶梢涡、叶片脱落涡、叶根涡和毂涡，叶片脱落涡向下游传输过程中能量迅速扩散；

3) 导管的存在直接影响螺旋桨梢涡的分布，梢涡形态发生变化，涡量更大面积的分布于导管内壁面，导管有助于螺旋桨尾流的涡能回收从而提高系统推进效率。