基于DEM-CFD的玉米气吸式排种器种盘设计与试验

丁 力 杨 丽 张东兴 崔 涛 高筱钧

(1.中国农业大学工学院， 北京 100083； 2.农业农村部土壤-机器-植物系统技术重点实验室， 北京 100083)

0 引言

玉米作为我国第一大粮食作物，种植面积和产量逐年增加[1-2]。目前，农村劳动力缺失，为满足发展现代农业需要，农村土地分散经营已向适度规模经营转变，高速精量播种技术是推动玉米产业可持续发展的重要途径之一[3-4]。

图1 排种器结构图Fig.1 Structure diagrams of metering device1.前外壳体 2.清种毛刷 3.清种刀 4.种盘 5.快拆机构 6.后外壳体 7.密封垫 8.排种轴 9.卸种机构 10.挡种刷 11.卸种挡板 12.种层高度调节板 13.吸孔 14.型孔凸台 Ⅰ.充种区 Ⅱ.携种区 Ⅲ.投种区 Ⅳ.过渡区

玉米高速播种过程中，由于种盘高速回转，充种时间急剧下降，加上玉米种子形状大小差异较大，导致漏播指数显著上升。种盘作为充种的关键部件，众多学者对此做了大量研究[5-8]。但未考虑不易充种的玉米种子形状，没作针对性分析，且因排种器实际工作中种子运动十分复杂，难以进行准确分析计算。

近些年，离散单元法与计算流体动力学的耦合方法已应用于气力输送、气流清选和流化床等方面的研究[9-11]。目前，此方法主要应用在气吹式排种器上[12-13]。而离散单元法(DEM)与计算流体动力学(CFD)耦合的分析方法在气吸式排种器上应用则鲜有报道。本文采用离散元和计算流体动力学耦合的方法模拟气吸式排种器工作过程，分析玉米种子的运动受力情况，选取不易充种的种子形状颗粒，有针对性地优化设计种盘结构参数；通过二次旋转正交组合试验，求得所设计排种器在工作速度范围内的最佳参数匹配，并进行试验验证。

1 排种器结构模型

排种器的结构如图1所示，种盘是核心工作部件，由型孔凸台、吸孔和种盘本体构成，型孔凸台均匀分布于吸孔周边，整个种盘分为充种区Ⅰ、携种区Ⅱ、投种区Ⅲ和过渡区Ⅳ 4个区域。

工作时，种子在重力的作用下落入排种器腔室底部，转动的种盘上的型孔凸台扰动充种区内的种子，并托持种子形成一定的初速度；同时，吸孔在气室负压的作用下吸附种子，吸附于吸孔上的种子在负压吸附力和型孔凸台托持力的作用下随种盘一起转动，转过携种区到达卸种区时，在卸种机构的外力作用下快速进入导种管并落入已开好的种沟中，完成播种作业。

2 仿真建模与分析

气吸式排种器工作过程中，种子处在气流场、颗粒场和重力场并存的复杂环境中[14]，传统的CFD模拟和EDEM模拟受其方法所限，无法准确模拟排种器内复杂受力情况下种子颗粒的运动情况，从而难以准确模拟其工作过程。因此，采用DEM-CFD气固耦合的方式，通过仿真分析种子颗粒的实际位置和运动情况，从而得到更加丰富的模拟结果[15-16]。

2.1 模型建立

图2 玉米种子仿真模型Fig.2 Simulation models of maize grains

以应用较为广泛的郑单958种子为模型，按照文献[17]玉米种子分级方式，将玉米种子分为大扁形、小扁形和类圆形3类。为解决气固双向耦合过程中要求仿真颗粒体积小于流场网格最小体积这一问题[18]，采用Bonding黏结模型填充玉米籽粒，即利用多个体积小于流场网格的圆球型颗粒对SolidWorks建立好的玉米模型填充，并将所有填充的小球利用Bonding力黏结在一起，将黏结在一起的模型作为玉米籽粒的仿真模型[17]。利用Bonding力粘结的小颗粒填充的模型如图2所示。图中从左至右依次为种子的实物图、三维模型图和颗粒黏结模型图。

在SolidWorks中建立排种器简化模型，保存为step格式导入EDEM中，并开始颗粒生成，模型中各部分如图3所示。流体域通过ICEM-CFD划分结构化网格，并定义流体入口、出口，吸孔接触面都设置为interface，将划分好的网格保存成网格文件，导入Fluent软件，设置滑移网格，其中：吸孔部分为动区域，其他区域为静区域，网格划分及相关设置如图4所示。

