我国西北地区由于气候干旱、降雨稀少等特点，导致风沙活动频繁，严重影响人们的生产活动，防护林依靠其防风阻沙等作用，已经成为该地区生态保护的最有效措施之一。农田防护林带主要功能在于降低风速改善局地小气候，然而不同配置结构的林带，对林带前后风速流场的影响也有所不同。因此研究当前不同配置农田防护林带的防风效益，对区域林带优化配置起到重要作用^[1]。目前，对于林带前后风速流场的研究主要集中在野外风速观测、风洞试验和风流场数值模拟这3方面^[2]。相比于传统方法，数值模拟因其成本低、预测精度高、结果可视化好等优势越来越受到学者们的关注。Kaiser^[3]早在20世纪60年代就将数值模拟技术应用于防护林带防风效应研究，之后随着在空气动力学研究上取得的突破，数值模拟在林带流场分析中的应用越来越广泛。林带防护效应的研究多关注于林带外流场风速变化，而无需考虑林带内部气流的运动情况，多将林带视为二维栅格结构，以简化数值模拟的计算量^[4]。Bitog等^[5]采用k-ε湍流模型对沿海防风栅栏周围风速进行模拟，探讨不同结构配置栅栏对气流的影响，并结合风洞试验验证二维模拟的准确性；Zhan等^[6]采用二维模型模拟固沙防风林周围不同位置的风速，结合林带前后高塔风速观测数据，更好地揭示防护林改变气流运动轨迹、降低风速的作用机理。笔者借鉴已有林带数值模拟的研究经验，在对大风绕流林带前后流场分布特征研究中，采用二维林带模型进行数值模拟。

当前针对乌兰布和沙漠绿洲区农田防护林带的研究中，多采用野外风速观测的方法对林带或整个林网进行风速流场分析^[7]，试验过程费时费力，且对不同配置林带前后的风速流场变化情况以及防护效应间的差异研究不足。笔者结合前人的研究成果，以乌兰布和沙漠绿洲区4种典型农田防护林带为研究对象，采用数值模拟方法模拟气流经过林带时前后流场的分布情况，对比分析2种不同宽度以及带后布灌对林带防风性能的影响，进一步揭示不同配置林带的防风机理，以期为该地区农田防护林带结构优化与配置提供参考。

研究区位于内蒙古磴口县中国林业科学院沙漠林业实验中心(简称沙林中心)实验林场，地处内蒙古自治区巴彦淖尔市西南部，属乌兰布和沙漠东北缘沙漠绿洲交错带(E 106°35′30″~106°59′48″，N 40°17′9″~40°29′15″)。属于温带大陆性季风气候，气候干燥，风大沙多，年均气温8.8 ℃，年均降水量143 mm，年平均蒸发量2 493 mm，年日照时间＞3 300 h，无霜期130 d。林带群落结构包括新疆杨(Populus alba var. pyramidalis)、箭杆杨(Populus nigra var.thevestina)、小美旱杨(Populus Simonii)、白榆(Ulmus pumila)等乔木，以及沙枣(Elaeagnus angustifolia)、梭梭(Haloxylon ammodendron)等灌木。

实验选取沙林中心实验林场内的4条不同配置主林带为研究对象，其中包含3条新疆杨纯林带和1条新疆杨与沙枣的混交林带。

在2021年7月对林带基本结构参数进行野外调查中，选的4条林带内分别设置3个长10 m，宽度与林带宽度相同的样方，采用测高仪分别测量每个样方内树木的树高与枝下高(灌木只测量冠层高度)，采用皮尺测量株行距、冠幅和林带宽度。林带疏透度使用照相法，在距离林带5倍树高处，利用数码相机拍摄林带迎风面照片，并在Photoshop软件中采用林木像素点所占比例的方法获取林带疏透度。根据新疆杨树木形态，笔者将新疆杨纯林带平均枝下高为分界线，分别计算上层冠层疏透度与下部疏透度。林带基本结构参数见表 1。

表 1(Tab. 1)

