大气边界层近地层是大气边界层中紧接地面的一层，是人类赖以生存和活动的地方。地球与大气的相互影响主要是通过近地层进行的[1]。风力机在利用风能的同时，改变了大气循环，使大气边界层内的风速降低、湍流强度增高。且大气边界层内湍流引起的各种物理量的交换和输送会改变局地能量分布和水汽交换，使空气温度、湿度等发生变化，进而可能会改变某些环境物理参数[2]。所以研究风力机对边界层近地层的影响，对风电场周边的环境变化研究有一定的意义。
胡菊[3]利用区域气候模式(Regional Climate Model 4.1，RegCM4.1)对酒泉千万千瓦基地的大型风电场长期气候效应进行了模拟研究，得出大型风电场建设30年边界层内温度、湿度等的变化特征。Roy等[4]应用区域大气模拟系统 (Regional Atmospheric Modeling System，RAMS)研究风电场对当地水文气象的潜在作用，结果证明风电场会对大气动力学等产生相应影响。Porté Agel[5-6]及Abkar[7-8]对整个风电场区域内的大气边界层流动进行了大涡模拟，分别采用不同的模型来研究风力机对大气边界层的影响。Wang等[9]利用全球气候模式 (Community Climate Model Version 3，CCM3)进行模拟实验，结果表明风电场增加了近地层动量的拖曳作用并且降低了风速，也加强了地表与大气之间的热湿交换。Adams[10]认为风电场阻碍风的流动，对风施加了额外的拖曳阻力，通过在其尾迹中产生湍流、降低风速等方式改变大气边界层，使地表通量、热湿对流等发生变化，进而对气候产生间接影响。Zhang等[11]通过风洞实验得出结论：尽管整体上风电场所引起的表面热通量变化较小，但是由此造成的不均匀的空间扰动却非常明显，这可能会影响到大气与地面之间动量、热量及水汽的传递运输。类似相关研究[12-14]也证明风电场可能会使其周边地区的环境气候发生变化。
计算流体力学(Computational Fluid Dynamics,CFD)是研究风力机较常用的方法，通常用于风力机水平方向的尾流流场分析及气动性能模拟等[15-16]。将CFD应用于大气边界层的模拟研究也较常见。但是通过CFD将两者结合，研究风力机运行对近地层影响的文献较少。本文采用FLUENT对单台风力机近地层的流场进行数值模拟，从不同的高度上分析两个基本参数即流速和湍动能的变化情况，为研究风电场引起的潜在环境影响提供指导。
1 计算模型及网格处理 1.1 几何模型利用叶素-动量理论在GAMBIT中建立1.2MW风力机模型，风力机参数如表 1所示。
| 风轮直径/m | 60 |
| 桨叶数量/个 | 3 |
| 轮毂中心高度/m | 60 |
| 轮毂直径/m | 2.8 |
| 风轮转速/(r·min-1) | 19.27 |
| 翼型系列 | NACA634/FX66S196 |
旋转流场是包含叶片在内的扁圆盘，由于主要考虑风力机下游的边界层近地层的流场变化情况，故忽略塔架的影响，整体流场区域是包含扁圆盘在内的整个流场。将整体区域简化为长方体状，流向沿z轴，风力机平面平行于xy平面，模型长宽高分别为1200m×300m×300m，风力机叶轮位于入口后3d处，d为风轮直径，以保证来流到达叶轮时已充分发展，模型如图 1所示。
|
| 图 1 整体流场区域尺寸 Fig. 1 Size of the whole flow field |
1.2 网格划分
针对上述建好的模型进行网格划分，采用四面体非结构化网格，应用尺寸函数对叶片、轮毂表面进行局部加密处理，流场部分采用相对较稀疏的非结构化网格。通过比较3种网格数量的模拟结果来验证网格无关性，以轮毂中心线上的湍动能(TKE)变化为参考，网格无关性验证如图 2所示。
|
| 图 2 网格无关性验证 Fig. 2 Grid independence calculation |
因此，考虑到计算时间及计算精度，选用418 万的网格模型进行计算，网格划分情况如图 3所示。
|
| 图 3 内部及外部流场区域的网格划分 Fig. 3 Grid distribution of the internal and external flow field |
1.3 边界条件
进口：边界层近地层速度分布遵循大气边界层指数方程。在我国，地貌按照地表粗糙度分为A、B、C、D四类[17]，不同的地貌类型对应着不同的地表粗糙度指数及风速廓线。本文针对A类地貌进行模拟，应用用户自定义函数(User-Defined Function,UDF)编译风速廓线。A类地貌对应的边界层指数方程即入口速度边界条件为：
|
(1) |
入口湍流边界条件：采用经常使用的湍动能k与比耗散率ω的一组湍流参数，湍动能取值为1，比耗散率取值为1。通过公式计算在风力机轮毂高度处的湍流强度为7%，湍流粘性系数为1.2。
出口：自由出流，边界条件为outflow。
叶片及轮毂：采用无滑移边界条件，设为旋转固体壁面wall。
流场底面：无滑移边界条件，固体壁面wall，A类地形所对应的地表粗糙高度为0.01m。
流场侧面及顶面：对称边界条件，模型顶部大气边界设置为零滑移壁面的对称边界，即symmetry。
选择旋转坐标系对风轮进行模拟，采用MRF(Multiple Reference Frame Model)模型，风力机匀速转动，转速为19.27r/min，风力机周边流体以相同的角速度旋转，其余流场静止。
由于风力机在低雷诺数下运行，马赫数也相对较低，故将模型简化为三维稳态不可压缩空气流动，流动过程与外界无换热(即不考虑能量方程)。采用Segregated隐式求解器，应用SST k-ω紊流模型，基于定常雷诺时均N-S方程(RANS)进行数值模拟，即用连续方程和动量方程来描述。速度与压力之间的耦合采用SIMPLIC算法实现，应用有限体积法对控制方程进行离散，对流项差分格式采用二阶迎风格式，计算收敛标准为连续方程的无量纲残差和所有变量降到1×10-3以下。
1.4 模拟方法验证为了验证数值模拟的正确性，将下游距离风力机6倍直径距离处的风力机尾流区速度剖面与风洞实验结果进行对比，如图 4所示。由图 4可知模拟结果在一定程度上与风洞实验相似，且Helmis[18]等学者曾指出风洞实验过高的预估了近尾流区的尾流效应，分析结果使风轮中心的速度亏损较大，所以在本文的模拟速度的计算结果要比实验结果略偏大。
|
| 图 4 CFD模拟的尾流区域速度剖面与风洞实验计算结果比较图 Fig. 4 Comparison between the velocity profiles from the CFD simulation and the wind tunnel experiment |
2 计算结果与分析
大气边界层近地层的主要特点是速度分布的规律性及运动的湍流特性，湍流交换过程决定了各种变量的分布及变化情况。所以从下游不同距离处速度和湍动能沿高度的变化情况来加以分析。
2.1 速度分布气流经过旋转的风力机后，由于风力机对来流的阻碍作用，使得下游流场发生很大变化。位于风轮后不同距离处的xy平面的速度分布云图如图 5所示。
|
| 图 5 不同截面上的速度云图 Fig. 5 Contours of velocity magnitude on different sections |
从图 5中可见，每个截面的尾流状态都与风轮的外形轮廓相似，尤其是近风轮区域，轮毂周围风速明显下降。随着空气继续向后流动，影响的面积在扩大，但是趋势逐渐变缓，不同距离处与风轮旋转平面相同高度处的速度变化速率逐渐变小，速度分布逐渐趋于来流速度分布情况，轮毂轴线高度处的速度逐渐增加，但直至流场出口即距离风力机叶轮17 d后，其速度仍低于相同高度处的来流速度。下游不同距离处速度随高度的变化的具体数值如表 2所示。
| 高度 | 20m | 40m | 60m | 80m | 100m |
| 2d | 10.702 | 10.273 | 10.325 | 12.222 | 13.202 |
| 4d | 10.586 | 10.642 | 10.735 | 12.128 | 13.173 |
| 6d | 10.522 | 10.805 | 11.050 | 12.145 | 13.182 |
| 8d | 10.572 | 10.947 | 11.179 | 12.132 | 13.108 |
| 10d | 10.584 | 11.022 | 11.248 | 12.047 | 13.052 |
| 12d | 10.603 | 11.026 | 11.337 | 12.132 | 13.020 |
| 14d | 10.629 | 11.037 | 11.428 | 12.212 | 13.101 |
| 17d | 10.802 | 11.720 | 12.453 | 12.914 | 13.267 |
从表 2中可以看出，在同一高度处，40m、60m处的速度均随距离的增加而增加，其他高度处的速度分布则不明显。在距离风力机同一位置处，除2d外，其它处的速度大体上随高度的增加而增加，具体中心线上速度变化情况如图 6所示。
|
| 图 6 下游不同距离处的风速的竖直分布图 Fig. 6 Diagram of the vertical profiles of wind speed at different height of different sections of downstream |
图 6表明，入口的速度分布情况符合采用的大气边界层函数公式。通过对比在风轮下游2d、6d、10d、14d、17d位置处中心线上不同高度的速度变化情况可以更加明显的看出：经过轮毂后，中心速度有很大的衰减；在轮毂以下，速度衰减得较快，所以此处也是最易出现漩涡的区域；轮毂以上高度虽然也有衰减，但是因为来流速度随高度的增加而增加，所以轮毂上部速度高于下部，因此轮毂上部分与下部分相比，衰减幅度略小。随着空气向下游的流动，低空低速气流与高空高速气流混合，风力机尾流区域与其周边区域的气流混合，使尾流的速度随距离的增大缓慢增加，但直至出口处仍然未达到来流速度。
2.2 湍动能分布情况湍动能是湍流强度的度量，是微气象学中的一个重要变量。湍流是始终占据边界层主导地位的流动，边界层内的气象要素的时空分布是湍流运动的直接结果。而由于风力机的运行搅动气流，在尾迹中产生湍流。模型的湍动能的模拟数据如表 3所示。
| 高度 | 20m | 40m | 60m | 80m | 100m |
| 2d | 0.194 | 0.223 | 0.376 | 0.3 | 0.292 |
| 4d | 0.138 | 0.138 | 0.123 | 0.190 | 0.239 |
| 6d | 0.115 | 0.104 | 0.118 | 0.151 | 0.180 |
| 8d | 0.106 | 0.089 | 0.101 | 0.127 | 0.155 |
| 10d | 0.108 | 0.081 | 0.093 | 0.114 | 0.136 |
| 12d | 0.098 | 0.078 | 0.087 | 0.107 | 0.121 |
| 14d | 0.097 | 0.077 | 0.086 | 0.099 | 0.111 |
| 17d | 0.095 | 0.076 | 0.082 | 0.090 | 0.098 |
从表 3中可以看出距离风力机不同距离处的湍动能沿竖直方向的分布情况。湍动能的变化情况可以分为两个区域分别进行讨论。在近尾流区域，湍动能变化较为剧烈。而在远尾流区域，湍动能由地面先减小直至叶轮底部，而后缓慢增加直至叶轮顶部偏上部位后几乎不再变化，具体变化情况如图 7所示。
|
| 图 7 下游不同距离处湍动能的竖直分布图 Fig. 7 Diagram of the vertical profiles of turbulent kinetic energy at different height on different sections downstream |
图 7为风力机叶轮下游不同位置处xy截面中心线上湍动能沿竖直方向的分布情况。从中图 7可以看出在与轮毂相同高度处，距离旋转叶轮2d处的湍动能在竖直向上的分布受风力机转动的影响比更远处要大的多，湍动能明显分布不均匀，流场扰动非常强烈。这说明离旋转叶轮越近，虽然轴向速度有很大的衰减，但是湍动能却有很大的增加，因为轮毂处速度为0，轮毂周围的气体混合要更加剧烈，所以湍动能较高。随着向下游的延伸，相比于靠近叶轮处，湍动能受到叶轮转动的影响逐渐降低，流动逐渐趋于缓和，所以在相同高度处，下游比上游湍动能要低。
在竖直方向上，由地面到旋转叶轮底部的过程中，由于地表粗糙度不为0，故近地面的速度由0突然增大，所以地面处湍动能较大，向上则脱离地面，与空气接触，湍动能逐渐减小。直至达到叶轮底部叶尖附近，此时速度突然增大，湍动能随之开始增强，而由底部叶尖到轮毂线速度则逐渐减小，但是由于与周围气流的混合搅动，相对于叶尖，轮毂处速度最低，混合越强，所以在轮毂高度处的湍动能也越大，混合也就越剧烈。同理，由轮毂中心至叶轮顶部叶尖附近，湍动能则逐渐降低。但是由于高空中的速度要远大于低空中叶轮底部的速度，所以高空中叶轮顶部的湍动能要高于低空叶轮底部的湍动能。进而由于叶尖处的扰动作用，轮毂顶端再往上部分的湍动能要相应的增大，直至一定高度后维持不变。
3 结 论利用数值模拟软件，结合大气边界层理论及风力机尾流理论，通过设置的表粗糙度，加载边界层内速度分布函数，针对单台风力机运行后对A类型地貌对应的大气边界层近地层的影响进行模拟计算，分析后得出以下结论：
1) 在靠近风力机叶轮附近，叶轮中心轴向速度衰减最快；随着尾流逐渐扩散，速度逐渐增加，但速度增加的速率逐渐变缓，尾迹范围逐渐扩大。
2) 由于风轮的扰动作用，边界层内的湍动能也发生变化。近尾迹区域，由地面至高空，湍动能呈现出先减小后增大，再减小再增大的趋势；远尾迹区域则先减后增；但随着向下游的延伸，湍动能逐渐降低且变化速率也逐渐降低。
3) 经过风力机后，原本速度分布较规律的大气边界层近地层受到扰动，速度和湍动能均发生剧烈变化，且在下游距离风轮17d的位置处仍未恢复到原来状态，而近地层的状态影响着动量、热量、水汽等的分布，决定着气象要素的变化，所以风力机对周边环境的影响可能不仅仅局限在风速和湍动能方面。
本文利用CFD模拟风力机下游的大气边界层近地层的流动变化过程，从最基本的两个参数进行分析，得到了风力机下游近地层的速度和湍动能在不同距离处沿高度的变化规律，可以用于优化风能工程设计，同时对预测风电场对当地环境的影响等方面有一定的指导意义。
| [1] |
Liu H Z, Feng J W, Wang L, et al. Overview of recent studies on atmospheric boundary layer physics at LAPC[J].
Chinese Journal of Atmospheric Sciences, 2013, 2(37):467–476.
(in Chinese) 刘辉志, 冯建武, 王雷, 等. 大气边界层物理研究进展[J]. 大气科学, 2013, 2(37) : 467–476. |
| [2] |
Zhao Z C, Luo Y, Jiang Y. Advances in assessment on impacts of wind farms upon climate change[J].
Advances in Climate Change Research, 2011, 7(6):400–406.
(in Chinese) 赵宗慈, 罗勇, 江滢. 风电场对气候变化影响研究进展[J]. 气候变化研究进展, 2011, 7(6) : 400–406. |
| [3] |
Hu J. Numerical simulation research on impact of large-scale wind farms on regional climate[D]. Lanzhou: Lanzhou University, 2012. (in Chinese) 胡菊. 大型风电场建设对区域气候影响的数值模拟研究[D]. 兰州: 兰州大学, 2012. |
| [4] | Roy S B, Traiteur J J. Impacts of wind farms on surface air temperatures[J]. Proceedings of the National Academy of the Sciences of the United States of America, 2010, 107(42):17899–17904. DOI:10.1073/pnas.1000493107 |
| [5] | Porte-Agel F, Wu Y T, Lu H, et al. Large-eddy simulation of atmospheric boundary layer flow through wind turbines and wind farms[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2011, 99(4):154–168. DOI:10.1016/j.jweia.2011.01.011 |
| [6] | Porte-Agel F, Wu Y T. Atmospheric turbulence effects on wind-turbine wakes: an LES study[J]. Energies, 2012, 6(5):5340–5362. |
| [7] | Abkar M, Porte-Agel F. The effect of free-atmosphere stratification on boundary-layer flow and power output from very large wind farms[J]. Energies, 2013, 6(5):2338–2361. DOI:10.3390/en6052338 |
| [8] | Abkar M, Porte-Agel F. Mean and turbulent kinetic energy budgets inside and above very large wind farms under conventionally-neutral condition[J]. Renewable Energy, 2014, 6(5):142–152. |
| [9] | Wang C, Prinn R G. Potential climatic impacts and reliability of very large-scale wind farms[J]. Atmospheric Chemistry and Physics, 2010, 10(4):2053–2061. DOI:10.5194/acp-10-2053-2010 |
| [10] | Adams A S, Keith D W. Wind energy and climate: modeling the atmospheric impacts of wind energy turbines[C]//American Geophysical Union Fall Meeting Abstracts. San Francisco, American, 2007. |
| [11] | Zhang W, Corey D M, Porte-Agel F. Experimental study of the impact of large-scale wind farms on land-atmosphere exchanges[J]. Environmental Research Letters, 2013(8):1–8. |
| [12] | Zhou L M, Tian Y H, Roy S B, et al. Impacts of wind farms on land surface temperature[J]. Nature Climate Change, 2012, 2(7):539–543. |
| [13] | Roy S B. Simulating impacts of wind farms on local hydrometeorology[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2011, 99:491–498. DOI:10.1016/j.jweia.2010.12.013 |
| [14] | Barrie D B, Kirk-Davidoff D B. Weather response to a large wind turbine array[J]. Atmospheric Chemistry Physics, 2010, 10:769–775. DOI:10.5194/acp-10-769-2010 |
| [15] |
Li S H, Yue W P, Kuang Q F, et al. Numerical simulation on flow field of a 1.2 MW wind turbine[J].
Journal of Chinese Society of Power Engineering, 2011, 31(5):101–107.
(in Chinese) 李少华, 岳巍澎, 匡青峰, 等. 双机组风力机尾流互扰及阵列的数值模拟[J]. 中国电机工程学报, 2011, 31(5) : 101–107. |
| [16] |
Wang Q F, Chen B J, An Y R. 3D numerical simulation for aerodynamic performance of large scale wind turbine[J].
Acta Aerodynamica Sinica, 2011, 29(6):810–814.
(in Chinese) 王琦峰, 陈伯君, 安亦然. 大型风力机三维空气动力学数值模拟[J]. 空气动力学学报, 2011, 29(6) : 810–814. |
| [17] |
Wang T T. Large eddy simulation of atmospheric boundary layer flow based on FLUENT[D]. Beijing: Beijing Jiaotong University, 2011. (in Chinese) 王婷婷. 基于FLUENT的大气边界层风场LES模拟[D]. 北京: 北京交通大学, 2011. |
| [18] | Helmis C G, Papadopoulos K H, Asimakopoulos D N, et al. An experimental study of the near wake structure of a wind turbine operating over complex terrain[J]. Solar Energy, 1995, 54(6):413–428. DOI:10.1016/0038-092X(95)00009-G |