基于水热效应数值模拟的旱地玉米全膜双垄沟农艺优化

魏万成 戴飞 张锋伟 张仕林 史瑞杰 刘元祥 沈鹏云

(甘肃农业大学 机电工程学院，兰州 730070)

全膜双垄沟播技术是我国西北旱作农业的一项创新突破性技术，该技术采用地膜全覆盖模式，不但可以显著提高农作物产量，而且具有保温、保湿、抑制杂草生长，防止病虫害发生的作用[1]。该技术对玉米增产效应的研究主要集中在热效应与水分运移方面。其中，垄体结构是影响玉米全膜双垄沟保温性和保水性的关键因素之一，因此，研究垄体结构对其水热效应的影响程度十分重要。

目前，对于玉米全膜双垄沟水热效应的研究大多以大田试验为主，王红丽等[2-3]通过大田试验研究了全膜双垄沟水热效应及其对产量的影响，通过垄体内温度和含水率对比得出全膜双垄沟技术在保温保水性能上有明显的优势。但大田试验有试验周期长、环境干扰大、成本高等缺陷，通过数值模拟可以弥补上述不足。近年来，有限元技术在土壤水热效应数值模拟领域已经取得了一系列成果，通过数值模拟研究土壤热平衡和降雨入渗的方式已经被广泛认可。其中ABAQUS软件在土壤研究方面提供了较为完善的土体本构模型，在其Soils分析步中可以解决土壤热平衡和渗流问题而被广泛使用[4-6]。依据现有的旱地全膜双垄沟种床构建农艺要求，起垄后的垄体结构仅对垄体宽度和高度2个参数进行了限定。对于起垄后垄体和垄沟的形状大小、坡度、沟宽比等参数并没有给出具体要求[7]，因此，有必要进一步优化全膜双垄沟种床构建农艺要求。

本研究拟采用数值模拟的方法，借助ABAQUS软件建立4种垄体模型，并对其热平衡和雨水入渗过程进行数值模拟，以寻求同时具备良好保温性和保水性的垄体结构。将现有的垄体结构和配套的起垄铺膜以及播种机具机进一步优化，以期为玉米全膜双垄沟种床构建农艺要求及其配套起垄覆膜机具研发标准化提供借鉴和参考。

1 模型与方法

1.1 垄体模型

参考现有的玉米全膜双垄沟种床构建农艺要求和相关文献研究，目前国内垄体的构建主要以梯形和拱形2种形状为主[7-10]。因此，在ABAQUS软件中分别建立4种全膜双垄沟垄体模型(图1)。

T1模型的垄体和垄沟均为梯形结构，大垄宽70 cm，垄高15 cm，小垄宽40 cm，垄高20 cm，垄沟底宽8～10 cm；

T2模型的垄体为梯形结构，大垄宽70 cm，垄高15 cm，小垄宽40 cm，垄高20 cm，垄沟为拱形结构，拱腹弧半径为10 cm，沟宽8～10 cm；

T3模型的垄体和垄沟均为拱形结构，大垄拱顶弧半径为48 cm，垄宽70 cm，垄高15 cm，小垄拱顶弧半径为22 cm，垄宽40 cm，垄高20 cm，垄沟拱腹弧半径为10 cm，沟宽8～10 cm；

T4模型的垄体为拱型结构，大垄拱顶弧半径为48 cm，垄宽70 cm，垄高15 cm，小垄拱顶弧半径为22 cm，垄宽40 cm，垄高20 cm，垄沟为梯形结构，沟底宽8～10 cm。

4种垄体模型土层深度值均为40 cm。

1.地膜；2.大垄覆土带；3.垄体 1.Plastic-film; 2.Large ridge surface covering the soil; 3.Body of ridgeH为小垄高度，cm；h为大垄高度，cm；B、b分别为大、小垄体宽度，cm；c为垄沟宽度，cm；R为垄沟拱腹弧半径，cmHis a large ridge height, cm;his a small ridge height, cm;Bandbare large and small ridge widths, cm;cis furrow widths;Ris furrow radius, cm.
图1 不同结构垄体模型
Fig.1 Model of different structure ridge body

1.2 计算模型与方法

土壤是一种典型的多孔介质，按照传热学理论，分析土壤时可以忽略空气的传热作用，导热规律基本符合常规固体内部导热规律[11]。在ABAQUS软件中只要建立准确的模型、添加与热性质相关的材料参数、选择合适的分析步和正确的温度边界条件，就能得出一定时间内模型内部热平衡规律。研究时可以将其抽象成半无限大平板的一维温度响应问题[12]。考虑土壤的水分、导热系数、比热和热物性变化时，对应的一维热扩散基本方程为：

