基于离散元法的椰糠颗粒建模及仿真参数标定

王刚 王倩 薛忠 郭昌进 宋刚 王槊

摘  要：為确定椰糠在离散元模拟过程中合理的仿真参数，本文基于离散元法Hertz-Mindlin（no slip）接触模型在仿真软件EDEM 2018中建立椰糠颗粒模型并生成颗粒工厂，通过测量比对物理试验和仿真试验椰糠堆积角的方法对其仿真参数进行标定。首先通过物理试验测得椰糠堆积角、堆积密度等本征参数，利用EDEM 2018内嵌GEMM数据库以及相关文献，综合分析得到椰糠待标定因素的高低水平;然后进行Plackett-Burman试验设计筛选出对椰糠堆积角影响显著的因素依次为：椰糠-椰糠静摩擦系数、椰糠-椰糠滚动摩擦系数和椰糠剪切模量;再进行Box-Behnken试验设计建立堆积角与3个显著性影响因素的回归模型，运用Design expert软件优化功能，以堆积角物理试验45.69°为目标，对回归方程寻优求解，得到显著性影响因素最佳组合：椰糠剪切模量为1.44 MPa，椰糠-椰糠静摩擦系数为1.1 椰糠-椰糠滚动摩擦系数为0.15;最后通过试验验证表明此组合参数可用于椰糠物料的离散元仿真，能为椰糠输送、混合等机械装备的设计研发提供理论参考。

关键词：椰糠;离散元;参数标定;堆积角;试验设计

中图分类号：S22      文献标识码：A

Modeling and Simulation Parameters Calibration of Coconut Coir Particles Based on DEM

WANG Gang1， WANG Qian1， XUE Zhong2， GUO Changjin2， SONG Gang1， WANG Shuo1

1. Tropical Agricultural Machinery Research Institute， Chinese Academy of Tropical Agriculture Sciences， Zhanjiang， Guangdong 524091， China; 2. Key Laboratory of Agricultural Equipment for Tropical Crops， Ministry of Agriculture and Rural Affairs， Zhanjiang， Guangdong 524091， China

Abstract： In order to determine the reasonable simulation parameters of coconut coir in the DEM process， based on the DEM Mindlin （no slip） contact model， a coconut coir model is established in the simulation software EDEM 2018 and the particle factory is generated. The simulation parameters were calibrated by comparing the physical test and virtual test of coconut coir repose angle. Firstly， the intrinsic parameters of coconut coir such as repose angle and bulk density were measured by physical tests， and the level of factors to be calibrated were obtained by comprehensive analysis using the GEMM database embedded in EDEM 2018 and relevant literature. A Plackett–Burman test design was carried out to screen out the factors that have significant influence on the repose angle of coconut coir： static friction coefficient of coconut coir with coconut coir， rolling friction coefficient of coconut coir with coconut coir， shear modulus of coconut coir. Through Box–Behnken experimental design， a regression model of repose angle and three significant influencing factors were established， using the optimization function of design expert software， taking the physical test of stacking angle of 45.69° as the objective， the regression equation was optimized and solved， and the best combination of significant influencing factors were obtained. Coconut coir shear modulus was 1.44 MPa， static friction coefficient of coconut coir with coconut coir was 1.11， and rolling friction coefficient of coconut coir with coconut coir was 0.15. The experimental results showed that the combined parameters could be used in the discrete element simulation of coconut coir materials， and provide theoretical references for the design of coconut coir conveying， mixing and other mechanical equipment.

Keywords： coconut coir; discrete element method; parameters calibration; repose angle; experimental design

DOI： 10.3969/j.issn.1000-2561.2021.09.036

椰糠是椰子加工过程中的副产物，是从椰子外壳纤维脱落下来的一种天然有机质媒介，椰糠因其保温、保湿、疏松、透气的优良理化性质，在设施农业中正逐步取代草炭成为应用最为广泛的栽培基质[1-4]。研究表明，单一依靠椰糠作为栽培基质是存在缺陷的，通常需要与其他物料如有机肥等按照一定比例搭配形成混合基质，才能保证栽培基质的碳氮比、养分含量、保水通气性等性状更适应作物生长[5-6]。