图3 仿真简化模型示意图Fig.3 Simplified model diagram of simulation1.清种刀 2.种盘 3.后外壳体 4.前外壳体 5.颗粒工厂 6.颗粒

图4 网格划分及相关设置示意图Fig.4 Diagram of grid partition and related settings1.出口(静区域) 2.吸孔(动区域) 3.入口(静区域)

2.2 接触模型

离散元法描述的是颗粒之间碰撞的过程，即接触的产生和发生作用的过程。离散元中根据接触方式的不同可分为硬颗粒接触和软颗粒接触，软颗粒接触方式允许颗粒之间接触点间出现重叠部分，并根据接触颗粒的物理属性和法向重叠量、切向位移计算接触力，而颗粒间内摩擦力即切向力与法向力有着密切的联系，因此本文采用软颗粒接触模型。考虑到玉米种子表面无粘附力，本文选取Hertz-Mindlin无滑移接触模型。

2.3 仿真参数确定

DEM-CFD气固耦合仿真需保证在EDEM中设置的种盘转速和Fluent中设置配合流体区域的型孔转速相同，考虑到所设计排种器为高速排种器，作业速度高达14 km/h，换算成种盘转速，设置种盘角速度为3.62 rad/s，在进种口设置颗粒工厂。为保证仿真真实性，同时节约仿真时间，经过多次尝试，选取160粒种子作为颗粒工厂生成的颗粒数，其中大扁形85粒、小扁形35粒、类圆形40粒。

排种器整体采用有机玻璃3D打印而成，根据所用材料确定相关参数如表1所示[19]。在EDEM中设置时间步长为1×10-5s，Fluent中时间步长应为EDEM时间步长的50～100倍，选取0.000 5 s作为Fluent中时间步长；设置Fluent步数为6 000步，即仿真时间为3 s，设置Max Iterations/Time Step 为80，即每个时间步最多迭代80次；每0.01 s保存一次数据。

表1 模拟所需物理和力学特性参数Tab.1 Physical and mechanical parameters required for simulation

2.4 仿真结果分析

通过仿真观察，标记吸附效果较好的3种类型种子Bonding模型，分析其受力状况和运动速度，如图5所示。提取每个Bonding模型小颗粒所受到的曳力总合，即整个种子所受到的曳力，如图5a所示。为了进一步观察种子的运动情况，提取标记Bonding模型每个颗粒的速度均值作为种子的整体速度，如图5b所示。

图5 种子受力情况及运动状态Fig.5 Stress and movement state of seeds

从图5a中可以看出，种子受到的曳力先急剧增大，然后缓慢减小，这是由于种子在吸孔周边受到气流的作用，离吸孔越近，所受到的曳力越大，被吸附于吸孔上时曳力达到最大值；随后，随着种盘的转动，吸附种子的吸孔离进气口越来越远，从而导致风压降低，因此，曳力不断减小，当吸附种子的吸孔转过气室范围，曳力消失。从图5a中也可明显看出3种种子由于颗粒迎风面积的影响，所受到的曳力大小顺序为大扁形、小扁形、类圆形，即种子越大在相同的气流场受到的力也越大，曳力最大值分别为0.015、0.013、0.006 N。

从图5b中可知，种子的速度先增加，然后趋于稳定，这是由于种子在气流的作用下迅速朝吸孔运动，当吸附稳定后随种盘一起转动，速度为种盘吸孔的线速度；当吸孔转过气室后，吸孔气流消失，种子在自身重力的作用下自由落体，直至掉落在模型外，无法检测到。进一步分析曲线变化可知，3种种子的速度增加快慢顺序为小扁形、类圆形、大扁形。对种子进一步分析有

FD=maa

(1)

vt=v1+aat

(2)

式中FD——曳力，Nm——种子质量，kg

aa——加速度，m/s2

vt——吸附型孔速度，m/s

v1——种子初速度，m/s

t——时间，s

由式(1)、(2)可以看出，种子速度的增加不仅与种子的曳力相关，还与种子自身质量相关。在充种区，速度能迅速提升到种盘线速度的种子容易被吸附[20]。因此，3种种子充种能力大小顺序为小扁形、类圆形、大扁形。为提高排种器作业质量，确保排种器对种子的适应性，需针对大扁形不易吸附的问题，优化设计种盘结构。

3 种盘优化设计

种盘作为与种子直接接触的部件，对充种性能起着至关重要的作用，优化种盘结构参数，可有效提高排种器作业质量[21-25]。

3.1 种子在流场中的受力

