水稻茎秆接触物理参数测定与离散元仿真标定

侯 杰，谢方平，2*，王修善，2，刘大为，2，陈志刚

（1.湖南农业大学 机电工程学院，湖南 长沙 410128；2.湖南省智能农机装备重点实验室，湖南 长沙 410128）

【研究意义】水稻是我国主要粮食作物之一，其种植面积广、产量大[1]，机械化收获水平逐年提高[2]，但是由于水稻生物学特性变化，水稻机械化收获清选过程中物料堵塞问题严重[3]，解决该问题的前提是建立较为精确的水稻模型，模拟其在联合收割机中的运动规律，深入研究堵塞原因，可为后续联合收割机优化设计提供理论基础。【前人研究进展】离散元法是专门用来解决不连续介质问题的数值模拟方法[3]，在离散元仿真模拟中，参数标定是最关键的问题之一，EDEM 模拟谷物脱粒与清选时[4-9]，茎秆等的参数设置直接影响仿真结果的可信度。近年来，离散元参数标定方法不断增多[10]，刘凡一等[11]采用三点弯曲试验对小麦短茎秆粘结参数进行离散元标定，结果表明：茎秆的弹性模量随着粘结半径、单位面积法向刚度和单位面积剪切刚度的增大而增大；张涛等[12]以径向堆积角相对误差值为评价指标，应用正交方法标定玉米秸秆离散元仿真需要输入的接触参数；王韦韦等[13]通过堆积角试验和模孔压缩试验标定玉米秸秆接触参数，试验结果表明：仿真试验堆积角与物理堆积角相对误差为0.68%，仿真试验与物理试验模孔压缩比相对误差为0.98%，表明标定的离散元参数可用作玉米秸秆离散元仿真研究；马彦华等[14]以休止角为目标值对苜蓿秸秆压缩过程中离散元仿真研究所用接触参数进行标定；王万章等[15]通过斜面法对比仿真与物理试验，对小麦茎秆接触参数进行标定；廖宜涛等[16]通过堆积角和响应面试验标定了油菜茎秆离散元接触和粘结参数。

【本研究切入点】较多研究者对茎秆离散元参数标定提供了可行的方法，但少有对水稻茎秆接触参数和粘结参数的全部标定，且试验设计所用方法未涉及两个显著性影响因素分析。【拟解决的关键问题】为减小脱粒清选中水稻茎秆EDEM 模拟仿真误差，笔者以带有叶鞘包裹的水稻茎秆为研究对象，以茎秆堆积角、茎秆弯曲弹性模量为评价指标，选取EDEM 软件中HBP（Hertz-Mindlin with bonding particle）模型，结合Plackett-Burman 试验、Central-Composite 试验和Box-Behnken 试验对水稻茎秆接触参数和粘结参数进行标定，对比分析物理试验与仿真试验结果进行验证，以期为水稻植株在脱粒清选等环节离散元建模及仿真研究提供参数依据。

1 材料与方法

1.1 水稻茎秆本征参数

选用湖南省常德市桃源县叶家冲村早稻，品种为‘中早35’。试验在湖南农业大学物料收获实验室完成，随机抽10 根水稻，选取整株水稻第一节（靠近穗头部分）由叶鞘包裹的茎秆做样本（图1），剪去叶片、穗头与剩余茎秆（水稻植株第一个结点之后部分），利用卷尺和数显游标卡尺测得茎秆本征参数均值，其长度为201.50 mm、上端（靠近穗头的一端）外径为2.69 mm、下端（远离穗头的一端）外径为3.94 mm、上端壁厚为0.44 mm、下端壁厚1.04 mm，试验茎秆真实密度为152.3 kg/m3，试验时茎秆含水率为79.04%。

图1 试验样本Fig.1 Test sample

1.2 堆积角试验

依据文献[12，15]，采用自制可抽钢板的水稻茎秆径向堆积角测定装置进行堆积角试验，钢板材料为不锈钢，测定装置长宽高分别为400，200，200 mm。试验时所有茎秆长度取200 mm，平均外径4 mm，平均壁厚0.7 mm，将一定数量的短茎秆放入堆积角测定装置一侧，堆积宽度为40 mm，竖直向上迅速抽出钢板，水稻茎秆整体向另一侧运动，与水平钢板形成堆积角，应用相机拍摄记录茎秆堆积后形成的照片，堆积角测量与处理过程如图2所示。