混合机械设备是椰糠栽培基质机械化生产的关键设备之一，其生产过程属于固体混合。因为固体混合区别于液体搅拌，不具有扩散特性，须施加外力才能强制流动，固体颗粒与混合设备（或其相关接触部件）的接触关系、颗粒运动特性以及微观作用机理等直接关系到混合设备的作业性能和工作效率，所以不同的固体混合设备没有统一的制造标准，往往根据不同的功能和使用要求单独设计，而设计专用设备，就要充分了解物料特性[7-8]。

目前，实际应用的椰糠混合基质生产设备大多为直接套用市面普通的饲料混合设备或粉粒体混合设备，生产过程中存在混合不均匀或功率不符等问题。现阶段仍缺乏对椰糠基质混合过程及混合性能的系统研究，对椰糠混合基质的研究也主要集中在基质不同配比对农作物生长特性的影响规律等方面，椰糠作为一种典型散体物料，科研院所对其开展物料离散元特性方面的研究较少，而近年来，针对土壤、颗粒肥料、农业物料等离散性物料特性方面的研究相对较多，并广泛运用了离散元法（discrete element method，简称DEM）进行分析，可对椰糠物料离散元特性的分析提供重要的借鉴与参考。在土壤方面，邢洁洁等[9]通过试验并利用离散元建模仿真确定砖红壤颗粒特性最佳组合，为研发适用于海南热区砖红壤地高性能耕作装备关键触土部件的设计优化提供支撑;Ucgul等[10]利用离散元软件中Hertze- Mindlin及Hysteretic Spring接触模型，研究分析了有黏结力和无黏结力的情况下土壤颗粒的塑性形变等问题。在颗粒肥料方面，罗帅等[11]使用离散元参数标定方法筛选蚯蚓粪基质颗粒参数，为研究蚯蚓堆肥处理及蚓粪利用相关设备奠定基础;Bangura等[12]采用离散元法对肥料进行建模通过分析排肥过程比较了螺旋槽轮和直槽轮的排肥速度与排肥均匀性。在农业物料方面，郑智旗[13]运用离散元分析粉碎后的玉米秸秆在输送装置内的运动规律，通过改变抛送板叶片倾角、旋转轴转速和喂入量等模拟研究参数改变对输送性能的影响;Ghodki等[14]建立了大豆離散元模型，通过与大豆形状、静止角等特征进行比较，获得了该模型的接触参数。因此，本文在对国内外散体物料离散元参数标定方法运用进行总结的基础上，针对椰糠物料具体实际，采用物理试验和建模仿真相结合的手段，基于离散元法分析确定椰糠的颗粒物料特性参数，为后续设计高效的椰糠栽培基质混合专用设备、优化混合结构关键部件和混合工艺参数、提高椰糠栽培基质生产效率与产品品质提供理论依据。

1  材料与方法

1.1  仿真原理及模型选取

椰糠离散元模拟中需要的物料特性参数主要分为三类：一是材料本征参数，包括颗粒大小、堆积密度、剪切模量和泊松比等，这是椰糠自身的特性参数，和外界无关;二是材料接触参数，包括碰撞恢复系数、静摩擦系数和滚动摩擦系数，这是椰糠颗粒相互接触才会发生作用的参数，与接触双方都有关系;三是接触模型参数，跟接触模型选取有关，它描述椰糠颗粒之间接触时颗粒的行为，一些特殊的接触模型需要特殊的模型参数[15-17]。本文选取“Hertz-Mindlin （no slip）”颗粒接触模型，该模型在力的计算方面高效准确[18-19]，其基本原理如图1所示。颗粒之间作用力包括法向弹性力 、法向阻尼力 、切向弹性力 、切向阻尼力 、滚动摩擦力 。

（1）

其中

（2）

（3）

式中： 为等效弹性模量， 为法向重叠量， 为等效接触半径， 、 为接触颗粒的泊松比， 、 为接触颗粒的弹性模量， 、 为接触颗粒的接触半径。

（4）

其中

（5）

（6）

（7）

式中： 为阻尼比， 为法向刚度， 为等效质量， 为相对速度的法向分量， 为恢复系数， 为接触颗粒的质量。

（8）

其中

（9）

式中： 为法向刚度， 为法向重叠量， 为等效剪切模量。

（10）

式中： 为相对速度的切向分量。

（11）

式中： 为滚动摩擦系数， 为接触点到颗粒质心的距离， 为颗粒在接触点处单位角速度矢量。

