基于离散元法(Discrete Element Method, DEM)的薄膜–土壤复合模型开发及验证

《Agriculture》：Development and Validation of a Film–Soil Composite Model Based on the Discrete Element Method

　　

残膜(Residual film)回收是缓解农业"白色污染(white pollution)"和保障土地可持续利用的关键途径。目前残膜回收机的研发主要依赖理论分析和田间性能试验，缺乏计算仿真模型支撑常导致机械性能欠佳，严重制约了回收设备的设计与优化。针对此，研究人员提出一种基于EDEM软件构建耕层残膜离散元(Discrete Element Method, DEM)模型并进行实验验证的方法。首先通过田间试验测定不同深度下的土壤紧实度和残膜分布，并测试残膜极限拉力为模型构建提供基础数据。采用EDEM中Hertz–Mindlin with bonding接触模型选取并优化残膜本征参数，结合Box–Behnken试验设计构建单位面积法向刚度(normal stiffness per unit area)、临界法向应力(critical normal stress)及粘结半径(bond radius)与薄膜极限拉力的二次回归模型。最优参数组合为：法向刚度=1.11×106N·m?3，临界法向应力=2.45×106Pa，粘结半径=0.03 mm。该参数下理论预测极限拉力为1.18 N，仿真值为1.16 N，相对误差仅1.69%，验证了单膜模型有效性。进一步基于田间实测数据，采用"降雨(rainfall)法"建立膜–土耦合模型并进行虚拟贯入试验，采用多元牛顿–拉夫逊迭代法(Multivariate Newton–Raphson iteration method)标定参数。土壤颗粒最优粘结参数为：单位面积法向刚度=9.6×105N/m2，单位面积切向刚度=9.6×105N/m2，临界法向应力=5.38×105Pa，临界切向应力=5.38×105Pa，粘结半径=4.3 mm。模拟平均贯入阻力为59.61 N，与田间实测值56.28 N相比相对误差为5.91%。结果表明所开发的膜–土耦合离散元模型可有效应用于耕层残膜回收机起膜与抛送过程仿真，为残膜回收机构设计与优化提供重要建模支撑。

论文解读：基于离散元法的薄膜–土壤复合模型开发及验证——以棉田耕层残膜为例

研究背景与意义

残留地膜(residual film / mulch film)长期滞留于耕层会导致土壤环境恶化，即农业"白色污染"。新疆棉田广泛采用地膜覆盖栽培，连续多年覆膜后耕层（0–200 mm）积累大量残膜，目前主流的滚筒式、齿梳式及铲式残膜回收机在实际作业中效果仍不理想。传统农机研发主要依赖理论力学分析与田间试制验证，难以精确获知作业过程中残膜与机具间复杂的接触力（大小、方向及作用点），且无法可视化膜–土–机具的动态交互过程，缺乏可靠的离散元(Discrete Element Method, DEM)仿真模型成为设计优化的瓶颈。已有少数研究尝试单独建立残膜DEM模型，但未实现与土壤颗粒的有效耦合。因此，研究人员以新疆棉田耕层残膜–土壤系统为对象，基于EDEM软件分别构建独立残膜DEM模型及膜–土耦合复合模型，并通过物理实验标定与验证，旨在为残膜回收装置的设计与结构优化提供微观尺度的数值仿真平台。该论文发表于《Agriculture》。

主要技术方法概述

研究人员于新疆巴州尉犁县典型棉田采用五点对角线取样法采集耕层0–200 mm分层（0–50、50–100、100–150、150–200 mm）土样，筛分清洗称重测定残膜含量；以万能试验机按GB/T 1040.3-2006标准对残膜条带进行拉伸试验获取极限拉力（均值1.18 N）；用数字式土壤紧实度仪测定耕层贯入阻力。残膜及土壤离散元模型均选用EDEM 2022内置Hertz–Mindlin with bonding接触模型（粘结接触模型，bonded contact model——颗粒间通过圆柱状粘结键连接以模拟材料内聚与抗拉/抗剪行为）。残膜模型以单层有序排列球状颗粒模拟，选取单位面积法向刚度、临界法向应力及粘结半径为因素，以极限拉力为评价指标，通过Box–Behnken设计(BBD)及二次回归拟合寻优，经拉伸仿真验证；膜–土耦合模型采用"降雨法"分层生成土床并将残膜模型随机嵌入耕层，以田间贯入阻力为参照，采用多元Newton–Raphson迭代法标定土壤颗粒粘结参数，通过虚拟锥尖贯入试验验证。

