CN109376389B - 一种基于2D_k Jensen模型的三维尾流数值模拟方法 - Google Patents

Abstract

本发明提供了一种基于2D_k Jensen模型的三维尾流数值模拟方法，可获得新型的三维尾流模型，其特征在于，包括：步骤1、计算入流风切变曲线u₀(z)以及湍流强度在垂直方向上的分布型I₀(z)；步骤2，计算尾流膨胀系数k_x,z；步骤3，基于原Jensen模型对尾流风速进行预测，得到初始尾流风速u^*(x,z)；步骤4，计算尾流膨胀半径r_x,r；步骤5，利用步骤4所得尾流膨胀半径r_x,r，基于2D_kJensen模型提出的cosine型速度分布，对步骤3所得的初始尾流风速u^*(x,z)进行三维修正，得到修正后的尾流风速分布u(x,y,z)。本发明获得的三维尾流模型继承了工程模型的优点，对流向、横风向和垂直向的尾流速度亏损均体现了良好的预测精度，且能够有效地反映垂直方向速度亏损的非对称特性，计算精度甚至优于基于CFD方法的数值模拟结果。

Description

一种基于2D_k Jensen模型的三维尾流数值模拟方法

技术领域

本发明属于风力发电技术领域。

背景技术

风力机尾流效应会引起下游风速降低、湍流强度增加、风切变加剧等现象。风速减小使得下游机组的输出功率降低，强湍流和附加的风剪切也会影响下游机组的疲劳载荷、结构性能和使用寿命。对于规模有限的风电场而言，其内的风力机不可避免地处于周围机组的尾流之中运行，尾流效应成为风电场机组排布和优化时应考虑的重要因素。由于机组布置优化过程中需要不断地评估排列方式、安装间距变化后的影响，不可能通过资源消耗大、耗时冗长的计算流体力学(CFD)数值模拟来实现，此时必须借助工程尾流模型。

尾流模型是基于风力机空气动力学原理发展得到的数学解析模型，主要用于计算尾流区的速度分布，其简便性适合于解决风电场发电量评估和微观选址等风工程问题。从1980年起，学者们先后提出了一系列的工程尾流模型，其中， Jensen尾流模型(也称为PARK模型)因具有形式简单、计算效率高、应用方便且精度可接受等特点，被认为是开创性尾流模型，在风工程领域得到了普遍认可及广泛应用。但是，Jensen模型属于一维模型，其理论相对简单，存在一些不足之处。之后，田琳琳、杨祥生等人对Jensen模型进一步修正改进，陆续提出了几种二维(2D)Jensen模型，用于计算尾流区顺风向和横风向的速度分布情况。

风电技术近年来迅速发展，主力机型从兆瓦发展到多兆瓦级，风轮直径达到了百米量级。随着风轮扫风面积和轮毂高度的增加，风切变和湍流强度的影响更加明显。对于风轮直径超过100m的多兆瓦级风电机组，其风轮范围内的风速差别可达30％。此时，采用塔架机舱处的单点风速来计算风力机性能的常规方法会造成设计功率与实际输出功率存在较大偏差，需要整体考虑叶轮平面内风速分布。为了解决该问题，2017年3月，风电行业标准IEC 61400-12-1： 2017引入了风轮等效风速的概念，将风轮扫掠面多个高度的风速大小进行加权平均，得到能够更准确反映风轮扫掠面动能通量的等效风速，提高了风力机功率预测精度。相应地，亟需一种除尾流区流向和横风向之外，还能够同时反映垂直向速度分布的三维尾流模型，用于计算尾流区的三维速度分布，进而计算尾流区机组的入流等效风速，最终提高整个风电场的发电功率预测精度。但据调研，当前国内外普遍使用一维尾流模型，少数使用二维模型，几乎没有相应的三维尾流模型供参考和使用。

发明内容

鉴于以上不足，本发明在田琳琳等人提出的2D_k Jensen模型的基础上，综合考虑风速和湍流强度在垂直高度的分布特性及其对尾流的影响，提出一种新型的三维尾流模型。

本发明提供的技术方案为：

一种基于2D_k Jensen模型的三维尾流数值模拟方法，其特征在于，包括：

步骤1，计算入流风切变曲线u₀(z)以及湍流强度在垂直方向上的分布型 I₀(z)，计算式如下：

