考虑粗糙度影响的单圆角柱体水动力特性影响研究*

2023-09-25于定勇房子昂杨子慧

中国海洋大学学报（自然科学版） 2023年10期

于定勇, 房子昂, 杨子慧

(中国海洋大学工程学院, 山东青岛 266100)

圆角柱体结构物在海洋采油平台、桥梁建设等工程中有大量应用,而水工结构物在投入使用后受海水冲刷侵蚀、生锈和生物附着等因素影响,表面粗糙度将发生改变,受力和水动力环境也发生相应的变化;而基于减阻减振的目的,人为改变柱体表面粗糙度也是被动控制的手段之一。

Achenbach[1]最早通过风洞试验对粗糙柱体表面压强进行测量,得出在临界区,表面粗糙的圆柱所受阻力低于表面光滑的圆柱,并且粗糙度的增加会导致流动分离点向柱体后方移动。然后Achenbach等[2]扩大了雷诺数(Reynolds number,Re)NRe范围,使用滚花机在柱体表面附加纹理结构,通过研究粗糙度对总阻力的影响[3],发现表面粗糙圆柱绕流具有更低的临界雷诺数。Guven等[4]使用砂纸模拟了均匀流条件下粗糙圆柱边界层特性和受力特征,发现粗糙度与柱体所受阻力及其表面压强分布具有高相关性,且会影响到边界层的厚度与分离特性。

根据柱体表面不同粗糙度对边界层有明显影响且具有减阻作用的结论,以往的学者们又对局部加糙的柱体模型进行了进一步的研究。Leung等[5]在柱体表面布置了局部粗糙带,结果表明非对称沟槽圆柱体的尾迹和边界层同光滑表面处于不同的流态,漩涡尾迹向下移动。Chang等[6]研究了高雷诺数3.0×104≤NRe≤1.2×105条件下,局部粗糙带的布设位置、粗糙高度和粗糙覆盖面积对圆柱涡激振动增强和抑制的影响。Wang等[7]在圆柱表面一定角度附加了两条沟槽状的对称局部粗糙带,当4.0×104≤NRe≤1.2×105时,利用数值模拟研究沟槽高度对涡激振动的影响。Michálek等[8]采用粒子图像测速(PIV)方法研究了Re在5 140～11 800之间,粗糙圆柱周围的流动及其尾迹特性。研究发现当Re较低时,表面粗糙度的增加会减小尾迹长度,增加尾迹宽度;在Re较高时,表面粗糙度的增加虽然没有减小尾迹长度,但是仍有尾迹宽度增加的现象,只是这种现象不明显。综上可以发现,有关粗糙度的研究已经从全覆盖加糙发展到局部加糙对水动力的影响,从绕流研究发展到涡激振动研究领域,但是研究对象均为圆柱体。

van Hinsberg等[9]发现虽然关于表面粗糙度单独作用于圆柱的影响得到广泛研究,但直到2018年,关于表面粗糙度和圆角半径共同作用下的柱体绕流研究仍尚无资料,于是在NRe=12×106的条件下,为圆角半径为0.16和0.29的圆角柱体增加一个值为1×10-3的粗糙度,并设置两个来流方向。结果表明,圆角半径越大,柱体阻力系数越低,涡脱频率越高,脱落强度越低。

综上所述,关于圆角柱体表面增加粗糙度对绕流特性影响的研究较少,且以超临界雷诺数的空气动力学为主,对亚临界条件下水动力特性研究较少,而NRe=3 900下柱体绕流问题是典型的亚临界雷诺数问题,因此本文选取NRe=3 900进行计算。早期的研究手法多为风洞试验,对粗糙度的模拟多借助砂纸、肋条、电缆等实物,易导致粗糙度不均匀,影响试验准确性。随着计算流体力学的快速发展,使用CFD软件来进行数值模拟可提高计算准确度和效率。圆角柱体结构特殊,柱体表面由直边和圆角两部分组成,因此有必要研究局部加糙位置和粗糙度大小对圆角柱体绕流的影响。在实际工程应用中,水流存在多角度来流的情况,对圆角柱体而言,迎流面和分离点粗糙情况将发生改变。因此,研究不同来流角度对加糙圆角柱体绕流的影响也十分必要。

