基于X射线断层扫描的小麦籽粒模型构建及离散元参数标定

朱 康， 惠延波， 周 颖， 王宏晓， 张云龙， 李 辉， 陈艳雷， 郑金锋， 白路遥

(河南工业大学先进制造研究所，郑州 450001)

离散元法又称离散单元法，是一种用于模拟并分析散体介质系统动力学行为的数值方法，广泛应用于地质工程、土木工程、农业工程[1]。在农业工程中，农业装备会涉及到大量的谷物颗粒的运动，运用离散元法研究两者之间的相互作用关系，可以优化农业装备的机械结构、提升装备工作效率，为现代农业装备的数值化设计提供理论基础[4]。

在进行离散元仿真时，模型的精度和仿真参数的设置会对仿真结果产生较大的影响[11]。刘凡一等[9]使用游标卡尺测量小麦籽粒的平均长轴、短轴，将小麦近似地看作规则的椭球型，在EDEM中采用5球组合的方式构建小麦离散元模型对小麦离散元参数进行了标定，但在实际中，小麦籽粒形状是不规则的，通过数量较少的球形颗粒的堆叠，无法准确地反映小麦籽粒的真实轮廓。赵武云[10]对玉米种子进行切片，每片厚度约为1 mm，采集玉米种子剖面图像，根据剖面图像提取种子外轮廓样条曲线，在SolidWorks中进行轮廓曲线放样绘制玉米种子三维模型，但切片数量过少，无法完整地表示玉米的轮廓特征。张荣芳等[11]借助三维扫描和逆向工程的方法获取水稻种子外形轮廓，构建不同半径的填充球离散元模型对种间静摩擦因数和动摩擦因数进行了标定，但三维激光扫描仪针对尺寸较小的物体扫描的精度较低，无法完整地获取物体三维轮廓。

因此，本实验以小麦为研究对象，基于X射线断层扫描技术获取1 024张不同灰度值的二维切片，采用图像滤波、图像分割等算法构建出小麦籽粒的三维外轮廓模型；在EDEM中导入外轮廓模型，分别填充4种不同数量的球(17、49、80、126球)来拟合小麦籽粒的真实形状，构建了更加精准的小麦籽粒离散元模型；通过休止角实际实验与EDEM仿真实验，运用Plackett-Burman实验、最陡爬坡实验、响应面优化实验标定了小麦离散元仿真参数；最后对标定的参数进行仿真分析并与实际实验进行对比，验证了本方法构建的离散元模型标定参数的准确性，为其他不规则谷物模型的构建提供了参考。

1 小麦籽粒模型构建

1.1 实验材料

实验材料选取丰德存麦5号，含水率12.5%。

1.2 图像采集

采用YXLON Y.Cheetah型微米X射线断层扫描系统(Micro-CT)采集小麦籽粒三维灰度图像，如图1所示。随机选取外形尺寸完好的小麦籽粒，通过X射线断层扫描系统对其进行扫描，单颗籽粒1次扫描约20 min，可以无损地获取小麦籽粒的1 024张二维灰度图像。

图1 二维灰度图像

1.3 图像重构

从二维灰度图像可以看出，小麦籽粒的不同组织对X射线的吸收不同，呈现的灰度值也就不同，因此可以有效地区分各组织的轮廓[12]。通过二维灰度图像的堆叠，可完整地构建出小麦籽粒轮廓。采用非局部均值滤波算法，去除图像噪声，提高图像质量；采用交互式阀值分割算法分离出小麦籽粒的外轮廓组织，如图2所示；最后基于体绘制及面绘制方法重构出小麦籽粒三维外轮廓模型。构建好的三维模型存在表面粗糙、孔洞等问题，将三维模型转化为点云格式导入Geomagic Design中进行孔洞的填补、表面光滑处理等操作，最终生成实体模型[13]，如图3所示。

图2 图像处理

图3 小麦籽粒模型构建

1.4 离散元模型构建

将X射线断层扫描方法构建的小麦籽粒外轮廓模型导入EDEM中，以此为基础拟合小麦籽粒的真实形状，构建出17、49、80、126球组合的小麦籽粒仿真模型，如图4所示。

图4 小麦籽粒离散元仿真模型

2 离散元仿真参数标定

2.1 休止角实际实验

休止角是物料由高处滑落至平面上堆积形成的锥体的斜面与水平面之间的夹角，是表征颗粒流动、摩擦特性等特性的宏观参数[14]。运用离散元法研究的谷物的休止角，可以掌握谷物堆积形成粮堆的机制，了解颗粒与颗粒、颗粒与装备的之间的力学特性。

本次离散元参数标定采用休止角实验与离散元仿真实验结合的方法，实验方法参照GB/T 16913—2008，实验装置如图5所示，由无底圆筒(长度120 mm，内径40 mm)、有机玻璃板(150 mm×150 mm)、实验台组成。无底圆筒与底面有机玻璃板重合，将小麦填满无底圆筒后，以0.02 m/s的速度向上匀速提升圆筒，小麦自然落下并在底面堆积，锥面与底面之间的夹角即为休止角，每组实验重复5次。采集籽粒堆积的8个方向的图像，使用Matable读取采集的图像，对图像进行二值化、边缘提取、边缘直线拟合来测量小麦堆积的休止角[15]，如图6所示。对采集的数据去除极值，取其平均值，最终测得实际休止角为20.326 5°。

