一种通过添加ε源项来校正标准k-ε模型的方法与流程

本发明属于风力发电领域，涉及一种通过添加ε源项来校正标准k-ε模型的方法。

背景技术：

风能是一种清洁的可再生资源，因此，风力发电技术近几年来在我国得到了快速发展，并推动了风电场的大量建设。早起的风电厂大多建设在风力资源丰富、地形较为平坦的地方，但是随着风电场的规模越来越大、数量越来越多、风电场开始建设在复杂地形区域。因此，采用标准线性模型模拟风电场条件，不能很好地再现复杂地形的风资源条件。而采用cfd技术进行模拟则变得越来越常见并且有良好的应用前景。综合考虑计算精度，运算时长，以及计算资源等条件，雷诺时均法仍是目前最经济最普遍的计算方法。考虑地表粗糙度的影响，采用带有壁面函数的标准k-ε湍流模型是目前工程计算中普遍采用的方法。

山丘地形流动最显著的特征是分离流动，与平板地形下的附着流动有很大区别。有些可以较好模拟平板地形大气边界层的方法，却不能很好地模拟分离区。因此，直接将平板地形的方法运用到山丘地形中是不合适的，

技术实现要素：

本发明的目的在于克服上述现有技术的缺点，提供了一种通过添加ε源项来校正标准k-ε模型的方法，该方法能够有效提高k-ε模型模型山丘地形风资源的准确性。

为达到上述目的，本发明所述的通过添加ε源项来校正标准k-ε模型的方法包括以下步骤：

1)获取山丘地形下的测风数据；

2)构建山丘地形模型，再进行网格划分，然后利用标准k-ε湍流模型计算；

3)分析步骤2)的计算结果与实际测风数据之间的误差，并根据分析结果设定ε源项，然后利用ε源项对k-ε湍流模型进行修正，最后利用修正后的k-ε湍流模型进行复杂地形风资源计算。

步骤1)的具体操作为：

确定待研究的山丘地形，再选取两个激光雷达的安装位置，然后再安装激光雷达，再利用激光雷达测量山丘地形的测风数据。

选取两个激光雷达的安装位置的具体过程为：

确定当地主风向方向，然后在山丘迎风面的山脚位置处放置第一个激光雷达来确定入口边界条件，在山丘背风面的山脚位置处放置第二个激光雷达。

步骤2)中还包括：确定入口条件中变量的取值范围。。

本发明具有以下有益效果：

本发明所述的通过添加ε源项来校正标准k-ε模型的方法在具体操作时，根据分析结果设定ε源项，然后利用ε源项对k-ε湍流模型进行修正，最后利用修正后的k-ε湍流模型进行复杂地形风资源计算，以弥补采用传统修正方法计算不准确的不足，以定性的观测不同来流条件和不同地形对流动的影响。

进一步，两个激光雷达即可完成对山丘地形风资源的监测，极大的节约人力物力资源。

附图说明

图1为风电场布机的点位图；

图2a为lidar387的风向玫瑰图；

图2b为lidar366的风向玫瑰图；

图3a为来流wnw扇区的计算域的选取图；

图3b为计算域边界的处理图；

图3c为计算域的网格划分图；

图4为280°进口两雷达所在平面的速度云图；

图5a为lidar387的风速廓线与实测数据对比图；

图5b为lidar366的风速廓线与实测数据对比图。

具体实施方式

下面结合附图对本发明做进一步详细描述：

本发明所述的通过添加ε源项来校正标准k-ε模型的方法包括以下步骤：

1)获取山丘地形下的测风数据；

2)构建山丘地形模型，再进行网格划分，然后利用标准k-ε湍流模型计算；

3)分析步骤2)的计算结果与实际测风数据之间的误差，并根据分析结果设定ε源项，然后利用ε源项对k-ε湍流模型进行修正，最后利用修正后的k-ε湍流模型进行复杂地形风资源计算。

步骤1)的具体操作为：

确定待研究的山丘地形，再选取两个激光雷达的安装位置，然后再安装激光雷达，再利用激光雷达测量山丘地形的测风数据。

选取两个激光雷达的安装位置的具体过程为：确定当地主风向方向，然后在山丘迎风面的山脚位置处放置第一个激光雷达来确定入口边界条件，在山丘背风面的山脚位置处放置第二个激光雷达。

步骤2)中还包括：确定入口条件中变量的取值范围。

实施例一

本实施例的具体操作过程为：

1)确定地形，选却两个激光雷达的安装位置，然后在激光雷达安装位置安装激光雷达，利用激光列大测量实际地形下的监测数据；

由于地形包含的种类繁多，不能一一枚举，本实施例以河北某复杂地形风场为例，说明本发明的适用性，其他种类的复杂地形方法类似。

怀来某风电场位于115°33’e～115°38’e、40°22’n～40°26’n之间，图1为风电场布机点位图。