(1)

式中：C为比热容，J/(m3·K)；k1为传热系数；t为时间，h；x为深度，cm。当不考虑热物性的变化时式(1)简化为：

(2)

对于单纯的雨水入渗模拟可忽略土壤中的气体和根系吸水对土壤水分的影响，根据渗流达西定律和能量守恒，用非饱和一维Richards渗流控制方程来描述垄体内土壤水分运动[12-15]：

(3)

式中：θ为单位体积土壤含水量，cm3/cm3；z为土壤深度，cm；K为土壤导水率，cm3/h；h为土壤压力水头，cm。

由于一定范围内农田土壤的属性基本一致，因此，为了减小计算成本和时间将3维模型转换为2维模型并给模型添加准确的材料参数。仿真时需要给模型添加的土壤材料参数见表1[15-17]。

表1 土壤材料参数Table 1 Parameters of soil

2 模拟过程与结果分析

2.1 雨水入渗模拟

2.1.1前处理及其边界条件

根据模拟要求建立垄体模型，并添加材料参数。由于土壤是典型的岩石分化颗粒在自然条件下堆积的产物，且在ABAQUS中Mohr-Coulomb土体本构模型主要适用于颗粒状材料，因此，在材料属性中选择Mohr-Coulomb土体本构模型[11,18]。添加土壤渗透系数、初始孔压、孔隙比、水的容重等材料参数，并指派截面给模型。划分网格类型为六面体结构化网格，单元类型为CPE4P，采用soils分析步。

在ABAQUS软件中降雨的施加是通过入渗强度即单位流通量来表示，即在相应的模型表面定义表面流体并制定入渗强度的大小和方向，同时排除了降雨造成的积水现象。依据全膜双垄沟膜面产流，沟垄集流的农艺效应，模型中设置雨水仅能通过垄沟膜孔下渗到深层土壤，垄体上表面为不透水地膜层。根据相关文献和案例，为了达到较好的仿真效果节约计算成本，在ABAQUS渗流仿真时入渗强度一般为恒定值，此处假设垄沟膜孔的入渗强度为0.02 m/h，降雨时间设为72 h[19]。

2.1.2垄体内渗流速度分布

模拟计算结束后得到各模型内雨水渗流速度的分布状态(图2)。在相同的降雨边界条件下，由于垄沟的集雨效应，垄体内土壤渗流速度大小和方向均发生了改变，垄沟附近土壤渗流速度明显大于其他位置，雨水在经过垄沟时产生了明显的侧渗现象，当下渗到一定深度时侧渗现象随之减弱。4种垄体模型中最大入渗速度均位于垄沟与垄体的交界处，造成这一现象的主要原因是降雨集中到垄沟内由于垄沟结构改变了雨水入渗的初速度使其优先产生侧渗。这一结果充分说明全膜双垄沟播技术能将有限的自然降水入渗到玉米根系附近，增大了根系附近土壤的湿润面积，起到了良好的集雨效果。

蓝色箭头代表降雨后土壤空隙中渗流速度的矢量分布。
The blue arrow represents the vector distribution of the effective velocity of rainwater in the soil voids after rainfall.
图2 不同垄体模型雨水渗流速度分布
Fig.2 Effective velocity distribution of rainwater seepage under different models

2.1.3垄沟内雨水渗流速度变化

为了进一步探明4种模型垄沟内雨水渗流速度的变化规律，在后处理中得到渗流速度在土层深度方向的变化和分布状态(图3)。垄体结构对垄沟内雨水渗流速度产生了一定影响，表层渗流速度明显大于深层，在20 cm的浅层土壤内各模型垄沟内渗流速度差异最大。T1模型垄沟内渗流速度从0.1 m/h迅速下降到0.01 m/h。在 3 cm 土层内就接近最小值，超过3 cm后渗流速度始终小于其他3种模型。T2模型垄沟内雨水渗流速度从0.08 m/h开始缓慢下降，超过3 cm后变化趋势基本和T4模型基本一致。在相同的降雨边界条件下，T2模型垄沟内平均渗流速度为0.23 m/s，均大于其他3种模型(图3(b))，说明T2模型的垄体结构在保证侧渗的同时还具有较强的下渗能力。

图3 不同模型垄沟内雨水渗流速度分布
Fig.3 Distribution of rainwater seepage velocity in different ridge models

2.1.4垄体含水率分布

含水率是土壤持水能力的关键指标，也是土壤水热效应研究中重要的观测对象[2]。ABAQUS软件渗流模拟功能可以通过后处理得到降雨后流体体积占土体总体积之比(图4)又称体积含水率[20]。在降雨结束后湿润土壤和干土之间交界面较为明显。垄沟附近土壤体积含水率明显高于深层土壤，水分富集区主要集中在垄沟附近，随着土层深度的增加体积含水率逐渐下降，且小垄垄体内平均含水率高于大垄。