表 1 林带结构参数统计 Tab. 1 Statistics of the shelterbelts structure parameters 林带编号 Shelterbelt No. 配置 Configuration 林带行数 Rows 平均枝下高 Mean HCB/m 平均高度 Mean height/m 宽度 Width/ m 下部疏透度 Lower porosity 冠层疏透度 Canopy porosity Ⅰ 新疆杨 Populus alba var. pyramidalis 4 2 25 20 0.8 0.3 Ⅱ 新疆杨 Populus alba var. pyramidalis 2 2 25 5 0.8 0.3 Ⅲ 新疆杨/沙枣 Populus alba var. pyramidalis/Elaeagnus angustifolia 2/1 2/- 20/3 6/4 0.8/- 0.4/0.1 Ⅳ 新疆杨 Populus alba var. pyramidalis 2 2 20 6 0.8 0.4 注：疏透度为林带观测的平均值。Notes: The porosity is the average of the shelterbelt. HCB: Height of clear branch.

表 1 林带结构参数统计 Tab. 1 Statistics of the shelterbelts structure parameters

以林带高度H作为划分单元，流场内其他相应距离(如林带宽度)均以H表示，Ⅱ和Ⅲ林带H值设为25 m，Ⅳ的混交林带以及无灌木林带H值为20 m，详见图 1。由表 1可知，乔木林带冠层疏透度和下部疏透度有很大不同，因此将林带平均枝下高的上下2个部分作为单独的求解区域^[8]。采用CFD ICEM绘制非结构化的四面体网格，全局网格尺寸大小设置为0.5 m。求解域下壁面采用标准壁面函数处理，要求下壁面与沿壁内节点的量纲为1距离＞11.225，最终确定第1层网格高度为0.04 m时满足计算要求。

图 1(Fig. 1)

H：林带高。 H: Shelterbelt height. The same below. 图 1 几何模型示意图 Fig. 1 Schematic diagram of geometric model

气流经过林带所形成的流场属于湍流场，选择标准k-ε模型^[9]。而在林带区域可近似为多孔渗透结构，通过改变来流气体的动量参数模拟林带对气流的阻碍作用。在二维模拟的源项设置中，根据J.D. Wilson^[10]提出的动量汇模型可将流场内的动量方程表示为：

在x方向上：

在y方向上：

式中：u_i为速度矢量；x_i为尺度矢量；s_i为林带区域不同方向的动量源项；k_r为x方向上林带压力损失系数，是与疏透度、雷诺数和马赫数有关的函数，k_r=A(1-β²)/β²，其中β为实测林带疏透度，A为雷诺数，A=0.52；u为局地时均速度。

将进口边界条件设为速度入口，用来定义进口速度以及气流湍动能k和耗散率ε的初始值。本研究将近地面边界层内的风速廓线用对数分布来表示，即x方向的进口初始风速

式中：u₀为动力速度(摩擦速度)；y₀为零风面位移(地面粗糙度)，根据本研究样地情况，u₀=0.65 m/s，y₀=0.12 m，K=0.42^[11]。在基于压力的稳态条件下，选择SIMPLE数值求解算法，压力向采用标准差分格式，动量项、湍动能项和湍动能扩散项均采用一阶迎风格式。各求解量的收敛精度设为10^-3，计算迭代步数设为1 000步。

数值模拟结束后，将结果导出至后处理软件CFD-Post中，并使用“export”功能导出林带前后距地面0.1、0.3、0.5、0.7、1和1.2H高度处的水平相对风速值和流场内整体水平风速值，同时将带后风速降为初始风速70%以下的范围作为林带的有效防护距离^[12]。相关数据采用Excel软件进行整理，并使用Origin 2021软件绘图。

为验证数值模拟方法中求解域网格划分和参数设置的精准度，采用野外实测旷野风速与数值模拟风速进行对比。旷野风速观测数据采集于2021年6—8月，利用研究区荒漠-绿洲过渡带风沙监测塔，该塔主要监测绿洲外围环境中气流与沙尘情况，塔身分别在距地面高1、2、4、8、12、24、36和48 m处装有风速风向采集仪(AV-30WS)，能够以30 d的采样频率每1 min自动记录1组数据。将采集9个高度上的实测风速与数值模拟中相同高度处的旷野风速分别进行拟合，绘制出各自的风速廓线，并对两者数据差异化进行相关性分析。

如图 2所示，数值模拟风速曲线与野外实测风速曲线基本吻合，模型拟合度R²均＞0.99，风速拟合结果与罗凤敏等^[13]在该地区风速研究结果一致；此外对模拟风速与实测风速的回归分析中，二者拟合度R²=0.997，均方误差(MSE)为0.008 16，证明该数值模拟方法满足后续实验的基本条件。