应用Origin2019软件对拍摄图像进行图2过程处理，由拟合直线斜率计算茎秆堆积角。试验图像处理时，为减少拍摄角度带来的误差，分别以水平钢板为X 轴，竖直固定钢板为Y 轴进行曲线拟合，试验重复5 次，进而得到堆积角均值为47.67°，拟合直线决定系数为0.937 9，接近于1，表明拟合方程可靠程度高，物理试验所测堆积角可用于后续仿真试验寻优及参数验证。

图2 堆积角测量与处理Fig.2 Measurement and processing of stacking angle

1.3 三点弯曲试验

笔者采用三点弯曲试验对茎秆弯曲特性进行分析测量[17-20]，基于弯曲特性中的弹性模量对茎秆粘结参数进行标定。试验采用CMT6104 电子万能试验台，其主要性能参数如表1 所示。试验时将茎秆两端固定在水平支座上，支座间距为30 mm，试验通过万能力学试验装置以0.06 m/s对直径为6 mm 的钢制圆柱压头从茎秆中心加载（图3所示），弯曲弹性模量Eb按下式计算可得[20]：

图3 三点弯曲试验Fig.3 Three point bending test

表1 CMT6104主要性能参数Tab.1 Main performance parameters of CMT6104

式中Eb是茎秆弯曲弹性模量，MPa；F是加载力，N；S是两支撑间的长度，mm；W是茎秆的弯曲挠度，mm；I是茎秆横截面相对中性轴的惯性矩[20]，mm4。

根据生物学特性，水稻茎秆横截面可近似为圆环[20]，其惯性矩I按下式计算：

式中d是茎秆外径，mm；t是茎秆壁厚，mm。

试验进行10次取均值，根据式 （1）和式（2）可得茎秆弯曲弹性模量为3 740.71 MPa。

2 仿真模型

2.1 茎秆模型创建

笔者采用多球圆环排列方式构建单层水稻茎秆模型[20]，由离散元软件的颗粒工厂和颗粒体力API颗粒替换构建茎秆模型，通过内置的HBP 模型沿轴向粘接。由于真实秆壁很薄，若采用真实模拟，则计算量很大，故在模拟中放大了秆壁厚度[20]，依据文献[20]，仿真放大了的秆壁厚度对试验结果无显著影响，因此，根据物理试验茎秆基本特征参数将仿真模型茎秆壁厚设为1 mm、外径为4 mm、长度为200 mm。

图4 水稻茎秆离散元模型Fig.4 Discrete element model of rice stalk

2.2 茎秆-钢板静摩擦系数测定

测定时将待测茎秆轴向放置（保证茎秆只滑动不滚动）在带有三量数显倾角仪的钢板上，三量数显倾角仪通过双面胶固定在钢板上，钢板一端与具有60°倾角的固定钢板平齐，钢板下放置一滑块，通过推动滑块以调节钢板的倾斜角度，当水稻茎秆与不锈钢板发生位移时停止推动，记录此时数显角度测试仪的读数α，按照式（3）计算摩擦系数，试验重复10次取均值，测定过程如图5所示，测得茎秆-钢板静摩擦系数为0.4。

图5 静摩擦系数测定试验Fig.5 Test for determination of static friction coefficient

式中fs为静摩擦系数；α为静摩擦临界角，（°）。

2.3 茎秆-钢板滚动摩擦系数测定

茎秆-不锈钢板滚动摩擦系数测定方法与静摩擦系数测定方法相似，测定时将茎秆横向放置在钢板上，拍摄滑块前进钢板缓慢抬起的视频过程，当茎秆发生滚动后停止拍摄，用Premiere 软件导出茎秆开始滚动那一帧图片，读出此时数显角度测试仪的读数β，按照式（7）计算滚动摩擦系数，试验重复10次取均值，测定过程如图6所示，测得茎秆-钢板滚动摩擦系数为0.3。

图6 滚动摩擦系数测定试验Fig.6 Rolling friction coefficient measurement test

式中M为滚动摩擦力偶矩，N*m；FN2为支撑面的正压力，N；fr为滚动摩擦系数；G为茎秆的重力，N；β为滚动摩擦临界角，（°）；r为茎秆的半径，mm。

2.4 仿真参数选取与模型建立

