全膜双垄沟覆膜土壤离散元接触参数仿真标定

戴 飞 宋学锋 赵武云 张锋伟 马海军 马明义

(1.甘肃农业大学机电工程学院， 兰州 730070； 2.甘肃洮河拖拉机制造有限公司， 定西 730500)

0 引言

全膜双垄沟播技术是我国西北地区旱作农业的一项突破性技术，要求先完成大小垄体种床的耕整、覆膜与覆土准备，后在小垄垄沟内进行播种作业[1-2]。其中，全膜双垄沟起垄、覆膜-覆土环节农艺要求较为特殊和复杂，首先需构建大垄垄体与小垄垄体相间的异形垄体覆膜种床；其次以小垄垄体为中心基准，铺设宽度为1 200 mm的白(黑)色地膜(厚度为0.01 mm)，即用整张地膜覆盖小垄整体及其两侧垄沟，且分别均等向两侧各1/2大垄垄体延伸进行覆膜作业，并同时完成两侧膜边、垄沟内的覆土镇压。同时，当下一组小垄垄体覆膜作业时，其左侧1/2大垄垄体覆土膜边正好与上一组右侧大垄垄体覆土膜边对接，并在中心覆膜土带的连接下完成大垄垄体的全膜覆盖，整张地膜所覆盖双垄底部总体宽度约为1 100 mm，覆盖边界经历了由垄到沟、沿沟至垄的双重跨越[3]。因此，在上述农艺要求下，覆膜种床覆土环节至关重要，研究全膜双垄沟机械化覆土作业过程对于其种床高效构建具有重要意义。

近年来，离散元法(DEM)及其仿真软件EDEM在农业工程领域得到了广泛研究与应用，主要在大田耕作、谷物排种、农业物料参数标定、作物收获等方面取得了一系列重要研究进展[4-8]。种床覆膜土壤作为全膜双垄沟覆土装置的主要作用对象，是一种典型的离散物质，可借助离散元法对其接触参数及机械化覆土作业过程进行研究。其中，覆膜土壤颗粒的本征参数包括形状尺寸、密度、泊松比和剪切模量，该类参数可通过查阅文献或试验测得；土壤颗粒的接触参数主要包括恢复系数、静摩擦因数和滚动摩擦因数等，较难直接测量，可借助离散元法采用虚拟试验进行标定与优化筛选。目前，国内外相关学者应用离散元法，通过物料堆积休止角及滑动摩擦角完成了对土壤、稻谷、饲料、马铃薯、玉米等物料颗粒间的静摩擦因数、滚动摩擦因数及恢复系数等接触参数的标定[7,9-12]；同时，依托具体研制样机开展了各类农机具触土部件与大田土壤关于上述接触参数的标定[13-16]，但由于两者互作方式、土壤类型、触土部件材料等多因素的影响，研究标定的土壤重要接触参数仍处于一个区间值，且基于全膜双垄沟覆土装置作用下有关土壤接触参数的研究鲜见报道。

因此，为进一步提升基于离散元法关于全膜双垄沟覆土作业过程研究的准确性，在参考已有土壤-土壤、土壤-农具重要接触参数(静摩擦因数、滚动摩擦因数及恢复系数)区间值的基础上，分别应用三因素三水平正交回归模拟试验，优化全膜双垄沟机械化覆土装置作用下覆膜土壤的离散元模型参数，建立优化参数下的机械化覆土作业仿真模型并进行试验验证，以期为西北旱区农田耕作土壤离散元仿真参数设置提供参考与借鉴。

1 试验标定参数确定

1.1 作业装置结构组成

双垄耕作施肥喷药覆膜机主要由悬挂架、喷药装置、旋耕刀组、取土铲、镇压装置、挂膜装置、刮板式升运器、覆土侧流槽、覆土直流槽和肥料箱等组成，如图1所示。

图1 双垄耕作施肥喷药覆膜机结构图Fig.1 Structure diagrams of operation machine for tillage-fertilization and spraying-filming on double ridges1.喷药装置 2.悬挂架 3.肥料箱 4.挂膜装置 5.镇压装置 6.覆土侧流槽 7.覆土直流槽 8.刮板式升运器

1.2 标定参数选取

全膜双垄沟机械化覆土装置作业过程如图2所示，在作业机旋耕刀组及其取土铲的共同作用下将切削土壤高速抛送至刮板式升运器并后输至覆土装置。覆土装置将土壤分为两条路径，并在覆土侧流槽、覆土直流槽及后置镇压轮的共同配合下完成全膜双垄沟覆膜种床两侧膜边、垄沟内的覆土镇压[1]。