1 计算模型及参数设置

1.1 控制方程

本文采用Fluent中的SSTk-ω模型进行模拟。连续方程和动量守恒方程Navier-Stokes(N-S) 方程为:

∇·U=0,

(1)

(2)

式中:U表示速度;p表示压力;ρ表示流体密度;υ表示流体动力黏度。

湍动能k和比耗散率ω方程为:

(3)

(4)

式中:k和ω分别表示湍动能和比耗散率;Gk表示湍动能产生项;Gω表示ω的生成项;Гk和Гω分别表示k和ω的有效扩散系数;Yk和Yω分别表示k和ω的耗散;Dω表示交叉扩散项。

本文采用有限体积法,压力与速度的耦合应用SIMPLE算法,对流差值应用二阶迎风格式。连续性方程与动量方程收敛残差为小于10-4,时间步长为0.001 s,每个时间步内迭代次数为30。

1.2 模型参数

如图1所示本文计算域相关参数:特征长度D=0.01 m;顺流向(X方向)长度为20D;垂直水流方向(Y方向)长度为10D;展向高度(Z方向)为πD。加密区域为距离圆柱1.5D的网格。雷诺数NRe=U0D/υ=3 900,流体介质的流速设置为υ=1.01×10-6m2/s,入口处速度设置为U0=0.391 m/s。边界条件设置如下:计算域入口设定为速度入口(Velocity-inlet);出口设定为压强出口(Pressure-outlet);计算域的左、右边界设定为对称边界(Symmetry);为真实模拟桩柱上方连接防波堤底面,下方接触海底的实际情况参照文献[10],将上、下边界设置为无滑移壁面(No-slip wall),柱体表面设置为无滑移壁面(No-slip wall)。

图1 计算模型(a)和全局网格(b)示意图

如图2所示,圆角半径柱体以R+=R/D来表示。其中:R为圆角半径;D为柱体边长;来流角度为θ。如图3所示,在圆角和直边两位置分别增加粗糙度,圆角加糙为在圆角柱体四个圆角位置增加粗糙度,直边加糙为在圆角柱体四个直边位置增加粗糙度。

图2 圆角半径示意图

图3 局部加糙位置布设图

2 模型验证

2.1 光滑单圆柱模型验证

为验证网格无关性分析及数值方法的准确性,给出了3套网格(A1、A2和A3)的设置情况(见表1)。选取在NRe=3 900条件下光滑单圆柱体绕流的工况,本文的数值模拟结果与前人研究结果的对比如表2所示。

表1 网格设置情况

表2 光滑单圆柱数值计算结果验证

根据乔永亮等[11]的研究,SSTk-ω模型平均阻力系数要大于 LES 的模拟结果,且SSTk-ω模拟结果的回流区长度只有试验值的一半左右。由表2所示,本文采用的SSTk-ω方法计算得到的平均阻力系数(Average resistance coefficient)Cd-ave和斯特鲁哈数(Strouhal number,St)均大于前人采用大涡模拟(LES)方法得到的计算结果。在同为SSTk-ω的计算方法下,计算得到的三套网格(A1、A2和A3)的Cd-ave和St同已有研究结果吻合度较高。在NRe=3 900条件下,由于缺少采用SSTk-ω模拟方法得到的圆柱中心顺流向时的均流速分布资料,因此将本文结果同前人采用PIV和LES方法的研究结果进行对比(见图4),SSTk-ω模拟结果的回流区长度约为PIV试验值的一半,符合乔永亮等[11]的研究结论。综上,说明本文选取的网格和数值模型准确性较高,考虑到计算效率,本文所有网格设置与A2网格相同。

图4 柱体中心顺流向时均速度分布图

2.2 加糙单圆柱模型验证

为验证本文选取的数值模型在模拟加糙柱体绕流时的可行性,本文模拟了雷诺数NRe=20 000,且粗糙度KS/D=0.003 25的条件下,单椭圆柱体的绕流现象,并计算了柱体的平均阻力系数和升力均方根系数。如表3所示,本文将计算结果同钱权等[20]的粗糙圆柱数模进行对比,各水动力参数相差均小于4%,吻合度较高,说明本文采用的数值方法具有较高的可靠性和准确性。