研究结果

2.1 残膜含量测定(Determination of Residual Film Content)

通过五点对角线法分层取样分析，耕层残膜含量在0–50、50–100、100–150、150–200 mm深度分别为6.57、12.37、11.63、5.30 g/m²，残膜主要富集于50–150 mm土层，占总量的67.17%。

2.2 耕层残膜物理性质测定(Determination of Physical Properties of Topsoil Residual Film)

测得残膜间动/静摩擦系数均值为0.422/0.438，残膜–土壤间为0.456/0.411；数显千分尺测得残膜厚度约0.01 mm，考虑田间褶皱取仿真等效厚度为0.02 mm（单层颗粒排列）。

2.3 残膜拉伸试验(Tensile Testing of Residual Film)

经100次有效拉伸试验，残膜极限拉力均值为1.18 N，以此作为DEM残膜粘结参数标定目标值。

2.4 耕层土壤贯入阻力测定(Determination of Soil Penetration Resistance in the Plow Layer)

0–100 mm土层平均紧实度0.256 kPa，100–200 mm增至0.856 kPa，呈随深度增加而增大趋势，取耕层0–200 mm平均贯入阻力56.28 N为耦合模型验证基准。

3.1 残膜建模与分析(Modeling and Analysis of the Residual Film)

以单位面积法向刚度(X₁)、临界法向应力(X₂)、粘结半径(X₃)三因素Box–Behnken试验建立二次回归模型（R²=0.9775，p<0.01，失拟项不显著），各因素及X₁X₃交互项、X₃²影响极显著(p<0.01)，影响程度排序X₃>X₃²>X₁X₃>X₁>X₂。优化得残膜最优粘结参数：单位面积法向刚度1.11×10⁶N·m^?3，临界法向应力2.45×10⁶Pa，粘结半径0.03 mm。理论预测极限拉力1.18 N，EDEM仿真均值1.16 N，相对误差1.69%；仿真能再现拉伸变形→启裂→完全断裂过程，力–时间曲线与实测高度吻合，证实残膜DEM模型可靠。

3.2 膜–土复合模型验证(Validation of the Soil–Film Composite Model)

设土粒法向与切向刚度相等、临界法向与切向应力相等，经多元Newton–Raphson迭代标定得土壤颗粒最优粘结参数：单位面积法向刚度=单位面积切向刚度=9.6×10⁵N/m²，临界法向应力=临界切向应力=5.38×10⁵Pa，粘结半径=4.3 mm。虚拟贯入试验平均贯入阻力59.61 N，与田间均值56.28 N相对误差5.91%，力–深度曲线趋势一致，膜–土耦合DEM模型满足工程仿真精度要求。

讨论与结论翻译总结

讨论指出：本研究首次构建了经田间数据标定、可反映真实耕层状态的残膜单独模型及膜–土耦合DEM模型，优于仅建单膜模型之既有研究；局限在于仅基于尉犁县单一土质校准，未引入膜–土颗粒界面粘结及非球形土团，且以极限拉力为唯一失效判据、土壤参数仅由静态贯入试验标定，后续将拓展至多区域土质、加入膜–土界面作用及动态作业工况验证。结论如下：

(1) 系统测定了耕层土壤紧实度、含水量、分层残膜含量及不同厚度残膜极限拉力，为膜–土耦合DEM模型建立奠定基础。

(2) 基于Hertz–Mindlin with bonding接触模型在EDEM中校准残膜DEM模型，通过Box–Behnken设计构建单位面积法向刚度、临界法向应力及粘结半径与极限拉力的二次回归模型，优化得最优参数组合：单位面积法向刚度1.11×10⁶N·m^?3、临界法向应力2.45×10⁶Pa、粘结半径0.03 mm，仿真极限拉力1.16 N，较理论预测误差1.69%，拉伸变形过程与实测一致，残膜模型可靠。

(3) 基于田间数据建立膜–土耦合模型并经贯入试验验证，土壤颗粒最优粘结参数为单位面积法向及切向刚度均为9.6×10⁵N/m²、临界法向及切向应力均为5.38×10⁵Pa、粘结半径4.3 mm；模拟平均贯入阻力59.61 N，与田间实测56.28 N相对误差5.91%，表明该膜–土复合模型可用于残膜回收机挖掘、抛送及输送过程仿真，为耕层残膜回收设备设计与优化提供稳健建模工具。