式中，z为距离地面的高度，z_hub为轮毂高度，z₀为地表粗糙度，u(z_hub)为风力机轮毂高度处的风速，I(z_hub)为风力机轮毂高度处的湍流强度；

步骤2，计算尾流膨胀系数k_x,z，计算式如下：

k_x,z＝k₀·I_wake,x,z/I₀(z)

式中，k₀为初始尾流膨胀系数，I_wake,x,z为风力机下游x位置，z高度处的湍流强度；

所述I_wake,x,z的计算式如下：

其中,

式中，x为风力机下游距离,D为风轮直径,C_T为风力机推力系数,I_add,x为尾流区x位置处的附加湍流强度，包含由于风力机的存在而产生的剪切湍流和附加的机械湍流；

步骤3，基于Jensen模型对尾流风速进行预测，计算初始尾流风速u^*(x,z)，计算式如下：

式中，ɑ为风力机轴流诱导因子，与风力机的推力系数C_T有关，计算式为

s为量纲化下游位置，计算式为s＝x/D；

步骤4，计算尾流膨胀半径r_x,r，计算式如下：

r_x,r＝k₀·x+r₁

式中，r₁为风轮后的特征尾流半径，表达式为

其中 r_d为风轮半径；

步骤5，建立三维尾流模型，利用步骤4所得尾流膨胀半径r_x,r，基于2D_k Jensen模型提出的cosine型速度分布，对步骤3所得的初始尾流风速 u^*(x,z)进行三维修正，得到修正后的尾流风速分布u(x,y,z)，其表达式为：

u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

f_x,r为风速分布函数，其表达式为f_x,r＝cos(2π/r_x,r·r+π)，r为以风轮中心线为圆心的辐射半径，表达式为：

式中，y为风力机尾流任一点的横向坐标，y_wake和z_wake分别为风力机下游x位置的尾流中心坐标。

在上述方案的基础上，进一步改进或优选的方案还包括：

所述步骤1中，风力机轮毂高度处的风速u(z_hub)和湍流强度I(z_hub)通过风速计测得。

步骤2中，初始尾流膨胀系数k₀通过考虑地表粗糙度计算获得，计算式为：

k₀＝0.5/ln(z_hub/z₀)。

有益效果：

经分别与风场实测数据和风洞试验结果比较发现，本发明三维尾流数值模拟方法所建立的三维尾流模型继承了工程模型的形式简单、容易编码、计算高效等优点，对流向、横风向和垂直向的尾流速度亏损均体现了良好的预测精度，且能够有效地反映垂直方向速度亏损的非对称特性，其计算精度甚至优于基于 CFD方法的数值模拟结果。

附图说明

图1是一维Jensen模型、二维2D_k Jensen模型和三维尾流模型的尾流速度型示意图；

图2是Nibe B风力机下游4.0D、7.5D处轮毂高度平面的径向风速廓线；

图3是Nibe B下游4.0D、7.5D处垂直方向风速廓线；

图4是陆上环境和海上环境的来流风速和湍流强度的垂直分布图；

图5a是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处轮毂平面的径向风速廓线(陆上环境)；

图5b是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处轮毂平面的径向风速廓线(海上环境)；

图6a是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处垂直向风速廓线(陆上环境)；

图6b是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处垂直向风速廓线(海上环境)；

图7a是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处垂直向风速亏损廓线(陆上环境)；

图7b是Vestas 2MW风力机下游5.0D、7.0D、10.0D和15.0D位置处垂直向风速亏损廓线(海上环境)。

具体实施方式

为了进一步阐明本发明的技术方案，下面结合附图与具体实施例对发明的技术方案进行详细说明。

一种基于2D_k Jensen模型的三维尾流数值模拟方法，建立三维尾流模型的过程如下：

1.给出入流条件

考虑风切变效应以及湍流强度在垂直方向上的变化，根据风力机的轮毂高度z_hub，当地地表粗糙度z₀，还有测风塔或机舱风速计测算得到的风力机轮毂高度处的风速u(z_hub)和湍流强度I(z_hub)，结合对数律规律，给出风力机的入流风速廓线u₀(z)和湍流强度分布型I₀(z)，

其中，z为距离地面的高度，下同。

2.计算尾流膨胀系数

首先根据风力机的轮毂高度z_hub以及当地地表粗糙度z₀，可以计算得到初始尾流膨胀系数k₀，