在风电场左下区域布置有两个测风激光雷达，分别为lidar387及lidar366，两激光雷达的直线距离为551.4m，二者连线与x轴正方向夹角约为160.2°。每个激光雷达都测量了从2019年4月1日至2019年5月6日为期36天的当地测风数据，包括风向、风速、风速最大值最小值及可利用率等，涵盖从距地面40m到230m间距不等的高度范围，根据两激光雷达的测风数据，制作16扇区的风向玫瑰图，如图2a及图2b所示。

由风向玫瑰图可得，该风电场lidar387处，90m高度主风向为ese,占总风向频率的13.72％，次主风向为wnw，占总风向频率的12.84％。lidar366处，90m高度主风向为wnw，占总风向频率的13.11％，次主风向ese，占总风向频率的12.08％。主风向和次主风向方向较明显。对风机的布置有利，同时风向的集中有利于风机的稳定运行。

入口速度廓线为：

其中，u表示来流速度，z表示对应高度，z0为地表粗糙长度，z0可根据风电场下垫面种类确定，κ为卡门常数，默认取0.4，u*为摩擦速度。

2)对地形建模并进行划分网格，利用cfd软件中的标准k-ε湍流模型进行计算，设计算域按照正南正北方向设置，计算域的选取原则为尽量使lidar387和lidar366的位置在所选计算域的中间，同时计算域的入口边界尽平整，方便设置来流边界条件。

参考图3a、图3b及图3c，wnw扇区下，计算域的尺寸为13800m×13500m×2000m，网格节点数为199×199×205，x、y方向均匀分布，z方向第一层网格高度为0.05m，延展比为1.04，网格总数为811万，计算域的西、北两面为进口边界条件，东、南两面为出口边界条件，底面为无滑移固壁，顶部为对称边界，数值计算采用ansys-fluent软件，采用标准k-ε湍流模型。

对两个主风向分别进行计算，发现在入口处给定对应扇区内的一个方向，气流在流场中经过起伏的地形均会发生一定偏转，等到流经lidar387和lidar366位置处时，速度方向都比入口处方向沿顺时针偏转了约20°左右，因此，将wnw扇区的来流方向定为280°，将ese扇区的来流方向定为100°。

当进口为280°时，气流流动到两测风雷达处的速度方向均为302°。根据现有的测风数据，以302°为基础选取振幅分别为2°和5°的两组范围，求得每个范围内不同高度下的平均速度，根据测风数据随高度变化的规律，选择摩擦速度u*＝0.23m/s，粗糙长度z0＝0.05m，图4为计算结果云图，截取的平面为两激光雷达所在平面，左侧黄线为lidar387位置，右侧黄线为lidar366位置，流动方向从左向右。

由云图可知，流动发展到两测风雷达处，基本保持均匀大气边界层的流动，流动状态基本相似，最大的区别在于lidar387位于坡前迎风面由于气流加速而边界层较薄，而lidar366位于坡后背风面气流减速而边界层偏厚，速度廓线对比如图5a及图5b所示。

3)对步骤2)的计算结果与实测数据之间的误差进行分析，通过添加源项来进行修正，以提高数值模拟的准确度，图中横坐标代表风速，纵坐标代表距地面高度，观察发现默认设置下，计算u*＝0.23m/s，z0＝0.05m在100m以上拟合相对较好，在100m以下比测风数据偏大，调整参数为u*＝0.22m/s，z0＝0.05m进行计算，发现在100m以下速度相对之前减小了，同时100m以上风速也减小了，偏离了测量数据，在调整入口参数无法优化的情况下，通过对算例添加源项sε，改进近壁面速度，经过计算确认，最合适的源项如下：

计算结果u*＝0.22m/s，z0＝0.05m，e如图5所示，在添加源项之后，增大了近壁面处的耗散率，使得近壁面初速度相对不加源项减小，而上部的速度基本不变，对lidar387和lidar366都有较好的改善作用。

表1及表2分别为lidar387和lidar366处模拟结果与测风数据的误差对比，其中，高度表示测风雷达测量的高度，从40m到230m不等距分布，实测风速表示对应高度处，302°±5°风向范围内的平均速度。标准风速为算例u*＝0.23m/s，z0＝0.05m，改进风速为算例u*＝0.23m/s，z0＝0.05m加源项，标准误差和改进误差分别为各自风速与实测风速间的误差，计算方法为标准风速(改进风速)减去实测风速，然后求绝对值，最后在除实测风速，由结果可以得到，无论lidar387还是lidar366，加源项改进后的方法整体上都要比不加源项的要好。其中，lidar387处标准风速各个高度的总误差为48.44％，而改进后总误差为36.21％，在轮毂高度90m处，改进误差也比标准误差小3％；lidar366处标准风速各高度的总误差为62.44％，而改进后总误差为46.44％，同样在轮毂高度90m处，改进误差比标准误差小4.23％，因此采用改进模型不论是在关键局部还是整体都提高了模拟的准确性，有助于更好地进行风电场风资源分析。

表1

表2

此实施例仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。在不脱离本发明主旨和范围的前提下，本发明还会有进一步的改进，所作改进仍在要求保护的本发明范围内，因此，本发明的保护范围应该以权利要求的保护范围为准。