基于致动盘方法分析潮汐能水轮机水动力性能

刘葳兴，刘 磊，陈雨龙，王 鑫，张之阳

（江苏海洋大学海洋工程学院，江苏 连云港 222005）

潮汐能对比其他可再生能源具有再生、可预测、环保等优势[1]，从潮汐能中提取能量是当下研究热点，而潮汐能水轮机是一种从潮汐能中提取能量转换为其他能量的装置。模型试验、数值仿真可用于水轮机的尾流场研究[2-3]。水轮机的数值仿真方法常分为全尺寸仿真、缩放尺寸仿真和致动盘方法。然而，从节约计算时间考虑，可用使致动盘方法捕捉水轮机工作时尾流。用致动盘方法替代真实叶片进行计算时，常用体积力作为附加源项，并规定体积力分布形式及其大小[4]。致动盘方法在捕捉叶片几何形状变化对流场造成的影响时有所欠缺。捕捉翼型对流场带来的影响可使用叶片元素动量理论（blade element momentum theory，BEM）方法[5-6]。张之阳等[7]采用BEM 分析水平轴潮流能涡轮机的水动力特性，其方法节省大量计算时间，但在尾流细节处研究不足。而应用致动盘模型（actua‐tor disc model），可在计算时间上、节省计算资源和捕捉尾流细节特征上相平衡。张玉全等[8]基于曼彻斯特大学循环水槽的实验数据，利用致动盘理论求解Navier-Stokes 方程获得尾流场数据，并将两者对比，发现致动盘方法更易获得尾流场情况。Myers等[9]用多孔圆盘代替潮流能水轮机叶片，对其性能开展研究，发现水轮机的推力、机架离拖曳池底面距离与水池壁面粗糙度对尾流都有不同影响。Sun等[10]使用Fluent软件建立多孔圆盘的二维和三维模型，发现水轮机运行中流场后方的液面高度有所下降与尾流的速度亏损有关。Edmunds等[11]使用非线性致动盘模型并对传统源项进行改进，分析叶片尖端损失及尾流场。致动盘替代实际叶片进行仿真，常用于风机尾流场分析，因风机与水平轴水轮机在结构上具有一定相似性，故可将致动盘模型应用于水轮机尾流场分析。本研究采用体积力附加源项方式对潮汐能水轮进行数值仿真，并与爱丁堡大学中的拖曳池试验作对比，评估致动盘模型在水轮机上的应用。

1 实验装置与方法

1.1 试验拖曳池

本实验数据参考爱丁堡大学Flowave 拖曳池中进行的水轮机尾流场实验[12]。实验拖曳池为直径25 m 的圆形结构，工作水深H=2 m。拖曳池中，沿整个水池壁面有28个叶轮单元，每个单元可控制水流流速而形成一个循环流动系统。通过设置叶轮转向在水池中央试验区的任一方向产生水流流动，可实现最大水流流速为1.6 m/s。水池上方有XY 行车方便安装水轮机。由于实验中没有设施可控制湍流强度，因此，本实验湍流强度与爱丁堡大学进行水轮机实验所选择的湍流强度相同（该湍流强度为欧洲海洋平均流速的7%[13]）。

1.2 涡轮机和仪器设备

固定在拖曳池底部的水平轴水轮机，叶片直径D=1.2 m，轮毂中心距水池底h=1 m。其翼型为NACA63812、叶片的几何形状可参考文献[13]。水轮机上有安装测量叶片弯矩R，转子扭矩Q和推力T的传感器。水轮机中的伺服电机可使用编码器来控制转子的转动角度θ，进而使叶片在不同尖端速度比（tip-speed ratio，TSR）范围内运行。

水轮机流场布置见图1（a），拖曳池水流方向沿x轴，水池纵向为y轴，深度方向为z轴，转子中心位置记（x,y,z）=（0,0,0）。用Vectrino Profiler声学多普勒测速仪（ADV）对流量采样，以转子中心为水平面（z=1 m 截面），采样点在x轴、x=0 m、x=2.4 m 上分布，见图1（b）。采样频率为100 Hz，采样时间256 s。测量的流量速度分别为水平速度、纵向速度和垂直速度。实验所用水流速度为0.8 m/s。