k₀＝0.5/ln(z_hub/z₀)

湍流强度对尾流恢复快慢有着重要影响，具体而言，较强的湍流有助于尾流与周围自由流进行对流交换，通过流动掺混使得尾流的速度亏损得到快速恢复；相反，当湍流强度较小时，尾流与自由流的动量交换作用偏小，从而尾流恢复较慢。鉴于上述湍流强度对尾流特性的影响机制，不得不考虑尾流区湍流强度对尾流膨胀系数的修正。事实表明，风力机尾流区不仅存在入流湍流的影响，还包含由于风力机的存在而产生的剪切湍流和附加的机械湍流等。为此，引入附加湍流强度I_add,x的计算公式：

式中，x为风力机下游距离，D为风轮直径，C_T为风力机推力系数，下同。

用附加湍流强度I_add,x，结合入流湍流强度I₀(z)，计算得到尾流区的湍流强度I_wake，x，z，

并以此修正之前求出的初始尾流膨胀系数k₀，得到新型三维尾流模型的尾流膨胀系数k_x,z，表达式如下：

k_x,z＝k₀·I_wake,x,z/I₀(z)

3.预测初始尾流风速

原Jensen尾流模型假设尾流是线性扩张的，尾流风速随风向x的距离变化而变化，在径向方向风速呈常数分布。这里考虑来流风切变u₀(z)，计算尾流区 x位置处的尾流风速u^*(x,z)，

u^*(x,z)＝u₀(z)[1-2a/(1+k_x,z·2s)²]

其中，ɑ为轴向诱导因子(下同)，k_x,z为修正后得到的尾流膨胀系数，s为量纲化下游位置，表达式为s＝x/D。

4.计算尾流膨胀半径

基于Jensen模型关于尾流是线性扩张的假设，风力机尾流膨胀半径r_x,r，可以通过如下线性表达式进行计算：

r_x,r＝k₀·x+r₁

式中，k₀为初始尾流膨胀系数，r₁为风轮后的特征尾流半径，由致动盘理论知其与风轮半径r_d相关，表达式为：

5.建立三维尾流模型

基于2D_k Jensen模型提出的cosine型速度分布，结合求得的尾流膨胀半径r_x,r，对初始尾流风速u^*(x,z)进行三维修正，得到修正后的尾流风速分布函数u(x,y,z)，

u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

其中，f_x,r为风速分布函数，f_x,r＝cos(2π/r_x,r·r+π)；r为以风轮中心线为圆心的辐射半径，表达式为：

式中，y为风力机尾流任一点的横向坐标，y_wake和z_wake分别为风力机下游x 位置的尾流中心坐标。

综上所述可得，新型三维尾流模型主要分为两步来实现：

预测步：u^*(x,z)＝u₀(z)[1-2a/(1+k_x,z·2s)²]

校正步：u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

一维Jensen模型、二维2D_k Jensen模型和新型三维尾流模型的尾流速度型示意图详见图1。

6.新模型验证

6.1用Nibe风电场实测数据验证模型有效性

Nibe风电场位于丹麦北部沿海地区，地表粗糙度估算为z₀＝0.070m，其内含两台风力机Nibe A和Nibe B，额定功率为630kW，风轮直径为 D＝40m，轮毂高度为z_hub＝45m。此外，Nibe B风力机下游2.5D、4.0D、6.0D、7.5D分别安装了四台测风设备，Nibe A则位于NibeB下游5.0D位置处。 Taylor等人对外场长期观测数据进行统计分析，在报告中给出了NibeB 风力机单独运行时各测风设备的一分钟平均尾流数据。由于2.5D位置为近尾流区，工程尾流模型不适用；6.0D位置会受到紧上游1.0D处Nibe A 风力机(处于停机状态)的影响，本文将剔除这两组数据，选取Nibe B 下游4.0D和7.5D位置为研究对象，验证新型3D尾流模型的有效性。轮毂高度处入流风速为u(z_hub)＝8.55m/s，入流湍流强度为I(z_hub)＝10％，相应条件下推力系数为C_T＝0.82。

步骤1，确定入流风速廓线u₀(z)和湍流强度分布型I₀(z)，根据对数分布律，得到

式中，z_hub为轮毂高度，z₀为地表粗糙度，u(z_hub)以及I(z_hub)分别为已知的轮毂高度处入流风速与湍流强度。