图 2(Fig. 2)

图 2 实测与模拟风速对比 Fig. 2 Comparison between measured and simulated wind speed

如图 3所示，4条林带对气流分布有相似的影响作用，并在林带前后形成不同的速度分区^[14]。气流由入口向林带运动的过程中，因受林带影响在带前出现指数型减速区。并在气流穿过林带时在林带冠层产生分流，结合图 4典型林带局布风速流线图可见，一部分气流从林带上方越过时与上方气流汇集，出现风速抬升区；而穿过林带冠层的气流受林带内部的动能消耗，使风速大大减弱，并在带后形成弱风区；另有小部分气流从林带下部穿过形成气流加速区，其原因在于气流在经过林带下部时产生汇流，形成狭管气流，从而致使流速激增。图 3林带后方蓝色椭圆区域表明4条林带在带后均出现反向气流，出现不同程度的速度回流区，并均随高度升高回流逐渐减弱。之后随着气流逐渐远离林带，风速逐渐恢复到初始水平，形成风速恢复区。

图 3(Fig. 3)

图 3 林带前后水平风速流场分布 Fig. 3 Horizontal wind speed field distributions around the shelterbelt

图 4(Fig. 4)

中间矩形框图为林带在流场中的简化图，框图内横线表示林带平均枝下高位置。灰色箭头表示气流，箭头稠密代表气流的汇聚情况，箭头方向表示气流流向，箭头长短表示流速大小。 The rectangular frame in the middle is a simplified diagram of the shelterbelt in the wind speed field, and the horizontal line in the frame represents the mean HCB of the shelterbelt. The gray arrows represent the airflow, the density of the arrow represents the convergence of the airflow, the direction of the arrow represents the direction of the airflow, and the length of the arrow represents the wind speed. 图 4 林带周围风速流场 Fig. 4 Wind speed field around the shelterbelt

图 5为0.8H宽度林带(林带Ⅰ)和0.2H宽度林带(林带Ⅱ)距地面6种不同高度处的相对风速变化曲线。2条林带在带前风速的变化趋势基本相同，其差异主要发生在带后，且以近地面0.1H和0.3H高度处差异最为明显。气流在经过林带下部时由于疏透度较大导致风速增大，以较窄的0.2H宽度林带为例，在带后1H范围内增至初始风速的73.4%，而0.8H宽度林带气流在运动过程中风速继续衰减。0.8H宽度林带在6个高度上的回流程度都要强于0.2H宽度林带，且均伴随着高度的增加逐渐减弱；6个高度中在更靠近地面的0.1H高度处带后回流最为剧烈，0.8H宽度林带在带后7H处达到回流峰值，最大风速为初始风速的110.6%；0.2H宽度林带在带后10H处达到回流峰值，其值超过初始风速的69.5%。结合表 2林带有效防护距离差异，可以看出，在林带高度以下的范围内0.2H和0.8H宽度林带的有效防护距离分别为17.8H和15.5H。

图 5(Fig. 5)

z：距地面高度。下同。 z：Height from ground. The same below. 图 5 不同宽度林带各高度水平相对风速对比 Fig. 5 Comparison of horizontal relative wind speed at each height in different width of shelterbelt

表 2(Tab. 2)

表 2 不同配置结构林带在不同高度处的有效防护距离 Tab. 2 Effective shelter distances at different heights from shelterbelt with different structure  H 高度 Height 0.8H宽度林带 0.8H width 0.2H宽度林带 0.2H width 有灌木林带 With shrub 无灌木林带 No shrub 0.1 16.2 18.5 22.1 24.6 0.3 15.8 18.1 22.5 23.8 0.5 15.6 17.6 20.8 23.6 0.7 14.6 17.2 20.6 23.1 1 13.3 15.8 19.1 22.1 1.2 12.2 14.6 17.8 20.6

表 2 不同配置结构林带在不同高度处的有效防护距离 Tab. 2 Effective shelter distances at different heights from shelterbelt with different structure