充种过程中，玉米种子处在气流场、颗粒场和重力场并存的复杂环境中，其中在气流场中，种子颗粒会同时受到来自气流场的曳力、浮力、压力梯度力、Basset力、Magnus力和Saffman力等多个力的作用[18]。由于实际工作中，种子通常处于低速运动状态，因此根据各力的定义，可以忽略其他力，重点考虑曳力对种子的作用。曳力的表达式为

(3)

式中A——物体在流动方向上的投影面积，m2

CD——绕流阻力系数

ρ——流体密度，g/cm3

U——来流速度，m/s

大扁形种子被吸附的姿态有3种可能性：竖躺、竖直和横躺，如图6所示。

图6 种子被吸附姿态示意图Fig.6 Schematics of seed adsorbed posture

表2 大扁形种子不同吸附姿态参数Tab.2 Different adsorption attitude parameter values of big flat seed

从表2中理论计算数据可以得出，大扁形种子不同吸附姿态受到的曳力均值为0.015 9 N，与气固耦合仿真数据较为接近，证明了仿真模拟的真实性。

3.2 型孔凸台高度设计

型孔凸台作为种盘上的关键结构，能够扰动种群、托持种子，起到辅助充种作用。为了计算型孔凸台高度，选取排种器设计的最高速度14 km/h，种盘型孔线速度为0.25 m/s进行计算。如果其间距满足充种要求，则低于14 km/h时都能满足。

为了便于分析，将大扁形种子受力简化，如图7所示。种子在吸孔周边将受到曳力的作用，如图7a所示，曳力对种子产生的加速度为a1，种子自身的加速度为g，如图7b所示。由牛顿第二定律公式有

FD=ma1

(4)

式中大扁形种子质量m利用精确到0.1 g的电子秤称量，每个重复放置100粒种子，重复3次，得出种子质量均值为m=0.42 g，FD取最小值(1.06×10-2N)进行计算，如果其能满足，则任何姿态都能满足。求得a1=25.2 m/s2。

由加速度合成定理，种子的加速度为

(5)

式中g——重力加速度，取9.81 m/s2

求得a=27.1 m/s2。

图7 大扁形种子受力和加速度分析Fig.7 Force and acceleration analyses of big flat seed

种子如果没有被吸附发生掉落，将从上一型孔凸台边缘掉落至下一吸孔区域，则需满足

(6)

(7)

式中L——玉米最大尺寸(14 mm)

n——型孔凸台厚度

v0——种盘吸孔处线速度

n取1.5 mm，即L+n=15.5 mm，a=27.1 m/s2、a1=25.2 m/s2。v0取0.25 m/s，求得t=0.025 s，h≤0.007 87 mm。

型孔凸台越高，对种子的托持作用越大，但会造成清种越困难，导致重播严重。因此，型孔凸台高度不易过高。由文献[27]玉米总体尺寸可知，大扁形种子厚度均值5.34 mm，选取1～5 mm高度区间，具体参数还需进一步试验分析。

3.3 型孔凸台角度设计

种盘型孔凸台的目的是能够托持种子，使种子具有一定的初速度，起到更容易被吸孔吸附的作用，因此，型孔凸台的角度也将影响对种子的托持效果。

对于种盘上吸孔附近的种子，设其质心在O′点，首先分析充种起始位置，即种子质心O′与种盘圆心O连线与0°位置，其夹角为α=30°。种子从上部充入吸孔时，所受到的力为：重力G=mg、随种盘转动产生的离心力J=mrω2(r为种子质心所在处半径，取70 mm)、沿型孔斜面的摩擦力f=μFN、斜面对种子的支持力FN，如图8所示。

图8 充种区种子受力分析Fig.8 Stress analysis of seed filling area

坐标系O′XY如图8所示，X方向沿着型孔凸台倾斜角的方向。为了使型孔凸台能够更好地托持种子，只要X方向上的加速度aX和Y方向上的加速度aY都为零，型孔凸台就能对种子起到托持作用。将各力投影到坐标系O′XY上，得到α=30°受力平衡方程

mgcos(α+θ)+mrω2=μFN

(8)

mgsin(α+θ)+mrω2=FN

(9)

联立式(8)、(9)得

gcos(30°+θ)+rω2cosθ-μ(gsin(30°+θ)+
rω2sinθ)=0

(10)

式中θ——型孔凸台边缘a1和a2两点间连线与种子离心力J反方向夹角，即型孔凸台角度，(°)