图4 不同模型垄体内体积含水率分布
Fig.4 Distribution of volumetric moisture content under different models

2.1.5垄沟含水率变化

根据目前有关农田土壤含水率的研究和土壤学理论，含水率一般都是通过烘干法测定，该方法得到的含水率为质量含水率[21]。因此，为了使模拟结果更具有科学性将体积含水率通过式(4)转换为质量含水率：

(4)

式中：ω为质量含水率；ρd为土壤干密度，g/cm3；ρw为水密度，g/cm3；θw为体积含水率。

垄沟内含水率变化和平均含水率分布见图5。随着土层深度的增加垄沟内含水率逐渐下降，在浅层土壤内下降速率略有不同。根据相关研究，全膜双垄沟播技术能提高15～20 cm内土壤含水率，可有效缓解早春干旱对水分的需求[3]。分析垄沟内平均含水率(图5(b))可知；降雨结束时T2模型内土壤平均含水率最高(含水率为18%)，T3模型最低(15%)，说明T2模型的垄体结构相对于其他3种模型在保水性能上有明显的优势。

图5 不同模型垄沟内含水率分布
Fig.5 Distribution of moisture content in different ridge models

2.2 热平衡模拟

2.2.1前处理及其边界条件

温度值直接反映了土壤热量的多少，根据玉米全膜全膜双垄沟播技术适用地区土壤特性和相关研究，本研究默认土壤空隙结构和温度不受含水率变化影响，土壤环境为理想状态，同时不考虑地热能对热平衡的影响，土壤热源仅为太阳热辐射[22]，且下层土壤的热量主要来源于表层热交换。因此，分别对各垄体模型添加表1中的材料热参数，划分四面体结构网格，单元类型为DC2D4单元，选择Heat transfer热传递分析步[11]，假设辐射时间为6 h，初始温度值为0 ℃。

本研究仅探讨垄体结构对土壤热量传递的影响，因此需要统一的温度边界条件。以玉米出苗期土壤温度为主要依据，根据已有的大田试验研究，出苗期内垄体表层温度与垄沟表层温度差为3.0～4.0 ℃，与覆土带温差为1.4～2.0 ℃。以此为依据分别对4种模型设置相同的温度边界条件[16,21]：垄沟表层温度为24.8 ℃，垄体地膜层温度为28.0 ℃，大垄垄体覆土带表层温度为26.0 ℃。

2.2.2垄体土壤温度分布

在后处理中得到4种垄体模型土壤温度分布状态(图6)，相同的温度边界条件下各模型垄体内热量分布截然不同，浅层土壤内温度的分布和垄体模型结构存在一定的关联，当深度达到一定值时温度分布差异越来越小。主要原因是随着土层深度的增加热辐射逐渐减弱，土壤能接收到的表层热量值也逐渐下降。说明地膜覆盖主要改变了表层土壤的传递能力和热能分布，因此有必要对垄沟内土壤温度进行进一步研究。

图6 不同模型垄体温度分布
Fig.6 Distribution of temperature under different models

2.2.3垄沟内温度变化

全膜双垄沟播技术要求先起垄覆膜，后在垄沟内进行穴播作业，播种深度为3～6 cm。提取垄沟内温度数据得到垄沟内温度变化和平均温度分布(图7)。随着垄沟内土层深度的增加温度值逐渐减小，耕层越深、温度值差异越小，当深度达到40 cm后各模型垄沟内温度均接近21±0.5 ℃。造成这一现象的主要原因是不同垄体模型结构土壤表层热辐射面积不相等，因而土壤表层所接收到的热量也不相等。随着土壤深度的增加土壤接收表层热辐射能力减弱，差异也逐渐消失。

根据不同模型结构下垄沟内土层平均温度(图7(b))对比情况可以得出：在相同的温度边界下T4模型垄沟内平均温度最高(温度为23.5 ℃)，T2模型最低(22.8 ℃)，说明T4模型的垄体结构相对于其他3种模型在保温性能上有明显的优势。

图7 不同模型垄沟内温度分布
Fig.7 The distribution of temperature under different ridge models

2.2.4垄体土壤热通量分布

土壤热通量直接影响着土壤温度变化的速度和时间，分析土壤热通量可以进一步了解土壤热量的收支和传递状况。因此，热通量也是土壤热平衡研究中重要的观测指标[23]。同样的方式通过后处理功能得到各模型土壤热通量分布状态(图8)。4种垄体模型中热通量值均为正值，说明在40 cm土层内温度始终在向下传递。土壤表层热通量大于底层，说明表层土壤单位时间内经过单位横截面积上的热量大于深层土壤。