步骤2，计算尾流膨胀系数k_x,z，由于湍流强度直接影响尾流恢复快慢，需要同时考虑入流大气湍流以及尾流区附加湍流来对初始尾流膨胀系数k₀进行修正，得到尾流膨胀系数k_x,z：

k_x,z＝k₀·I_wake,x,z/I₀(z)

式中，初始尾流膨胀系数k₀可以通过考虑地表粗糙度计算获得，

k₀＝0.5/ln(z_hub/z₀)

I_wake,x_,z为风力机下游x位置，z高度处的湍流强度，计算式如下：

其中,

式中，I_add,x为尾流区x位置处的附加湍流强度，包含由于风力机的存在产生的剪切湍流和附加的机械湍流等；D为风轮直径；C_T为风力机推力系数。

步骤3，基于原Jensen模型对尾流风速进行预测，得到初始尾流风速u^*(x,z)，

u^*(x,z)＝u₀(z)[1-2a/(1+k_x,z·2s)²]

式中，ɑ为风力机轴流诱导因子，与风力机的推力系数C_T有关，计算公式为

s为量纲化下游位置s＝x/D。

步骤4，计算尾流膨胀半径r_x,r，

r_x,r＝k₀·x+r₁

式中,r₁为风轮后的特征尾流半径，表达式为

其中 r_d为风轮半径。

步骤5，利用步骤4所得尾流膨胀半径r_x,r，基于2D_k Jensen模型提出的cosine型速度分布，对步骤3所得的初始尾流风速u^*(x,z)进行三维修正，得到修正后的尾流风速分布u(x,y,z)，

u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

其中，f_x,r为风速分布函数，f_x,r＝cos(2π/r_x,r·r+π)；r为以风轮中心线为圆心的辐射半径，表达式为：

式中，y_wake和z_wake分别为风力机下游x位置的尾流中心坐标。

经过分析图2和图3可知：

经典Jensen模型为一维高帽(top-hat)型速度分布，在尾流中心高估了尾流速度，而在尾流区两侧严重低估了尾流速度，整体上与实测结果偏差较大。新型3D尾流模型给出的速度型形状与大小均与实测数据吻合较好，显著优于一维Jensen模型结果，且优于采用AD/RANS-RSM方法得到的数值模拟结果，显示了优良的预测性能。

6.2用Vestas V80-2MW风力机实测数据验证模型有效性

为了研究大气湍流强度对风力机尾流的影响，Porté-Agel等人采用旋转致动盘结合大涡模拟(AD-R/LES)的方法开展了不同粗糙度(反映为大气湍流强度)条件下单台风力机的尾流变化研究。这两种地表粗糙度和湍流强度，分别近似对应于陆上环境和海上环境，来流大气条件及风力机气动特性见表1。

表1测试算例及相应工况参数

步骤1，确定入流风速廓线u₀(z)和湍流强度分布型I₀(z)，根据对数分布律，得到

式中，z_hub为轮毂高度，z₀为地表粗糙度，u(z_hub)以及I(z_hub)分别为已知的轮毂高度处入流风速与湍流强度。

步骤2，计算尾流膨胀系数k_x,z，由于湍流强度直接影响尾流恢复快慢，需要同时考虑入流大气湍流以及尾流区附加湍流来对初始尾流膨胀系数k₀进行修正，得到尾流膨胀系数k_x,z：

k_x,z＝k₀·I_wake,x,z/I₀(z)

式中，初始尾流膨胀系数k₀可以通过考虑地表粗糙度计算获得，

k₀＝0.5/ln(z_hub/z₀)

I_wake,x,z为风力机下游x位置，z高度处的湍流强度，计算式如下：

其中,

式中，I_add,x为尾流区x位置处的附加湍流强度，包含由于风力机的存在产生的剪切湍流和附加的机械湍流等；D为风轮直径；C_T为风力机推力系数。

步骤3，基于原Jensen模型对尾流风速进行预测，得到初始尾流风速u^*(x,z)，

u^*(x,z)＝u₀(z)[1-2a/(1+k_x,z·2s)²]

式中，ɑ为风力机轴流诱导因子，与风力机的推力系数C_T有关，计算公式为

s为量纲化下游位置s＝x/D。

步骤4，计算尾流膨胀半径r_x,r，