图2 全膜双垄沟机械化覆土作业过程Fig.2 Operation process of mechanized covering soil on whole plastic film mulching on double ridges

由机械化覆土装置作业过程可以看出，覆膜土壤间、土壤与覆土侧(直)流槽间主要呈密相堆积与倾斜相对滑动，并形成了覆膜土壤颗粒流，以4条覆土带完成对种床垄沟及膜边的定位覆盖。因此，应用离散元法对覆土过程进行仿真分析时，土壤与土壤静摩擦因数、土壤与土壤滚动摩擦因数、土壤与土壤恢复系数、土壤与钢(溜土槽)静摩擦因数、土壤与钢(溜土槽)滚动摩擦因数、土壤与钢(溜土槽)恢复系数均是影响其仿真结果准确性的关键参数，本研究将标定上述6个主要接触参数。

1.3 标定参数范围确定

根据已有文献研究，相关标定参数的参考范围为：土壤与土壤静摩擦因数在0.30～0.70之间[5,9,13,17-20]，土壤与土壤滚动摩擦因数在0.14～0.40之间[5,9,13-14,16-17,20]，土壤与土壤恢复系数在0.20～0.60之间[5,15,19,21]；土壤与钢静摩擦因数在0.30～0.60之间[5,13,21-22]，土壤与钢滚动摩擦因数在0.04～0.20之间[5,14,16,21,23]，土壤与钢恢复系数在0.28～0.60之间[15-16,20]。

2 试验材料与方法

2.1 试验方法与仪器

依托全膜双垄沟机械化覆土作业过程，分别通过土壤休止角及土壤与钢(溜土槽)滑动摩擦角的实测试验与仿真试验，标定土壤与土壤静摩擦因数、土壤与土壤滚动摩擦因数、土壤与土壤恢复系数、土壤与钢(溜土槽)静摩擦因数、土壤与钢(溜土槽)滚动摩擦因数及土壤与钢(溜土槽)恢复系数6个关键参数。

图3a所示为自制土壤休止角测定装置，主要由漏斗、支架、土壤堆积台和抽片等组成。其中，抽片和漏斗嘴上的插槽配合使用，装置支架可调节漏斗与土壤堆积台之间的垂直距离。如图3b所示，将田间耕作土壤装入土壤休止角测定装置，待土壤在漏斗内静止后瞬时拉出装置下端部抽片，并等土壤坡面稳定后，垂直拍照，导入CAD 软件中标注土壤休止角(Yα)。

图3 土壤休止角测定Fig.3 Determination of repose angle of soil-soil1.漏斗 2.抽片 3.土壤堆积台 4.支架

图4a所示为自制土壤-钢(溜土槽)滑动摩擦角测定装置，主要由拉伸装置、力学试验机、翻转合页、不锈钢溜土槽和绳索等组成。其中，翻转合页随着力学实验机的拉伸可绕轴向旋转90°，且为保证测定参数的准确性，装置选用钢材与实际溜土槽材料一致。如图4b所示，将田间耕作土壤装入土壤-钢(溜土槽)滑动摩擦角测定装置，通过设置万能力学试验机向上拉伸速度(20 mm/min)，在绳索牵引下将平行连接有溜土槽的平板匀速提起，并使其缓慢倾斜。待溜土槽内土壤开始持续滑动时，平板倾斜角度即为土壤与钢(溜土槽)的滑动摩擦角，待土壤颗粒持续滑动稳定后，垂直拍照，导入CAD 软件中标注土壤与钢(溜土槽)的滑动摩擦角(Yβ)。

图4 土壤与钢滑动摩擦角测定Fig.4 Determination of sliding friction angle of soil-steel1.拉伸装置 2.绳索 3.翻转合页 4.不锈钢溜土槽 5.力学试验机

2.2 试验设计

依据全膜双垄沟机械化覆土过程中土壤与土壤、土壤与钢(溜土槽)的互作特性及其运动规律，选取土壤与土壤静摩擦因数x1、土壤与土壤滚动摩擦因数x2和土壤与土壤恢复系数x3作为试验因素，以仿真土壤休止角Yα作为评价指标；选取土壤与钢(溜土槽)静摩擦因数x4、土壤与钢(溜土槽)滚动摩擦因数x5和土壤与钢(溜土槽)恢复系数x6作为试验因素，以仿真土壤与钢滑动摩擦角Yβ作为评价指标，分别进行三因素三水平正交回归模拟试验，选取的各试验因素编码如表1、2所示，各实施17组响应面分析试验，每组仿真试验重复进行3次，取3次测试结果的平均值作为试验结果，通过应用Design-Expert 8.06进行数据处理分析[24]。