表3 加糙单圆柱数值计算结果验证

3 模拟结果及分析

3.1 受力特性

图5—9给出了不同加糙位置圆角柱体升、阻力系数和斯特鲁哈尔数随粗糙度大小和来流角度变化的规律。

((a) 平均阻力系数 Drag coefficient; (b) 升力均方根系数 Lift coefficient; (c) 斯特鲁哈数 St.)

在圆角柱体圆角位置增加粗糙度时,图5(a)给出了不同来流角度θ下Cd-ave随KS/D变化的曲线。研究发现所有来流角度的Cd-ave均随粗糙度的增加而逐渐降低,其中正向来流条件下减阻率最高,KS/D从0增大至0.02时的减阻率约为36.7%,其余来流角度下减阻率则均在20%左右。而图6(a)对比了不同KS/D下Cd-ave随θ变化的曲线。当θ从0°到30°变化时,各种大小KS/D下圆角柱体的Cd-ave均随θ的增大而增大,当θ>30°时,Cd-ave由增大的趋势转为减小趋势。其中,θ从0°到15°变化时Cd-ave曲线上升幅度较大。该现象出现的原因可能是在θ增大的过程中,垂直方向的柱体宽度逐渐随之增大,迎流面积随之增加,导致平均阻力系数上升。当θ=45°时,虽然迎流面积达到最大,但此时边界层分离位置均处于圆角加糙处,两加糙圆角对边界层分离的影响较θ=30°时增大,减阻作用加强,因此θ=45°的Cd-ave数值反而要小于θ=30°。

((a) 平均阻力系数 Drag coefficient; (b) 升力均方根系数 Lift coefficient; (c) 斯特鲁哈数 St.)

图5(b) 给出了不同θ下Cl-rms随KS/D变化的曲线。可以发现柱体的升力均方根系数(Root mean square coefficient of lift)Cl-rms在θ从7.5°到45°之间的条件下保持较平稳的过渡,整体趋势为随着KS/D的增加略有下降。当θ=0°且KS/D=0.01时,Cl-rms发生了突增现象。为研究该突增现象,图7给出了柱体在正向来流下,当KS/D为0.005或0.01或0.02时的升、阻力时程曲线,可以发现KS/D=0.005时,升、阻力系数变化较为稳定;而KS/D=0.01时,升、阻力系数发生了很明显的“拍”现象,升阻力系数波动强弱变化幅度较大,且存在周期性;KS/D=0.02时,升、阻力系数变化再次稳定,且波动幅度更小。而图6(b)对比了不同KS/D下Cd-ave随θ变化的曲线。随着θ的增加,整体而言各粗糙条件下圆角柱体的Cl-rms均随θ的增加呈先增大后减小的趋势。其中θ从0°到15°时呈增大趋势,θ从15°到45°时呈减小趋势,在该变化过程中,Cl-rms增大的幅度远大于减小的幅度。

图7 圆角粗糙柱体升阻力系数时程曲线随粗糙度大小的变化

图5(c) 和图6(c)分别给出了圆角加糙柱体的St随KS/D和θ的变化曲线。θ在7.5°～45°之间时,柱体的St曲线均随着KS/D的增加而小幅上升且均随θ的增加而下降。当θ=0°时,St的变化规律不明显,表现为:St曲线在KS/D=0.01有一处突增,但随着θ的增大,St曲线上升、下降的情况同时存在。

在圆角柱体直边位置增加粗糙度时,由Cd-ave、Cl-rms和St三者随KS/D(见图8)和θ(见图9)的变化曲线可知:随着KS/D的增加,所有θ下的Cd-ave、Cl-rms和St曲线波动幅度均较微小,且曲线的微弱变化也无明显规律(见图8);而Cd-ave、Cl-rms和St数曲线均随θ的增加,有着较明显的波动,均呈先上升后下降的变化趋势(见图9)。其中,Cd-ave曲线在θ=30°处发生转折,Cl-rms曲线在θ=15°处发生转折,这同圆角加糙柱体力系数的变化规律相似。在不同KS/D下,Cd-ave、Cl-rms和St曲线几乎重合。同圆角加糙柱体相比,直边加糙柱体下St曲线变化幅度较小,但是变化规律较为明显。说明直边加糙对圆角柱体受力特性的影响较小,θ的影响占主导地位。

