1 引言
为了应对全球气候变化和满足清洁能源的需求, 风力发电作为一种清洁且相对成熟的能源技术, 在过去的几十年中得到了快速的发展。根据中国能源局发布的《风电发展“十三五”规划》可知, 风电已成为中国继火电和水电之后的第三大电力能源（国家能源局, 2016）。截止2017年, 中国累计装机容量已达到188 GW, 居世界第一。在中国北方地区, 已有9个千万千瓦级的风电基地, 包括新疆哈密、 甘肃酒泉、 河北张家口和承德以及内蒙古等大型风电基地。根据中国目前的能源规划, 截止2020年, 风电累计装机量将超过210 GW, 预计占我国总发电量的6%左右。
风轮机作为一种风力发电设备, 能将大气中的空气动能转换为机械能并进一步转换为电能（李国庆等, 2016；李照荣等, 2018; 刘丽珺等, 2018）。一方面大气中气象条件（风速、 风向）会对风力发电机的运行存在影响, 另外一方面风力发电机能够影响大气边界层中质量、 动量、 能量以及水汽的传输（赵宗慈等, 2011；张双益等, 2017）。目前对风电场的研究主要是通过模型模拟和现场观测（雷达、 卫星等）等方法（胡菊, 2012）。大量的模拟研究表明, 风轮机群能够有效地减少轮毂高度处的风速并产生湍动动能（Turbulent Kinetic Energy, TKE）（Roy, 2004；Christiansen et al, 2005；Fitch et al, 2012；Jimenez et al, 2015）。Roy et al（2004, 2011）和Zhou et al（2012）发现, 风轮机产生的湍动动能增强了空气的垂直混合, 从而促进地面边界层物质和动量的垂直交换, 影响地表温度的变化。另外, 部分模拟实验还表明风轮机还会造成大气边界层高度的增加（Frandsen et al, 2006；Fitch et al, 2013；Porté-Agel et al, 2014）。Fiedler et al（2011）利用WRF模式, 结合1948 -2009年的气象数据, 对美国中部地区的风电场进行模拟, 结果指出大型风电场会造成区域内及其周边地区1%的降水增加。Li et al（2018）在对撒哈拉地区大型风电场模拟的过程中发现, 风能发电场会形成地表摩擦增大—降水增加—植被增加的正反馈机制, 造成局部降水的增加。
观测数据也表明风电场对气象条件存在着明显的扰动和影响。Christiansen et al （2005）使用SAR雷达对丹麦的两个大型海上风电场的进行观测, 发现在轮毂高度处存在着8%~9%的风速衰减。Zhou et al（2012）使用MODIS卫星遥感数据发现美国德克萨斯州的一个大型风电场内部和临近区域出现0.72 ℃的变暖效应。Xia et al（2016）研究指出, 风电场会对地表温度产生一定的增温效果, 且夏季夜间风电场所造成的地表温度升高效应略微高于冬季, 且地表温度变化的大小主要取决于风轮机所造成的湍动动能和背景湍动动能之比。风力发电场除了影响其临近区域, 还会对其下游地区风速、 地表温度和降水等产生明显的影响（Jimenez et al, 2015；Fiedler et al, 2011；Li et al, 2018；Wang et al, 2010；Xia et al, 2017；Sun et al, 2018；吴正人等, 2014）, 如Walsh-Thomas et al (2012)利用卫星遥感数据发现在风电场下游存在变暖效应。模拟和观测研究均发现, 对于单个风电场而言, 风电场产生的尾流能够延伸至风机下游超过20 km处, 且尾流的量级和距离取决于风电场的大小以及大气边界层的稳定性（Christiansen et al, 2005；Fitch et al, 2012；Roy, 2011；Barrie et al, 2010）。Mo et al (2017) 发现大型风电场对大气污染物的分布也有着明显的影响, 具体表现为在风电场上下游存在“边缘效应”: 在风机上游区, 由于地表粗糙度的增加所强迫出的内边界层, 使得空气运动减缓, 有利于大气污染物的聚集。其结果表明在风能场上游二氧化氮（NO2）出现增加, 风机下游区则相反。
为了探究中国北方大型风电场对区域内气象要素及其邻近和下游地区天气预报的影响, 我们使用了大气数值模式WRF（Weather Research and Forecasting Model）模拟了我国北方大型风电场对华北平原地区气象要素的影响, 并通过相关的敏感性实验分析了风电场的具体作用。我们设计了四个情景实验, 具体实验模拟方案如表1所示, 其中模型的相关参数化方案将在第2.2节介绍。
表1 模拟情景Table 1 Modeled scenarios |
| 情景 | 模拟方案 |
|---|---|
| CON | 无风电场 |
| S1 | 地表粗糙度参数化方案 |
| S2 | 风机曳力参数化方案 |
| S3 | 与S2一致但风电场面积扩大一倍 |
2 模式配置和参数化方案
2.1 模式配置
本研究采用WRF模式（3.8.1版本）, 模拟区域如图1所示（图1和文中其他涉及中国国界和行政区划的地图均基于中国科学院资源环境数据云平台2015年中国省级行政边界shape矢量数据直接绘制, 数据下载链接为http: //www.resdc.cn/data.aspx?DATAID=200）。模拟的中心网格点经纬度为42.9°N, 110.0°E, 采用一层网格嵌套, 其中外层网格包括90×50个网格点, 分辨率为50 km, 内层网格包括91×91个网格点, 分辨率为10 km。模式垂直方向上分为29层, 最顶层为50 hPa。为了更好地模拟风轮机和大气边界层之间的相互作用, 本研究中提高了模式低层的垂直分辨率（大约10, 60, 90, 120, 150, 180, 210 m等）。主要的物理机制参数化方案如表2所示。模型使用NCEP FNL (Final) 再分析资料作为气象驱动场, 水平分辨率为1°×1°, 模拟时间包括了2016年1月（在下文中用冬季指代）和7月（在下文中用夏季指代）。
表2 主要物理参数化方案Table 2 Major physical parameterization schemes |
| 物理过程 | 参数化方案 |
|---|---|
| 微物理过程 | Morrison double-moment scheme |
| 长波辐射 | RRTMG scheme |
| 短波辐射 | RRTMG shortwave |
| 近地面过程 | MYNN surface layer |
| 路面过程 | Noah Land Surface Model |
| 边界层方案 | MYNN 2.5 level TKE scheme |
| 积云过程 | Grell 3D ensemble scheme |
2.2 风电场参数化方案
本研究中模拟的风能发电场位于中国北方地区, 和风电场有关的信息主要来自中国自愿减排交易信息平台（China Certified Emission Reduction Exchange Info-Platform, CCER）和谷歌地图等, 其中轮毂平均高度为96 m, 叶片直径为113 m, 平均风轮机间距为678 m（6倍叶片直径）。在模拟研究中我们使用了两种风机参数化方案。第一种方案主要考虑到风机增加了地表粗糙度（Barrie et al, 2010；Mo et al, 2017；Keith et al, 2004；Kirk-Davidoff et al, 2008）。利用Lettau（1969）提出的粗糙度方案计算风轮机的地表粗糙度:
,
式中: h *为风机轮毂高度(单位: m)；s为风力涡轮机叶片的旋转区域(单位: km2)；S为风机密度（S=A/n, A为风电场所占区域, n为风机数量）。风机间距设置为风机转子直径的6倍, s和S经计算分别为0.01 km2和0.46 km2, 从而得到风电场的粗糙度z 0为1.04 m。这种方案优点主要在于其表征风轮机的简洁性和准确性（Wieringa, 1993；Petersen, 1997）。
第二种参数化方案将风轮机考虑为大气动能的汇和湍动动能的源（Roy, 2004；Fitch et al, 2012；Jimenez et al, 2015；Fiedler et al, 2011；Xia et al, 2016, 2017 Barrie et al, 2010； Mo et al, 2017；Cervarich et al, 2013）。这里所使用的风机参数化方案为风机曳力参数化方案（Fitch et al, 2012）。这个参数化方案将大气动能通过推力系数CT转换为机械动能, CT的值取决于风轮机的类型和风速。风轮机分别通过电力转换系数CP和湍流转换系数C TKE将机械动能转换为电能和湍动动能, 且C TKE = CT - CP, 该参数化方案忽略了风机自身的机械损耗能量。系数CT和CP可从风轮机制造商处得到, 本研究将CT设为0.16, 额定功率为1.5 MW（Mo et al, 2017）。假设风轮机的方向与风向垂直, 则风轮机产生的阻力 如下:
,
式中: CT为风轮机推力系数； ρ为空气密度； A=(π/4)D 2为转子横截面积（D为涡轮叶片直径）； V（u, v）为水平风速。则风力涡轮机造成的大气动能损失率 为:
,
风轮机产生的电能 为:
,
风轮机产生的湍动动能 为:
.
3 模拟结果
3.1 冬季月平均
图2为2016年1月风电场对轮毂高度处月平均风速的影响, 在三个风机情景模拟中, 均可以看到在风电场内部及其临近地区都出现了显著的的风速衰减, 其最大值可以达到6 m·s-1左右。由于风轮机将大气运动动能转换为机械能, 使风力发电机轮毂高度处的风速出现明显的衰减, 尤其是对于大型的风电基地而言更为明显（Roy, 2004；Christiansen et al, 2005；Fitch et al, 2012）。通过对比可以发现, 粗糙度参数化方案造成的风速衰减比风机曳力参数化方案造成的风速衰减更明显[图2(a), (b)], 且当风能场面积扩大时, 其所造成的风速衰减也会随之增大[图2(b), (c)]。
图3和图4分别为风电场在轮毂高度处（大约100 m）和地面位置对温度的扰动, 三个风能发电场的模拟情景（S1, S2和S3）均表现为温度的升高。地表和轮毂高度处温度扰动变化较为一致, 温度升高的区域主要位于风电场内部和其下风区, 其温度变化最大值接近了1.2 °C。地表空气温度变暖效应与之前的模拟和观测到的结果较为一致（Zhou et al, 2012；Li et al, 2018；Xia et al, 2017；Roy et al, 2010）, 但变化幅度更大。推断这主要是由于风电场规模的差异所造成的（Roy et al, 2010）。当风电场面积扩大时, 温度扰动也随之扩大, 最大可达到1.5 °C[图3(c), 图4 (c)]。除了地表温度的增加, 还在轮毂高度处发现了明显的变暖效应。由于风轮机在轮毂高度上空会产生湍动动能, 推断轮毂高度处的变暖效应主要是由于湍流逆温效应造成的: 假设初始气层温度递减率 小于干绝热直减率 , 初始位温垂直分布 , 由于风轮机轮毂高度处产生大量湍动动能, 增强空气的垂直混合, 混合后气层温度递减率 接近于 , 上层热量向下输送（盛斐轩等, 2003）。
图5为在风电场内部[网格点（22, 52）]、 风电场下游50 km [网格点（20, 53）]、 风电场下游300 km [网格点（15, 53）]处CON（无风电场）和S2（有风电场）情景中风速、 TKE和温度的垂直廓线, 可以看出, 在风电场内部, 与情景CON相比, 情景S2中地面到400 m处的风速形成了显著的衰减, 而更高层的风速则无明显的变化[图5(a)]。风电场内湍动动能[图5(b)]也出现了明显的增加, 变化的高度从地面一直延伸至1000 m处, 其中变化最显著的高度约为150 m, 位于轮毂高度上空。在风电场内部气温的扰动随高度发生明显的变化: 从地面到180 m, 扰动为正值, 即风电场使低层气温升高；从180~1000 m, 扰动为负值, 即风电场使高层气温降低。低层温度的垂直变化与图2和图3基本保持一致。而对于风电场的下游区, 风速衰减的程度逐渐减小[图5(d), (g)], 但风速变化的高度却有所升高, 如在图5(g)中风速发生衰减的高度延伸到了1000 m左右。在风电场下游50 km处, 湍动动能[图5(e)]在不同高度处表现出相反的趋势: S2情景中, 从地面至180 m高度湍动动能小于CON情景, 而高于180 m高度时则大于CON情景。在风电场下游300 km处, 有风电场的S2情景中地面到300 m高度湍动动能出现了一定程度的下降, 却并不显著[图5(h)]。推断这种变化主要是由于风机叶片和临近区域（60~400 m）风速衰减使其下方风切变减小, 进而导致湍动动能的减少, 而对于上层的湍动动能则相反。在风电场下游50 km处温度[图5(f)]的垂直变化趋势类似于风电场内部[图5(c)], 但扰动的大小出现一定程度的下降, 而在风电场下游300 km处[图5(i)]从地面到600 m出现了较小的增温, 但增温的量级和变化高度进一步减小。
图5 冬季风速、 TKE和温度在风电场内部（a~c）、 风电场下游50 km（d~f）、 风电场下游300 km（g~i）上CON情景和S2情景的垂直廓线Fig.5 The vertical profiles of wind speed, TKE and temperature for CON scenario and S2 scenario within wind farms (a~c), at 50 km downstream of wind farms (d~f ) and at 300 km downstream of wind farms (g~i ) in winter |
3.2 夏季月平均
从风电场对轮毂高度处风速的影响（图6）可以看出, 冬季在风电场内部及其临近地区存在明显的风速衰减, 其空间分布以及不同模拟情景之间的差异变化与冬季在一定程度上保持一致, 但风速衰减的量级偏小, 最大的风速衰减量约3 m·s-1。这主要由于冬季盛行西北风, 风向较为一致, 而夏季风场的变化则与局地环流有关, 风向和风速变化较大, 而风向和风速的频繁变化会对风力发电机的运行状态和效率造成较大影响（张慧宁, 2010）。图7和图8分别显示了风电场对轮毂高度处和地表位置温度的影响。夏季（7月）风电场所引起的温度变化在空间分布和变化程度上都与冬季（1月）的情形存在显著差别。受季风的影响, 我国华北地区夏季盛行偏东南风, 图7和图8中温度升高的区域位于风电场的临近地区, 即河北省中部地区（王奕丹等, 2019）。同时, 夏季风电场所造成的温度扰动也明显弱于冬季（1月）, 最大值约为0.4 °C[图7(b), 图8(b)]。对于不同的情景模拟（S1, S2和S3）, 温度的变化也存在较大的差别, 风机曳力参数化方案（S2）的影响大于地表粗糙度参数化方案（S1）, 而且当风能场面积扩大一倍后, 风机对温度造成的扰动反而随之减小[图7(c), 图8(c)]。需要指出的是, 在内蒙古东部, 温度在轮毂高度处和地表位置都出现了一定程度的下降, 其最大值约为-0.3 °C, 且对于风电场面积扩大的模拟情景而言更为显著。
图6 2016年7月风电场对轮毂高度（约100 m）风速的影响（单位: m·s-1）Fig.6 The impact of wind farms on the wind speed at hub-height (about 100 m) in July 2016.Unit: m·s-1 |
图7 2016年7月风电场对轮毂高度（约100 m）温度的影响（单位: °C）Fig.7 The impact of wind farms on the temperature at hub-height (about 100 m) in July 2016.Unit: °C |
4 结论和讨论
本研究使用耦合的大气模式WRF模拟了中国北方大型风电场对于大气的影响。在本次模拟实验中地表粗糙度参数化方案和风机曳力参数化方案表现出较为明显的一致性。通过一系列敏感性实验得到如下主要结论:
（1） 中国北方大型风电基地显著地减少了局地的大气动能, 在风电场内部和其临近区域造成了风速的明显衰减, 其风速衰减效应会随着风电场距离的增加而减少。同时受不同季节环流的影响, 在夏季风速的衰减现象比冬季偏弱。
（2） 由于风轮机在轮毂高度及更高处产生了湍动动能, 增加了空气的垂直混合形成湍流逆温效应, 最终使地面和轮毂高度（约100 m）处的温度增加, 这种现象在冬季（1月）表现的更加明显。未来风电场（本研究中将风电场面积加倍）对大气造成的影响在不同季节存在一定程度的差异, 在冬季（1月）会促进轮毂高度处和地表位置温度的上升, 而在夏季（7月）导致风电场上游区增温趋势的减弱以及下游地区降温趋势的增强。
风轮机不仅增加了相应区域内的地表粗糙度, 且在将大气动能转换为机械能时使大气边界层中动能减少, 在风电场内部和临近地区存在明显的风速衰减。另外风轮机在将风能转换为电能时, 还会生成大气湍动动能, 增加动量、 热量以及其他物质的垂直混合, 导致气象要素（如温度）的变化。考虑到风电场的地表粗糙度增加效应和对大气动能的转换效应, 在模型中分别使用地表粗糙度参数化方案和风机曳力参数化方案来表征风电场的影响。另外基于风电场的持续发展, 在模型中通过将风电场面积增加一倍来研究未来风电场对大气的影响。需要强调的是, 风电场对于局地气候的影响还存在很多的不确定性, 包括数值模式模拟, 风机参数配置（如风机性能、 数量、 高度等）以及相关参数化方案的不确定性, 因此对于大规模的风电场在局地以及更大尺度上的影响还需要进一步研究。