表1 仿真土壤休止角试验因素编码Tab.1 Coding of simulation test factors of repose angle of soil-soil

表2 仿真土壤与钢滑动摩擦角试验因素编码Tab.2 Coding of simulation test factors of sliding friction angle of soil-steel

2.3 仿真模型建立

由于全膜双垄沟覆膜土壤是在作业机旋耕刀组的高速转动切削过程中产生，其土壤颗粒的基本结构主要以球形块状颗粒为主，可由EDEM软件自带的球形颗粒单元进行土壤模型建立。为保证仿真与实际土壤的一致性，设置EDEM球形填充单元的半径为3 mm，覆膜土壤颗粒及几何体的离散元模拟参数设置见表3[5,9,13-14,15-23]。

表3 覆膜土壤颗粒及几何体的离散元模拟参数Tab.3 DEM simulation parameters of covering soil on plastic film and geometry

采用EDEM 2018软件进行土壤休止角测定(图5)及土壤与钢滑动摩擦角测定(图6)的建模仿真，仿真模型选用Hertz-Mindlin (no slip)接触模型，rayleigh时间步长设置为20%，计算时间步长为0.01 s。待仿真结束后，将模型调整到正视视角，采用屏幕量角器Screen Protractor 4.0工具分别对土壤休止角、土壤与钢滑动摩擦角进行测量。

图5 土壤休止角仿真模拟试验Fig.5 Simulation test of repose angle of soil-soil

图6 土壤与钢滑动摩擦角仿真模拟试验Fig.6 Simulation test of sliding friction angle of soil-steel

3 试验结果与分析

3.1 实测试验结果

应用如图3、4所示的仪器与方法分别进行土壤休止角、土壤与钢滑动摩擦角的测定。每类试验重复5次，得到土壤休止角均值为31.2°，土壤与钢滑动摩擦角均值为22.6°。

3.2 仿真试验结果

3.2.1土壤休止角

以各影响因素编码值为自变量，以仿真结果测得的土壤休止角(Yα)为评价指标，构建不同试验组的几何体模型导入到EDEM中进行仿真试验，试验设计及结果如表4所示，X1、X2、X3为各变量编码值。

表4 土壤休止角仿真试验设计与结果Tab.4 Simulation test design and result of repose angle of soil-soil

借助Design-Expert 8.06软件对试验结果进行分析，获得编码值表示的土壤休止角Yα二次回归模型为

(1)

对上述二次回归模型进行方差分析和回归系数显著性检验，结果如表5所示。

表5 土壤休止角回归方程方差分析Tab.5 Variance analysis of soil-soil repose angle regression equation

注：*表示差异显著(P<0.05)，** 表示差异极显著(P<0.01)。下同。

3.2.2土壤与钢滑动摩擦角

以各影响因素编码值为自变量，以仿真结果测得的土壤与钢滑动摩擦角(Yβ)为评价指标，构建不同试验组的几何体模型导入到EDEM中进行仿真试验，试验设计及结果如表6所示，X4、X5、X6为各变量编码值。

借助Design-Expert 8.06软件对试验结果进行分析，获得编码值表示的土壤与钢滑动摩擦角Yβ二次回归模型为

表6 土壤与钢滑动摩擦角仿真试验设计与结果Tab.6 Simulation test design and result of sliding friction angle of soil-steel

(2)

对上述二次回归模型进行方差分析和回归系数显著性检验，结果如表7所示。

表7 土壤与钢滑动摩擦角回归方程方差分析Tab.7 Variance analysis of soil-steel sliding friction angle regression equation

3.3 最优标定参数确定

应用Design-Expert 8.06软件Optimization-Numerical模块，将3.1节实测试验目标参数土壤休止角Yα=31.2°、土壤与钢滑动摩擦角Yβ=22.6°代入回归方程(1)、(2)，分别对回归模型进行所对应目标的优化求解，获得土壤休止角最优标定参数组合为：x1=0.68，x2=0.27，x3=0.21；土壤与钢滑动摩擦角最优标定参数组合为：x4=0.31，x5=0.13，x6=0.54。因此，全膜双垄沟机械化覆土仿真模型中土壤离散元所需参数如表8所示。

表8 覆土作业过程仿真离散元模拟参数Tab.8 Discrete element simulation parameters of covering soil process

图7 全膜双垄沟覆土过程仿真模拟结果Fig.7 Simulation operation process of covering soil on whole plastic film mulching on double ridges