1.2  试验材料与物理试验

1.2.1  试验材料选取  试验中所用的材料采用椰糠砖充分泡发后的椰糠物料，椰糠砖尺寸为200 mm×100 mm×50 mm、质量为（550±10）g、泡发量为12 L/kg，相关术语取自由中国农业机械化协会设施农业分会发布的团体标准T/CAMA 12—2019《无土栽培椰糠》[20]。

1.2.2  含水率测量  在待测的椰糠物料堆中采用五点法取样，每个取样点称取3～5 g样品，分别放入MA45水分分析测试仪的样品盘中，并确保样品在样品盘中均匀分布。根据5次样品测量的平均值，计算出待测椰糠的含水率为86.57%。

1.2.3  堆积密度测量  称取一定质量的椰糠将其沿喂料漏斗自然下落到量筒中，通过用质量除以量筒内椰糠的体积读数计算出待测椰糠的堆积密度，重复试验5次取平均值为370 kg/m3。

1.2.4  真实密度测量  称取一定质量的椰糠装入盛有一定体积水的100 mL量筒中，量筒增加的体积读数就是椰糠的真实体积，用椰糠质量除以增加的体积，重复5次试验，计算取平均值得到椰糠的真实密度为1 074 kg/m3。

1.2.5  颗粒尺寸测定  称取100 g椰糠样品，使用分级筛进行颗粒分级，用电子天平称取不同粒径范围对应的质量，重复5次试验，获得椰糠粒径质量分数分布，为仿真时椰糠颗粒建模和颗粒工厂生成提供依据，结果如表1所示。

1.2.6  堆积角的测定  堆积角也称安息角或休止角，是表征颗粒物料流动、摩擦等特性的宏观参数，与颗粒表面积、颗粒密度、颗粒摩擦系数及含水率均有关系[21-22]。当大量颗粒物料被倾倒于水平面上堆积为锥体，锥体表面与水平面所成内角即为堆积角，堆积角物理试验是颗粒物料离散元参数标定的常用手段。采用注入法将椰糠从漏斗上方慢慢加入，漏出的物料在水平面上形成圆锥状堆积体的倾斜角即为椰糠的堆积角，重复试验5次，取平均值为45.69°。

1.3  颗粒建模与参数选取

1.3.1  建立椰糠仿真模型  由上述表1可以看出，在选取的椰糠样本中，颗粒粒径1.25～5 mm的质量占总质量的93.18%，其中粒径1.25～2 mm的质量分数为60.38%，占大多数。在离散元参数标定中，若完全按照颗粒真实尺寸进行建模，仿真效率太低，研究表明：相比真实的颗粒尺寸与形状，颗粒之间以及颗粒与几何体之间的相互作用更为关键[23-24]。因此，在建立椰糠颗粒模型时，以粒径1.25～2 mm的颗粒为基础，考虑缩短数值模拟时间和简化模型因素，结合椰糠实际外观形态创建其离散元颗粒模型，如图2所示。

1.3.2  标定参数范围确定  椰糠物料特性的研究目前尚不完善，其离散元仿真的参数尤其匮乏。EDEM软件内置颗粒材料数据库（Generic EDEM Material Model Database）中包含了数千种具有代表性的颗粒物料，如岩石、土壤和矿石等。通过选择仿真规模和堆积密度范围、输入物料堆积角，获得椰糠颗粒之间的接触参数参考值范围：恢复系数为0.15～0.75，静摩擦系数为0.20～1.16，滚动摩擦系数为0.05～0.15。考虑到椰糠作为农业物料及其作为栽培基质与土壤接近，主要参照文献[25–30]中各类物料参数取值或范围，综合考虑设定椰糠泊松比为0.1～0.5，剪切模量为1.0～10.0 MPa，椰糠与钢材恢复系数为0.1～0.7、静摩擦系数为0.20～0.90、滚动摩擦系数为0.05～0.30。通过查阅文献取钢材泊松比为0. 剪切模量为7.9×1010 Pa，密度为7 850 kg/m 重力加速度为9.81 m/s2。各参数标定试验因素及水平如表2所示。

1.4  数据处理

采用Design Expert 10软件设计试验方案并对试验数据进行方差分析和显著性分析，因试验