r_x,r＝k₀·x+r₁

式中,r₁为风轮后的特征尾流半径，表达式为

其中 r_d为风轮半径。

步骤5，利用步骤4所得尾流膨胀半径r_x,r，基于2D_k Jensen模型提出的cosine型速度分布，对步骤3所得的初始尾流风速u^*(x,z)进行三维修正，得到修正后的尾流风速分布u(x,y,z)，

u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

其中，f_x,r为风速分布函数，f_x,r＝cos(2π/r_x,r·r+π)；r为以风轮中心线为圆心的辐射半径，表达式为：

式中，y_wake和z_wake分别为风力机下游x位置的尾流中心坐标。

经过分析图4-7b 可知：

新型3D模型对流向、横风向和垂直向的尾流速度亏损均体现了良好的预测精度，且能够有效地反映垂直方向速度亏损的非对称特性，其计算精度甚至优于基于CFD方法的数值模拟结果。另外，通过高湍流强度和低湍流强度两个算例，该模型很好地体现了环境湍流强度越大，尾流速度恢复得越快的规律，这与风力机空气动力学理论和观测到的尾流发展规律相符。

本发明经与多个风场实测数据对比，发现新提出的三维尾流模型可以很好地模拟尾流区的风速，不仅在精度上与实验结果接近一致，而且在三维分布上也更加符合真实流场，同时继承了工程模型的形式简单、编码容易、计算高效等优点。

以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，本发明要求保护范围由所附的权利要求书、说明书及其等效物界定。

Claims (3)

1.一种基于2D_k Jensen模型的三维尾流数值模拟方法，其特征在于，包括：

步骤1，计算入流风切变曲线u₀(z)以及湍流强度在垂直方向上的分布型I₀(z)，计算式如下：

式中，z为距离地面的高度，z_hub为轮毂高度，z₀为地表粗糙度，u(z_hub)为风力机轮毂高度处的风速，I(z_hub)为风力机轮毂高度处的湍流强度；

步骤2，计算尾流膨胀系数k_x,z，计算式如下：

k_x,z＝k₀·I_wake,x,z/I₀(z)

式中，k₀为初始尾流膨胀系数，I_wake,x,z为风力机下游x位置，z高度处的湍流强度；

所述I_wake,x,z的计算式如下：

其中,

式中,x为风力机下游距离，D为风轮直径,C_T为风力机推力系数,I_add,x为尾流区x位置处的附加湍流强度，包含由于风力机的存在而产生的剪切湍流和附加的机械湍流；

步骤3，基于Jensen模型对尾流风速进行预测，计算初始尾流风速u^*(x,z)，计算式如下：

u^*(x,z)＝u₀(z)[1-2a/(1+k_x,z·2s)²]

式中，ɑ为风力机轴流诱导因子，与风力机的推力系数C_T有关，计算式为

s为量纲化下游位置，计算式为s＝x/D；

步骤4，计算尾流膨胀半径r_x,r，计算式如下：

r_x,r＝k₀·x+r₁

式中，r₁为风轮后的特征尾流半径，表达式为

其中r_d为风轮半径；

步骤5，建立三维尾流模型，利用步骤4所得尾流膨胀半径r_x,r，基于2D_k Jensen模型提出的cosine型速度分布，对步骤3所得的初始尾流风速u^*(x,z)进行三维修正，得到修正后的尾流风速分布函数u(x,y,z)，其表达式为：

u(x,y,z)＝(u₀(z)-u^*(x,z))·f_x,r+u^*(x,z)

f_x,r为风速分布函数，其表达式为f_x,r＝cos(2π/r_x,r·r+π)，r为以风轮中心线为圆心的辐射半径，表达式为：

式中，y为风力机尾流任一点的横向坐标，y_wake和z_wake分别为风力机下游x位置的尾流中心坐标。

2.根据权利要求1所述的一种基于2D_k Jensen模型的三维尾流数值模拟方法，其特征在于：

步骤1中，风力机轮毂高度处的风速u(z_hub)和湍流强度I(z_hub)采用风速计测得。

3.根据权利要求1所述的一种基于2D_k Jensen模型的三维尾流数值模拟方法，其特征在于：

步骤2中，初始尾流膨胀系数k₀通过考虑地表粗糙度计算获得，计算式为：

k₀＝0.5/ln(z_hub/z₀)。

Images