笔者以带后布灌林带(林带Ⅲ)和无灌木林带(林带Ⅳ)为研究对象，分析在其他条件相同情况下带后灌木对风速的影响差异。如图 6所示，林带迎风侧气流变化受带后灌木带的影响较小，风速差异主要出现在带后近地面处。在0.1H高度处，有灌木林带在带后0.7H范围内风速降为0，而无灌木林带则出现明显的气流加速，并在带后1H处的风速达到最大，约为初始风速的70.5%；之后随着距地面高度的增加，气流受灌木带的影响逐渐减弱，风速变化趋势也趋于一致。在带后8H范围内，有灌木林带风速总体上都要低于无灌木林带，带后回流也要强于无灌木林带，且由于带后回流的位置靠前，因而导致风速恢复的更早。从表 2可见，在6个高度上的有效防护距离有灌木林带均小于无灌木的林带，在较为关注的林带高度以下区域，有灌木林带的防护距离约为21.5H，无灌木林带的防护距离约为23.7H。

图 6(Fig. 6)

图 6 带后有无灌木林带各高度水平相对风速对比 Fig. 6 Comparison of horizontal relative wind speed at each height in shelterbelt with and without shrub

笔者通过对林带风场数值模拟，绘制出林带前后水平风速分布云图，结果显示大风气流经过林带时会产生不同的速度分区，主要形成林带上方风速抬升区、林带下方气流加速区以及带后回流区^[11]。研究结果与前人林带流场研究结果类似。李雪琳等^[15]通过风洞试验对不同宽度配置结构林带迎风面、带中和背风面的风速进行测定，结果显示不同配置林带上方均出现高风速区，以及带后近地面处小面积加速区。杨文斌等^[16]在对灌木带的风洞试验中也得到相同的结果，在冠幅间的通道中会加快气流运动形成短暂加速区，而越过上方的气流会加速形成风速抬升区。杨学军等^[17]通过数值模拟对防护林带进行研究，模拟在固定宽度下林带对风场的影响，发现由于受到林带的阻碍作用，气流在翻过林带时会形成风速增加的区域，且在下部近地面处出现气流加速区。

带后回流也称带后湍流，是林带风速流场的研究重点，诸多研究表明其强度大小不仅改变带后的气流结构，也会影响着带后气流恢复^[18]。对于林带后回流的产生机理，Wilson^[19]首先从林带疏透度的角度进行解释，得出疏透度小的林带对气流的阻碍作用更大，加快带后的动能交换，从而导致回流。文中所选的4条林带疏透度分别为0.3、0.3、0.4、0.4，研究结果发现即使是疏透度相同的林带，带后回流程度仍存在不同，可见在林带其他条件相同的情况下林带宽度与带后布灌均会对回流产生造成影响，周军莉等^[20]通过对相同疏透度下不同宽高比的林带进行数值模拟，得出宽高比越大的林带，气流进入回流会加快，且气流恢复也越快，降低林地有效防护距离。此外根据以往研究，带后布灌的配置形式则在整体上减小了林带下部的疏透性，增大林带下方的阻力，从而增强带后湍流作用^[21]。

防护林因其防风作用，一般布置在受保护区边缘^[22]。笔者发现，0.8H宽度林带和带后布灌林带能有效削减带后近地面风速，原因在于2种配置均不同程度上增大了林带对来流风的阻力，这有助于强风条件下对带后农作物的防护。此外，王元等^[21]对灌草与林带搭配下的风场数值模拟研究中发现，带后灌草带的存在有效阻挡因前方窄林带下部的加速气流，增强林带后缘近地面防护效果。结合前文分析，宽林带和带后布灌也会增大带后回流强度，从而降低林带有效防护距离，可见不同林带配置结构在防风效应上各有利弊，因此笔者认为在农田防护林规划设置中，应综合考虑不同配置结构林带的防护特征进行合理布置。本研究实验结果揭示4种典型林带的防护特征，证明二维数值模拟在林带流场分析中的可靠性；然而笔者仅从林带光学疏透度出发，并未探讨真实情况下林带三维结构以及气流在林带内部的运动情况，这也将成为今后的研究重点。

采用数值模拟方法能够很好描述林带周围流场的分布特征，且较为详细地阐明流场不同区域的风速变化规律：1)不同配置结构林带水平风速总趋势为先减小后增大，并在带后流场中出现明显的速度分区，包括林带上方风速抬升区、带后弱风区、下部气流加速区、带后回流区和气流恢复区；2)相比于0.2H宽度窄林带与带后无灌木带的配置形式，较宽的0.8H宽度林带与带后布灌的乔灌混交林带对来流风速的削减更为显著，有效降低带后近地面处风速，但带后回流也更为剧烈，这导致气流恢复速率加快，降低林带有效防护距离。