注：1.实验台 2.小麦 3.无底圆筒 4.有机玻璃板。图5 小麦休止角实验

图6 Matable图像处理

2.2 休止角仿真实验

2.2.1 仿真参数

本次实验所需的仿真参数如表1所示[16]。

表1 小麦离散元仿真参数

2.2.2 仿真实验

为缩短仿真时间，休止角仿真实验选用17球小麦模型，在EDEM中构建长120 mm、内径40 mm、上下开口的圆筒。在圆筒顶端以5 000颗/s的速度随机生成小麦颗粒直至填满圆筒，等待颗粒稳定，此过程时间为1 s；随后以0.02 m/s的速度向上提升圆筒，直至颗粒全部落下并形成稳定的堆积体，此过程过程仿真时间为3 s，如图7所示。仿真模型选用Hertz-Mindlin (no slip)接触模型，时间步长为20%，仿真总时长为5 s。

图7 小麦籽粒休止角离散元仿真

2.3 参数标定

2.3.1 Plackett-Burman实验

在Plackett-Burman实验设计的基础上，以小麦籽粒的休止角为响应值，选取影响小麦籽粒休止角的八个参数进行评价，筛选出影响小麦休止角的显著因素[19]。每个因素取低水平(-1)和高水平(+1)2个水平，共12组实验，如表2所示。

表2 Plackett-Burman实验设计的因素水平及编码

Plackett-Burman实验设计及结果如表3所示。运用Design-Expert软件进行求解得到的结果如表4所示。根据各因素的P值大小判断出影响小麦休止角的关键因素依次为：E(小麦-小麦静摩擦系数)>G(小麦-小麦滚动摩擦系数)>F(小麦-有机板静摩擦系数)。其中小麦-小麦静摩擦系数、小麦-小麦滚动摩擦系数及小麦-有机板静摩擦系数均达到极其显著水平。在后续的最陡爬坡实验中只考虑这3个参数的影响，结合其他仿真参数进行最陡爬坡实验和响应面优化实验。

表3 Plackett-Burman实验设计及结果

表4 Plackett-Burman实验设计各因素效应评价

2.3.2 最陡爬坡实验

基于Plackett-Burman实验筛选出较显著的参数，设计最陡爬坡实验来快速逼近各显著参数的最佳区域。根据小麦-小麦静摩擦系数、小麦-小麦滚动摩擦系数和小麦-有机板静摩擦系数3个因素的效应大小设计相应的步长，实验设计及结果如表5所示。在第二组时，休止角的相对误差最小，因此选择第二组为水平中心点，第一组为低水平，第三组为高水平进行响应面优化实验。

表5 最陡爬坡实验设计及结果

2.3.3 响应面优化实验

根据Plackett-Burman和最陡爬坡实验结果，以小麦-小麦静摩擦系数、小麦-小麦滚动摩擦系数和小麦-有机板静摩擦系数为自变量，休止角θ为因变量建立三因素三水平实验，如表6所示。

表6 Box-Behnken实验设计及结果

通过Design-Expert求解得到休止角与3个显著参数的二次多项式回归模型为：

θ=26.51+4.13A+1.14B+3.17C-0.19B-1.18AC-1.05BC-1.32+0.37B2-4.57C2

式中：θ为休止角；A为小麦-小麦静摩擦系数；B为小麦-小麦滚动摩擦系数；C为小麦-有机板静摩擦系数。

回归模型方差分析结果如表7。其中，回归模型P<0.000 1(极其显著)，失拟项P=0.52>0.05(不显著)，说明该回归模型可靠，拟合程度高。由A、B、C的P值可判断出：三者对休止角的影响都极其显著。此外，模型的决定系数R2=0.983 8，说明回归方程相关性较好；校正决定系数adjR2=0.962 9，说明预测值与实际值具有较高的相关性；变异系数CV=3.72%，说明模型合理。此模型可用于对休止角的影响因素进行分析。

表7 回归模型方差分析结果

运用Design-Expert中的Optimization功能，以实际休止角为目标值，求解得到最优参数为小麦-小麦静摩擦系数为0.185、小麦-小麦滚动摩擦系数为0.046、小麦-有机板静摩擦系数数为0.203。

3 验证实验

分别用构建的4种小麦籽粒模型对标定的离散元参数进行仿真，并于实际实验进行对比分析，得到的仿真结果结果及相对误差如表8所示。从表8中可以看出，4种模型的仿真结果相近，与实际休止角的相对误差小于1.30%，表明本文基于X射线断层扫描方法构建的模型、标定的参数可用于离散元仿真。

表8 实验验证结果

4 结论

基于X射线断层扫描技术，运用图像滤波、图像分割等算法重构了小麦籽粒的外轮廓模型；在EDEM中，分别填充4种不同数量的球(17、49、80、126球)来拟合小麦籽粒的真实形状，构建了更加精准的小麦离散元模型。

基于X射线断层扫描技术构建的离散元仿真模型，通过休止角实验与离散元仿真实验结合，运用Plackett-Burman实验、最陡爬坡实验、响应面优化实验，标定了小麦离散元仿真所需要的参数；用构建的4种小麦籽粒模型对标定的离散元参数进行仿真分析并与实际实验进行验证，结果表明本文方法构建的离散元模型及标定的参数可用于离散元仿真。