4 试验验证

4.1 休止角、滑动摩擦角验证

为验证覆膜土壤颗粒及几何体最优参数组合的准确性，采用表3、8所示参数值进行离散元仿真试验，模拟土壤休止角及土壤与钢滑动摩擦角堆积情况。每类试验重复5次，得到离散元模拟参数标定后土壤休止角、土壤与钢滑动摩擦角均值为30.4°和21.9°，两者相对误差分别为2.6%和3.1%，表明经过模型参数标定与优化，全膜双垄沟覆膜土壤离散元模型较为准确。

4.2 覆土作业过程仿真与验证

为进一步验证全膜双垄沟覆土装置采用覆膜土壤颗粒优化标定接触参数后作业模型的准确性，采用离散单元法对装置覆土作业过程进行数值模拟。其中，覆土土壤颗粒选取球体颗粒建模，其直径设定为3 mm，土壤与土壤、土壤与钢(溜土槽)接触模型选择Hertz-Mindlin(no-slip)，仿真试验参数设置如表3、8所示。仿真时间步长1.405×10-5s，是瑞利时间步的40%，仿真共进行2 s。根据全膜双垄沟覆土装置结构设计，覆土侧流槽角为75°，覆土直流槽角为48°，刮板升运带式提土器单侧刮板为10个，两刮板间距离设置为130 mm，运动通过EDEM中coupling server面板，由动力学耦合来控制。依据覆土装置作业参数优化值，其前进速度为0.70 m/s，刮板升运带式提土器线速度设置为0.67 m/s[1]。颗粒工厂是160 mm×180 mm的长方形平面，单侧颗粒工厂每秒生成土壤颗粒数为1.447×107个。

图7所示为全膜双垄沟覆土装置在覆膜土壤颗粒及几何体最优标定参数组合设置下的种床覆土仿真过程。

在全膜双垄沟覆土模拟仿真前0.15 s内为土壤颗粒生成与刮板升运带式提土器填充时段(图7a)。当仿真时间为0.44 s 时(图7b)，刮板升运带上方土壤已经足量填充并开始逐渐进入覆土侧、直流槽；0.73 s时(图7c)覆土装置各溜土槽处均出现少量、不连续的土壤颗粒流；0.95～1.56 s时(图7d、7e、7f)覆土侧、直流槽内土壤逐步形成连贯的覆土颗粒流，且土壤颗粒流量逐渐增大；1.73～2.00 s(图7g、7h)全膜双垄沟覆土装置趋于稳定的覆土输送状态，所形成4条覆土带的土壤颗粒流厚度相对均匀一致，达到机械化覆土的作业要求。全膜双垄沟覆土过程仿真模拟(图8a)与田间覆土验证试验(图8b)基本一致，表明相关接触参数的优化标定可靠有效。

图8 覆土作业效果验证Fig.8 Verification test of covering soil effort1.覆土侧流槽 2.种床垄沟覆土带 3.种床小垄 4.种床大垄覆土带 5.种床大垄 6.覆土直流槽

5 结论

(1)结合EDEM软件进行土壤休止角测定及土壤与钢滑动摩擦角测定的建模仿真。通过三因素三水平正交组合试验，得出各因素对土壤休止角的影

响主次顺序依次为：土壤与土壤静摩擦因数、土壤与土壤滚动摩擦因数和土壤与土壤恢复系数；得出各因素对土壤与钢滑动摩擦角的影响主次顺序依次为：土壤与钢静摩擦因数、土壤与钢滚动摩擦因数和土壤与钢恢复系数。

(2)分别建立了各关键接触参数与土壤休止角、土壤与钢滑动摩擦角的二次多项式回归模型；以自制土壤休止角、土壤与钢滑动摩擦角试验装置测定结果作为优化的目标值，获得全膜双垄沟覆膜土壤离散元最优接触参数组合为：土壤与土壤静摩擦因数0.68、土壤与土壤滚动摩擦因数0.27、土壤与土壤恢复系数0.21、土壤与钢静摩擦因数0.31、土壤与钢滚动摩擦因数0.13和土壤与钢恢复系数0.54。

(3)为验证所标定全膜双垄沟覆膜土壤接触参数的可靠性，对模拟仿真与实际试验的土壤休止角、土壤与钢滑动摩擦角进行了对比，两者相对误差分别为2.6%和3.1%；应用离散单元法进行全膜双垄沟覆土装置在覆膜土壤颗粒最优标定参数组合设置下的种床覆土过程仿真模拟，通过与实际作业效果对比，验证了仿真试验与回归模型的有效性。