参数较多，故先设计Plackett-Burman试验筛选出对堆积角影响显著的试验因素，然后根据试验结果对显著性影响因素高（1）、中（0）、低（-1）3个水平设计Box-Behnken试验进行响应面分析，得到堆积角和显著性影响因素之间的回归模型，最后通过参数优化寻找最优解和试验验证完成参数标定。

2  结果与分析

2.1  Plackett-Burman 试验设计与结果

Plackett-Burman设计法是一种两水平的试验设计方法，通过比较各个因素两水平的差异与整体的差异来确定因素的显著性。以椰糠堆积角为响应值，延用表2中8个因素高低水平，并设置3个虚拟参数用于误差分析，共进行12次试验，试验设计及结果如表3所示。

通过方差分析得到8个因素对堆积角的影响效果和显著性，如表4所示。分析可得：A、B、D、E、G对椰糠颗粒的堆积角有正效应，即堆积角随该5个因素的增大而增大;C、F、H则有负效应，即堆积角随该3个因素的增大而减小;其中对堆积角影响极显著（P<0.01）的因素是椰糠-椰糠静摩擦系数（D）和椰糠-椰糠滚动摩擦系数（E），影响显著（P<0.05）的因素是椰糠剪切模量（B），而其余因素影响不显著;模型的P<0.0 说明回归模型极显著，预测值 =0.994 0，说明模型能很好地反映各因素的变化情况。

2.2  Box-Behnken试验设计与结果

2.2.1  试验设计与结果  Box-Behnken设计是一种响应曲面设计类型，根据Plackett-Burman试验结果，取椰糠泊松比为0.3、椰糠-椰糠恢复系数为0.45、椰糠和钢材恢复系数为0.4、椰糠和钢材静摩擦系数为0.55、椰糠和钢材滚动摩擦系数为0.18，以椰糠-椰糠靜摩擦系数（D）、椰糠-椰糠滚动摩擦系数（E）、椰糠剪切模量（B）为试验因素，堆积角为响应指标，延用表2中3个因素的高中低水平，进行Box-Behnken试验设计共进行17模拟仿真，结果如表5所示。

2.2.2  堆积角回归模型分析  应用Design expert软件建立堆积角与3个显著性影响因素的二阶模型：

（12）

该回归方程的P=0.00 决定系数 = 0.9492，校正决定系数 =0.8838，信噪比为13.47，说明该堆积角回归方程极显著，相关性好，可以对目标堆积角进行预测;变异系数CV=3.43%，失拟项P=0.4007>0.05，表明该方程拟合程度好，对该模型试验结果进行回归方差分析，如表6所示。结果表明，椰糠-椰糠静摩擦系数（D）、椰糠-椰糠滚动摩擦系数（E）对堆积角影响极显著;椰糠-椰糠静摩擦系数的二次方项（ ）、椰糠-椰糠滚动摩擦系数的二次方项（ ）对堆积角影响显著。

2.2.3  参数优化与试验验证  运用Design expert软件优化功能，以椰糠堆积角试验值45.69°为目标值，在设置的试验因素水平范围内对堆积角回归模型进行寻优，所得优化解若干组。以优化解作为因素水平值进行堆积角仿真试验，并与物理试验堆积角进行对比验证，得到试验形状相似的一组最优解，即椰糠剪切模量（B）为1.44 MPa，椰糠-椰糠静摩擦系数为（D）为1.1 椰糠-椰糠滚动摩擦系数（E）为0.15，其余非显著性因素取中间水平。此优化解下的堆积角仿真试验结果45.25°与物理试验结果45.69°对比如图3所示，仿真结果与真实情况一致，说明该优化参数可为样品椰糠与混合设备的进一步离散元仿真提供依据。

3  讨论

（1）利用椰糠+有机肥制成混合基质替代草炭，不仅可以“变废为宝”，避免生态资源的枯竭，还可以减少环境的二次污染。采用离散元法研究机械设备与物料相互作用的研究虽然不少，在混合设备领域也有相关报道，但主要集中在饲料、颗粒肥料或者混凝土的研究方面，缺乏针对椰糠基质综合利用的混合搅拌研究。现有文献中针对椰糠物料特性的研究主要分布在结构及比例、分解度、灰分含量、孔隙度、pH、矿质营养等方面，如前所述椰糠离散元模拟中所需的本征参数、接触参数等必备特性参数缺乏数据支撑，而这些参数以试验方法测量往往难以获得。本文选择利用EDEM 2018离散元软件建立椰糠颗粒模型，预设相关参数进行椰糠堆积角的模拟仿真，通过与物理试验相结合进行参数标定，确定椰糠离散元模拟所需的准确参数，研究结果在国内外椰糠加工设备基础研究方面尚未见报道。