根据文献[21]，试验中不锈钢板密度为7 860 kg/m3、泊松比为0.3、剪切模量为7.9×1010Pa。在茎秆堆积角仿真试验时，茎秆生成方式为Dynamic，生成速率设为5 000 个/s，生成总数量为100 个，时间步设为20%，时间间隔为0.01 s，网格尺寸为3.0R，茎秆生成时间为0.2 s，待茎秆模型生成后取出隔板进行堆积，堆积时间为2 s。在茎秆三点弯曲仿真试验时，仿真几何模型与实测模型相同，时间步设为10%，时间间隔为0.01 s，网格尺寸为3.0R，茎秆弯曲仿真时间为2 s，建立仿真模型如图7所示。

图7 仿真试验模型Fig.7 Simulation test model

3 试验设计与结果分析

3.1 Plackett-Burman试验及显著性分析

Plackett-Burman 试验是基于目标响应与各因素间关系，比较每个因素2 水平间的差异来确定因素显著性[14]。试验以水稻茎秆堆积角θ和弯曲弹性模量Eb为响应值，根据国内外文献对影响参数设计如表2 所示[5-10,13-16]，各参数取高低2 个水平，采用1 个中心点，对仿真接触参数和粘结参数的显著性进行评估试验。

表2 Plackett-Burman仿真试验参数列表Table2 Plackett-Burman simulation test parameter list

由Plackett-Burman 仿真试验得到的响应值结果如表3 所示，其中X1～X8为参数编码值，利用Design Expert 软件对响应值θ和Eb进行方差分析，得到各影响参数的显著性如表4 所示。由表4 可知，茎秆-茎秆静摩擦系数us与茎秆-茎秆动摩擦系数u（tP＜0.01）对堆积角θ影响极显著；茎秆-茎秆滚动摩擦系数ut与茎秆切向粘结刚度kt（b0.01＜P＜0.05）对茎秆弹性模量Eb影响显著，茎秆粘结半径Rb与茎秆法向粘结刚度kn（bP＜0.01）对茎秆弹性模量Eb影响极显著。

表3 Plackett-Burman试验结果Tab.3 Plackett-Burman test scheme and results

表4 Plackett-Burman试验方差分析Tab.4 Analysis of variance of Plackett-Burman test

3.2 Central-Composite试验及回归模型

为对影响实测茎秆接触参数指标θ的us和ut进行标定，依据Plackett-Burman试验结果，以us和ut为影响因素，θ为评价指标，选取3个中心点进行Central-Composite试验，试验因素设计与结果如表5。仿真试验中对θ影响不显著的其它参数均取中间值。根据表5得到关于堆积角θ的多元回归拟合方程为：

表5 Central-Composite试验设计与结果Tab.5 Central-composite test design and results

Central-Composite 试验回归模型方差分析如表6 所示，从表6 可知，该拟合回归模型P＜0.01，说明堆积角与所得回归方程有密切关系；其中各个参数的P值均小于0.05、表明了回归模型的有效性。试验精确度Adeq Precision=8.991 8＞4，说明模型可行且具有良好的精确度；决定系数R2=0.925 3 趋近于1，表明物理试验与仿真试验拟合度好，试验误差小。在Desirability 值为1（可靠度最高）的一百组回归解中寻找与实测堆积角相近的一组，可得茎秆-茎秆静摩擦系数为0.35、茎秆-茎秆滚动摩擦系数0.016。

表6 Central-Composite试验θ方差分析Tab.6 Analysis of variance of Central-composite test θ

3.3 Box-Behnken试验及回归模型

为对影响实测茎秆粘结参数指标Eb的Gb、Rb和knb进行标定，基于Plackett-Burman 试验因素范围以及Central-Composite 试验确定的us和ut，依据Box-Behnken 试验原理，以Rb、knb和ktb为影响因素，Eb为评价指标，选取3 个中心点进行3 因素3 水平Box-Behnken 试验，试验设计与结果如表7。仿真试验中除已标定参数外，对Eb影响不显著的其它参数均取中间值。根据表7，应用Design Expert 软件得到关于弹性模量Eb的多元回归拟合方程为：

表7 Box-Behnken试验设计与结果Tab.7 Box-Behnken test design and results

Box-Behnken 试验回归模型方差分析如表8所示，分析可知，该拟合线性回归模型P＜0.01，说明弹性模量与所得回归方程是极显著的；其中Rb、knb和ktb的P值均小于0.01、说明这3个参数对于弹性模量的影响极显著，表明了回归模型的有效性。

