船舶操纵性与其航行的安全性和经济性密不可分．近年来，随着船舶大型化发展，人们在船舶设计的过程中越来越注重评估船舶的操纵性能．相比于单桨船，双桨船具有更好的操纵性能，但其在操纵运动的过程中船尾流场复杂，左右舷的螺旋桨载荷和舵效有着明显的差异，因此须要考虑操纵运动下螺旋桨的水动力性能，从而提高双桨船操纵预报的准确性，同时可以为操纵运动下双桨船螺旋桨载荷的分配提供依据．近年来CFD(计算流体力学)被广泛运用于船舶操纵运动预报[1]，基于CFD的船舶操纵运动试验数值模拟主要分为约束模操纵运动试验数值模拟和自航船舶操纵运动数值模拟[2]．文献[3]对DTMB 5415M进行约束模试验数值模拟，并对MR(多次工况拟合)法和SR(单次工况拟合)法的计算结果进行了对比；文献[4]采用约束模试验数值模拟对ONRT在操纵状态下的船、桨和舵各部分水动力性能进行了系统研究；文献[5]采用实桨法在网格较为粗糙的情况下进行了DTMB5415M的回转运动和Z形操纵运动数值模拟；文献[6]基于HUST-Ship对DTMB5415M和KCS分别进行了回转运动和Z形运动的数值模拟．值得注意的是，以上预报方法均要模拟螺旋桨的旋转，虽然采用实桨法能真实反映船-桨-舵之间的相互作用，但须要耗费很大的计算资源，而采用螺旋桨体积力法可以在满足一定精度的条件下大幅度提高计算效率．在船舶操纵运动过程中，船尾伴流的不均匀性对桨推力、侧向力有着较大的影响，须要充分考虑．文献[7]通过螺旋桨敞水性能计算建立数据库体积力模型，将其用于某集装箱船的约束模试验数值模拟；文献[8-9]采用三种螺旋桨体积力法对某双桨单舵船进行了回转运动预报，分析了船桨舵系统各部分水动力的大小．总体来看，通过自航船舶操纵运动数值模拟对船舶操纵运动直接进行预报更准确，但一次仅能求取一个结果；而采用约束模试验数值模拟结合操纵运动方程预报船舶的操纵性，在求取水动力导数后可以求得所有工况的结果，但求取完整的水动力导数耗费大量的时间，远高于一次自航船舶操纵运动数值模拟用时．本研究采用商业软件STAR CCM+，以DTMB5415M为研究对象，首先根据约束模试验数值模拟建立船舶操纵运动下的螺旋桨修正体积力模型，接着采用螺旋桨修正体积力模型和原螺旋桨体积力模型分别对其进行操纵运动数值模拟研究，最后对两种方法的预报结果进行了对比．1 数学模型和数值求解1.1　目标船型研究对象为DTMB5415M，其垂线间长Lpp为4.002 m，船宽B为0.537 m，吃水T为0.173 m，排水体积∇ 为0.189 m3，方型系数CB为0.507，设计航速U为2.59 m/s．图1给出了几何模型示意图及本研究采用的大地坐标系o0-x0 y0 z0和随船坐标系o-xyz．10.13245/j.hust.230410.F001图1DTMB5415M几何模型和坐标系1.2　螺旋桨体积力模型1.2.1　原螺旋桨体积力模型STAR CCM+中提供了一种螺旋桨体积力法，以桨前方某个截面处的轴向平均流入速度作为进速带入螺旋桨敞水性能曲线以求得螺旋桨水动力，但其未考虑船尾伴流的不均性．该螺旋桨体积力模型在桨盘面处推力和扭矩的分布形式可参考文献[10]．1.2.2　螺旋桨修正体积力模型对于船舶操纵运动，由于漂角β和转首角速度r'的存在，螺旋桨的入流速度也会发生复杂的变化，其水动力性能和在均匀来流中会有着明显的差异，仅按照原体积力中的螺旋桨水动力求取方法并不能反映伴流不均匀性的影响．本研究通过DTMB5415M的约束模试验建立船舶操纵运动下的螺旋桨修正体积力模型．1.3　数值理论方法本研究采用工程中应用最广泛的RANS(雷诺平均方程)方法对船体绕流场进行求解，湍流模型采用SST k-ω湍流模型，自由面捕捉求解采用VOF(流体体积函数)法．船体的六自由度运动通过DFBI(动态流体物体相互作用)模块计算获得，采用重叠网格法处理船-舵之间的相对运动及船体在自由液面附近的大幅度运动．采用分离流计算模型和SIMPLE(压力耦合方程组的半隐式方法) 算法对速度项和压力项进行分离求解．其中时间瞬态项采用一阶离散，对流项采用二阶迎风离散，扩散项采用中心差分格式．当模拟实桨的旋转运动时，采用滑移网格法，为保证求解精度，时间步长对应于螺旋桨旋转2°所需时间；采用螺旋桨体积力法时所取时间步长参考ITTC(国际拖曳水池会议)给出的建议，取为0.005Lpp/U．1.4　计算域设置计算域共分为背景域、船体域、舵域和螺旋桨域4个部分，计算域的具体大小和边界条件设置如图2所示．若采用螺旋桨体积力法，则计算域中仅有背景域、船体域和舵域．10.13245/j.hust.230410.F002图2计算域示意图1.5　网格划分采用切割体网格对计算域进行网格划分，在船体壁面附近生成棱柱层网格．为比较准确地捕捉船体周围的绕流场，须要对船体附近和开尔文波附近的区域进行加密，最终全附体DTMB5415M的网格数量为5.91×106，船体表面第一层网格使得无因次壁面距离y+值在60~80，船体周围的网格如图3所示．若采用螺旋桨体积力法，则仅须在螺旋桨力源处进行适当网格加密，该工况的网格数量约为4.65×106．另外本研究还进行了敞水桨的水动力计算，具体设置和上述类似，在此不再赘述．10.13245/j.hust.230410.F003图3船体周围网格2 不同工况下螺旋桨水动力性能验证进行不同工况下螺旋桨水动力性能的验证有两个目的，其一是验证数值模拟方法的有效性，其二是得到一些基本规律，为操纵运动下的螺旋桨修正体积力模型建模提供依据．2.1　敞水桨水动力性能为方便谈论来流方向对螺旋桨水动力的影响，定义进速系数Ji、无因次水动力KTi和力矩KQi及螺旋桨敞水效率η0，有Ji=UinD;KTi=Tiρn2D4;KQi=Qiρn2D5;  η0=KTxKQxJx2π, (1)式中：i表示轴方向，可取值为x，y和z，分别对应于轴、横和垂向；Ti和Qi为对应方向上的力和力矩；Ui为对应方向上的速度；n为螺旋桨转速．计算中n=20 r/s，改变Ux使得 Jx从0.2~0.8变化．求得敞水性能曲线如图4所示 ．10.13245/j.hust.230410.F004图4螺旋桨敞水性能曲线从图4中可以看到：计算值和试验值基本符合，KTx和KQx的误差在5%以内，其中当Jx在0.6附近时计算值和试验值最接近．考虑到实际操纵运动下螺旋桨Jx在0.5~0.7范围内，故认为该方法可以满足操纵运动下螺旋桨水动力性能的计算精度要求．2.2　斜流下桨水动力性能螺旋桨转速保持为n=20 r/s，Jx取为0.3，0.5和0.7，Jy取为0.0，0.1，0.2，0.3和0.4，通过二者两两组合共15个工况．由于对船舶操纵运动影响较大的主要是螺旋桨的推力Tx、扭矩Qx和侧向力Ty，故不对其他方向的力和力矩进行具体分析．斜流下螺旋桨的水动力性能计算结果如图5和图6所示．从图5中可见：随着Jy增加，Tx和Qx增加，两者的变化基本一致．在低进速系数重载荷的情况下，Tx和Qx随着Jy增大呈现出非线性增长的趋势，而在中高进速、轻载荷的情况下基本为线性增长．对于本研究的螺旋桨，其在船舶操纵运动下的螺旋桨实际进速系数范围约为0.5~0.7．根据图5的计算结果可知：当Jx=0.5和10.13245/j.hust.230410.F005图5螺旋桨推力、扭矩系数和横向进速系数关系1—Jx=0.3；2—Jx=0.5；3—Jx=0.7(下同)．10.13245/j.hust.230410.F006图6螺旋桨侧向力系数和横向进速系数关系Jx=0.7时，Tx和Qx随Jy变化的曲线基本平行，即Tx和Qx的变化量仅和Jy的大小相关．从图6中可知：Ty可以认为和Jx无关，而和Jy近似呈线性关系，即和Uy近似呈线性关系．2.3　操纵运动下螺旋桨水动力性能为方便和试验数据[11]进行对比，对β=10°和r'=0.3这两个典型操纵运动工况下螺旋桨水动力进行计算验证．U=1.53 m/s，n=10.83 r/s，具体计算结果如表1和表2所示．其中直航CFD计算值表示设计航速下β=0°,r'=0工况下螺旋桨CFD计算结果，该结果可以作为参照反应两个典型操纵运动工况下螺旋桨水动力相对于直航工况下的变化趋势．10.13245/j.hust.230410.T001表1操纵运动下的螺旋桨推力系数KTxβ=10°r'=0.3左舷桨右舷桨左舷桨右舷桨直航CFD计算值0.097 20.097 20.097 20.097 2CFD计算值0.096 40.085 10.085 40.111 9试验值0.097 50.088 90.093 20.117 110.13245/j.hust.230410.T002表2操纵运动下的螺旋桨扭矩系数KQxβ=10°r'=0.3左舷桨右舷桨左舷桨右舷桨直航CFD计算值0.019 50.019 50.019 50.019 5CFD计算值0.020 10.018 20.018 80.021 7试验值0.017 30.016 10.019 10.021 8从表1和2可知：计算值在变化趋势和计算精度上和试验值相比符合良好，KTx和KQx的变化趋势基本一致，这种规律和敞水、斜流条件下的螺旋桨水动力性能一致．对于β=10°，船体左舷为迎风侧而右舷为背风侧，左舷桨的载荷相比于直航时变化较小，右舷桨的载荷比直航时减小10%；对于r'=0.3，船体首部左舷为背风侧而船体尾部左舷为迎风面，船体右舷相反，此时左舷桨的载荷相比于直航时减小15%，而右舷桨的载荷比直航时增大15%．虽然对于两种工况下船体尾部左舷均为迎风面，但左右两舷螺旋桨的载荷相比于直航时的变化完全相反，可能是船体周围流场的复杂性造成．3 螺旋桨修正体积力模型建模船舶在操纵运动过程通常处于回转和漂角耦合的状态，螺旋桨的水动力性能更为复杂．采用2.3节方法，对不同β和r'组合工况下螺旋桨水动力性能进行了模拟，包含船舶实际操纵运动下所遇到的所有操纵工况，另外同时采用了STAR CCM+中提供的原螺旋桨体积力模型和实桨法的计算结果进行了对比．U=2.59 m/s，n=18.85 r/s，β取值为0.0°，2.5°，5.0°，7.5°和10°，r'取值为0.0，0.1，0.2，0.3和0.4，两两组合共计20个工况．从图5和表1-2中可知：Qx和Tx的变化趋势基本一致，故仅给出Tx的计算结果．图7给出了操纵运动下不同计算方法求得的Tx，图中：Tx1为左舷桨的Tx，Tx2为右舷桨的Tx(后文相同符号代表同样含义)．从图7中可以看到：在操纵工况下采用原螺旋桨体积力模型所得到的螺旋桨水动力性能和实际情况有较大区别，对于左舷桨，当r'不变时，Tx1随着β的增大而增大，而原螺旋桨体积力模型求得的Tx1则是先减小而缓慢增大；对于右舷桨，当r'不变时，总体上Tx2随着β的增大而减小，而原螺旋桨体积力求得的Tx2则和r'的取值大小有关．以工况β=10°,r'=0.4为例，实际情况下Tx1＞Tx2，但原螺旋桨体积力法求得的结果正相反，左、右舷桨的推力相差分别达15%和42%．总的来说，迎风面桨(左舷桨)推力变化较小，在12%以内，而背风面桨(右舷桨)变化很大，变化幅度达50%．10.13245/j.hust.230410.F007图7操纵运动下不同计算方法的桨推力图8给出了操纵运动下采用真实桨计算求得的螺旋桨侧向力Ty计算结果，图中：Ty1为左舷桨的，Ty2为右舷桨的．由图8可知：Ty1随r'和β的增加呈线性增长的规律，而Ty2虽然呈现出非线性增长的规律，但总体变化趋势基本接近于线性增长．另外从侧向力的大小看，当β=10°,r'=0.4时，左右舷桨Ty的大小分别约为左右舷桨Tx的1/3和1/4，因此其大小不能忽略．10.13245/j.hust.230410.F008图8操纵运动下桨侧向力总而言之，原体积力模型仅以螺旋桨前方的某个截面处的轴向平均流入速度作为进速来求得螺旋桨的水动力性能，并不能较好地得到船舶操纵运动下的螺旋桨水动力性能，因此须要构建更精确的体积力模型．虽然舵角的变化也会对桨的水动力产生影响，但根据文献[12]中的计算结果，舵角对螺旋桨的水动力性能影响较小，当进行螺旋桨水动力建模时可以忽略．为方便对操纵运动下的螺旋桨水动力进行建模，本研究将Tx和Ty表示为和船舶运动状态U，r'及β有关的多项式形式，推力和扭矩在桨盘面内的分布形式采用Hough和Ordway[10]提出的H-O体积力分布方式，侧向力在桨盘面内则采用均匀分布．根据2.2节中的计算结论可知：Tx决定于该航速U下直航状态下的螺旋桨推力T0(U)及Jy，而Jy决定于船舶的运动状态U，β和r'，故可以将Tx近似表达成Tx=T0(U)+Uf1(β,r')的形式，其中f1(β,r')为关于β和r'的多项式．根据2.2节的计算结论可知：Ty和Jy近似成正比而与Jx无关，而Jy决定于船舶运动状态U，β和r'，故可以将Ty近似表达成Ty=Uf2(β,r')的形式，其中f2(β,r')为关于β和r'的多项式．4 自航船舶操纵运动数值模拟进行DTMB5415M回转运动和Z型运动数值模拟，其中基于原螺旋桨体积力模型的操纵运动数值模拟用BF表示，基于螺旋桨修正体积力模型的操纵运动数值模拟用新BF表示．因文献[5]中采用实桨法进行了以上两种运动的数值模拟，因此本研究将计算结果与之进行对比．4.1　35°舵角回转运动数值模拟从图9中可以看到，相比于BF，新BF的计算结果明显和文献[5]的结果符合较好：对于运动轨迹中的回转半径，新BF和文献[5]的结果相差为2.5%，BF和文献[5]的结果相差为10.1%；对比回转运动定常阶段的计算结果，对于U，虽然BF的结果和文献[5]的结果更加符合，但两种方法的结果和文献[5]的结果相差都在3.5%以内；而对于横摇角ϕ，文献[5]的结果为-2.8°，BF和新BF的结果分别为-1.2°和-2.6°，BF的结果仅为文献[5]的43%；对于r'，新BF，BF和文献[5]的结果相差分别为3.0%和9.5%．10.13245/j.hust.230410.F009图9回转运动参数图10给出了回转运动下Tx，Ty和舵侧向力Ry时历曲线，其中Ry1为左舵的，Ry2为右舵的．由图10(a)和(b)可知新BF的左、右舷桨推力结果和文献[5]的结果基本一致，相差都在5%以内，而BF的桨推力结果和文献[5]的结果相差较大，特别是回转运动定常阶段右舷桨推力和文献[5]的结果相差达到了47%．值得说明的是实际桨载荷是高频振荡的，因此在本研究中对文献[5]的结果进行了光滑处理．由于桨推力计算的准确性会影响桨后舵的舵效，因此可以看到新BF的左舵侧向力结果明显和文献[5]10.13245/j.hust.230410.F010图10回转运动螺旋桨载荷、舵侧向力时历曲线的符合更好，但右舵侧向力的精度并没有提高，这可能是因为右舵处于背流面其周围流场复杂，不仅受桨的影响，船体及附体对其舵效影响也不可忽视．观察图10(c)和(d)可知：左、右舷桨侧向力的合力在回转运动定常阶段大约是舵侧向力合力的1/3，说明桨侧向力的大小不可忽略，会产生对船桨舵系统影响较大的横摇力矩和转首力矩：在回转阶段，桨侧向力和舵侧向力相反，因此其产生的转首力矩会减弱舵的实际效果，使得船体转首角速度减小，同样由于桨侧向力产生的横摇力矩会使得船舶回转定常阶段的横摇角增大．以上分析和从图9中的得到的结果基本一致．图11给出了回转运动定常阶段螺旋桨后0.25D10.13245/j.hust.230410.F011图11回转运动定常阶段螺旋桨后0.25D处伴流场处流场分布图，图中：Z为垂向位置；Y为横向位置．可以看到左舷桨后0.25D的轴向速度，新BF的计算结果相比于BF的轴向速度较大，而右舷桨则相反．该截面处的流速的大小影响着舵的实际进流速度，因此对于左舷舵的进速，新BF的计算结果相比于BF较大，而对于右舷舵则相反，这和图10(c)和(d)中舵力大小的结果一致．4.2　20°/20°Z型运动数值模拟图12给出了Z型运动参数时历曲线，图中Ψ为转首角．由图可知：相比于BF的各项参数计算结果，新BF的各项参数和文献[5]相比在时间上稍有滞后，但该变化并不明显，而各项参数峰值的结果和文献[5]的结果更符合．10.13245/j.hust.230410.F012图12Z型操纵运动参数图13给出Z型运动下Tx，Ty和Ry时历曲线．由图可知：新BF的左右舷桨推力结果和文献[5]结果变化趋势基本一致，但在峰值上相差稍大达到了8.9%，可能是因为Z型运动是一个频繁操舵的非定常运动，桨周围的流场时刻变化，而新BF仅是根据定常操纵运动状态下建立的螺旋桨体积力模型，不能准确考虑这种动态变化过程，另外舵角的转动可能对于桨的水动力也有一定影响．BF的左右舷桨推力结果和文献[5]结果变化趋势完全相反，峰值上最大相差达到22%．对于舵侧向力，新BF在计算精度上有所提高，但和文献[5]的结果相比仍偏小，特别是5~10 s这样的不转舵阶段．10.13245/j.hust.230410.F013图13Z型运动螺旋桨载荷、舵侧向力时历曲线相比于回转运动，桨侧向力和舵侧向力的方向只是在部分时刻相反，但总体上桨侧向力仍旧减弱了舵效，只是效果不如回转运动中明显，这也和新BF各项参数的时历曲线相比BF稍微滞后的结果一致．总之，虽然新BF能够较准确地求得Z型运动下桨及舵的水动力，对于Z型运动特征参数的预报也更准确，但对于运动参数的预报在时间上会稍有滞后．一方面可能是新BF不能考虑Z型运动下船体周围流场的非定常性和复杂性，另一方面可能是力源分布的简化及体积力法不考虑螺旋桨实际几何模型造成．5 结论a．均匀来流下，Tx和Qx随着Jy的增大而增大，且其变化量仅与Jy有关；Ty则和Jx几乎无关，仅和Jy呈正线性相关；在操纵运动下，β和r'对桨载荷的影响情况不同，Tx和Qx变化比较复杂，Ty则与β和r'的变化呈线性正相关．b．通过操纵工况下螺旋桨的水动力性能计算可建立螺旋桨修正体积力模型，该模型能较好描述实际操纵运动下螺旋桨的水动力性能．运用该螺旋桨修正体积力模型可以较准确快速地对双桨船的操纵运动进行预报，相比原螺旋桨体积力模型考虑了伴流不均匀性对Tx和Ty的影响，且计算效率相比实桨法大幅提高，各种运动参数及桨、舵的受力情况均和实际情况符合较好．c．Ty的大小对于双桨船操纵运动的影响不可忽略．对于回转运动定常阶段，Ty会降低舵的实际效果，使得r'减小，而对于ϕ，其影响达到了100%以上；对于Z型运动，Ty同样会影响实际舵效，但对各项运动参数峰值影响较小，其影响体现在各项参数的时间滞后上．d．螺旋桨修正体积力模型是根据定常操纵运动状态下所建立的，采用力源的分布也比较简单，因此对于Z型运动这样频繁操舵的非定常操纵运动，可能不能显著提高预报的精度，须要进一步采用更精确的螺旋桨体积力模型进行预报．