（2）本文通过分析文献以及结合前期试验测量的堆积角，根据EDEM 2018材料库确定椰糠泊松比、剪切模量，椰糠-椰糠恢复系数、静摩擦系数、滚动摩擦系数，椰糠-钢材恢复系数、静摩擦系数、滚动摩擦系数这8个因素的取值范围。然后进行了Plackett-Burman试验设计筛选出对椰糠堆积角影响显著的接触参数和本征参数包括：椰糠-椰糠静摩擦系数、椰糠-椰糠滚动摩擦系数和椰糠剪切模量。通过Box-Behnken试验设计建立了堆积角与显著性影响因素的回归模型，方差分析表明了椰糠-椰糠静摩擦系数和椰糠-椰糠滚动摩擦系数对堆积角影响极显著，椰糠-椰糠静摩擦系数的二次方项、椰糠-椰糠滚动摩擦系数的二次方项对堆积角影响显著。以椰糠物理试验堆积角为目标，对回归方程寻优求解，得到显著性影响因素最佳组合为椰糠剪切模量1.44 MPa，椰糠-椰糠静摩擦系数1.1 椰糠-椰糠滚动摩擦系数0.15。最后通过最优组合参数试验验证，证实了标定的参数准确、可靠，建立的模型正确、可行。

（3）离散元法可以准确高效地分析颗粒的运动情况，并具有可视化特点，本文研究结果对椰糠输送、混合等机械装备的研发提供了基础数据。但农业物料的各种特性往往与含水率有很大关系，本文取样以椰糠砖泡发量为依据，虽尽可能缩小了与实际应用中的差距，但受不同栽培基质配方特性的要求，需进一步摸清不同含水率对椰糠堆积角的影响规律，从而确定不同含水率下椰糠物料的特性参数，实现工业化生产。另外，上述试验得出：对椰糠堆积角影响显著的是椰糠-椰糠静摩擦系数、椰糠-椰糠滚动摩擦系数和椰糠剪切模量，它们都是颗粒物料特性，与钢材特性关联不大。而椰糠混合基质生产肯定至少包含2种物料，如椰糠和有机肥等，那么有必要进一步研究椰糠-有机肥的颗粒接触特性，这样才能使整个离散元仿真模拟系统更加完善。

参考文献

[1]  谭海燕， 张振荣， 林峰莹， 等. 椰糠基质限根栽培下樱桃番茄肥料利用率研究[J]. 湖北农业科学， 2020， 59（17）： 46-48.

[2] 邱志豪， 汤柔颖， 韩莹琰， 等. 椰糠与蛭石混配基质理化性状及其对生菜生长的影响[J]. 北京农学院学报， 2019， 34（4）： 62-65.

[3] 肖守华， 赵  西， 肖真真， 等. 以椰糠为基质的设施甜瓜无土栽培基质配方筛选[J]. 山东农业科学， 2019， 51（1）： 61-64.

[4]  Barbara M， Sofia M， Sabrina R， et al. Coconut coir as a sustainable nursery growing media for seedling production of the ecologically diverse quercus species[J]. Forests， 2020， 11（5）： 522.

[5]  劉  佳， 季延海， 王宝驹， 等. 椰糠复合基质对温室番茄生长及品质的影响[J]. 江苏农业科学， 2019， 47（17）： 150-154.

[6]  王跃华， 张明科， 惠麦侠， 等. 不同椰糠配比基质对白菜幼苗生长的影响[J]. 南方农业学报， 2019， 50（12）： 2749-2754.

[7]  李萍萍， 吴  硕， 张西良， 等. 番茄秸秆与发酵菌种立式螺旋混合设备参数优化[J]. 农业机械学报， 2016， 47（11）： 114-120.

[8]  吴  硕. 番茄秸秆混料立式螺旋带式混合方法及性能试验研究[D]. 镇江： 江苏大学， 2016.