图8 不同垄体模型热通量分布
Fig.8 The distribution of heat flux under different models

2.2.5垄沟内热通量变化

提取各模型垄沟内热通量数据，得到垄沟内土壤热通量变化和分布状态(图9)。垄体结构对垄沟内热通量变化影响较大，在有限的深度内(<40 cm)，4种垄体模型垄沟内的热通量没有出现负值，说明热量始终在向下传递只是传递的程度有所差异。其中T3和T4模型在浅层土壤的热通量值较大且下降速率较为明显，分别从119和87 W/m2开始急剧下降，当深度超过10 cm时趋于稳定，此后T4模型热通量值始终大于T3模型。而对于T1模型和T2模型来说表层土壤热通量值都较小，随着土层深度的增加热通量下降趋势较为缓慢，分别从45和68 W/m2开始缓慢下降，直到超过深度30 cm才趋于稳定。T4模型垄沟内热通量平均值最大(热通量为30 W/m2)，T1模型最小(19 W/m2)。显然T4模型有利于大量的热量被传递到深层土壤，防止热量散失。

图9 垄沟内热通量分布
Fig.9 The distribution of heat flux in different ridge models

3 优化方案

根据雨水入渗和热平衡数值模拟分析,全膜双垄沟垄体结构对垄体内土壤温度和含水率产生了影响。以T4模型结构为基础，对玉米全膜双垄沟垄体结构、配套起垄铺膜机具、配套播种机具提出了优化方案，具体如下：

1)垄体结构优化。以现有的玉米全膜双垄沟种床构建农艺要求为基础，结合数值模拟结果以及目前垄沟集雨种植模式中垄沟大小和沟宽比的确定方法[24-25]。以T4模型为基础构建垄体，起垄后的垄体呈拱型结构，大垄拱顶弧半径为48 cm，垄宽70 cm，垄高15 cm；小垄拱顶弧半径为22 cm，垄宽40 cm，垄高20 cm；垄沟为梯形结构，沟底宽8～10 cm，沟底与大垄宽度比约1∶7，与小垄宽度比约1∶4。起垄后要求垄体结构规整、表面光洁、膜土耦合效果良好。

2)配套起垄铺膜机具的优化。在传统的机械化起垄铺膜作业基础上提出一种垄体精整化作业模式，既在起垄铺膜机上安装振动式垄体整形装置，在其振动挤压的作用下将垄体土坯固定成形并使表面光洁平整。具体实施方式为；提取T4模型结构数据，根据仿生学原理设计出一种仿形整形器，其形状结构与T4模型一致。起垄作业时，通过整形器的持续性击打和夯实将垄体土坯定形。

3)配套播种机具的优化。在播种作业时，对于目前主流的自走式玉米全膜双垄沟直插式穴播机而言，行走轮宽度应该小于垄沟宽度(8～10 cm)，以防止播种作业时挤压和破坏垄体结构。优化鸭嘴结构使其在播种后膜孔半径达到4～5 cm，以便降雨能及时通过垄沟入渗防止径流。

4 结 论

本研究通过数值模拟的方式研究了不同垄体结构对垄体内含水量及温度的影响，得出如下结论：

1)全膜双垄沟垄体结构直接影响了垄体内土壤的水分运移状况和含水率分布，通过数值模拟可以直观的观察到雨水的侧渗现象，4种模型中T2模型的测渗能力最强。T2模型垄沟内平均含水率达到18.0%，相对与其他模型含水率增加了1.1%～3.0%，其次T4模型垄沟内含水率也达到了16.9%。说明这2种模型的垄体结构更有利于降雨在垄沟内及时入渗，保水效果明显。

2)全膜双垄沟垄体结构对垄体内土壤热平衡过程也产生了影响，主要体现在温度和热通量方面。40 cm土层内各模型垄体内土壤热通量都为正值，说明在相同的温度边界条件下，热量的传递方向是一致的。垄体结构仅改变了土壤热量的传递程度，其中T4模型垄沟内平均热通量为30 W/m2高于其他3种模型，因此T4模型垄沟内平均温度也最高(温度值为23.5 ℃)，相对于其他3种模型T4模型垄沟内平均温度提高了0.3～0.5 ℃。说明T4模型的垄体结构能将更多热量储存在垄沟内，显著提高了保温性能。

3)T4模型的垄体结构能同时提升垄沟内土壤温度和含水率，在全膜双垄沟种床构建农艺优化过程中应该以T4模型为基础构建统一的垄体结构，并对其配套的起垄铺膜机具和播种机具进行优化。