μ——滑动摩擦因数，玉米种子与3D打印材料型孔间的滑动摩擦因数在0.135～0.287之间，取中间值0.21[19]

α——充种开始角度，α=30°

θ的取值还由种盘角速度ω决定，因此需先求得排种盘角速度的变化范围，选取高速作业速度14 km/h时种盘转速作为计算所用角速度，ω=3.62 rad/s。以上各式在Matlab中进行计算，代入各值，求得θ=50°。

为保证在充种区内，不同充种角度下型孔凸台角度都能满足对种子的托持作用，选取充种区另一角度α=60° 进行分析，采用相同的分析方法，得到α=60°受力平衡方程

mgcos(θ-30°)+mrω2sinθ=FN

(11)

mgsin(θ-30°)-mrω2cosθ=μFN

(12)

联立式(11)、(12)得

rω2cosθ-gsin(θ-30°)+μ(gcos(θ-30°)+
rω2sinθ)=0

(13)

代入上面数据，求得θ=45°。

最后，选取最大充种角度α=90°进行分析，得到受力平衡方程

mgsinθ-mrω2cosθ=μFN

(14)

mgcosθ+mrω2sinθ=FN

(15)

联立式(14)、(15)得

μ(gcosθ+rω2sinθ)+rω2cosθ-gsinθ=0

(16)

在Matlab中代入各值，得出θ=17°。

通过对不同充种角度下型孔凸台角度的分析，得出种子在充种区范围内型孔凸台角度最佳范围为17°～50°。但由于在排种盘转动的情况下，种群间的相互作用和种群与排种器之间的相互作用所造成的复杂的运动情况，还需通过试验进一步分析确定最佳角度。

4 试验

4.1 试验条件

选用仿真模型所用的郑单958未分级的种子，自制 9 种不同规格的种盘，通过更换不同的种盘实现排种器结构参数的改变；排种检测装置选用中国农业大学自主研发的排种器性能检测仪，风压测定选取RE-1211型风压计，试验装置如图9所示。

图9 排种试验装置Fig.9 Seed metering test device

4.2 试验方法

根据GB/T 6973—2005《单粒(精密)播种机试验方法》，每组试验采集251粒种子进行统计，每组重复3次，取平均值作为试验结果记录分析，以重播指数、漏播指数、合格指数为排种性能评价指标，设置理论株距为25 cm，设定排种器作业速度为14 km/h，吸种负压为-3 kPa，选取影响气吸式排种器充种关键因素种层高度[28-29]，结合前文分析的型孔凸台角度、型孔凸台高度因素，设置型孔凸台角度为17°～50°，型孔凸台高度1～5 mm，种层高度36～72 mm，进行三因素二次旋转正交组合试验。试验因素编码如表3所示。试验结果如表4所示，X1、X2、X3为因素编码值。采用方差分析法进一步分析，结果如表5所示。

表3 试验因素编码Tab.3 Factors and codes of test

表4 试验设计及结果Tab.4 Experiment design and result

4.3 试验结果分析

根据表4数据，运用数据处理软件Design-Expert 8.0.6对试验数据进行多元回归拟合，建立合格指数、漏播指数、重播指数与三因素的二次多项式回归模型，回归方程为

表5 方差分析Tab.5 Variance analysis

