附件3.科技报告正文填写模板.doc
附件3.科技报告正文填写模板.doc
目录 1 研究概述 ...............................................................1 1.1 研究背景和意义 ......................................................1 1.2 总体目标 ............................................................3 1.3 考核指标 ............................................................4 2 边界建模方法研究 .......................................................4 3 颗粒建模方法研究 .......................................................7 4 接触作用力学模型研究 ..................................................13 5 接触力学模型中参数的研究...............................................14 6 三维离散元法算法研究和软件设计.........................................26 6.1 颗粒入口的研究 .....................................................26 6.2 接触叠合量求解算法研究 .............................................27 6.3 接触相对位移和转角的求解算法研究 ...................................31 6.4 接触作用力求解算法研究 .............................................35 6.5 颗粒运动方程求解 ...................................................36 6.6 机械部件性能评价算法研究 ...........................................37 6.7 动态显示等其它算法研究 .............................................38 7 新型 CAE 软件的应用研究 ................................................39 8 结论 ..................................................................54 附录 A 发表论文情况 ......................................................57 I 插图清单 图 1 三种离散元法边界建模方法简图 .........................................2 图 2 基于离散元法的新型 CAE 软件体系结构 ...................................3 图 3 基于 CAD 模型的三维离散元法边界建模方法 ...............................5 图 4 球颗粒与带有虚圆形平面边界的接触作用仿真 .............................6 图 5 由组合内窝孔精密排种器的二维 CAD 模型建立的该排种器三维离散元法分析模型7 图 6 由组合内窝孔精密排种器的三维 CAD 模型建立的该排种器三维离散元法分析模型7 图 7 由超二次方程建立的颗粒模型 ...........................................7 图 8 由球颗粒组合方法建立的颗粒模型 .......................................8 图 9 自制的静摩擦系数测试装置 ............................................20 图 10 自制的碰撞恢复系数测试系统 .........................................21 图 11 散粒物料由矩形进料口进入过程简图 ...................................27 图 12 求解两颗粒接触点的几何势方法示意图 .................................27 图 13 坐标系建立及球颗粒空间位置简图 .....................................29 图 14 球颗粒与圆形平面接触判定简图 .......................................31 图 15 球体相对位移合转角计算简图 .........................................31 图 16 球体与平面边界的接触 ...............................................32 图 17 接触边界速度计算简图 ...............................................33 图 18 中心差分格式时步之间的关系 .........................................37 图 19 单颗粒作用力变化曲线 ...............................................37 图 20 单颗粒速度变化曲线 .................................................38 图 21 颗粒群体流量变化曲线 ...............................................38 图 22 机械部件工作阻力变化曲线 ...........................................38 图 23 大豆种子跌落碰撞过程及弹起高度的仿真和试验结果对比 .................39 图 24 多种属性颗粒和多种属性边界的分析计算 ...............................39 图 25 具有弹性边界的排种器工作过程仿真分析 ...............................39 图 26 由 CAD 模型得到的组合内窝孔精密排种器播种大豆种子过程仿真分析与实验对比 .........................................................................40 II 图 27 由 CAD 模型得到的组合内窝孔精密排种器播种玉米种子过程仿真分析与实验对比 .........................................................................40 图 28 由 CAD 模型得到的组合内窝孔精密排种器播种小麦种子过程仿真分析与实验对比 .........................................................................40 图 29 由 CAD 模型得到的垂直圆盘型孔轮式精密排种器播种大豆种子过程仿真分析与实 验对比 ...................................................................40 图 30 由 CAD 模型得到的垂直圆盘型孔轮式精密排种器播种玉米种子过程仿真分析与实 验对比 ...................................................................41 图 31 由 CAD 模型得到的窝眼轮式精密排种器播种大豆种子过程仿真分析与实验对比41 图 32 由 CAD 模型得到的窝眼轮式精密排种器播种玉米种子过程仿真分析与实验对比41 图 33 由 CAD 模型得到的槽轮式精密排种器播种小麦种子过程仿真分析与实验对比 .41 图 34 由 CAD 模型得到的第一种外槽轮式排肥器工作过程仿真分析与实验对比 .....41 图 35 由 CAD 模型得到的第二种外槽轮式排肥器工作过程仿真分析与实验对比 .....42 图 36 由 CAD 模型得到的第三种外槽轮式排肥器工作过程仿真分析与实验对比 .....42 图 37 由 CAD 模型得到的施肥开沟器工作过程仿真分析与实验对比 ...............42 图 38 由 CAD 模型得到的船铲开沟器工作过程仿真分析与实验对比 ...............42 图 39 由 CAD 模型得到的芯铧开沟器工作过程仿真分析 .........................42 图 40 组合内窝孔精密排种器充种角、清种角和投种角的示意图 .................45 图 41 组合内窝孔排种器和型孔轮式排种器的种子面倾角示意图 .................45 图 42 组合内窝孔精密排种器内特殊点种子位置示意图 .........................46 图 43 窝眼轮式排种器内特殊点种子位置示意图 ...............................47 III 附表清单 表 1 三种大豆种子的三维尺寸和千粒重 .......................................8 表 2 三种大豆种子颗粒的密度(G/CM3) .........................................8 表 3 三种大豆原种的含水率 .................................................9 表 4 三种玉米种子的三维尺寸 ...............................................9 表 5 三种玉米种子的千粒重(G)...............................................9 表 6 三种玉米原种的含水率(%) .............................................10 表 7 三种玉米种子颗粒的密度(G/CM3) ........................................10 表 8 三种小麦种子的三维尺寸、等效直径和球形率 ............................10 表 9 三种小麦种子原种的含水率(%) .........................................10 表 10 三种小麦种子的千粒重(G).............................................11 表 11 三种小麦种子颗粒的密度(KG/M3) .......................................11 表 12 三种肥料原始样品含水率的测量结果 ...................................11 表 13 三种肥料的三维尺寸和千粒重 .........................................12 表 14 三种肥料的颗粒密度(G/CM3) ...........................................12 表 15 实验用壤土的粒级分布(%) ............................................12 表 16 实验用壤土的含水率 .................................................12 表 17 实验用壤土的密度 ...................................................13 表 18 不同加载速度时吉科豆原种的刚度系数和破碎力 .........................14 表 19 加载速度为 1MM/MIN 时三种大豆种子的刚度系数和破碎力...................14 表 20 不同含水率时三种大豆种子的弹性模量(MPA).............................15 表 21 三种玉米种子的刚度系数(N/MM) ........................................15 表 22 同一含水率时三种玉米种子的弹性模量(MPA).............................16 表 23 不同含水率时三种小麦种子的弹性模量(MPA).............................16 表 24 加载速度为 1MM/MIN 时三种小麦种子的刚度系数(N/MM) .....................16 表 25 加载速度为 6MM/MIN 时三种肥料颗粒的刚度系数和破碎力...................16 表 26 不同含水率时三种肥料颗粒的弹性模量(MPA).............................17 表 27 离散元模拟双轴试验得到的三种土壤微观参数 ...........................17 IV 表 28 三种土壤样本的弹性模量 .............................................17 表 29 大豆原种与不同摩擦表面之间的动摩擦系数 .............................18 表 30 压力和速度一定时大豆种子与不同摩擦表面间的动摩擦系数 ...............18 表 31 压力和速度一定时不同含水率下玉米种子的动摩擦系数 ...................19 表 32 压力和速度一定时在不同含水率小麦种子与不同材料表面间的动摩擦系数 ...19 表 33 原含水率时三种肥料颗粒与不同摩擦表面之间的动摩擦系数 ...............20 表 34 不同含水率时大豆种子与不同材料表面间的静摩擦系数 ...................21 表 35 不同含水率下玉米种子的静摩擦系数 ...................................21 表 36 不同含水率时小麦种子与不同材料表面间的静摩擦系数 ...................22 表 37 在不同含水率时三种肥料与不同材料表面间的静摩擦系数 .................22 表 38 吉科豆在不同含水率和不同高度时的碰撞恢复系数 .......................23 表 39 吉新豆在不同含水率和不同高度时的碰撞恢复系数 .......................24 表 40 吉豆在不同含水率和不同高度时的碰撞恢复系数 .........................24 表 41 三种玉米种子与铸铁的碰撞恢复系数 ...................................25 表 42 三种小麦种子与铸铁的碰撞恢复系数 ...................................25 表 43 北京生物肥料在不同含水率、不同高度、不同碰撞表面时的碰撞恢复系数 ...26 表 44 土壤颗粒之间的摩擦系数及颗粒与边界之间的摩擦系数 ...................26 表 45 三维仿真模拟时大豆种子的参数选取 ...................................43 表 46 开沟过程仿真模拟时三种土壤微观参数选取 .............................43 表 47 土壤颗粒之间的摩擦系数及颗粒与边界之间的摩擦系数 ...................43 表 48 三维仿真模拟时玉米种子的参数选取 ...................................43 表 49 三维仿真模拟时小麦种子的参数选取 ...................................44 表 50 三维仿真模拟时肥料颗粒的参数选取 ...................................44 表 51 组合内窝孔精密排种器播大豆种子时粒数的仿真与试验对比 ...............47 表 52 组合内窝孔精密排种器播玉米种子时粒数的仿真与试验对比 ...............47 表 53 组合内窝孔精密排种器播小麦种子时粒数的仿真与试验对比 ...............48 表 54 组合内窝孔精密排种器播大豆种子时清种角、充种角和投种角的仿真与试验对比 .........................................................................48 表 55 组合内窝孔精密排种器播玉米种子时清种角的仿真与试验对比 .............48 V 表 56 组合内窝孔精密排种器播小麦种子时清种角、投种角和种子面倾角仿真与试验对 比 .......................................................................48 表 57 组合内窝孔精密排种器播大豆时单粒种子合速度模拟值与试验值对比(M/S) ...49 表 58 组合内窝孔精密排种器播玉米时单粒种子合速度模拟值与试验值对比(M/S) ...49 表 59 组合内窝孔精密排种器播小麦时单粒种子合速度模拟值与试验值对比(M/S) ...49 表 60 型孔轮式精密排种器播大豆种子时粒数的仿真与试验对比 .................49 表 61 型孔轮式精密排种器播玉米种子时粒数的仿真与试验对比 .................50 表 62 型孔轮式精密排种器播大豆种子时种子面倾斜角的仿真与试验对比 .........50 表 63 型孔轮式精密排种器播玉米种子时种子面倾斜角的仿真与试验对比 .........50 表 64 窝眼轮式精密排种器播玉米种子时粒数的仿真与试验对比 .................50 表 65 窝眼轮式精密排种器播大豆时单粒种子合速度仿真与试验对比 .............51 表 66 槽轮排种器播小麦种子时种子流量的仿真值与试验值比较 .................51 表 67 第一种槽轮排肥器仿真值与试验值比较 .................................51 表 68 第二种槽轮排肥器仿真值与试验值比较 .................................51 表 69 第三种槽轮排肥器仿真值与试验值比较 .................................52 表 70 含水率为 22%、开沟深度为 30MM 时芯铧开沟器工作阻力的二维模拟值与试验值对 比 .......................................................................52 表 71 含水率为 22%、工作速度为 0.2M/S 芯铧开沟器工作阻力的二维模拟值与试验值对 比 .......................................................................52 表 72 开沟深度 50MM、工作速度 0.2M/S 时施肥开沟器工作阻力的三维模拟与试验相比52 表 73 刚度系数变化时组合内窝孔精密排种器播小麦的模拟值与试验值对比(M/S) ...53 表 74 刚度系数变化时槽轮排种器播小麦的模拟值与试验值对比(M/S) .............53 表 75 已毕业硕士研究生名细 ...............................................59 VI 1 研究概述 1.1 研究背景和意义 在自然界和工农业生产领域，大量存在着散粒物料与机械部件的接触作用和散粒物 料的流动过程。在农业生产领域，例如耕地、开沟、播种、施肥、脱粒、分离、清选、 粉碎、干燥、输送、仓储、分级、加工和包装等过程中，始终存在着散粒物料与机械部 件的接触作用和散粒物料的流动过程。在众多工业生产领域，例如制药、食品、化工、 水泥、材料、矿山、能源和岩土工程等领域，也大量存在着散粒物料与机械部件的接触 作用和散粒物料的流动过程。 散粒物料的流动除与颗粒本身的形状、尺寸和物理力学性质有关外，还与相关机械 部件的结构形式、尺寸参数和运动参数等有关。一个好的机械设计，应使散粒物料按照 预期的方式流动，减少流动过程中不必要的损伤，节省机械动力消耗，此时必须考虑散 粒物料与机械部件的接触作用和颗粒群体的动力学问题。 散粒物料的性质介于固体与流体之间，虽然颗粒结构简单，却有复杂的物理力学性 质。二十世纪 70 年代后，许多物理学家、力学家和应用数学家开始对颗粒运动的物理 机制发生兴趣并进行了深入研究，建立了两类颗粒动力学理论：①基于连续介质力学的 理论，如颗粒动理论、摩擦塑性模型和光滑粒子法等；②基于离散介质力学的理论，如 软颗粒模型（简称离散元法） 、硬颗粒模型和 Monte Carlo 方法等。 以往采用基于连续介质力学的理论，研究散粒物料与机械部件的接触作用和颗粒群 体动力学问题，只能把颗粒群体作为一个整体来考虑，不能考虑颗粒群体中颗粒间的接 触作用及其对颗粒运动的影响，也无法分析颗粒群体中颗粒之间时而接触时而又分离的 不连续性，因此不能很好地解决该问题。目前进行相关机械设计时，大都依靠设计者的 经验和试验方法，既费时费力又得不到理想的效果。据估计仅由散粒物料输送所造成的 相关设备利用损失就达 40%，远未达到优化设计和节省能源的要求。 随着计算机技术的发展，基于离散介质力学的理论，愈来愈引起人们的重视。其中 Cundall 提出的离散元法，是计算密相颗粒群体力学行为的一种数值方法。该方法把散 粒物料简化成具有一定形状和质量颗粒的集合，赋予接触颗粒间及颗粒与接触边界（机 械部件）间某种接触力学模型和模型中的参数，以考虑颗粒之间及颗粒与边界之间的接 触作用和颗粒与边界的不同物理力学性质。离散元法采用动态松弛法、牛顿第二定律和 时步迭代求解每个颗粒的运动速度和位移，因而特别适合于求解非线性问题；当采用不 1 同的接触模型时，还可以分析颗粒结块、颗粒群聚合体的破坏过程、多相流动甚至可以 包括化学反应和传热等问题。正是由于诸多优点，使得离散元法已成为研究散粒物料与 边界相互作用和颗粒群体动力学问题的一种通用方法，并在岩土工程、采矿工程、化工 过程、制药工程、食品工程和农业工程等领域得到较多应用。 在采用离散元法分析散粒物料与机械部件的接触作用时，需要建立机械部件（边界） 的离散元法分析模型。目前已报道三种离散元法边界建模方法，参见图 1。一是函数建 模，即对于比较简单的边界如一个圆柱面，可以采用圆柱面方程建立边界的离散元法分 析模型。二是排列颗粒的方法，即把一定尺寸的球颗粒按照边界形状排列建立边界模型， 一般也用于简单边界的建模。第三种称“有限壁（Finite Wall Methods）”方法，该方法 把与颗粒接触的边界表面离散成小三角形平面的组合，是复杂边界建模常用的方法。但 “有限壁”方法是一种近似方法：①求解颗粒与边界的接触点（接触力作用点）和接触 叠合量（用于计算接触作用力）是近似的；②在小三角形平面的连接处一阶导数不连续， 致使所求颗粒与边界的法向和切向接触作用力具有突变性；③在采用赫芝模型求解接触 作用力时，对于曲面边界也只能取无穷大的曲率半径，这些均与实际的边界情况差别较 大。为了提高建模精度，往往还需要较多数量的三角形平面，由此增加了离散元法的计 算时间。国外 2 个著名的离散元法商品软件 PFC 和 EDEM，即采用函数建模和“有限 壁”的边界建模方法。而建立一种高效、精确和适应性广的离散元法边界建模方法，已 成为离散元法实际应用需要解决的关键问题之一。 图 1 三种离散元法边界建模方法简图 为此，本研究组提出由机械部件的 CAD 模型（CAD 软件设计图）建立机械部件的 离散元法分析模型的方法和实边界与虚边界的方法，在此基础上实现 CAD 软件与自主 研制的离散元法分析计算软件集成，从而构建一种集设计与性能分析评价为一体的新型 CAE 软件，该软件的体系结构见图 2 所示。从已经掌握的文献资料看，本研究提出的边 界建模方法和构建的 CAE 软件尚未见他人报道。 本研究提出的边界建模方法和构建的 CAE 软件，把边界设计与其离散元法性能分析 2 有机结合，因而其特点是：①在设计阶段通过修改机械部件的 CAD 模型，能对不同结 构和尺寸的机械部件进行性能分析和评价，由此实现机械部件结构方案和尺寸参数的优 化；②由机械部件的 CAD 模型能进行机械部件工作过程的动态仿真，由此分析机械部 件的工作原理和工作过程，从而可以发明新原理和新结构的机械部件，这是现有与散粒 物料接触作用的机械部件的设计方法和设计技术不能实现的。 本研究不仅可以丰富离散元法，具有科学意义，而且对于建立相关机械部件优化设 计和虚拟设计的新方法、减少相关机械部件试验环节和试验次数、缩短相关机械部件新 产品的研发周期、降低相关机械部件工作过程中的损失和节省能源消耗等，均有较大意 义。例如，我国每年生产粮食 5 亿吨，每年耕翻处理土壤上千亿吨，约四分之三的化工 原料和二分之一的化工产品都是散粒物料等等。如果能提高 1%的产出率或降低 1%的能 源消耗，其经济效益和社会效益也是巨大的。 农业是我国的基础产业，农业生产好坏关系到国计民生和国家安全。作为农业生产 装备的农业机械，对于我国农业的生产状况、农业现代化进程和增加农民收入等都是至 关重要的。但我国农机产品的研究和设计水平还相对落后，导致许多农业生产中急需的 农机产品或没有或性能较差。利用信息技术提升我国农机产品的研究和设计水平，使其 实现跨越式发展，是解决该问题的好时机，这也是本研究的意义所在。 图 2 基于离散元法的新型 CAE 软件体系结构 1.2 总体目标 各选典型颗粒农业物料玉米种子、大豆种子、小麦种子、颗粒肥料和土壤（每种各 选三种）及其相关农机工作部件排种器、排肥器和开沟器（每种各选三种）为对象，研 3 究由农机工作部件的 CAD 模型生成其三维离散元法分析模型的方法，研究散粒物料的 颗粒之间及颗粒与农机工作部件之间接触作用的力学模型及模型中参数的选用方法，在 此基础上实现 CAD 软件与自主研制的三维离散元法分析软件集成，从而构建一种农机 工作部件的通用数字化设计方法及其软件，为农业机械工作部件的优化设计和虚拟设计 提供一种新方法。 1.3 考核指标 ①数字化设计方法及其软件在微机环境下运行，具有工作部件分析模型自动生成及 颗粒模型、力学模型、分析计算参数、性能评价指标等自动选择添加和工作部件工作过 程动态显示等功能； ②采用该方法设计出 2~3 种新型农机工作部件，在 1~2 家研究和设计单位应用该方法 和软件； ③新方法及其软件的相关技术资料。 2 边界建模方法研究 为构建该新型 CAE 软件，如何由农机工作部件的 CAD 模型建立其三维离散元法分 析模型，是需要研究解决的首要问题。为此，研究组对散粒物料与农机工作部件的接触 作用进行了分析。一般来说农机工作部件的结构是比较复杂的，但与散粒物料接触作用 的是农机工作部件中某些零部件的表面，这些表面又是由一些基本面元素组成，如平面、 球面、圆柱面、圆锥面、圆环面、螺旋面等完整图元和扇形平面、弓形平面、球冠面、 球台面、圆台面等非完整图元。通过读取农机工作部件 CAD 模型中与散粒物料接触作 用的基本面元素（图 2 中称图元），来建立农机工作部件的离散元法分析模型，是一种 可行和通用的方法。 边界建模方案①：考虑到目前二维 CAD 软件应用广泛和使用简便，本研究采用由 农机工作部件的二维 CAD 模型，通过交互式模型引导的方法，重建农机工作部件中与 散粒物料接触作用的图元及参数，然后通过添加运动属性和材料特性参数实现图元的几 何信息、运动信息和材料信息绑定，最后把图元信息（图元类型+几何信息+运动属性+ 材料特性）存储到数据库中，以供离散元法计算和计算结果仿真显示使用，建模方法和 步骤参见图 3。其中设置运动属性和材料特性参数，是将不同运动方式和不同物理机械 性质的图元给定不同特征值，以便在离散元法计算和计算结果仿真显示时采用不同的处 理方式。 4 图 3 基于 CAD 模型的三维离散元法边界建模方法 边界建模方案②：由三维 CAD 软件 PRO/E，建立农机工作部件的三维 CAD 模型， 采用 CAD 软件本身提供的接口和二次开发，通过人机交互读取 CAD 模型中与散粒物料 接触作用的图元及参数，通过人机交互添加图元的运动属性和材料特性，最后把图元信 息（图元类型+几何信息+运动属性+材料特性）存储到数据库中。 在上述边界建模方案中，存在的一个问题是有些零部件表面不是完整的图元，例如 一个平面上钻有圆孔或一个柱面上开有沟槽等。为了在建模时除去这些缺失部分（图 2 中称冗余） ，本研究组提出实边界和虚边界的方法。实边界表示农机工作部件中实际存 在的面，这些面与散粒物料接触时有接触力产生。虚边界表示实边界上缺失部分，散粒 物料与虚边界接触时没有力产生，因而可以穿过虚边界（参见图 4）。分析可知，由于制 造工艺的限制，农机工作部件表面的缺失部分，往往也可以由图元表示。分析还可知， 定义虚实二种边界的方法，还可以解决不同图元间的相贯问题。 5 图 4 球颗粒与带有虚圆形平面边界的接触作用仿真 采用上述方法建立边界模型并存储到数据库后，当用户分析该机械部件的工作情况 时，只需采用离散元法软件打开该机械部件的分析模型即可。 根据上述边界建模方法，本研究组初步建立了平面（包括圆形平面、椭圆形平面、 扇形平面、弓形平面、多边形平面）、球面（包括完整球面、球冠面、球台面） 、圆柱面 （包括完整和非完整圆柱面） 、圆锥面（包括完整和非完整圆锥面） 、圆环面（包括完整 和非完整圆环面）等，十余种图元的读取和离散元法计算方法。采用这些图元，已能建 立较多机械部件的离散元法分析模型。当一个机械部件包含其它种类图元时，只要按照 上述方法建立新图元的读取方法和计算方法，即可由该机械部件的 CAD 模型建立其离 散元法分析模型。当然对于复杂曲面边界，也可以采用上述图元组合的方法，来建立该 边界的分析模型。因此本研究提出的边界建模方法具有通用性，与“有限壁”边界建模 方法相比，具有更大的灵活性、更高的计算精度和更短的计算时间。 在对二维 CAD 软件 AutoCAD 和三维 CAD 软件 PRO/E 进行二次开发的基础上，通 过人机交互实现了与散粒物料接触作用的图元识别、参数读取、运动属性及材料特性参 数绑定，还通过数据库共享实现了 CAD 软件与自主研制的三维离散元法分析计算 软件集成，从而研制出基于 CAD 模型的三维离散元法边界建模软件。采用上述方法和 软件分别由机械部件的二维 CAD 模型和三维 CAD 模型，建立的机械部件三维离散元法 分析模型见图 5 和图 6 所示。 6 图 5 由组合内窝孔精密排种器的二维 CAD 模型建立的该排种器三维离散元法分析模型 图 6 由组合内窝孔精密排种器的三维 CAD 模型建立的该排种器三维离散元法分析模型 3 颗粒建模方法研究 目前国外大多数三维离散元法软件均采用球颗粒模型，因该颗粒模型接触检测算法 简单。对于非球形的复杂形状颗粒，一般采用 2 种建模方法。一是由超二次方程 （|x/a|n1+|y/b|n2+|z/c|n3-1=0）建立颗粒模型；另一种是采用球颗粒组合的方法。 为能适用于多种散粒物料的分析，本研究分别采用超二次方程和球颗粒组合的方法， 建立了不同颗粒模型见图 7 和图 8 所示。其中玉米种子和小麦种子的离散元法分析模型 的研究，尚未见他人报道。 图 7 由超二次方程建立的颗粒模型 7 图 8 由球颗粒组合方法建立的颗粒模型 本研究分别采用烘箱干燥法、排水法和游标卡尺及筛分等方法，测试分析了 3 种大 豆种子（吉科豆 1 号、吉新豆 1 号和吉豆）、3 种玉米种子（吉峰 218、长单 228 和平安 18） 、3 种小麦种子（吉 12、丰强 11 和小冰麦 33）、3 种颗粒肥料（北京生物、大庆 3 和山东 4）和 3 种土壤（含水率为 3%、12%和 22%）颗粒的几何形状、尺寸分布、千粒 重和颗粒的密度，测试结果见表 1~表 17 所示。通过加湿的办法还得到了不同的含水率。 表 1 三种大豆种子的三维尺寸和千粒重 表 2 三种大豆种子颗粒的密度(g/cm3) 8 表 3 三种大豆原种的含水率 表 4 三种玉米种子的三维尺寸 三维尺寸 长(mm) 宽(mm) 厚(mm) 等效直径(mm) 球形率 吉峰 218 9.57 7.55 4.87 7.04 0.74 长单 228 9.58 8.74 5.95 7.89 0.83 平安 18 9.84 7.99 5.40 7.48 0.77 吉峰 218 0.865 0.617 0.90 0.42 0.08 长单 228 0.688 0.59 0.70 0.48 0.05 平安 18 0.687 0.700 1.060 0.58 0.07 平均值 (mm) 标准 差 表 5 三种玉米种子的千粒重(g) 第一次平均值 第二次平均值 第三次平均值 三次平均值 吉峰 218 214.57 217.76 213.61 215.31 长单 228 296.55 298.64 299.24 298.14 平安 18 266.17 265.92 267.89 266.66 9 表 6 三种玉米原种的含水率(%) 编号 吉峰 218 长单 228 平安 18 1 12.94 12.77 13.36 2 12.90 12.96 13.31 3 12.82 12.92 13.45 4 12.74 12.89 13.49 平均含水率 12.85 12.96 13.40 表 7 三种玉米种子颗粒的密度(g/cm3) 编号 吉峰 218 长单 228 平安 18 1 1.174 1.260 1.270 2 1.166 1.255 1.269 3 1.180 1.256 1.276 4 1.163 1.254 1.282 5 1.176 1.258 1.284 平均密度 1.172 1.257 1.276 表 8 三种小麦种子的三维尺寸、等效直径和球形率 长(mm) 宽(mm) 厚(mm) 等效直径(mm) 球形率 平均值 7.04 3.24 2.96 4.07 0.58 标准差 0.34 0.18 0.17 0.17 0.01 平均值 7.0 3.19 2.86 3.99 0.57 标准差 0.36 0.21 0.19 0.21 0.02 平均值 7.06 3.39 2.96 4.14 0.59 标准差 0.36 0.27 0.19 0.23 0.02 小冰麦 33 丰强 11 吉 12 表 9 三种小麦种子原种的含水率(%) 试样 1 的 试样 2 的 试样 3 的 试样 4 的 含水率 含水率 含水率 含水率 含水率 均值 吉 12 12.41 12.39 12.37 12.38 12.39 丰强 11 12.27 12.24 12.38 12.37 12.31 10 12.32 小冰麦 33 12.42 12.32 12.39 12.36 表 10 三种小麦种子的千粒重(g) 第一次 第二次 第三次 第四次 第五次 均值 44.0525 45.2751 46.8912 47.105 47.8659 46.238 丰强 11 38.6912 39.431 38.7972 40.315 39.0483 39.257 吉 12 43.6236 43.0656 43.1064 43.103 43.5603 43.292 小冰麦 33 表 11 三种小麦种子颗粒的密度(kg/m3) 第一次 第二次 第三次 第四次 第五次 均值 小冰麦 33 1385.299 1363.708 1359.165 1389.528 1367.597 1373.059 丰强 11 1338.796 1355.017 1342.464 1348.328 1355.844 1348.09 吉 12 1325.945 1341.607 1334.563 1342.773 1344.454 1337.869 表 12 三种肥料原始样品含水率的测量结果 肥料品种 样品 + 样品盒重量(g) 样品盒重量(g) 二次测量误差 含水率 64.7310 第二次称重(g) 63.5468 0.10% 2.99% 0.12% 2.96% 63.5011 0.11% 2.97% 烘干前称重(g) 烘干冷却后 烘干冷却后 1 27.1448 65.9427 北京 2 27.1322 64.5418 64.7429 第一次称重(g) 63.5645 生物 3 27.1556 64.5152 63.5121 2.97% 平均 1 26.5801 73.7983 73.6346 73.6245 0.10% 2.64% 大庆 2 26.3430 72.8008 71.7650 71.7521 0.11% 2.62% 3号 3 26.5702 72.4772 72.2317 72.1183 0.11% 2.63% 2.63% 平均 山东 4号 1 28.0972 66.3624 65.3569 65.2317 0.12% 1.93% 2 28.1326 65.3819 64.2834 64.1356 0.15% 1.89% 3 27.9624 65.2378 64.3582 64.2168 0.14% 1.91% 1.91% 平均 11 表 13 三种肥料的三维尺寸和千粒重 长(mm) 宽(mm) 厚(mm) 等效直径(mm) 球形率 千粒重(g) 平均 标准 平均 标准 平均 标准 平均 标准 平均 标准 平均 标准 值 差 值 差 值 差 值 差 值 差 值 差 北京生物 4.21 0.56 4.35 0.53 4.25 0.52 4.11 0.45 0.97 0.036 182.19 1.29 大庆 3 号 4.41 0.54 4.45 0.51 4.36 0.48 4.36 0.36 0.98 0.037 203.64 1.21 山东 4 号 4.45 0.57 4.53 0.37 4.55 0.39 4.31 0.37 0.97 0.035 171.96 1.20 品种 表 14 三种肥料的颗粒密度(g/cm3) 密度 1 2 3 平均 标准差 北京生物 1.970 1.978 1.969 1.973 0.049 大庆 3 号 1.775 1.774 1.771 1.772 0. 021 山东 4 号 1.288 1.368 1.361 1.363 0. 044 品种 表 15 实验用壤土的粒级分布(%) ≥2 2~0.9 0.9~0.3 0.3~0.05 0.05~0.074 ≤0.074 第一次 20.56 67.22 389.53 205.36 25.36 69.67 第二次 18.73 66.38 391.88 203.87 24.37 72.47 第三次 21.46 65.68 391.74 206.35 27.64 64.83 平均数 20.25 66.43 391.05 205.20 25.79 68.99 百分比(%) 2.6 8.54 50.28 26.39 3.16 8.87 粒径(mm) 重量(g) 表 16 实验用壤土的含水率 编号 湿土重 (g) 干土重 (g) 水分重 (g) 含水率 (%) 1 80.80 78.74 2.06 2.616 2 68.42 66.67 1.75 2.625 3 69.70 67.88 1.82 2.681 4 72.46 70.62 1.84 2.605 5 68.09 66.41 1.68 2.530 6 76.85 74.99 1.86 2.480 2.590 平均 12 表 17 实验用壤土的密度 编号 土样重 (g) 体积 (cm3) 密度 (g/cm3) 1 79.88 59.692 1.338 2 78.87 59.692 1.321 3 80.20 59.692 1.344 4 80.86 59.692 1.355 5 79.91 59.692 1.339 1.339 平均 由测试结果可知，种子和肥料颗粒的三维尺寸分布，近似于正态分布。通过分析还 可知，大豆种子颗粒近似于椭球形，但球形率较高，可用球或椭球颗粒模型模拟；玉米 种子是非对称的多面体形状，可用球颗粒组合的方法建立其颗粒模型；小麦种子近似于 椭球形，可用椭球颗粒模型模拟；肥料颗粒近似于椭球形，但球形率较高，可用球或椭 球颗粒模型模拟；土壤颗粒形状比较复杂，但为简化通常也可用球或椭球颗粒模型模拟。 本研究所建立的几种典型颗粒模型如图 7 和图 8 所示。 4 接触作用力学模型研究 研究组研究了散粒物料颗粒之间及颗粒与农机工作部件之间接触作用的形式，完成 了常用 7 种力学模型（线性粘弹性模型、非线性粘弹性模型、赫芝模型、双线性模型、 表面粘附力模型、湿颗粒模型和液桥模型）的算法研究和软件设计，以使该软件能适应 多种情况下农机工作部件的分析计算。还研究了分别采用线性模型和赫芝模型时，对离 散元法分析结果的影响。由研究可知，分别选取线性模型和赫芝模型计算接触力时，对 离散元法的分析结果影响不大。 本研究还改进了湿颗粒模型用于分析土壤，改进后的湿颗粒模型为：当法向叠合量δ n>0 时，法向接触力由线性粘弹性模型计算如下 Fn (t) Fnk (t) Fnd(t) Fna(t) 式中 Fnk (t) knn ，当两颗粒接触时 kn 取颗粒间接触的刚度系数，当颗粒与边界接 触时 kn 取颗粒与边界接触的刚度系数； Fnd(t) cnvn ；当两颗粒接触时 cn 取颗粒间接 触的阻尼系数，当颗粒与边界接触时 cn 取颗粒与边界接触的阻尼系数，vn 为接触处的法 向相对速度； Fna(t) ka (Ca n ) 。 当 – Ca ≤δn ≤0 时，法向接触力（拉力）为 13 Fn (t) kbCb 当 δn< - Ca 时，Fn(t)= 0。 5 接触力学模型中参数的研究 散粒物料的刚度系数、弹性模量、泊松比、碰撞恢复系数、破碎力和动静摩擦系数 等已有研究报道。为了给接触力学模型中参数的选取提供依据，本研究深入系统地研究 了散粒物料这些参数的测试方法及参数的影响因素和变化规律。 本研究采用电子万能试验机测试系统测试了三种大豆种子、三种玉米种子、三种小 麦种子和三种颗粒肥料的刚度系数、破碎力和弹性模量，见表 18~表 26 所示。采用三轴 试验及直剪试验和离散元法模拟双轴试验及直剪试验的方法，确定 3 种含水率土壤颗粒 的刚度系数和弹性模量，见表 27 和表 28 所示。 表 18 不同加载速度时吉科豆原种的刚度系数和破碎力 平放 加载速度 立放 刚度系数(N/mm) (mm/min) 破碎力(N) 刚度系数(N/mm) 破碎力(N) 平均值 标准差 平均值 标准差 平均值 标准差 平均值 标准差 0.5 125.152 16.397 127.904 17.2 79.9 13.795 79.812 8.4 1 126.978 15.001 132.400 10.2 77.0 8.726 80.236 9.1 2 127.901 13.791 139.967 16.1 84.6 14.886 83.426 10.3 表 19 加载速度为 1mm/min 时三种大豆种子的刚度系数和破碎力 平放 品种 吉 科 豆 吉 新 豆 吉豆 含水率 刚度系数(N/mm) 立放 破碎力(N) 刚度系数(N/mm) 破碎力(N) 平均值 标准差 平均值 标准差 平均值 标准差 平均值 标准差 9.2% 126.978 15.006 132.400 10.2 80.236 10.189 77.017 8.7 12.4% 50.484 7.890 83.547 12.0 37.836 12.269 47.170 9.9 16.3% 43.311 7.643 61.347 9.6 33.120 8.322 40.377 7.0 20.0% 28.481 8.500 43.557 8.8 23.902 7.697 29.023 6.7 23.9% 10.515 3.961 20.800 5.3 12.314 5.232 22.203 7.2 9.1% 144.483 14.215 161.423 12.8 106.751 13.423 68.917 10.1 12.1% 58.412 11.721 77.727 13.8 49.148 1.489 46.323 8.7 16.0% 39.231 9.213 42.333 7.7 24.823 1.988 34.768 9.2 20.2% 18.268 7.508 35.473 8.8 21.286 3.417 27.847 7.9 24.0% 12.562 5.569 24.243 8.1 14.371 4.707 23.183 7.3 9.4% 142.00 14.201 161.233 16.8 101.393 14.713 75.393 12.2 12.3% 47.582 11.177 65.493 12.8 39.304 8.775 42.663 10.1 14 16.2% 22.903 7.889 44.010 10.9 26.958 8.408 39.423 8.6 20.3% 13.110 5.415 26.697 8.9 14.777 5.375 21.947 8.7 24.1% 11.384 4.802 21.450 7.4 14.184 4.037 21.373 7.9 表 20 不同含水率时三种大豆种子的弹性模量(MPa) 品种 吉科豆 吉新豆 吉豆 含水率 9% 平均值 153.7 标准差 18.5 平均值 163.7 标准差 13.7 平均值 165.4 标准差 15.4 13% 59.6 9.7 54.8 10.1 45.8 9.3 17% 29.1 7.5 26.8 8.3 18.5 6.1 21% 18.2 5.9 11.4 4.3 9.8 3.0 25% 6.3 2.2 7.9 2.7 7.2 2.7 表 21 三种玉米种子的刚度系数(N/mm) 品种 吉峰 218 长单 228 平安 18 含水率(%) 加载速度 种子放置方向 平放 侧放 立放 12.85 (mm/min) 0.5 1 184.92 194.39 152.53 148.96 88.24 95.43 3 186.71 162.59 86.81 15.48 17.66 1 1 167.57 155.765 129.58 102.52 82.66 67.16 19.32 1 145.56 87.17 58.70 21.09 1 135.7 79.03 54.15 12.96 0.5 1 209.83 208.23 163.31 169.24 104.19 102.23 3 196.26 177.39 95.14 15.82 17.71 1 1 183.06 173.97 148.09 126.03 92.52 85.05 19.4 1 157.29 102.31 66.45 20.92 1 148.21 98.01 62.20 13.40 0.5 1 197.19 202.20 161.32 165.5 98.94 95.78 3 203.93 172.11 93.34 15.32 17.14 1 1 180.43 162.34 137.5 122.1 79.59 75.89 18.93 1 156.62 100.3 63.67 20.65 1 148.56 93.30 58.90 表 22 同一含水率时三种玉米种子的弹性模量(MPa) 15 含水率(%) 吉峰 218 长单 228 平安 18 13 358 409 382 表 23 不同含水率时三种小麦种子的弹性模量(MPa) 含水率(%) 小冰麦 33 丰强 11 吉 12 4 93.56 149.6 139.7 12 87.43 127.05 108.62 16 61.65 95.32 86.47 18 26.82 41.32 34.98 24 10.75 20.18 18.21 表 24 加载速度为 1mm/min 时三种小麦种子的刚度系数(N/mm) 品种 立放 放置方式 平放 横放 3.82 40.05 56.35 52.36 12.36 29.75 48.06 40.22 15.63 23.46 29.02 27.68 18.92 8.08 13.04 10.77 23.67 1.87 4.13 3.35 4.64 35.15 54.87 50.85 12.31 26.34 45.85 39.83 15.70 20.13 27.23 24.98 18.08 6.73 11.11 8.68 23.69 1.47 4.07 3.03 4.37 33.32 47.99 46.63 12.39 25.04 39.43 35.94 15.86 19.73 25.66 22.04 17.99 5.67 10.54 6.32 24.30 1.38 2.81 2.81 含水率(%) 小冰麦 33 丰强 11 吉 12 表 25 加载速度为 6mm/min 时三种肥料颗粒的刚度系数和破碎力 颗粒直径为 3mm 时 品种 含水率 2.97% 北京 生物 4% 8% 12% 大庆 3号 2.62% 4% 刚度系数(N/mm) 颗粒直径为 4mm 时 破碎力(N) 刚度系数(N/mm) 破碎力(N) 平均值 标准差 平均值 标准差 平均值 标准差 平均值 标准差 27.117 22.484 2.631 2.890 32.663 28.547 2.099 1.920 31.971 27.836 3.285 3.269 36.163 37.170 2.082 2.915 19.311 2.643 23.347 2.635 23.120 3.322 30.377 2.037 18.481 2.500 21.557 3.281 21.902 3.697 29.023 3.057 39.304 33.412 3.774 3.721 44.336 37.727 3.368 3.128 45.637 43.148 3.152 2.489 52.558 46.323 2.381 2.735 16 8% 12% 1.91% 山东 4号 4% 8% 12% 29.231 3.213 42.333 3.537 34.823 2.988 38.768 3.208 25.268 3.508 35.473 3.418 31.286 3.417 37.847 3.192 94.637 87.582 3.412 3.177 99.339 95.493 2.938 2.834 98. 767 89.304 2.403 3.775 108.266 96.663 4.136 4.012 75.903 2.889 84.010 2.929 76.958 3.408 89.423 3.986 68.110 3.415 76.697 2..936 72.777 3.375 81.947 3.787 表 26 不同含水率时三种肥料颗粒的弹性模量(MPa) 品种 大庆 3 号 北京生物 山东 4 号 平均值 标准差 平均值 标准差 平均值 标准差 含水率 原湿度 132.7 9.5 123.6 6.7 135.4 5.4 4% 86.6 9.7 74.8 9.1 85.8 9.3 8% 53.1 7.5 46.8 8.3 58.5 6.1 12% 28.2 5.9 22.4 4.3 21.8 5.0 表 27 离散元模拟双轴试验得到的三种土壤微观参数 样本 含水率 3%土壤样本 含水率 12%土壤样本 含水率 22%土壤样本 1 2 1 2 1 2 法向刚度系数 微观参数 3000~5000 5000~8000 5000~8000 8000~10000 15000~30000 20000~40000 切向刚度系数 2500~4200 4200~7500 4500~7500 7500~9000 10000~25000 15000~35000 法向阻尼系数 0.10~0.18 0.18~0.25 0.18~0.25 0.25~0.32 0.39~0.58 0.45~0.62 切向阻尼系数 0.08~0.14 0.14~0.23 0.16~0.23 0.20~0.28 0.32~0.48 0.39~0.55 注：1 指土壤颗粒之间的微观参数 2 指土壤颗粒与开沟器之间的刚度系数单位(N/m) 阻尼系数单位(N s/m) 表 28 三种土壤样本的弹性模量 土壤样本 弹性模量（MPa） 3% 1% 3.05 12% 1% 2.25 22% 1% 1.18 本研究还采用微观摩擦学试验机测试系统，测试了三种大豆种子、三种玉米种子、 三种小麦种子和三种颗粒肥料与不同表面的动摩擦系数，见表 29~表 33 所示；采用自制 的静摩擦系数测试装置和碰撞恢复系数测试装置（参见图 9 和图 10），测试了三种大豆 种子、三种玉米种子、三种小麦种子和三种颗粒肥料与不同表面的静摩擦系数和碰撞恢 复系数，见表 34~表 43 所示；采用三轴试验及直剪试验和离散元法模拟堆积的方法，确 定 3 种含水率土壤颗粒的摩擦系数，见表 44 所示。 17 表 29 大豆原种与不同摩擦表面之间的动摩擦系数 品 压 速度 种 力 (mm (N) /min) 500 1 2 3 平均 1 2 3 平均 1 2 3 平均 0.244 0.241 0.229 0.238 0.169 0.153 0.160 0.161 0.115 0.112 0.123 0.117 1000 0.231 0.245 0.241 0.239 0.153 0.156 0.163 0.157 0.112 0.111 0.118 0.114 1500 0.241 0.247 0.239 0.242 0.154 0.160 0.159 0.158 0.107 0.115 0.119 0.114 10 500 0.237 0.223 0.229 0.230 0.141 0.149 0.143 0.144 0.108 0.109 0.114 0.110 15 500 0.222 0.210 0.227 0.220 0.135 0.137 0.131 0.134 0.111 0.108 0.103 0.107 500 0.221 0.218 0.229 0.226 0.153 0.154 0.158 0.155 0.117 0.111 0.121 0.116 1000 0.215 0.227 0.228 0.223 0.152 0.148 0.157 0.152 0.124 0.114 0.118 0.119 1500 0.223 0.227 0.231 0.227 0.147 0.158 0.156 0.154 0.118 0.121 0.113 0.117 10 500 0.213 0.205 0.208 0.210 0.144 0.139 0.143 0.142 0.111 0.108 0.113 0.111 15 500 0.198 0.203 0.195 0.199 0.135 0.141 0.138 0.138 0.109 0.112 0.104 0.108 500 0.234 0.231 0.228 0.231 0.157 0.155 0.162 0.158 0.131 0.116 0.124 0.124 1000 0.223 0.231 0.228 0.227 0.149 0.156 0.158 0.154 0.120 0.129 0.125 0.125 1500 0.227 0.227 0.236 0.230 0.153 0.161 0.155 0.156 0.121 0.125 0.123 0.122 10 500 0.226 0.218 0.216 0.220 0.142 0.157 0.143 0.147 0.111 0.108 0.115 0.111 15 500 0.213 0.213 0.217 0.214 0.146 0.134 0.139 0.140 0.108 0.107 0.103 0.106 吉 5 科 豆 吉 5 新 豆 吉 豆 5 摩擦表面 铸铁 塑料 有机玻璃 表 30 压力和速度一定时大豆种子与不同摩擦表面间的动摩擦系数 摩擦表面 铸铁 品种 含水率 1 1 2 3 平均 1 2 3 平均 0.244 0.244 0.226 0.238 0.169 0.153 0.160 0.161 0.109 0.112 0.129 0.117 12.4% 0.228 0.269 0.261 0.253 0.211 0.204 0.195 0.203 0.147 0.136 0.123 0.135 16.3% 0.293 0.281 0.279 0.284 0.221 0.223 0.218 0.221 0.166 0.151 0.157 0.158 20.0% 0.313 0.302 0.291 0.302 0.229 0.251 0.241 0.240 0.184 0.197 0.194 0.192 13.9% 0.348 0.343 0.325 0.334 0.257 0.267 0.259 0.261 0.249 0.230 0.216 0.232 9.1% 0.210 0.239 0.229 0.226 0.153 0.154 0.158 0.155 0.117 0.111 0.121 0.116 12.1% 0.256 0.225 0.282 0.254 0.182 0.179 0.207 0.189 0.135 0.129 0.130 0.131 16.0% 0.277 0.276 0.284 0.279 0.205 0.208 0.240 0.216 0.184 0.181 0.165 0.177 20.2% 0.311 0.308 0.299 0.306 0.228 0.223 0.243 0.231 0.151 0.196 0.219 0.189 24.0% 0.349 0.338 0.341 0.343 0.248 0.264 0.242 0.251 0.300 0.205 0.178 0.228 9.4% 0.247 0.223 0.223 0.231 0.157 0.155 0.162 0.158 0.131 0.116 0.124 0.124 12.3% 0.271 0.273 0.257 0.267 0.183 0.183 0.218 0.195 0.144 0.138 0.135 0.139 吉豆 16.2% 0.304 0.285 0.287 0.292 0.219 0.203 0.225 0.216 0.178 0.183 0.160 0.174 20.3% 0.319 0.325 0.307 0.317 0.224 0.218 0.249 0.230 0.180 0.225 0.219 0.208 24.1% 0.341 0.335 0.327 0.334 0.254 0.243 0.256 0.251 0.239 0.251 0.211 0.234 吉 科 豆 吉 新 豆 3 有机玻璃 平均 9.2% 2 塑料 表 31 压力和速度一定时不同含水率下玉米种子的动摩擦系数 18 品种 摩擦表面 含水率(%) 吉峰 218 长单 228 平安 18 铸铁 塑料 有机玻璃 12.85 0.274 0.187 0.135 15.48 0.308 0.201 0.178 17.66 0.345 0.256 0.206 19.32 0.356 0.275 0.236 21.09 0.377 0.296 0.269 12.96 0.325 0.196 0.150 15.82 0.334 0.261 0.196 17.71 0.363 0.287 0.228 19.4 0.394 0.293 0.251 20.92 0.412 0.332 0.287 13.40 0.305 0.212 0.155 15.32 0.311 0.235 0.186 17.14 0.353 0.264 0.217 18.93 0.367 0.283 0.242 20.65 0.387 0.315 0.267 表 32 压力和速度一定时在不同含水率小麦种子与不同材料表面间的动摩擦系数 品种 小冰麦 33 丰强 11 吉 12 塑料 摩擦表面 有机玻璃 铸铁 3.82 12.36 0.129 0.172 0.095 0.131 0.191 0.224 15.63 0.187 0.155 0.254 18.92 0.209 0.192 0.260 23.67 0.263 0.238 0.305 4.64 12.31 0.149 0.172 0.102 0.136 0.215 0.251 15.70 0.202 0.156 0.246 18.08 0.225 0.195 0.272 23.69 0.252 0.245 0.313 4.37 12.39 0.139 0.172 0.109 0.146 0.211 0.237 15.86 0.200 0.172 0.262 17.99 0.237 0.209 0.288 24.30 0.275 0.261 0.331 含水率(%) 表 33 原含水率时三种肥料颗粒与不同摩擦表面之间的动摩擦系数 19 品 种 压 力 (N) 速度 (mm /min) 大 庆 3 号 山 东 4 号 铸铁 塑料 有机玻璃 1 2 3 平均 1 2 3 平均 1 2 3 平均 0.293 0.295 0.293 0.137 0.145 0.142 0.142 0.091 0.086 0.084 0.089 8000 0.287 0.303 0.301 0.308 0.306 0.139 0.146 0.138 0.14 0.088 0.089 0.085 0.086 10000 0.307 0.297 0.306 0.303 0.153 0.161 0.155 0.153 0.091 0.097 0.099 0.094 12000 0.315 0.311 0.316 0.154 0.151 0.152 0.094 0.097 8000 0.348 0.343 0.344 0.178 0.173 0.175 0.152 0.175 0.099 1.0 0.313 0.335 0.106 0.108 0.103 0.096 0.105 1.5 8000 0.359 0.363 0.355 0.358 0.201 0.204 0.205 0.202 0.110 0.115 0.111 0.113 6000 0.307 0.316 0.313 0.313 0.161 0.163 0.162 0.162 0.091 0.086 0.088 0.090 8000 0315 0.317 0.326 0.323 0.149 0.156 0.152 0.153 0.098 0.099 0.095 0.098 10000 0.323 0.333 0.328 0.326 0.161 0.163 0.165 0.163 0.101 0.097 0.099 0.100 12000 0.343 0.315 0.328 0.151 0.152 0.153 0.103 0.102 8000 0.352 0.355 0.356 0.135 0.141 0.138 0.152 0.138 0.104 1.0 0.330 0.354 0.109 0.112 0.104 0.103 0.108 1.5 8000 0.383 0.379 0.378 0.380 0.181 0.186 0.182 0.184 0.105 0.102 0.103 0.103 6000 0.373 0.374 0.373 0.374 0.197 0.203 0.199 0.199 0.101 0.096 0.099 0.100 8000 0.393 0.391 0.398 0.203 0.206 0.208 0.108 0.109 0.115 10000 0.387 0.383 0.386 0.191 0.193 0.195 0.111 0.107 0.109 12000 0.405 0.406 0.401 0.206 0.203 0.155 0.155 8000 0.424 0.428 0.425 0.232 0.231 0.234 0.203 0.232 0.116 1.0 0.403 0.426 0.204 0.152 0.149 0.148 0.115 0.150 1.5 8000 0.452 0.451 0.448 0.450 0.267 0.265 0.263 0.265 0.181 0.176 0.179 0.180 6000 北 京 生 物 摩擦表面 0.5 0.5 0.5 0.394 0.385 图 9 自制的静摩擦系数测试装置 20 0.206 0.193 0.112 0.109 图 10 自制的碰撞恢复系数测试系统 表 34 不同含水率时大豆种子与不同材料表面间的静摩擦系数 铸铁 品种 含水率 1 9.2% 摩擦表面 塑料 2 3 有机玻璃 平均 1 2 3 平均 1 2 3 平均 0.355 0.341 0.358 0.351 0.198 0.194 0.211 0.201 0.147 0.152 0.155 0.151 12.4% 0.376 0.374 0.385 0.378 0.225 0.231 0.239 0.232 0.176 0.179 0.187 0.181 16.3% 0.411 0.395 0.399 0.402 0.259 0.250 0.244 0.251 0.211 0.208 0.198 0.206 20.0% 0.417 0.421 0.429 0.422 0.273 0.275 0.281 0.278 0.223 0.229 0.224 0.225 13.9% 0.444 0.431 0.433 0.436 0.308 0.298 0.296 0.301 0.271 0.265 0.264 0.267 9.1% 0.344 0.335 0.353 0.344 0.188 0.183 0.191 0.187 0.169 0.175 0.171 0.172 12.1% 0.373 0.354 0.362 0.363 0.227 0.209 0.216 0.217 0.188 0.191 0.192 0.190 16.0% 0.386 0.395 0.393 0.391 0.251 0.258 0.261 0.257 0.194 0.209 0.211 0.205 20.2% 0.415 0.411 0.409 0.412 0.271 0.273 0.281 0.275 0.234 0.237 0.226 0.232 24.0% 0.449 0.447 0.458 0.451 0.312 0.317 0.305 0.311 0.277 0.268 0.259 0.264 9.4% 0.334 0.317 0.321 0.324 0.18 0.187 0.199 0.191 0.165 0.169 0.157 0.164 12.3% 0.367 0.351 0.363 0.360 0.222 0.212 0.209 0.214 0.211 0.199 0.197 0.202 吉豆 16.2% 0.381 0.396 0.378 0.385 0.233 0.229 0.246 0.236 0.234 0.225 0.224 0.228 20.3% 0.409 0.399 0.411 0.406 0.266 0.258 0.254 0.259 0.240 0.239 0.251 0.243 24.1% 0.431 0.436 0.429 0.432 0.296 0.287 0.299 0.294 0.277 0.269 0.274 0.273 吉 科 豆 吉 新 豆 表 35 不同含水率下玉米种子的静摩擦系数 品种 吉峰 218 长单 228 摩擦表面 含水率(%) 铸铁 塑料 有机玻璃 12.85 15.48 0.362 0.389 0.270 0.324 0.237 0.266 17.66 0.436 0.335 0.303 19.32 0.441 0.372 0.334 21.09 0.467 0.393 0.355 12.96 15.82 0.394 0.434 0.302 0.343 0.232 0.265 17.71 0.449 0.354 0.294 19.4 0.467 0.381 0.337 21 平安 18 20.92 0.492 0.412 0.376 13.4 15.32 0.396 0.412 0.302 0.347 0.241 0.275 17.14 0.429 0.368 0.299 18.93 0.457 0.378 0.323 20.65 0.467 0.395 0.364 表 36 不同含水率时小麦种子与不同材料表面间的静摩擦系数 品种 塑料 摩擦表面 有机玻璃 铸铁 3.82 0.18 0.15 0.24 12.36 0.23 0.20 0.29 15.63 0.25 0.22 0.31 18.92 0.28 0.26 0.34 23.67 0.36 0.34 0.41 4.64 0.18 0.15 0.24 12.31 0.22 0.19 0.27 15.70 0.25 0.23 0.30 18.08 0.27 0.26 0.34 23.69 0.34 0.33 0.37 4.37 0.18 0.15 0.24 12.39 0.23 0.20 0.29 15.86 0.25 0.22 0.31 17.99 0.28 0.26 0.34 含水率(%) 小冰麦 33 丰强 11 吉 12 表 37 在不同含水率时三种肥料与不同材料表面间的静摩擦系数 品种 北京 生物 大庆 3号 含水 率 摩擦表面 铸铁 塑料 有机玻璃 1 2 3 平均 1 2 3 平均 1 2 3 平均 2.97% 0.287 0.293 0.295 0.293 0.167 0.163 0.161 0.165 0.144 0.148 0.147 0.145 4% 0.321 0.322 0.323 0.322 0.171 0.179 0.175 0.176 0.148 0.142 0.146 0.147 8% 0.332 0.335 0.336 0.334 0.181 0.183 0184 0.183 0.166 0.151 0.157 0.159 12% 0.349 0.345 0.347 0.347 0.199 0.193 0.191 0.196 0.175 0.179 0.176 0.178 2.62% 0.373 0.354 0.362 0.363 0.223 0224 0.218 0.223 0.117 0.121 0.121 0.119 4% 0.385 0.381 0.389 0.382 0.242 0.239 0.237 0.242 0.145 0.139 0.140 0.142 8% 0.407 0.406 0.402 0.404 0.265 0.268 0.264 0.266 0.134 0.131 0.138 0.136 22 山东 4号 吉农 3号 12% 0.411 0.418 0.419 0.415 0.278 0.273 0.273 0.275 0.151 0.146 0.149 0.148 1.91% 0.389 0.386 0.388 0.388 0.244 0.250 0.259 0.251 0.144 0.148 0.145 0.146 4% 0.418 0.409 0.411 0.416 0.273 0.254 0.262 0.263 0.153 0.155 0.161 0.158 8% 0.417 0.421 0.429 0.423 0.283 0.289 0.284 0.285 0.172 0.176 0.179 0.177 12% 0.433 0.439 0.434 0.436 0.291 0.287 0.296 0.292 0.183 0.183 0.189 0.187 1.26% 0.397 0.393 0.393 0.395 0.217 0.215 0.212 0.213 0.161 0.156 0.154 0.156 0.421 0.413 0.417 0.418 0.223 0.229 0.228 0.226 0.174 0.168 0.165 0.168 0.424 0.425 0.417 0.423 0.251 0.255 0.257 0.255 0.178 0.173 0.168 0.172 0.432 0.435 0.437 0.436 0.254 0.261 0.259 0.267 0.185 0.180 0.189 0.187 4% 8% 12% 表 38 吉科豆在不同含水率和不同高度时的碰撞恢复系数 碰撞表面 铸铁 塑料 有机玻璃 大豆 含水率 跌落高度(cm) 30 15 50 9.2% 平均值 0.579 标准差 0.025 平均值 0.538 标准差 0.021 平均值 0.506 标准差 0.022 12.4% 0.497 0.036 0.486 0.026 0.481 0.010 16.3% 0.481 0.029 0.473 0.019 0.465 0.014 20.0% 0.467 0.032 0.463 0.012 0.456 0.013 23.9% 0.458 0.021 0.445 0.012 0.442 0.014 9.2% 0.508 0.032 0.494 0.025 0.489 0.014 12.4% 0.469 0.029 0.449 0.026 0.441 0.013 16.3% 0.448 0.030 0.439 0.021 0.436 0.013 20.0% 0.442 0.033 0.430 0.027 0.414 0.013 23.9% 0.427 0.031 0.420 0.016 0.402 0.011 9.2% 0.565 0.033 0.522 0.018 0.501 0.015 12.4% 0.488 0.017 0.481 0.014 0.474 0.011 16.3% 0.473 0.028 0.467 0.023 0.456 0.013 20.0% 0.461 0.031 0.446 0.022 0.453 0.017 23.9% 0.455 0.022 0.444 0.026 0.440 0.012 9.2% 0.472 0.030 0.455 0.014 0.440 0.010 12.4% 0.440 0.043 0.440 0.020 0.432 0.016 16.3% 0.438 0.031 0.433 0.019 0.423 0.017 20.0% 0.432 0.027 0.423 0.022 0.419 0.025 23.9% 0.431 0.024 0.419 0.022 0.415 0.011 23 表 39 吉新豆在不同含水率和不同高度时的碰撞恢复系数 碰撞表面 铸铁 塑料 有机玻璃 大豆 含水率 跌落高度(cm) 30 15 50 9.2% 平均值 0.593 标准差 0.042 平均值 0.547 标准差 0.043 平均值 0.545 标准差 0.015 12.4% 0.512 0.039 0.500 0.015 0.492 0.023 16.3% 0.494 0.027 0.488 0.016 0.485 0.008 20.0% 0.481 0.039 0.478 0.015 0.464 0.008 23.9% 0.477 0.029 0.465 0.019 0.457 0.016 9.2% 0.545 0.047 0.533 0.021 0.521 0.014 12.4% 0.471 0.020 0.469 0.026 0.452 0.018 16.3% 0.440 0.031 0.436 0.020 0.429 0.012 20.0% 0.434 0.023 0.426 0.027 0.414 0.015 23.9% 0.420 0.029 0.414 0.016 0.411 0.018 9.2% 0.584 0.024 0.537 0.022 0.538 0.017 12.4% 0.503 0.031 0.492 0.019 0.478 0.013 16.3% 0.485 0.038 0.484 0.026 0.472 0.018 20.0% 0.474 0.021 0.468 0.024 0.465 0.019 23.9% 0.463 0.027 0.462 0.021 0.451 0.012 9.2% 0.467 0.018 0.453 0.015 0.443 0.025 12.4% 0.463 0.021 0.448 0.009 0.431 0.006 16.3% 0.448 0.021 0.443 0.007 0.427 0.008 20.0% 0.440 0.028 0.440 0.006 0.426 0.006 23.9% 0.434 0.028 0.432 0.010 0.422 0.005 表 40 吉豆在不同含水率和不同高度时的碰撞恢复系数 碰撞表面 铸铁 塑料 有机玻璃 大豆 含水率 跌落高度(cm) 30 15 50 9.2% 平均值 0.635 标准差 0.028 平均值 0.577 标准差 0.029 平均值 0.554 标准差 0.008 12.4% 0.541 0.041 0.513 0.017 0.495 0.009 16.3% 0.523 0.016 0.490 0.025 0.488 0.017 20.0% 0.486 0.016 0.478 0.012 0.477 0.016 23.9% 0.481 0.029 0.471 0.016 0.463 0.015 9.2% 0.568 0.024 0.556 0.023 0.510 0.015 12.4% 0.467 0.037 0.460 0.021 0.441 0.013 16.3% 0.432 0.026 0.431 0.023 0.420 0.014 20.0% 0.416 0.031 0.405 0.019 0.401 0.017 23.9% 0.405 0.028 0.402 0.015 0.400 0.012 9.2% 0.623 0.023 0.567 0.022 0.547 0.016 12.4% 0.537 0.021 0.492 0.025 0.476 0.019 16.3% 0.511 0.034 0.476 0.018 0.471 0.022 20.0% 0.477 0.033 0.471 0.015 0.466 0.013 23.9% 0.475 0.028 0.468 0.020 0.461 0.015 9.2% 0.507 0.028 0.491 0.009 0.457 0.013 12.4% 0.492 0.018 0.468 0.018 0.451 0.011 24 16.3% 20.0% 0.473 0.462 0.022 0.019 0.455 0.454 0.012 0.011 0.451 0.445 0.011 0.011 23.9% 0.451 0.017 0.447 0.007 0.444 0.013 表 41 三种玉米种子与铸铁的碰撞恢复系数 种子品种 吉峰 218 长单 228 平安 18 跌落高度(mm) 含水率(%) 100 200 300 12.85 15.48 0.732 0.692 0.698 0.676 0.686 0.658 17.66 0.644 0.638 0.622 19.32 0.630 0.616 0.604 21.09 0.615 0.602 0.596 12.96 15.82 0.703 0.664 0.691 0.634 0.672 0.601 17.71 0.612 0.615 0.58 19.4 0.592 0.572 0.561 20.92 0.579 0.563 0.543 13.4 15.32 0.698 0.677 0.671 0.639 0.655 0.63 17.14 0.619 0.621 0.599 18.93 0.603 0.57 0.562 20.65 0.587 0.553 0.548 表 42 三种小麦种子与铸铁的碰撞恢复系数 种子品种 小冰麦 33 丰强 11 吉 12 含水率(%) 跌落高度 10(cm) 跌落高度 20(cm) 跌落高度 30(cm) 3.82 0.571 0.569 0.535 12.36 0.532 0.504 0.462 15.63 0.512 0.490 0.434 18.92 0.444 0.423 0.387 23.67 0.378 0.354 0.338 4.64 12.31 0.565 0.528 0.550 0.500 0.520 0.477 15.70 0.506 0.479 0.447 18.08 0.457 0.450 0.412 23.69 0.418 0.371 0.325 4.37 0.580 0.547 0.517 12.39 0.543 0.528 0.509 15.86 0.483 0.460 0.427 17.99 0.438 0.4299 0.389 24.30 0.383 0.334 0.311 25 表 43 北京生物肥料在不同含水率、不同高度、不同碰撞表面时的碰撞恢复系数 碰撞表面 含水率 10 平均值 2.97% 4% 铸铁 8% 12% 2.97% 4% 塑料 8% 12% 2.97% 4% 有机玻璃 8% 12% 2.97% 肥料颗粒 4% 8% 12% 跌落高度(cm) 20 40 标准差 0.025 平均值 0.538 标准差 0.021 平均值 0.506 标准差 0.022 0.036 0.486 0.026 0.481 0.010 0.228 0.029 0.473 0.019 0.465 0.014 0.212 0.032 0.463 0.012 0.456 0.013 0.508 0.032 0.494 0.025 0.489 0.014 0.469 0.029 0.449 0.026 0.441 0.013 0.448 0.030 0.439 0.021 0.436 0.013 0.442 0.033 0.430 0.027 0.414 0.013 0.565 0.033 0.522 0.018 0.501 0.015 0.488 0.017 0.481 0.014 0.474 0.011 0.473 0.028 0.467 0.023 0.456 0.013 0.461 0.031 0.446 0.022 0.453 0.017 0.472 0.030 0.455 0.014 0.440 0.010 0.440 0.043 0.440 0.020 0.432 0.016 0.438 0.031 0.433 0.019 0.423 0.017 0.432 0.027 0.423 0.022 0.419 0.025 0.248 0.237 表 44 土壤颗粒之间的摩擦系数及颗粒与边界之间的摩擦系数 土壤样本 颗粒与颗粒之间摩擦系数 颗粒与边界之间的摩擦系数 3% 1% 0.3076 0.3~0.5 12% 1% 0.2568 0.3~0.5 22% 1% 0.2235 0.3~0.5 还分析了 15 种散粒物料的物理力学参数的影响因素，通过离散元法计算结果与工作 部件试验结果的对比，确定了参数的选取方法。 6 三维离散元法算法研究和软件设计 6.1 颗粒入口的研究 为与农机工作部件的实际工作情况相同，本研究提出了基于进料口的颗粒样本生成 方法，其思想是：根据机械部件的实际工作情况，在计算区域内确定一个或几个具有一 定形状和尺寸的颗粒生成区域——进料口（参见图 11） ，也可由 CAD 模型读取进料口； 然后在进料口内每一个计算时步生成一个颗粒，颗粒中心位置随机确定，尺寸按一定分 布规律随机产生；如果生成的颗粒不与其它已生成的颗粒或边界重叠则生成成功，随即 26 将该颗粒加入到存储颗粒的数组中；下一时步该颗粒在重力作用下下落，并按离散元法 计算该颗粒的运动速度和位移。 图 11 散粒物料由矩形进料口进入过程简图 图 12 求解两颗粒接触点的几何势方法示意图 6.2 接触叠合量求解算法研究 本研究采用“几何势方法(Geometric Potential Algorithms)” ，建立了颗粒之间及颗粒 与农机工作部件图元之间接触判断算法和接触叠合量求解算法。对于椭球颗粒“几何势 方法”的思想是（见图 12）：分别求出两椭球颗粒进入对方的最远点，如图中的 A、B 两点，如果 A、B 两点分别在两个椭球内，则两个椭球颗粒接触，A、B 两点的距离为 接触叠合量，否则不接触。接触时求解 A、B 点的数学模型分别为 FB min xn1 yn2 zn3 s.t. F A an1 bn2 cn3 1 0 xn1 yn2 zn3 F 1 A an1 bn2 cn3 min s.t. F 0 B 式中 FA 为在颗粒 A 的局部坐标系下，表示颗粒 A 的超二次方程；FB 为在颗粒 A 的 局部坐标系下，表示颗粒 B 的超二次方程。 27 当求解颗粒与边界图元的接触作用点和叠合量时，以圆形平面边界图元为例， 求解 A、 B 两点的数学模型分别为 xn1 yn2 zn3 l x m y n z t min 1 0 0 0 0 n n n a 1 b 2 c 3 min n n n x 1 y 2 z 3 s . t . 1 s . t . l x m y n z t 0 n1 0 0 0 0 n2 n3 0 a b c (x xo0 )2 ( y yo0 )2 (z zo0 )2 R0 (x xo0 )2 ( y yo0 )2 (z zo0 )2 R0 式中(xo0, yo0, zo0)、R0、 l 0 , m0 , n0分别为在颗粒的局部坐标系下，圆形平面图元的圆 心点坐标、半径和法向量。求解其它种类图元时，只需更换上述的图元方程即可。 下面以球颗粒为例，介绍颗粒之间接处判定算法和接触叠合量求解算法。首先应建 立如图 13 所示的坐标系。全局坐标系在整个计算和演示过程中是固定不变的，坐标系 的原点由用户选取，为了计算方便并获得更好的视图效果，推荐将整个分析区域的中心 作为坐标原点，如图 13(a)中的 xyz 坐标系。局部坐标系一般取被分析球体 i 的球心作为 坐标原点，取由球体 i 的球心到球体 j 的球心连线为 X 轴；在与 X 轴垂直的平面内，过 球心 i 取一平行于 x-y 平面的直线为 Y 轴， Z 轴由右手螺旋法则确定，局部坐标系如图 13(b) 所示。此时，局部坐标系和全局坐标系中各物理量间有如下对应关系。 lx X x Y T y； T m l x2 my2 y 2 2 Z z l xnz l x my 0 l x2 my2 my lx nz l x2 my2 l x2 my2 mynz 式中{ l x , my , nz }是 X 轴在全局坐标系中的方向余弦，且 l x xj xi Lij ， my yj yi Lij ， nz z j zi Lij ； Lij x x y y z z ；ri、rj、 2 i j 2 i (xi, yi, zi)和(xj, yj, zj)分别为球体 i 和 j 的半径及中心点坐标。 28 j 2 i j 图 13 坐标系建立及球颗粒空间位置简图 对于两球体 i 和 j，当球心距离 Lij 满足 Lij x x y y z z r r 2 i 2 j 式时， 则认为两球接触，且接触叠合量为： ij ri r j i 2 j i j i j x x y y z z 。 2 i j 2 i j 2 i j 当球颗粒与边界接触时，以圆形平面边界图元为例，参见图 14 所示。设圆形平面的 圆心为 O0(xo0, yo0, zo0)，半径为 R0，法向量为 l 0 , m0 , n0。颗粒球心点 Q(xq, yq, zq)，半径 为 r，球与圆形平面边界接触判断如下。 （1） 边界运动否，如果边界运动首先更新边界参数，然后进入（2） 。如果边界不 运动，不用更新边界参数，直接进入（2）。 （2）计算球心点 Q 到圆形平面 O 所在平面的距离，该圆形平面的点法式方程为 l 0 x xo0 m0 y yo0 n0 z zo0 0 化为一般式为 l 0 x m0 y n0 z t0 0 式中 t0 l 0 xo0 m0 yo0 n0 zo0 。求出点 Q 到平面的距离 dqo0 l 0 xq m0 yq n0 zq t0 l 02 m02 n02 如果 dqo0 r 则可能接触，进行下一步判断；如果 dqo0 r 则无接触，判断结束退出。 （3）求出 Q 在圆面上的投影点 Q0(xq0, yq0, zq0)为 xq0 xq l 0 k0 2 yq0 yq m0 k0 2 z z n k 2 q 0 0 q0 29 式中 k0 2l 0 xq m0 yq n0 zq t0 l 0 m0 n0 。由此求出 Q0 点 O0 点间距离为 dq0o0 2 2 2 x x y y z z 2 qo o0 2 qo o0 2 qo o0 如果 dq0o0 R0 r 则可能接触，进行下一步判断；如果 dq0o0 R0 r 则无接触，判 断结束退出。 （4） 判断该圆形平面对应的虚边界数组是否为空，是空则进入（6） ；不为空则调 用虚边界的接触判断。 （5）如果与虚边界没有接触点则返回进入（6）；如果与虚边界有接触点，且接触点 在虚边界内，则颗粒与实边界无接触，判断结束退出；如果与虚边界有接触点，且接触 点在虚实边界的交线上，设接触点为 Cv(xcv, ycv, zcv)，求得 xcv xo0 2 ycv yo0 2 zcv zo0 2 do0cv 如果 do0cv R0 则接触，局部坐标系 X 轴为由 Q 指向 Cv，叠合量 do0cv ，判断结 束退出。否则 do0cv R0 无接触，判断结束退出。 （6）实圆形平面边界的判断和叠合量计算： 如果 dq0o0 R0 说明 Q0 在圆面内，可能发生面接触，接触点为 C0 Q0 ，叠合量为 r dqo0 ，局部坐标系 X 轴为由 Q 指向 Q0，判断结束退出。 如果 R0 dq0o0 R0 r ，说明可能发生圆面圆弧边界与球颗粒的接触，接触点 C0 为 xc0 xo0 xq0 xo0 R0 dq0o0 yc0 yo0 yq0 yo0 R0 dq0o0 z z z z R d o0 q0 o0 0 q0o0 c0 叠合量 r QC0 ，将 Q、C0 代入求得两点间距离为 QC0 x x y y z z 2 c0 2 q c0 q 局部坐标系 X 轴为由 Q 指向 C0，判断结束退出。 30 2 c0 q 图 14 球颗粒与圆形平面接触判定简图 图 15 球体相对位移合转角计算简图 限于篇幅，上面仅介绍了球颗粒间及球颗粒与圆形平面边界之间的接触处判定算法 和接触叠合量求解算法，椭球和超球颗粒及其它种类边界图元的求解算法丛略。 6.3 接触相对位移和转角的求解算法研究 以球颗粒为例，当时间增加一个时步 Δt 后，设两球心沿 x、y、z 轴方向的平动位 移增量分别为 Δuix、Δuiy、Δuiz、Δujx、Δujy、Δujz，绕 x、y、z 轴的转角增量分别为 Δφix、Δφiy、Δφiz、Δφjx、Δφjy、Δφjz，则有 ix Δ t Δ uix x iyt uiy y 31 uiz ziz Δ t ix ixt jxt u jx x jx jxt iy iyt jyt u jy y jy jyt iz iz t u jz zjz t jz jz t 由此可得局部坐标系中单元 i 相对于单元 j 的相对平动位移 ΔuX、ΔuY、ΔuZ 和转角增 量 ΔφX、ΔφY、ΔφZ 为 Δuix Δu jx 0 ΔuX Δu T Δu Δu Δ jy iZ iy Y ΔuZ Δuiz Δu jz Δ iY Δ X Δ ix Δ jx Δ jZ ； Δ Y T Δ iy Δ jy Δ iz Δ jz Δ jY Δ Z 0 并且规定：相对位移和转角的方向均以坐标系正方向和沿坐标轴逆时针方向为正。 其中 iZ iZ t iY iY t jZ jZ t jZ jZ t 局部坐标系下两球体的转动角速度 iX 、 iY 、 iZ 、 jX 、 jY 、 jZ 分别为 jx iX ix jX T ； T jy iY iy jY jz iZ iz jZ t t t t t 当球颗粒与边界 j 接触时，以平面边界为例， 求解 u jx ， u jy ， u jz ， jx ， jy ， t jz 的过程如下。参图 16 和 17 所示，vsx ，vsy ，v sz 表示边界的平动速度；ωzz 为边界 转动的角速度，其转轴为 Lzz，dzz 为转轴到接触点的距离；vzz 为旋转的线速度。 图 16 球体与平面边界的接触 32 图 17 接触边界速度计算简图 过已知点 P 的转轴 Lzz 的方程为 x xp y yp z zp lp mp np 式中 l p , mp , np 是 Lzz 在全局坐标系中的方向余弦，已知点 P 在全局坐标系下的坐 标为(xP, yP, zP)。转轴 Lzz 的方向规定为：以转轴 Lzz 的向量方向为右手拇指方向，角速度 方向与四指方向相同为正，相反则为负。接触点 C 到转轴 Lzz 的距离为 n y y m z z l z z n x x m x x l y y 2 dzz p c p p c 2 p p c p p c 2 p p c p p c p l p2 mp2 np2 过点 C 向转轴做垂线，垂点 D(xD, yD, zD）为 xD xp l p s yD yp mp s z z n s p p D 其中 s xc xp l p yc yp mp zc zp np 点的向量 a 在三个坐标轴上方向余弦为 33 l m n 。由垂点 D 指向接触 2 p 2 p 2 p la ma na xC xD xC xD 2 yC yD 2 zC zD 2 yC yD xC xD 2 yC yD 2 zC zD 2 zC zD xC xD 2 yC yD 2 zC zD 2 由此产生的线速度向量在全局下为转轴向量与向量 a 的差积。设转轴 Lzz 的方向向 量为 b，则线速度的方向向量为 i j k b a l p mp np mpna npma i npl a l pna j l pma mpl a k l a ma na 线速度的大小为 vzz ωzz dzz 则线速度在全局坐标中为 vzx mpna npma v n l l n v zy p a p a zz vzz l pma mpl a 由此求得边界在全局坐标中的速度为 vjx vzx vsx vjy vzy vsy vjz vsz vzz 进而得到边界在全局坐标中的位移增量为 t u vjx jx t u jy tvjy t u vjz jz 将角速度ωzz 变换到全局坐标下得 jx l p jy mp zz jz np 进而得到 t jx jx t jy t jy t jz jz 34 限于篇幅，上面仅介绍了球颗粒间及球颗粒与圆形平面边界之间的接触相对位移和 转角的求解算法，椭球和超球颗粒及其它种类边界图元的求解算法丛略。 6.4 接触作用力求解算法研究 以线性粘弹性模型接触力学模型为例，颗粒上所受的接触作用力，其法向力为 (t) (t) FX(t) F XK F XD (t) 式中 FX(t)为当前时步接触两体间的法向作用力； F XK 为接触的法向弹性力， (t) F XK K X X ，KX 为接触的法向刚度系数， K X K 1K 2 ，K1、K2 分别为接触两体的 K1 K 2 (t) 法向刚度系数； X 为接触的法向叠合量； F XD C X uX t ；CX 为法向粘性阻尼系数， CX 2lne m0 K 1,2 C1C2 ，C1、C2 分别为接触两体的法向粘性阻尼系数； C1,2 ， C1 C2 2 ln2 e m0 m1m2 ，m1、m2 分别为接触两体的质量，当颗粒与边界接触时，m0 等于颗粒的 m1 m2 质量，e 为碰撞恢复系数， uX 为两体接触处的法向相对位移，Δt 为时步。 接触颗粒上所受的切向作用力为 FY((tZ) ) FY((tZ) )e FY((tZ) )d 式中 FY((tZ) ) 为 t 时刻接触两体间切向作用力； FY((tZ) )e FY((tZ)et) FY((tZ) )e ； (t) t t FY((tZ) )e K Y ( Z ) uY ( Z ) ； FY ( Z )d CY ( Z ) uY ( Z ) t ；KY(Z)为接触的切向刚度系数 K Y(Z) k1k2 ，k1、k2 分别为接触两体的切向刚度系数； uY (Z) 为接触处的切向相对位 k1 k2 移；CY(Z)为切向粘性阻尼系数， CY (Z) 2lne m0 K Y (Z) 2 ln2 e 。 接触颗粒上所受的回转力矩为 (t ) M X(t) M Xe Md t (t) (tt) t t 式中 M Xe 为接触的回转 M Xe M Xe； M Xe K R X ； M d C R X t ； K R 刚度系数； X 为接触处的回转相对角位移。 35 2 2 接触颗粒上所受的作用力还应满足库仑—莫尔准则，即当 FY(t) FZ(t) F X(t) 时， 2 2 FY((tZ) ) 由上式计算；当 FY(t) FZ(t) F X(t) 时， FY((tZ) ) FX(t) FY((tZ) ) 2 2 FY(t) FZ(t) ， FY((tZ) ) 符号取为修正前切向力的符号。 6.5 颗粒运动方程求解 设颗粒 i 在 x，y，z 方向所受的作用力分别为 Fix，Fiy，Fiz；转动力矩(逆时针为正) 分别为 Mix，Miy，Miz，则有 t t Fix FiX 1 t t Fiy T FiY t t F F i z iZ t t M ix M iX 0 1 1 t t M iy T 0 ri T FiY t t M 0 FiZ iz 式中，Σ表示对所有与 i 单元相接触的边界求和。由此得到单元 i 在当前时步的合 力、合力矩，则得 t +Δt/2 时刻单元 i 的速度和角速为 t x Fix i tt 2 x i t-t 2 tt 2 t-t 2 t t i i y y m Fiy t i ztt 2 zt-t 2 i i Fiz t ixtt 2 ixt-t 2 M ix tt 2 t-t 2 t t iy iy I M iy t t 2 t t 2 t i M iz iz iz 式中 mi、Ii 分别为单元 i 的质量和转动惯量。由 t t 2时刻的速度可得 t t 时刻 单元 i 的位移和转角为 t xi tt xi x i tt 2 tt tt 2 t i y y t i i y t z tt z ztt 2 i i i t ixtt ix ixtt 2 tt tt 2 t iy iy tiy t tt tt 2 iz iz iz 根据单元的位移和转角可得单元的新位置，然后重新进行新一轮的接触判断、接触 力计算及位移和转角计算，如此重复直至计算终止。离散元法计算循环是以时步 t 向前 36 差分进行的，具体求解过程参见图 18 所示。 图 18 中心差分格式时步之间的关系 6.6 机械部件性能评价算法研究 在给定时间内根据每个颗粒的受力和位移，即可对农机工作部件性能进行分析和评 价。如统计颗粒运动通过某一截面的数量，可判断颗粒通过的质量或生产率；当颗粒受 力超过极限值为破碎；与农机工作部件接触的全部颗粒作用的合力为工作部件的工作阻 力等。根据上述思想本研究组还建立并实现了单颗粒运动速度、位移和作用力分析、颗 粒群体的速度场和作用力场分析、颗粒流量和机械部件工作阻力分析算法和求解模块， 以评价农机工作部件的性能，见图 19~图 22 所示。 图 19 单颗粒作用力变化曲线 37 图 20 单颗粒速度变化曲线 图 21 颗粒群体流量变化曲线 图 22 机械部件工作阻力变化曲线 6.7 动态显示等其它算法研究 本研究以 C++为开发工具，利用 OpenGL 软件及其函数绘制并显示颗粒和农机工作 部件图元的三维图形，采用 OpenGL 消隐功能去除虚边界，从而实现散粒物料与农机工 作部件图元运动过程仿真显示，由此可直观看到工作部件的工作原理、结构和尺寸参数 对颗粒运动的影响。 由于在短时间内，对软件做较大改动很难实现。因此探讨了基于 OpenMP 和 MPI 的 并行算法，并据此对程序进行了改进，结果表明改进效果明显。 本研究采用“白盒测试”与“黑盒测试”相结合的方法，进行软件的单元测试、子 系统测试和和系统综合测试，综合测试的初步结果见图 23 所示，利用数据库共享实现 了异构环境系统集成。 38 图 23 大豆种子跌落碰撞过程及弹起高度的仿真和试验结果对比 本研究建立了多种属性颗粒和多种属性边界的计算方法、变速运动边界和弹性边界 的计算方法等，以同时分析不同属性颗粒（如种子、肥料、土壤）与具有不同属性边界 （如金属、橡胶）机械部件的接触作用，见图 24 和图 25 所示。 图 24 多种属性颗粒和多种属性边界的分析计算 图 25 具有弹性边界的排种器工作过程仿真分析 7 新型 CAE 软件的应用研究 采用新型 CAE 软件对 15 种散粒物料（种子、颗粒肥料和土壤）及 11 种农机工作 部件（5 种排种器、3 种排肥器和 3 种开沟器）的工作过程进行了仿真分析，并将仿真 结果与农机工作部件工作过程的实验结果进行了对比，见图 26~图 39 所示。 39 图 26 由 CAD 模型得到的组合内窝孔精密排种器播种大豆种子过程仿真分析与实验对比 图 27 由 CAD 模型得到的组合内窝孔精密排种器播种玉米种子过程仿真分析与实验对比 图 28 由 CAD 模型得到的组合内窝孔精密排种器播种小麦种子过程仿真分析与实验对比 图 29 由 CAD 模型得到的垂直圆盘型孔轮式精密排种器播种大豆种子过程仿真分析与实验对比 40 图 30 由 CAD 模型得到的垂直圆盘型孔轮式精密排种器播种玉米种子过程仿真分析与实验对比 图 31 由 CAD 模型得到的窝眼轮式精密排种器播种大豆种子过程仿真分析与实验对比 图 32 由 CAD 模型得到的窝眼轮式精密排种器播种玉米种子过程仿真分析与实验对比 图 33 由 CAD 模型得到的槽轮式精密排种器播种小麦种子过程仿真分析与实验对比 图 34 由 CAD 模型得到的第一种外槽轮式排肥器工作过程仿真分析与实验对比 41 图 35 由 CAD 模型得到的第二种外槽轮式排肥器工作过程仿真分析与实验对比 图 36 由 CAD 模型得到的第三种外槽轮式排肥器工作过程仿真分析与实验对比 图 37 由 CAD 模型得到的施肥开沟器工作过程仿真分析与实验对比 图 38 由 CAD 模型得到的船铲开沟器工作过程仿真分析与实验对比 图 39 由 CAD 模型得到的芯铧开沟器工作过程仿真分析 42 仿真过程中，选择的散粒物料物理力学参数见表 45~表 50 所示。 表 45 三维仿真模拟时大豆种子的参数选取 颗粒与颗粒 颗粒与边界 5200 3500 0.35 0.28 0.23 0.11 10400 6900 0.43 0.35 0.27 0.45 法向刚度系数 (N/m) 切向刚度系数 (N/m) 法向阻尼系数 (N s/m) 切向阻尼系数 (N s/m) 动摩擦系数 静摩擦系数 颗粒密度 (kg/m3) 球颗粒直径(mm) 椭球颗粒长、宽、高(mm) 计算时步(s) 颗粒数 1.21 x 103 6.76 7.6，6.9，6.3 5 x 10-5 570 表 46 开沟过程仿真模拟时三种土壤微观参数选取 样本 含水率 3%土壤样本 1 2 含水率 12%土壤样本 1 2 含水率 22%土壤样本 1 2 微观参数 法向刚度系数 3000~5000 5000~8000 5000~8000 8000~10000 15000~30000 20000~40000 切向刚度系数 2500~4200 4200~7500 4500~7500 7500~9000 10000~25000 15000~35000 法向阻尼系数 0.10~0.18 0.18~0.25 0.18~0.25 0.25~0.32 0.39~0.58 0.45~0.62 切向阻尼系数 0.08~0.14 0.14~0.23 0.16~0.23 0.20~0.28 0.32~0.48 0.39~0.55 表 47 土壤颗粒之间的摩擦系数及颗粒与边界之间的摩擦系数 土壤样本 颗粒与颗粒之间摩擦系数 颗粒与边界之间的摩擦系数 3% 1% 0.3076 0.3~0.5 12% 1% 0.2568 0.3~0.5 22% 1% 0.2235 0.3~0.5 表 48 三维仿真模拟时玉米种子的参数选取 颗粒与颗粒 颗粒与边界 法向刚度系数(N/m) 4860 9720 切向刚度系数(N/m) 3888 7776 法向阻尼系数(N/m) 0.19 0.67 43 切向阻尼系数(N/m) 0.17 0.60 动摩擦系数 0.285 0.135 静摩擦系数 0.387 0.237 颗粒密度 (kg/m3) 1.172x 103 球颗粒直径(mm) 7.04 颗粒上底长、下底长、高度、厚 度(mm) 4.9， 7.55， 9.57，4.87 计算时步(s) 5 x 10-5 颗粒数 950 表 49 三维仿真模拟时小麦种子的参数选取 法向刚度系数 (N/m) 切向刚度系数 (N/m) 法向阻尼系数 (N s/m) 切向阻尼系数 (N s/m) 动摩擦系数 静摩擦系数 颗粒与颗粒 颗粒与边界 1967 1312 0.092 0.075 0.23 0.11 3934 2623 0.15 0.12 0.27 0.45 颗粒密度 (kg/m3) 球颗粒直径(mm) 椭球颗粒长、宽、高(mm) 计算时步(s) 颗粒数 1.373 x 103 4.16 7.18，3.32，3.04 10-4 3500 表 50 三维仿真模拟时肥料颗粒的参数选取 法向刚度系数 (N/m) 切向刚度系数 (N/m) 法向阻尼系数 (N s/m) 切向阻尼系数 (N s/m) 动摩擦系数 静摩擦系数 颗粒与颗粒 颗粒与边界 22000 18000 0.22 0.18 0.09 0.14 25000 21000 0.24 0.2 0.18 0.24 颗粒密度 (kg/m3) 1.973 x 103 球颗粒直径(mm) 椭球颗粒长、宽、高(mm) 计算时步(s) 颗粒数 3.79 4.01，3.67，3.46 10-4 5000 44 图 40 组合内窝孔精密排种器充种角、清种角和投种角的示意图 图 41 组合内窝孔排种器和型孔轮式排种器的种子面倾角示意图 在实验过程中，采用高速摄像系统拍摄了排种器的工作过程及颗粒的运动过程，采 用爆破分析软件 BLASTER’SMAS（BLASTER’S MOTION ANALYSIS SOFTWARE）分 析了颗粒运动速度和位移，并将实验分析结果与软件仿真计算结果进行了比较，见表 51~ 45 表 72 所示。其中各特征点的定义见图 40~图 43 所示。 图 42 组合内窝孔精密排种器内特殊点种子位置示意图 46 图 43 窝眼轮式排种器内特殊点种子位置示意图 表 51 组合内窝孔精密排种器播大豆种子时粒数的仿真与试验对比 转速(r/min) 线速度(m/s) 单粒率(%) 26.82 39.33 52.0 65.43 0.19 0.28 0.38 0.47 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 播吉科豆 100 100 0 100 100 0 100 100 0 93.7 100 播吉新豆 100 100 0 100 100 0 100 100 0 90.9 100 10.01% 播吉豆 100 100 0 100 100 0 100 100 0 88.2 100 13.38% 表 52 组合内窝孔精密排种器播玉米种子时粒数的仿真与试验对比 26.82 39.33 51.96 试验值(%) 80 70.97 60.47 仿真值(%) 65.4 78.25 75 相对误差(%) 18.3 10.25 24.2 试验值(%) 20 25.71 39.53 仿真值(%) 19.2 16.33 25 转速(r/min) 单粒率 双粒率 47 6.72% 相对误差(%) 4 36.48 36.32 试验值(%) 0 0 0 仿真值(%) 15.4 5.42 0 空穴率 相对误差(%) 表 53 组合内窝孔精密排种器播小麦种子时粒数的仿真与试验对比 转速(r/min) 仿真值(粒) 试验值(粒) 相对误差(%) 26.82 12 22 45 39.33 18 30 40 51.96 33 41 20 65.43 68 57 19 表 54 组合内窝孔精密排种器播大豆种子时清种角、充种角和投种角的仿真与试验对比 转速(r/min) 线速度(m/s) (°) 26.82 39.33 52.0 65.43 0.19 0.28 0.38 0.47 排种性能(°) 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 清种起始角 101.8 清种终止角 142.2 135.3 4.92% 145.3 149 2.55% 169.3 161 4.90% 185.3 178 充种角 18.4 18.7 1.63% 19.2 19 1.04% 19.2 21 9.38% 20.6 22.5 9.22% 投种角 30.8 29.3 4.87% 30.2 28.7 4.97% 27.4 28 2.19% 26.8 26.3 1.87% 99 2.65% 105.0 101.8 3.05% 116.0 108.8 6.21% 120.8 111.8 7.45% 3.94% 表 55 组合内窝孔精密排种器播玉米种子时清种角的仿真与试验对比 清种开始角 排种性能 转速(r/min) 清种终止角 试验值(°) 仿真值(°) 误差(%) 试验值(°) 仿真值(°) 误差(%) 13.85 92.3 83.9 9.1% 119.5 107.2 10.3% 26.82 95.0 90.4 4.8% 137.4 124.3 9.5% 39.33 98.4 94.2 4.3% 150.3 139.1 7.5% 51.96 101.2 97.3 2.9% 167.6 156.1 6.9% 表 56 组合内窝孔精密排种器播小麦种子时清种角、投种角和种子面倾角仿真与试验对比 (°) 26.82 转速(r/min) 排种性能(°) 仿真值 试验值 清种起始角 86.5 91.5 清种终止角 134.6 140 投种角 30.3 31 39.33 误差 5.46% 3.86% 2.26% 仿真值 试验值 93.2 95 145.7 157 28 29.5 51.96 误差 1.89% 48 7.20% 5.08% 仿真值 试验值 97.3 100 169 171 26.2 272 65.43 误差 2.70% 1.17% 90.37% 仿真值 试验值 误差 100 104 186.7 196 24.6 25.8 3.85% 4.74% 4.65% 种子面倾角 30.7 31.4 2.23% 32.5 32.3 32.8 0.62% 32.6 34 0.61% 33.5 1.49% 表 57 组合内窝孔精密排种器播大豆时单粒种子合速度模拟值与试验值对比(m/s) 转速(r/min) 线速度(m/s) 26.82 39.33 52.0 65.43 0.19 0.28 0.38 0.47 单粒速度 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 充种点 0.208 0.261 25.48% 0.312 0.334 0.327 0.376 14.98% 0.329 0.413 25.53% 群体点 0.117 0.068 41.88% 0.103 0.091 11.65% 0.109 0.075 31.19% 0.118 0.105 11.02% 护种点 0.207 0.198 4.35% 0.328 0.311 5.18% 0.634 0.582 8.20% 0.756 0.801 5.95% 投种点 0.232 0.238 2.59% 0.332 0.359 8.13% 0.427 0.453 6.09% 0.52 0.547 5.19% 7.05% 表 58 组合内窝孔精密排种器播玉米时单粒种子合速度模拟值与试验值对比(m/s) 13.85 26.82 39.33 51.96 群体点 清种点 0.089 0.261 0.083 0.379 0.097 0.521 0.095 0.629 护种点 0.250 0.386 0.517 0.686 投种点 0.148 0.332 0.487 0.671 群体点 清种点 0.065 0.260 0.072 0.391 0.074 0.583 0.069 0.692 护种点 0.235 0.396 0.518 0.704 投种点 0.136 0.309 0.451 0.623 群体点 清种点 36.9 0.4 15.3 3.1 31.1 10.6 37.7 9.1 护种点 6.4 2.5 0.2 2.6 投种点 4.4 7.4 8.0 7.7 转速(r/min) 仿真值(m/s) 试验值(m/s) 相对误差(%) 表 59 组合内窝孔精密排种器播小麦时单粒种子合速度模拟值与试验值对比(m/s) 26.82 转速(r/min) 39.33 51.96 65.43 分析种子 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 误差 充种点 0.24 0.25 5.27% 0.28 0.29 3.11% 0.33 0.31 7.28% 0.36 0.33 10.3% 群体点 0.24 0.21 9.46% 0.35 0.32 9.05% 0.39 0.34 14.7% 0.42 0.32 34.4% 护种点 0.25 0.28 9.75% 0.34 0.33 1.88% 0.44 0.45 2.89% 0.53 0.5 5.41% 投种点 0.26 0.26 1.57% 0.39 0.4 3.07% 0.51 0.52 1.18% 0.65 0.61 6.15% 表 60 型孔轮式精密排种器播大豆种子时粒数的仿真与试验对比 性 能 转速 指 (r/min) 标 13.85 单粒率(%) 双粒率(%) 空穴率(%) 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 误差 92.86 91.67 1.27 0 0 0 7.14 8.33 14.87 49 20.42 85 88.9 4.3 0 0 0 15 11.1 26 26.82 88.48 95.45 7.9 0 0 0 11.52 4.55 60 表 61 型孔轮式精密排种器播玉米种子时粒数的仿真与试验对比 性 能 转速 指 (r/min) 标 单粒率(%) 双粒率(%) 仿真值 试验值 误差 仿真值 试验值 13.85 80.7 87.18 7.3 0 20.42 76.7 85.08 9.8 26.82 71.43 81.29 12 空穴率(%) 误差 仿真值 试验值 误差 7.69 19.3 5.13 72.5 0 6.82 23.3 8.10 66.1 0 5.81 28.57 12.90 55.9 表 62 型孔轮式精密排种器播大豆种子时种子面倾斜角的仿真与试验对比 转速(r/min) 13.85 20.42 26.82 试验值(°) 26.5 25.8 20.8 仿真值(°) 23.8 18.8 15.0 相对误差(%) 10.2 27.1 27.9 表 63 型孔轮式精密排种器播玉米种子时种子面倾斜角的仿真与试验对比 转速(r/min) 13.85 20.42 26.82 试验值(°) 30.1 28.2 26.5 仿真值(°) 28.1 25.9 24.1 相对误差(%) 6.6 8.2 9.1 表 64 窝眼轮式精密排种器播玉米种子时粒数的仿真与试验对比 性 能 转速 指 (r/min) 标 单粒率(%) 双粒率(%) 空穴率(%) 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 13.85 2.2 6.39 59.8 97.8 72.37 25.9 0 21.24 26.82 2.6 11.5 74.1 97.4 72.13 25.7 0 16.37 50 误差 39.33 20 16.13 18.5 80 71.17 11 0 12.10 51.96 12.5 21.43 40.5 87.5 71.43 18.26 0 7.14 表 65 窝眼轮式精密排种器播大豆时单粒种子合速度仿真与试验对比 转速(r/min) 线速度(m/s) 13.85 26.82 39.33 52.0 0.06 0.11 0.16 0.22 单粒速度 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 试验值 模拟值 误差 投种(m/s) 0.307 0.293 4.56% 0..391 0.363 7.16% 0.444 0.430 3.15% 0.476 0.510 7.14% 拖带层 1(m/s) 0.025 0.030 20.00% 0.045 0.049 31.53% 0.038 0.056 47.37% 0.088 0.105 19.32% 拖带层 2(m/s) 0.025 0.039 56.00% 0.035 0.054 54.29% 0.064 0.087 35.94% 0.125 0.120 4.00% 拖带层 3(m/s) 0.009 0.013 44.44% 0.029 0.039 34.48% 0.033 0.052 57.58% 0.044 0.069 56.82% 表 66 槽轮排种器播小麦种子时种子流量的仿真值与试验值比较 转速(r/min) 仿真值(g/s) 试验值(g/s) 相对误差(%) 26.82 3.067 5.231 31 39.33 5.149 7.781 34 51.96 6.637 10.214 35 65.43 8.158 12.736 36 表 67 第一种槽轮排肥器仿真值与试验值比较 转速(r/min) 线速度(m/s) 13.42 26.35 32.46 44.36 0.15 0.24 0.36 0.43 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 1.516 1.286 15.13% 1.639 1.396 14.83% 1.780 1.489 16.32% 1.817 1.557 14.31% 进肥点单粒速度(m/s) 0.224 0.281 25.02% 0.242 0.298 23.29% 0.277 0.345 24.45% 0.311 0.392 26.08% 槽轮点单粒速度(m/s) 0.424 0.453 6.80% 0.503 0.632 25.59% 0.509 0.653 28.35% 0.638 0.799 25.25% 群体点单粒速度(m/s) 1.207 1.671 38.38% 1.288 1.637 27.13% 1.294 1.362 1.296 1.686 30.10% 排肥点单粒速度(m/s) 2.079 2.359 13.48% 2.217 2.534 14.32% 2.229 2.306 14.26% 2.256 2.564 13.67% 质量流率(g/s) 5.26% 表 68 第二种槽轮排肥器仿真值与试验值比较 转速(r/min) 线速度(m/s) 13.42 26.35 32.46 44.36 0.15 0.24 0.36 0.43 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 质量流率(g/s) 1.436 进肥点单粒速度(m/s) 0.237 槽轮点单粒速度(m/s) 0.420 1.427 4.16% 1.549 0.254 5.49% 0.267 1.513 2.83% 1.6370 1.587 5.24% 1.747 0.279 4.94% 0.321 0.284 4.22% 0.336 1.707 3.31% 0.310 4.46% 0.441 3.50% 0.504 0.521 3.91% 0.516 0.534 3.98% 0.610 0.634 3.82% 群体点单粒速度(m/s) 1.213 排肥点单粒速度(m/s) 2.228 1.253 7.98% 1.318 1.387 9.75% 1.648 1.711 9.72% 1.777 1.844 8.62% 2.232 2.75% 2.341 2.324 2.99% 2.354 2.362 2.64% 2.422 2.484 3.11% 51 表 69 第三种槽轮排肥器仿真值与试验值比较 转速(r/min) 线速度(m/s) 13.42 26.35 32.46 44.36 0.15 0.24 0.36 0.43 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 试验值 仿真值 误差 1.627 1.528 4.13% 1.687 1.643 4.82% 1.730 1.635 8.32% 1.767 1.732 4.26% 进肥点单粒速度(m/s) 0.246 0.261 4.57% 0.271 0.281 4.08% 0.340 0.354 4.65% 0.351 0.361 4.87% 槽轮点单粒速度(m/s) 0.437 0.449 5.22% 0.511 0.533 4.81% 0.527 0.538 5.63% 0.616 0.635 5.37% 群体点单粒速度(m/s) 1.238 1.309 8.39% 1.322 1.380 8.04% 1.655 1.705 8.63% 1.761 1.784 8.02% 排肥点单粒速度(m/s) 2.240 2.246 3.50% 2.326 2.338 3.68% 2.444 2.453 3.03% 2.535 2.552 3.18% 质量流率(g/s) 表 70 含水率为 22%、开沟深度为 30mm 时芯铧开沟器工作阻力的二维模拟值与试验值对比 水平作用力 工作速度(mm) 试验值(N) 模拟值(N) 差值(N) 误差(%) 0.2 51.5 40.95 10.55 25.76 0.5 55.4 44.1 11.3 25.62 0.8 56.9 45.9 11.0 23.97 垂直作用力 0.2 13.6 10.8 2.8 25.92 0.5 14.3 11.7 2.6 22.22 0.8 15.6 12.6 3.0 23.81 表 71 含水率为 22%、工作速度为 0.2m/s 芯铧开沟器工作阻力的二维模拟值与试验值对比 水平作用力 开沟深度(mm) 试验值(N) 模拟值(N) 差值(N) 误差(%) 30 51.5 40.95 10.55 25.76 60 76.6 64.35 12.25 19.04 90 127.7 106.65 21.05 19.73 垂直作用力 30 13.6 10.8 2.8 25.92 60 21.5 17.1 4.4 25.73 90 33.3 28.8 4.5 15.63 表 72 开沟深度 50mm、工作速度 0.2m/s 时施肥开沟器工作阻力的三维模拟与试验相比 水平作用力 含水率（%） 试验值(N) 模拟值(N) 差值(N) 误差(%) 3 22.1 20.54 1.56 7.6 12 23.4 22.43 0.97 4.32 22 38.8 33.7 5.1 15.13 52 垂直作用力 3 5.6 6.13 0.53 9.46 12 6.3 6.99 0.69 10.95 22 9.1 10.65 1.55 17.1 由比较可知：对于排种器来说，其排种量当种子较规整时如大豆种子，误差较小小 于 20%，当种子不规整时误差较大 10~50%；排种器充种、清种和投种性能仿真分析及 某一确定单粒种子的速度和位移仿真分析，与试验结果的误差在 5~20%之间，而群体内 某一种子的速度和位移仿真分析，与试验结果的误差较大在 10~50%之间；对于排肥器 来说，其排肥量的仿真分析，与试验结果的误差小于 20%，而群体内某一肥料颗粒的速 度和位移仿真分析，与试验结果的误差较大在 10~50%之间；对于开沟器来说，其工作 阻力的二维仿真分析，与试验结果的误差小于 30%。由于开沟器仿真时，所需计算的土 壤颗粒数较多，采用微机计算时间较长，因此本研究大量采用二维仿真，只进行了较少 数量的三维仿真参见表 72。上述工作初步证明了本研究所建立方法和软件的可行性。 本研究还研究了分别选取球颗粒模型、椭球颗粒模型和多面体颗粒模型及不同尺寸 参数时，对离散元法分析结果的影响。由研究可知，选取不同颗粒模型时，对离散元法 的分析结果影响较大。 本研究还研究了离散元法计算时，接触力学模型中的刚度系数、弹性模量、泊松比、 碰撞恢复系数和动静摩擦系数等，选取不同数值时对离散元法分析结果的影响，见表 73 和表 74 所示。由研究可知，在本研究的范围内，参数变化对离散元法的分析结果影响 不大。 表 73 刚度系数变化时组合内窝孔精密排种器播小麦的模拟值与试验值对比(m/s) 法向刚度系数 (N/m) 25000 1967 仿真值 试验值 相对误差 仿真值 试验值 相对误差 清种起始角(°) 86.3 95 9.16% 85.67 95 9.82% 清种终止角(°) 151 157 3.82% 146.67 157 6.58% 投种角(°) 27.3 29.5 7.46% 28 29.5 5.08% 排种量(粒) 18 30 40% 19 30 36.7% 表 74 刚度系数变化时槽轮排种器播小麦的模拟值与试验值对比(m/s) 刚度系数(N/m) 种子流量(g/s) 5x104 11.16x104 22.31x104 39.34x104 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 误差 仿真值 试验值 误差 5.72 29.44% 7.78 26.48% 5.58 7.78 28.28% 53 5.65 7.781 27.39% 5.49 7.781 本研究还采用新方法和软件进行了 2 种新型排种器（组合内窝孔精密排种器和垂直 圆盘型孔轮式排种器）的设计和工作过程仿真分析，见图 28~图 30 所示。本研究还采用 该方法和软件，通过改变力学模型，还进行了气力排种器、分离清选和玉米小麦脱粒等 其它种类农机工作部件的应用研究。 8 结论 （1）本研究组提出由机械部件的 CAD 模型（CAD 软件设计图）建立其三维离散元 法分析模型的方法和实边界与虚边界的方法，在对二维 CAD 软件 AutoCAD 和三维 CAD 软件 PRO/E 进行二次开发的基础上，通过人机交互实现了与散粒物料接触作用的图元 识别、参数读取、运动属性及材料特性参数绑定，在此基础上还实现了二维和三维 CAD 软件与自主研制的三维离散元法分析计算软件集成，从而构建了一种集设计与性能分析 评价于一体的新型 CAE 软件，其体系结构见图 2； （2）该新型 CAE 软件可在微机环境下运行，可由机械部件的二维或三维 CAD 模 型生成其三维离散元法分析模型，可生成多种几何模型的散粒物料，可采用多种力学模 型进行接触力的计算，通过人机交互方式输入计算参数，还具有颗粒运动速度场和作用 力场分析、单个颗粒运动速度和位移分析、颗粒群体流量和机械部件工作阻力分析、机 械部件工作过程动态仿真显示等功能； （3）以 3 种玉米种子（吉峰 218、长单 228 和平安 18） 、3 种大豆种子（吉科豆 1 号、吉新豆 1 号和吉豆）、3 种小麦种子（吉 12、丰强 11 和小冰麦 33）、3 种颗粒肥料 （北京生物、大庆 3 和山东 4） 、3 种土壤（含水率为 3%、12%和 22%的壤土）和 5 种 排种器（组合内窝孔精密排种器、垂直圆盘型孔轮式排种器、窝眼轮式玉米排种器、窝 眼轮式大豆排种器、槽轮式排种器） 、3 种排肥器（3 种结构外槽轮式排肥器）和 3 种开 沟器（芯铧式开沟器、施肥开沟器和船铲式开沟器）为对象，对上述 15 种散粒物料与 11 种农机工作部件的接触作用过程和工作过程进行了仿真分析，并将仿真结果与农机工作 部件工作过程的实验结果进行了对比，初步建立了农机工作部件的数字化设计方法； （4）采用该新方法设计并研制出 2 种新型排种器（组合内窝孔精密排种器和垂直圆 盘型孔轮式排种器），已有 4 个单位（①吉林农业大学工程技术学院；②吉林大学生物 与农业工程学院；③华南农业大学工程技术学院；④吉林省恒源农业机械有限公司）使 用该方法和软件； 主要技术成果及创新点 54 （1）提出由机械部件的三维 CAD 模型（CAD 软件设计图）建立机械部件的三维离 散元法分析模型的方法和实边界与虚边界的方法； （2）提出由机械部件的二维 CAD 模型（CAD 软件设计图）建立机械部件的三维离 散元法分析模型的方法； （3）采用数据库共享实现了 CAD 软件与自主研制的三维离散元法分析计算软件集 成，从而构建一种集设计和性能分析评价于一体的新型 CAE 软件； （4）分别采用超二次方程（|x/a|n1+|y/b|n2+|z/c|n3-1=0）和球颗粒组合的方法，建立了 玉米种子和小麦种子的离散元法分析模型，这方面的工作尚未见报道； （5）提出基于进料口的颗粒样本生成方法。 （6）提出离散元法软件中基于图元的边界类划分； （7）采用 OpenGL 和图形的伸缩变换，实现了椭球颗粒的运动过程动态显示； （8）系统地研究了三种散粒物料玉米种子、大豆种子和小麦种子的物理力学参数及 这些参数对离散元法分析结果的影响； （9）改进了湿颗粒间接触力学模型，以用于土壤颗粒的分析； （10）采用力和扭矩控制边界的方法，以分析机械部件中存在的弹性边界； （11）采用力控制边界的方法，以模拟土壤的坚实度； （12）建立了多种属性颗粒和多种属性边界的计算方法； （13）首次采用离散元法分析排种器的工作过程； （14）首次采用离散元法分析农业土壤——壤土与实际开沟器的工作过程； （15）初步建立了一种农机工作部件（排种器、排肥器和开沟器等）优化设计和虚 拟设计的新方法和及其软件等。 未来前景 目前进行与散粒物料接触作用的机械部件设计时，通常是根据设计要求和设计者的 经验拟定结构方案并进行详细设计，设计完成后再进行样机试制和试验研究，接着再根 据试制和试验中发现的问题修改设计及进行新一轮的样机试制和试验。一般一个新产品 从设计到最后定型，往往需要数十次甚至上百次的试验、修改才能完成，不仅费时、费 力，而且开发成本较高，还往往得不到理想的效果。 改变设计过程中这种被动局面，是设计者多年的梦想。科学技术的发展为我们提供 了一个从根本上改变这一局面的良机，那就是建立基于动力学分析的数字化设计方法， 55 这已被成功开发和应用的实例所证明。 目前国外有 2 个公司出售离散元法分析软件，美国的 ITASCA 和英国的 DEM Solutions 公司，价格在 50~100 万人民币不等，国内有多家单位已经购买。 本研究研制的新型 CAE 软件，不仅可以用于新型机械产品的开发，还可用于已有产 品的性能分析，不仅可以替代进口，而且还有望出口国外，其应用前景、经济效益和社 会效益是不可估量的。 正像基于连续介质力学的有限元法及其应用软件，已广泛应用于科学研究和工程设 计的各个领域一样。基于颗粒动力学的离散元法及其应用软件，也将使与有关的科学研 究和工程设计领域产生前所未有的变化，其应用前景是非常光明的。 研究中的主要问题、经验和建议 （1）由于本研究研制的是一个原型系统，因此在软件的结构方面还有许多问题需要 深入研究解决； （2）由于本研究研制的是一个原型系统，因此在软件的算法方面还有许多问题需要 深入研究解决； （3）当计算颗粒数较多时，计算时间还较长，如何提高计算速度还需深入研究。 56 附录 A 发表论文情况 发表研究论文 18 篇，其中 14 篇已被 EI 收录，分别为： （1）于建群，王刚，心男，付宏．型孔轮式排种器工作过程及性能的仿真分析．农 业机械学报，2011，42（12）：83~87(EI20120214672702)； （2）李政权，于建群，张尉林，付宏．内充式排种器工作过程和性能的离散元法仿 真分析．农业工程学报，2011，27（11）：32~36(EI20114714541922)； （3）于亚军，于建群，陈仲，付宏．三维离散元法边界建模软件设计．农业机械学 报，2011，42（8）：99~103(EI20113614302307)； （4）于建群，钱立彬，于文静等．开沟器工作阻力的离散元法仿真分析．农业机械 学报，2009，40(6)：53-57(EI20092912201619)； （5）冯占荣，于建群，付宏，卢洪涛．一种二维贴体网格改进算法的研究．西安交 通大学学报，2009，43(3)：105-109(EI20091512029456)； （6）付宏，乌兰，黄万风，于建群．基于图元的三维离散元法边界建模方法．计算 机集成制造系统，2008，14(12)：2328-2333(EI 090611901757)； （7）付宏，贾慧敏，张晓旭，于建群．基于超圆颗粒模型的二维离散元法计算方法 研究．吉林大学学报(工学版)，2008，38(6)：1383-1388(EI 085011777038)； （8）于建群，申燕芳，牛序堂，付宏．组合内窝孔精密排种器清种过程的离散元法 仿真分析．农业工程学报，2008，24(5)：105-109(EI082811368674)； （9）于建群，付宏．基于离散元法的数字化设计软件及其应用．中国农业机械学会 2008 年学术年会，山东济南，2008 年 9 月； （10）Yu Jian-qun，Fu Hong，Sun Yu-jing．Development and its Application of A Digital Design Software Based on Discrete Element Method．15th International Conference on Mechatronics and Machine Vision in Practice(M2VIP08)，2-4 Dec. 2008，Auckland， New-Zealand(EI 20091211965415)； （11）贾慧敏，付宏，张磊，于建群．基于 CAD 模型的边界离散元法仿真算法研 究．系统仿真学报，2007，19(20)：4607-4611(EI 074610916604)； （12）于建群，付宏．基于离散元法的数字化设计软件及其应用．第三届中国 CAE 工程分析技术年会，大连，2007 年 7 月；该文在本届年会被评为优秀论文； 57 （13）颜辉，于建群，付宏，寇幸幸．离散元法软件与 CAD 软件集成开发研究．现 代情报，2007，9：209-211； （14）于建群，付宏，李红，申燕芳．离散元法及其在农业机械工作部件研究与设 计中的应用．农业工程学报，2005，21(5)：1-6(EI 05259172437)．该文被评为农业工程 学报 2004-2005 年度优秀论文； （15）付宏，商慧，于建群．三维离散元法软件开发研究．华中科技大学学报(自然 科学版)，2005，33(增刊)：342-345(EI 06109747622)； （16）付宏，董劲男，于建群．基于 CAD 模型的离散元法边界建模方法．吉林大 学学报(工学版)，2005，35(6)：626-631(EI 06029637182)； （17）付宏，董劲男，刘振宇，危夏，于建群．精密排种器的集成分析设计系统开 发研究．计算机集成制造系统，2005，11(5)：744-750(EI 05289208169)； （18）滕安翔，付宏，于建群．OpenGL 下椭球的近似表示算法．吉林大学学报(信 息科学版)，2005，23(增)：154-156。 6）获国家发明专利 2 项，分别是： （1）组合内窝孔精密排种器，专利号 200410011100.0，授权时间 2007.10.17； （2）基于 CAD 模型的离散元法边界建模方法，专利号 200510016835.7，授权时间 2009.4.29。 7）申报国家发明专利 3 项，分别是： （1）一种离散元法边界建模的通用方法，专利申请号 201110303054.1； （2）基于离散元法的玉米脱粒过程分析方法，专利申请号 201110293826.8； （3）内充式玉米排种器，专利申请号 201110276013.8。 8）获软件版权登记 3 项，分别是： （1）基于离散元法的三维 CAE 软件，登记号 2008SR21750； （2）基于离散元法的二维 CAE 软件，登记号 2009SR02179； （3）基于二维 CAD 模型的三维离散元法边界建模软件，登记号 2009SR03204。 9）培养毕业硕士生 28 名见表 75。目前在读硕士生 16 名、在读博士生 6 名，全部 完成了课题任务规定的研究内容和考核目标。 58 表 75 已毕业硕士研究生名细 序 研究生姓 论文题目 号 1 名 张列南 2 专业 毕业时间 基于离散元法的大豆精密排种器 农 机 2009 年 6 月 杨明芳 的数字化设计方法研究 基于离散元法的玉米精密排种器 农 机 2009 年 6 月 3 李政权 的数字化设计方法研究 基于离散元法的小麦精密排种器 农 机 2009 年 6 月 4 陈 斌 的数字化设计方法研究 基于椭圆颗粒模型的离散元法算法 计算机 2009 年 6 月 5 武广涛 研究 组合颗粒模型的离散元法算法研究 计算机 2009 年 6 月 6 潘克峰 与软件设计 几种三维离散元法边界的建模方法 计算机 2009 年 6 月 7 钱立彬 研究 基于离散元法的开沟器的数字化 农 机 2008 年 6 月 8 吕 昊 与软件设计 设计方法研究 基于离散元法的排肥器的数字化 农 机 2008 年 6 月 9 刘 睿 设计方法研究 基于 AutoCAD2007 软件的三维离散 农 机 2008 年 6 月 10 李志勇 元法边界建模方法研究 基于椭球颗粒模型的离散元法基本 计算机 2008 年 6 月 11 李百青 理论及算法研究 基于非完整图元的三维离散元法边 计算机 2008 年 6 月 12 宋 鹗 界的计算方法研究 基于椭圆颗粒模型的离散元法基本 计算机 2008 年 6 月 13 卞小兵 理论及算法研究 基于二维 CAD 模型的边界离散元法 计算机 2008 年 6 月 14 寇幸幸 建模方法研究 基于 Pro/ENGINEER 软件的三维离 农 机 2007 年 6 月 15 杨霄鹏 散元法边界建模研究 基于动态八叉树的三维离散元法邻 计算机 2007 年 6 月 16 贾慧敏 居搜索方法研究 二维超圆形颗粒离散元法计算 计算机 2007 年 6 月 17 张 磊 仿真软件研制 圆盘颗粒模型的二维离散元仿真系 计算机 2007 年 6 月 18 么 鑫 统 模型的三维离散 基于边界二维 CAD 计算机 2007 年 6 月 19 乌 兰 的研究与开发 元法边界建模方法 三维离散元法计算仿真软件的研究 计算机 2007 年 6 月 20 许志宝 基于离散元法的大豆碰撞过程仿真 农 2006 年 6 月 分析 59 机 21 滕安翔 基于 AutoCAD 软件的三维离散元法 计算机 2006 年 6 月 22 商 慧 边界建模方法研究 三维离散元法计算仿真软件开发研 计算机 2006 年 6 月 23 申燕芳 究 基于离散元法的精密排种器设计 农 机 2005 年 6 月 24 董劲男 与工作过程的仿真分析 基于面向对象技术的离散元法 计算机 2005 年 6 月 25 王秋雁 分析设计软件开发研究 采用 OpenGL 实现机械部件的 计算机 2005 年 6 月 26 李 红 离散元法仿真分析 精密排种器的数字化设计与 农 机 2004 年 6 月 27 刘振宇 工作过程仿真分析 基于离散元法的精密排种器 计算机 2004 年 6 月 28 危 分析设计软件开发研究 基于 AutoCAD 二次开发的精密排种 计算机 2004 年 6 月 夏 器 分析设计软件开发研究 60