[9]  邢洁洁， 张  锐， 吴  鹏， 等. 海南热区砖红壤颗粒离散元仿真模型参数标定[J]. 农业工程学报， 2020， 36（5）： 158-166.

[10]      Ucgul M， Fielke J M， Saunders C. Three dimensional discrete element modeling DEM of tillage accounting for soil cohesion and adhesion[J]. Biosystems Engineering， 2015（129）： 298-306.

[11]     罗  帅， 袁巧霞， Gouda Shaban， 等. 基于JKR粘结模型的蚯蚓粪基质离散元法参数标定[J]. 农业机械学报， 2018， 49（4）： 343-350.

[12]      Bangura K， Gong H， Deng R L， et al. Simulation analysis of fertilizer discharge process using the Discrete Element Method （DEM）[J]. PLoS One， 2020， 15（7）： e0235872 .

[13]      郑智旗. 玉米秸秆捡拾粉碎沟埋还田机研究[D]. 北京： 中国农业大学， 2017.

[14]      Ghodki B M， Patel M， Namdeo R， et al. Calibration of discrete element model parameters： soybeans[J]. Computational Particle Mechanics， 2018， 6（1）： 3-10.

[15]      彭才望， 许道军， 贺  喜， 等. 黑水虻处理的猪粪有机肥离散元仿真模型参数标定[J]. 农业工程学报， 2020， 36（17）： 212-218.

[16]      张  珂， 于天行， 于文达， 等. 基于JKR黏结模型的混凝土离散元参数标定[J]. 混凝土， 2020（8）： 46-50， 55.

[17]      文恩杨， 李玉华， 牛子孺， 等. 蒜种颗粒离散元模型参数标定[J]. 农机化研究， 2021， 43（5）： 160-167.

[18]      武  涛， 黄伟凤， 陈学深， 等. 考虑颗粒间黏结力的黏性土壤离散元模型参数标定[J]. 华南农业大学学报， 2017， 38（3）： 93-98.

[19]     Cundall P A， Strack O D L. A discrete numerical model for granular assembles[J]. Geotechnique， 1979， 29（1）： 47-65.

[20]      团体标准T/CAMA 12—2019无土栽培椰糠[J]. 农业工程技术， 2019， 39（31）： 69-73.

[21]      吴  鹏， 张喜瑞， 李  粤， 等. 基于离散元的砖红土壤模型堆积特性研究[J]. 农机化研究， 2020， 42（9）： 145-150.

[22]      刘羽平， 张  拓， 刘  妤. 稻谷颗粒模型离散元接触参数标定与试验[J]. 中国农业科技导报， 2019， 21（11）： 70-76.

[23]      张  锐， 韩佃雷， 吉巧丽， 等. 离散元模拟中沙土参数标定方法研究[J]. 农业机械学报， 2017， 48（3）： 49-56.

[24]      Thomas S， Hassanizadeh S M， Hamed A， et al. The effect of particle shape on porosity of swelling granular materials： Discrete element method and the multi-sphere approximation[J]. Powder Technology， 2020， 360（15）： 1295-1304.

[25]      马彦华， 宋春东， 宣传忠， 等. 苜蓿秸秆压缩仿真离散元模型参数标定[J]. 农业工程学报， 2020， 36（11）： 22-30.

[26]      廖宜涛， 廖庆喜， 周  宇， 等. 饲料油菜薹期收获茎秆破碎离散元仿真参数标定[J]. 农业机械学报， 2020， 51（6）： 73-82.

[27]     吴文江， 郭  斌， 高占凤， 等. 基于DEM的沙土颗粒建模及参数标定研究[J]. 中国农机化学报， 2019， 40（8）： 182-187.

[28]      向  伟， 吴明亮， 吕江南， 等. 基于堆积试验的黏壤土仿真物理参数标定[J]. 农业工程学报， 2019， 35（12）： 116-123.

[29]      Mehrdad P， Colin H， Mojtaba G， et al. Effect of particle shape on flow in discrete element method simulation of a rotary batch seed coater[J]. Powder Technology， 2016， 296： 29-36.

[30]      Józef H， Marek M. Parameters and contact models for DEM simulations of agricultural granular materials： A review[J]. Biosystems Engineering， 2016， 147： 206-225.

責任编辑：崔丽虹