注：*表示差异显著(0.01

(17)

(18)

(19)

回归模型方差分析和显著性检验结果如表5所示。由表5可知，合格指数、漏播指数、重播指数模型的拟合度极显著(P<0.01)。且回归方程失拟不显著，与实际情况拟合较好。

对于合格指数回归方程，回归项X1X3的P值大于0.05，影响不显著，其它各项的影响均显著或极显著。失拟项P=0.584 5，说明不存在其它影响响应指标的主要因素。对于漏播指数回归方程，回归项X1、X3、X1X2、X1X3和X2X3的P值大于0.05，影响不显著，其他各项的影响都是极显著。失拟项P=0.884 3，说明不存在其它影响响应指标的主要因素。对于重播指数回归方程，回归项X1和X1X2的P值大于0.05，影响均不显著，其它各项的影响均显著或极显著。失拟项P=0.132 7，说明不存在其它影响响应指标的主要因素。剔除不显著的回归项，在保证回归模型显著、失拟项不显著的前提下，对合格指数、漏播指数和重播指数的回归方程进行重新拟合，回归方程为

(20)

(21)

(22)

通过对式(20)～(22)回归系数的检验得出，影响合格指数和重播指数的主次因素为:型孔凸台高度、种层高度、型孔凸台角度。影响漏播指数的主次因素为：型孔凸台高度、型孔凸台角度、种层高度。

4.4 各因素对排种合格指数的影响

通过Design-Expert 8.0.6对数据进行处理，可得到型孔凸台角度、型孔凸台高度和种层高度对合格指数的影响，其响应曲面如图10所示。分别固定种层高度、型孔凸台高度和型孔凸台角度因素为零水平，根据响应曲面图，分析其余2个因素间的交互作用对排种合格指数的影响。

图10 因素交互作用对合格指数的影响Fig.10 Impacts of interaction on qualified index

图10a为种层高度54 mm时，型孔凸台高度和型孔凸台角度对排种合格指数交互作用的响应曲面图。由图10a可知，型孔凸台高度一定时，随着型孔凸台角度的增大，排种合格指数先上升后下降。型孔凸台角度一定时，随着型孔凸台高度的增大，排种合格指数同样先上升后下降。型孔凸台高度为3.3 mm，型孔凸台角度为36°时，排种合格指数最高。

图10b是型孔凸台高度为3 mm时，型孔凸台角度和种层高度对排种合格指数交互作用的响应曲面图。由图10b可知，在型孔凸台角度为30°～40°，种层高度为50～60 mm时，排种合格指数较高。型孔凸台角度一定时，随着种层高度的增大，排种合格指数先上升后下降。种层高度一定时，随着型孔凸台角度的增大，排种合格指数同样先上升后下降。

图10c为型孔凸台角度为33.5°时，型孔凸台高度和种层高度对排种合格指数交互作用的响应曲面图。由图10c可知，在型孔凸台高度为2.6～3.4 mm，种层高度为50～60 mm时，排种合格指数较高。型孔凸台高度一定时，随着种层高度的增大，排种合格指数先上升后下降。种层高度一定时，随着型孔凸台高度的增大，排种合格指数同样先上升后下降。

由上述内容可知，各因素以及各因素之间的交互作用对排种器工作性能影响较大，为了获得最佳排种作业参数，优化型孔凸台高度、型孔凸台角度和种层高度，实现因素间参数的合理匹配是提高合格指数的关键[30-31]。

以合格指数为最终优化目标，根据JB/T 10293—2001《单粒(精密)播种机技术条件》中的要求，播种合格指数大于等于80.0%，重播指数小于等于15.0%，漏播指数小于等于8.0%。结合各因素边界条件，建立参数化数学模型，进行优化求解，其目标函数和约束条件为

Fmax=Y1-Y2-Y3

(23)

(24)

Fmax为最终优化的目标；运用Design-Expert 8.0.6软件的优化模块，对约束目标优化求解，得到型孔凸台角度为35.76°，型孔凸台高度为3.11 mm，种层高度为55.61 mm时，排种合格指数最高。

根据优化得到的最优参数，进行 5 次重复台架试验，并选取文献[27]设计的排种器作对比，试验结果平均值见表6。可见在最优参数作用下，实际结果与理论结果近似，合格指数和漏播指数均优于原排种器，可以实现有效播种。

表6 试验验证结果Tab.6 Test verification result %

5 结论

(1)针对气吸式排种器充种过程中种子运动十分复杂，难以进行准确分析计算的问题，通过建立玉米种子Bonding模型和划分流场网格，实现了DEM-CFD气固耦合仿真，并通过提取仿真过程中种子所受曳力和吸附速度运动参数，分析得出3种种子充种能力大小依次为小扁形、类圆形、大扁形。

(2)采用理论计算的方式分析了不同姿态大扁形种子在流场中的受力，验证了气固耦合仿真的真实性，建立了大扁形种子吸附充种过程数学模型，得到种盘型孔凸台高度和型孔凸台角度最佳参数范围。

(3)选取型孔凸台角度、型孔凸台高度和种层高度为试验因素，结合理论计算的各因素的取值变化范围，以排种合格指数、重播指数和漏播指数为评价指标，采用二次旋转正交组合试验方法进行台架试验，通过建立排种性能指标与各影响因素之间的回归模型，得出了影响各指标的因素主次顺序，并运用响应面法分析了试验因素对响应指标的影响。采用多目标优化方法，确定了排种器最佳参数组合: 型孔凸台角度为35.76°，型孔凸台高度为3.11 mm，种层高度为55.61 mm，排种合格指数最高。此时，排种器性能指标为: 合格指数 91.60%，漏播指数3.90%，重播指数4.50%。对优化结果进行验证试验，并与原排种器进行对比，验证结果与优化结果基本一致，且合格指数和漏播指数均优于原排种器。