表8 Box-Behnken试验Eb方差分析Tab.8 Analysis of variance of Box-Behnken test Eb

试验精确度Adeq Precision=28.912 0＞4，说明模型可行且具有良好的精确度；决定系数R2=0.993 4趋近于1，表明拟合方程可靠度高。在Desirability值为1得一百组回归解中寻找与实测弹性模量相近的一组，可得茎秆粘结半径Rb为2.331 mm、茎秆法向粘结刚度knb为4.7×109N/m3、茎秆切向粘结刚度ktb为3.7×109N/m3。

3.4 标定结果验证

根据茎秆接触参数和粘结参数的标定结果，应用EDEM 软件建立最佳参数离散元仿真模型，以验证标定试验结果的可靠度。试验时茎秆-茎秆静摩擦系数为0.35、茎秆-茎秆滚动摩擦系数为0.016、茎秆粘结半径为2.331 mm、茎秆法向粘结刚度为4.71×109N/m3、茎秆切向粘结刚度为3.67×109N/m3，茎秆泊松比为0.3，茎秆剪切模量为5.05 GPa，茎秆-茎秆恢复系数为0.3，试验重复5 次，模拟测得茎秆堆积角分别为47.79°、51.59°、43.87°、46.77°、50.86°，平均值为48.27°，与实际物理堆积角（47.67°）的相对误差为1.26%，应用SPSS 软件进行T检验，得P1=0.737；模拟测得茎秆弹性模量分别为3 672.21，3 819.74，3 595.82，3 661.68，3 819.74 MPa，平均值为3 713.84 MPa，与实际物理弹性模量（3 740.71 MPa）得相对误差为0.72%，应用SPSS 软件进行T检验，得P2=0.584。两次T检验所得P值均大于0.05，表明仿真试验与物理试验无显著性差异。综上可知，笔者所标定的水稻茎秆离散元仿真参数准确可靠，可用于水稻离散元仿真建模研究。

4 结论

笔者基于EDEM 离散元软件，选用HBP（Hertz-Mindlin with bonding）模型，通过Plackett-Burman 试验、Central-Composite 试验和Box-Behnken 试验对水稻短茎秆样本的接触参数和粘结参数进行标定，主要结论如下：

（1）根据Plackett-Burman试验方差分析可知，茎秆-茎秆静摩擦系数与茎秆-茎秆滚动摩擦系数对堆积角影响极显著，茎秆-茎秆滚动摩擦系数与茎秆切向粘结刚度对茎秆弹性模量影响显著，茎秆粘结半径与茎秆法向粘结刚度对茎秆弹性模量影响极显著。

（2）根据Central-Composite 试验方差分析和回归方程可得，对茎秆影响显著的2 个接触参数的一次项、交互项以及二次项对水稻茎秆堆积角影响显著；根据Box-Behnken 试验方差分析和回归方程可得，对茎秆影响显著的3 个粘结参数的一次项与交互项对水稻茎秆弹性模量影响极显著，茎秆粘结半径的二次项对水稻茎秆弹性模量影响显著。

（3）以水稻茎秆实际堆积角与弹性模量为目标，对回归方程进行寻优求解，得到显著性接触参数的最优值为茎秆-茎秆静摩擦系数0.35、茎秆-茎秆滚动摩擦系数0.016，显著性粘结参数的最优值为茎秆粘结半径2.33 mm、茎秆法向粘结刚度4.71×109N/m3、茎秆切向粘结刚度3.67×109N/m3。

（4）最优接触参数得出的仿真试验堆积角与物理试验堆积角的相对误差为1.26%，基于SPSS对仿真堆积角进行T试验，测得P1=0.737＞0.05，证明仿真堆积角与物理堆积角无显著性差异；最优粘结参数得出的仿真试验弹性模量与物理试验弹性模量的相对误差为0.72%，基于SPSS对仿真弹性模量进行T试验，测得P2=0.584＞0.05，证明仿真弹性模量与物理弹性无显著性差异。标定结果验证试验表明，基于HBP 接触模型所标定的接触参数与粘结参数可用于水稻茎秆离散元仿真，可为后续水稻整株仿真研究以及联合收割机物料优化设计提供参数依据。

致谢：湖南省高新技术产业科技创新引领计划项目（2020NK2002）同时对本研究给予了资助，谨致谢意！