图1 水轮机布置示意Fig.1 Schematic showing the turbine arrangement

水轮机实验时还需验证其运行过程中是否产生阻塞效应。水轮机阻塞率为叶片旋转一周面积S1比拖曳池面积S2的值，S1/S2=0.002 3。垂直阻塞率为叶片直径D比水池水深h的值，D/h=0.6。水平阻塞率为叶片直径D比水池直径w的值，D/w=0.048。可见，水轮机阻塞率、垂直阻塞率、水平阻塞率均很小，表明在潮汐能水轮机测试期间不会出现阻塞效应进而影响流场。

2 致动盘模型

致动盘方法基于一维动量理论，将水轮机叶轮视为一个可穿透的等效面积圆盘区域代替[14]。同时，将水轮机致动盘周围流场视为稳态、一维、不可压流场，并以体积力的形式实现致动盘对流场的作用[15]。

由于致动盘模型常用于风机尾流场研究，而在水轮机上应用较少，因此，本研究基于文献[16]，将OpenFOAM 的自带“simpleFoam”求解器进行修改自定义成新的致动盘求解器，并将fan 边界条件在OpenFOAM 中应用。假设在不可逆流体下，方程组可写成以下形式：

其中，xi代表笛卡尔坐标分量（x,y,z），Ui是笛卡尔平均速度分量，fi包括圆盘转子特征的附加源，即叶片对流体作用，本研究所用体积力作为附加源项。雷诺应力为，需应用湍流模型来建模得到控制方程。

将实验得到的推力T以体积力源项形式均匀分布于致动盘模型上。体积力与来流速度和推力系数有关。假设来流速度为U，对应功率曲线的空气密度和推力系数分别为ρ,CT，则水轮机叶轮单位面积的轴向力为

Δx为致动盘厚度，V为致动盘体积。其动量方程中，致动盘体积力源项表示为[17]

本研究采用RNGk-ε湍流模型，k代表湍流动能，ε代表这个能量的耗散率。

3 数值模拟

3.1 网格规划

仿真建模时对拖曳池和水轮机按1∶1 建模，以保证数值仿真的准确性。用OpenFOAM 6.0版本软件中“blockMesh”和“snappyHexMesh”功能划分网格。“blockMesh”工具用来生成一个初始块（网格域），并将其细分为不同的尺寸区域。块初始尺寸的底部面积为30 m2，高4 m 的圆柱体。接着“simpleGrading”参数设置为1。生成块状网格后，使用“snappyHexMesh”工具迭代细化“blockMesh”，将产生的拆分的“hex”网格变形到几何体上，从而近似地符合几何体。

将尾流区域规划为一个半径为0.7 m 的圆柱体，尾流长度以转子中心起到尾流下游9 m 处结束。对流体计算域不同区域的网格进行加密（图2），水轮机周围边界层网格密度，细化水平最高；机架上方网格密度次之；为节省计算资源，流体计算外域网格密度最低。计算域网格规划完成后进行网格敏感性分析，网格数量分4 个等级（表1）。对上述4 种数目的网格进行计算，计算结果以y=0 m，z=1 m截面速度分析，计算结果见图3。由图3 可见，当选用网格数量等级3后，流场中速度差异变化较小，故为节省计算时间后续仿真计算选择网格数量等级3进行求解。

图2 涡轮机网格细化示意Fig.2 Turbine mesh refinement diagram

表1 网格数量Table 1 Total number of meshes

图3 网格敏感度Fig.3 Mesh sensitivity

3.2 初始条件和边界条件