((a) 平均阻力系数 Drag coefficient; (b) 升力均方根系数 Lift coefficient; (c) 斯特鲁哈数 St.)

(a) 平均阻力系数 Drag coefficient; (b) 升力均方根系数 Lift coefficient;(c) 斯特鲁哈数 St

图10给出了不同KS/D和θ下圆角加糙柱体的平均压强云图。

图10 不同来流角度、不同大小粗糙度下圆角加糙柱体平均压强云图Fig.10 Mean pressure cloud of fillet roughed cylinder under different flow angles and different roughness

由图10可见,θ相同的情况下,圆角加糙柱体后方时负压值随粗糙度的增加而逐渐降低,且θ在0°～15°之间变化时,这种变化较为明显。在KS/D大小相同的情况下,随着θ的增加,两种局部加糙位置的圆角柱体均出现负压最大值点后延的现象,且θ越大,柱体后方平均负压区影响范围也越长,相应的负压值也越大。另外如图10(a) 所示,在θ=0°,KS/D=0.01时,柱体后方和两侧有负压值突变现象,同图5对比可见,平均压强云图的变化规律与升阻力系数的变化规律是一致的。

3.2 瞬时涡量

图11给出了圆角加糙柱体的瞬时涡量图。在相同来流角度下,随着KS/D增大,剪切层总体呈变长的趋势,柱体后方涡街逐渐变得混杂,且当θ<15°时,这种变化较为明显。其中,观察图11(a3)可以发现,KS/D=0.01时,柱体后方脱落涡旋之间的间隔变小, 即涡旋脱落频率明显加快。从前倒角处分离的剪切层在柱体后方相互作用,卷起并形成漩涡。相比其他KS/D下的柱体,漩涡脱落距离柱体更近,柱体两侧扰动更加剧烈,导致升力系数显著增大。

图11 不同来流角度、不同大小粗糙度下圆角加糙柱体瞬时涡量云图

在相同KS/D下,随着θ的增大,迎流面积的增大导致圆角柱体对流场的影响范围扩大,柱体后方尾流宽度逐渐变宽。另外,可以观察到,当非正向来流时,柱体两侧涡脱落并非同时进行,柱体下侧的边界层分离要早于上侧,这可能是由于上表面因存在倾斜角度而具有导流作用,延缓了涡脱落的时间。柱体两侧涡脱落时间不同步导致了涡街混杂程度的加剧。

图12给出了同一时刻(工况第8 s),θ=0°时四种KS/D柱体的瞬时流线图。可以发现,当KS/D=0、0.005和0.01时,柱体后方均形成了2个漩涡。而当KS/D=0.01时,2个漩涡大小相差较大,且两漩涡中心距离柱体较近,对柱体产生剧烈扰动,导致升力系数增大。继续增加粗糙度至KS/D=0.02,此时涡脱落位置距离柱体较远,柱体受到的扰动减弱,导致升力系数减小。图12中流线的变化可以进一步解释图5(b)中θ=0°且KS/D=0.01时,Cl-rms曲线出现突增现象的原因。

图12 正向来流下,不同大小粗糙度下圆角加糙柱体流线图

3.3 时均流场

为清晰地观察到柱体后方流场时均速度的分布情况,于计算域展向高度Z=πD/2处设置一截面,在圆角加糙柱体后方x/D=1.06位置布设垂直于流场方向的监测线。图13给出了该位置处x方向时均速度分布图。各KS/D和θ条件下,顺流向时均速度整体均以柱体中心为轴线呈对称分布。但在θ为7.5°或15°或30°的条件下,时均速度最小值向y轴负方向存在微弱的偏移,来流角度越小,偏移幅度越大。这种现象产生的原因可能是θ改变后,柱体两侧边界层分离不同步,从而导致只有θ为0°或45°时后方顺流向时均速度曲线是对称分布。另外,当KS/D逐渐增大时,时均速度曲线因θ变化而产生的偏移现象有所减弱。