水池壁面各部分按不同位置设置为质量流入、质量流出或无流流动。将图4 拖曳池中的28 个叶轮单元配置为入口、出口及壁面，质量流量大小根据表2 设定。表2 中质量流量（X轴流向）正值为流入，负值为流出。流入之和与流出之和为相同数值，则保持质量守恒。水池壁面速度设置为0 m/s，壁面函数用于k，ε。流体域顶部设置为无滑移壁面。除速度外，其余初始条件皆设置为边界条件，初始速度设置为0 m/s，动力粘度ν设定为1.666 7 ×10-6m2/s。控制方程求解选择以压力为基础的求解器，采用Simple 算法，空间离散采用二阶迎风格式，收敛标准为1 × 10-4。

图4 拖曳池质量流量配置示意Fig.4 Mass flow configuration diagram for towing pool

表2 质量流量设置Table 2 Mass flow rate setting

4 仿真结果与分析

本研究仿真结果分为两部分，一是拖曳池数值仿真，即在无水轮机情况下水流在拖曳池中自由发展；二是在致动盘模型下进行仿真。将两者仿真结果对比，分析水轮机尾流亏损情况，并为后续研究做基础。

4.1 自由发展区仿真

以水轮机转子中心z=1 为截面拖曳池数值仿真见图5。当水流速度以0.8 m/s 从壁面入口加速流向拖曳池中心时，水流会与壁面发生相互作用，引起湍流强度的加剧，从而导致速度亏损。水流在拖曳池中心充分发展时，流动速度更加均匀。水流经下游壁面流出时，部分水流相互作用时会产生出口回流，导致速度在壁面附近发散，又因水池上下壁面不是出入口，回流在此加剧湍流强度。

图5 拖曳池仿真Fig.5 Simulation diagram of fluid flow rate in the towing pool

以图1（b）采样点对拖曳池不同位置的流速和湍流数据进行采样评估（图6）。由图6（a）可见，上游水流速度沿x轴逐渐增大，到x=0 附近水流增速变缓，x=0 之后速度开始下降。这是由于上下游水流与壁面作用引起速度亏损，以及拖曳池上下湍流强度大引起水流流速在中心处增大和水流经充分发展在拖曳池中心处流速均匀。图6（a）中，实验与仿真流速误差在靠近x=0位置时减少，并随着尾流长度的增加而大，湍流强度两者误差较小。实验和仿真计算的水流速度和湍流强度之间的平均误差分别为0.3%和1.2%。

图6（b，c）显示以x轴向上的水流发展和扩散，在仿真模拟结果中可见速度和湍流强度曲线的对称性，然而在实验中存在不对称性。这可能是在实验过程中很难精确控制拖曳池中各个单元所需的水流流量，从而引起整个区域的平均流量和湍流强度变化。在x=0 m 和x=2.4 m处，实验和仿真数值的水流速度的差值分别为6.7%和8.2%。两截面还显示出外侧水流速度均大于中心水流速度，其原因是拖曳池外侧叶轮控制的质量流入大于中心处的质量流入。

图6 不同截面速度与湍流强度Fig.6 Different cross-sectional velocities and turbulence intensities

4.2 致动盘模型仿真

控制水轮机叶片转速，进而改变尖端速度比（TSR），实验中TSR 设置为4～7。由于推力T以体积力源项形式均匀分布于致动盘模型上，因此需要得到推力系数CT。由前文实验，得到推力T，带入公式（5），便可得到推力系数CT。不同的TSR 可以评估功率系数CP，进而得到水轮机最大功率。TSR 是由叶片旋转速度比来流速度计算得来，见公式（6）。CP由水轮机实际功率比上理论功率得来，见公式（7）。由图7 可见，TSR=6 时，功率系数CP最大，且只比TSR=5.5 时的CP略大，又因TSR=5.5 时实验和仿真的推力系数CT高度吻合，故综合考虑后续仿真实验采取TSR=5.5进行。

图7 不同TSR下CP与CT值的比较Fig.7 CP and CT results under a range of TSRs

其中，T为推力，ρ流体密度,U为来流速度。

其中，ω为旋转角速度。

其中，P为水轮机实际功率，A为致动盘面积。

图8（a）为z=1 m 截面水轮机速度云图，左图为速度细节，右图为流场区整体速度。图8（b）为湍流强度云图，左图为水轮机尾流区湍流强度细节，右图为流场区整体湍流强度。图10为y=0 m，z=1 m截面水轮机速度和速度亏损曲线。由图9 可见，速度和湍流数据的实验值与数值模拟高度吻合，两者误差分别为3%到5%。结合图8（a）和图9 可看出，水流速度由上游到致动盘前1.5 m（x=-1.5 m）处增加，接着由于转动盘及机架对水流有阻碍作用引起湍流而使水流速度下降。在致动盘上存在不连续的压强下降，致使尾流向两侧扩张。水流在经过旋转的转盘后，由于转盘吸收来流的部分能量并且来流与机架相互作用，使尾流动能小于来流动能，该区域内产生较大的湍流，同时降低流速并且尾流恢复严重不足。图8（a）显示，下游尾流的速度云图中涡流呈不对称状态，其原因可能是水轮机建模细节处网格划分不一致。图9 流量亏损显示，上游流量到制动盘前，水轮机与自由发展区要之间约有-0.018～0.13 m/s的流量亏损。

图8 水轮机仿真Fig.8 Tidal turbine simulation contour

图9 y=0 m,z=1 m截面速度Fig.9 y=0 m,z=1 m section velocity

观察尾流流场不同剖面图的速度变化是衡量尾流水动力特性的重要依据。图10 为x=-0.5、0、0.5、2.4、4.0 m，z=1 m 处截面图。其中x=2.4 m 处的实验数据和仿真结果高度一致，表明致动盘方法在捕捉尾流特性效果较好。从图10可知，水轮机从水流中提取动能，通过水轮机的部分流体会有减速效应，并产生一个对称“W 形尾流”。当水流在尾流中进一步扩散时，“W 形尾流”剖面消失。随着尾流的扩大，尾流剖面宽度将大于水轮机直径，并且由于湍流和机架作用尾流中心处的速度逐渐下降，而两侧尾流速度逐渐恢复。这说明尾流场两侧与周围环境流场的动量交换促使两侧速度恢复，而中心区域随着与湍流区越接近使其耗散的能量也就越大。

图10 TSR=5.5下不同截面速度Fig.10 TSR=5.5,different section velocities

5 结论

本研究应用实验和数值仿真的方法对比分析潮汐能水轮机尾流场特性。为保证仿真结果准确性，建模时将水轮机机架结构一并考虑进去。在实验设计过程中，通过拖曳池中有无水轮机对其流场影响分析。得到以下结论：

1）拖曳池数值仿真结果得到的速度与湍流强度与实验数据相吻合，但在纵向截面（z=1 m，不同y截面）上仿真模拟和实验之间有轻微的差异。这是因为在实验过程中，叶轮机械存在摩擦、阻力等因素，无法精确控制水流的输出与仿真过程一致，进而导致整个水池区域的平均流量和湍流强度出现一些变化。

2）对于致动盘水轮机模型，在TSR=5.5 时，y=0 m，z=1 m 截面上速度与湍流强度实验值与仿真值两者误差分别为3%到5%。随着尾流的扩大，x=0 m 截面时中心速度最大，x=4.0 m 截面时，其中心速度最小。且由于湍流和机架作用尾流中心处的速度逐渐下降后，两侧尾流速度逐渐恢复。

应用推力T作为致动盘模型附加源项，使仿真模拟具有计算时间少所需计算资源少优点，并且也能精准捕捉尾流场变化。圆形拖曳池相较于方形拖曳池在造波性能上有范围广、方向多、符合真实环境等优点，但有占用面积大、建设费用较高、不利于大型实验等缺点。由于圆形拖曳池中心区较小，不利于开展大规模阵列水轮机实验，因此后续将利用方形拖曳池对水平轴水轮阵列间距机进行实验研究。