基于光滑粒子流体动力学的旋耕切土数值模拟

2023-01-11钟嘉雄任述光吴明亮

湖南农业大学学报（自然科学版） 2022年6期

钟嘉雄，任述光,2*，吴明亮,2

钟嘉雄1，任述光1,2*，吴明亮1,2

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

为明确ⅠT245型旋耕刀(刀座式)切割重黏土过程，基于光滑粒子流体动力学(SPH)方法对旋耕刀–土壤进行建模，运用ANSYS/LS–DYNA软件进行仿真，分析单个旋耕刀在刀盘转速200 r/min、机具前进速度0.6 m/s、耕深0.12 m工况下的切土扭矩情况，仿真平均扭矩为5.84 N·m，以库伦定理和朗肯土压力理论为依据计算的理论平均扭矩为6.42 N·m，二者误差为8.98%。分别使用理论方法和有限元方法，以机具前进速度、刀盘转速和耕深为影响因素，以比功为评价指标，进行正交旋转组合试验，得到旋耕比功与机具前进速度、刀盘转速及耕深的回归方程，表明在符合农艺要求的情况下，应尽可能选择较大的机具前进速度和耕深，较小的刀盘转速。

刀座式旋耕刀；土壤切削；数值模拟；光滑粒子流体动力学

旋耕作业具有较强的切土和碎土能力，一次作业即可完成翻、耙、耪、平等耕地平整全过程[1–3]。旋耕刀是旋耕机主要工作部件，刀片的形状和参数，直接影响旋耕机的工作质量和功率消耗[4–5]。旋耕刀种类较多，目前国内以刀座式旋耕刀为主，研究旋耕刀切土的工作过程，对优化旋耕刀的结构参数，降低耗能有着重要意义[6–8]。采用田间试验方法对刀具–土壤相互作用的测定费时费力，试验结果受试验条件限制，而计算机数值模拟则可以克服试验条件差异，同时可以定量分析作用过程中各种动力参数，诸如作用力、作用力矩、功耗等随时间的变化，动态演示切土过程中土壤变形碎裂的进程。

数值模拟方法主要有离散元法、有限单元法、有限体积法与差分法等。离散元法在土壤切削应用的最大困难是如何在计算机中制造出符合真实土壤体力学特性的离散元堆积体[9]，而有限单元法在模拟高速切割破坏的过程中网格在大变形状况下会发生畸变，导致计算失败。

BOJANOWSKI[10]用有限元拉格朗日网格法(FE)、无网格伽辽金法(EFG)、光滑粒子流体动力学法(SPH)[11]以及多材料任意拉格朗日欧拉法(MM–ALE)进行了带砂粉质黏土的大型钢板贯入试验以及砂土的标准圆锥贯入试验，结果表明SPH法在模拟土体大变形中的表现稳健可靠。高建民等[12]研发了基于SPH的土壤高速切削仿真系统，并模拟了R245旋耕刀潜土逆旋切土的过程，展示了模拟结果，为高速切削数值模拟提供了新思路。LIMIDO等[13]用SPH方法在二维空间内模拟金属高速切削加工过程，切向和法向切削力的模拟结果与试验测量误差分别在15%和35%以内。韩煜杰等[14]基于SPH算法对立式旋耕刀进行土壤切削模拟，仿真平均扭矩为46.97 N·m，与试验测得数据的误差(34.4%)在允许范围内。SPH方法被认为是土壤大变形和破碎研究的理想方式[15–16]。笔者对刀座式ⅠT245型旋耕刀和重黏土进行建模，采用SPH方法进行旋耕切土的显式动力计算，探究旋耕切土过程中刀片与土壤的相互作用，以旋耕刀运动参数为试验因素，比功为评价指标，进行正交旋转组合试验，得出旋耕切土运动参数与比功的关系。

1　旋耕刀和土壤模型的建立

1.1　旋耕刀模型

ⅠT245型标准旋耕刀根据GB/T 5669—2017《旋耕机械刀和刀座设计》，采用SOLIDWORKS绘制三维图形，对旋耕刀的部分地方作了圆角处理，使过渡处更加平滑。为简化计算，并有更好的网格质量，去除了螺栓、螺母、垫片等复杂零部件。将设计的旋耕刀保存为STEP格式文件，导入前处理器LS–PREPOST中，使用实体单元，对刀片进行网格划分，刀片部分划分为3672个单元。由于旋耕刀弹性模量远大于土壤的弹性模量，为了减少计算量，刀的材料模型设置成*MAT_RIGID，即刚体材料，密度为7850 kg/m3，弹性模量为1.96×1011 Pa，泊松比为0.3。

1.2　土壤模型

土壤几何模型采用LS–PREPOST[17]建模，运用SPHGEN建立900×120×180 mm的立方体，使用SPH默认单元类型，一共720 000颗粒子。通过对LS–DYNA中可用土壤模型的比较，选用最适合模拟耕作土壤的*MAT_FHWA_SOIL模型[18]。使用土壤模型的主要参数包括土壤的密度、容重、体积模量、剪切模量、摩擦角、内聚力和水的密度等。用环刀法测得土壤容重，比重瓶法测得密度，燃烧法测得土壤含水率，通过三轴试验测得体积模量、剪切模量、摩擦角和内聚力。经测定，取自湖南农业大学工程实训中心土槽土壤的容重为1680 kg/m3、密度为2550 kg/m3、体积模量为4000 kPa、剪切模量为857 kPa、内聚力为38.55 kPa、内摩擦角为9.97°、含水率为22.4%。

1.3　旋耕刀–土壤接触模型的确定

设定土壤是一组具有密度、容重、体积模量、剪切模量等物理性质的粒子，旋耕刀以面的形式与土壤接触。考虑粒子与面之间的接触时，通常优先使用NODES_TO_ SURFACE的接触形式，再考虑接触类型。当应用侵蚀接触类型(ERODING)时，模型外部单元失效删除后，剩下的单元依旧能够形成新的表面继续处理接触，但该接触会消耗极大的算力；而在与自动接触类型(AUTOMATIC)进行旋耕切土对照试验时，发现在相同算力的服务器中计算时间大大缩短，耗时约为侵蚀接触类型的1/12，而计算结果基本一致。因而采用接触搜索算法。综合算力和效果，使用*CONTACT_AUTOMATIC_NODES_TO_ SURFACE接触模型，其中静摩擦系数与动摩擦系数通过摩擦系数试验测得，分别为0.44和0.24。

2　旋耕刀切土数值模拟

2.1　运动条件

以刀盘旋转中心为坐标原点，机具前进方向为轴方向，轴铅锤向上，刀尖运动方程为：

式中：为旋耕刀回转半径；为刀轴角速度；为切土时间；为旋耕机前进速度。

旋耕刀的运动可分解为旋转和移动2个部分。旋耕刀被设置成刚体，使用*BOUNDARY_PR ESCRIBED_MOTION_RIGID对其运动进行控制。

2.2　其余边界条件及求解设置

调整土壤模型与旋耕刀的距离，设置耕深。由于数值模拟中土壤建模区域远小于实际大小，为使土壤模型更符合实际，在土壤两侧添加无反射边界条件，在土壤底部添加固定约束。设置求解时间0.4 s、*DATABASE_ASCII_option中需要输出的值和频率后，得到完整的仿真模型(图1)，保存为K格式文件，递交ANSYS/LS–DYNA Slover求解。

图1　旋耕刀–土壤仿真模型

2.3　算例分析

以刀盘转速200 r/min、机具前进速度0.6 m/s、耕深0.12 m为例，图2为SPH数值模拟旋耕刀切土扭矩，可以看到切土过程共有3个波形，对应于旋耕刀3次切土过程，第1、2个波形的持续时间明显比第3个波形的持续时间长，这是由于前2次切土需要切割截面为弓形的土壤，第3次切割才是正常工况下的切土过程。

图2　SPH数值模拟旋耕刀切土扭矩

分析第3次切土过程，在LS–PREPOST中截取部分土壤进行观察，取宽度0.015 m的土层，以旋耕刀即将与土壤接触时刻作为分析的初始时刻(图3–a)，旋耕刀开始运动，逐渐与土壤接触，最大等效应力迅速增大，在0.02 s左右(图3–b)时达到最大值，随后最大等效应力逐渐减小(图3–c、图3–d、图3–e)，至0.08 s左右(图3–f)时完成切土过程。

图3　旋耕刀切土过程应力云图

根据旋耕刀的转速、机具前进速度、耕深以及土壤参数等，结合库伦定理和朗肯土压力理论[19]，设计Matlab程序，计算旋耕切土的扭矩，并与仿真数据作对比，结果如图4所示，SPH数值模拟的扭矩值从旋耕刀与土壤接触开始迅速增大，在0.01 s～0.03 s趋于平稳，后逐渐减小，与切土过程的分析结果一致。对比仿真曲线与理论曲线，在0.04 s以后，数值模拟曲线扭矩明显小于理论曲线的扭矩，这是因为在模拟计算中，土壤在切土后期已经发生了变形，土壤结构发生变化，所以扭矩下降速度较快；而理论分析中仅考虑了旋耕刀与土壤接触的面积，不能很好地反映土壤内部结构的变化。分别求得数值模拟与理论旋耕切土扭矩平均值为5.84 N·m和6.42 N·m，数值模拟数据比理论计算数据低8.98%。

图4　旋耕刀切土扭矩数值模拟与理论结果对比

3　旋耕切土运动参数的优化

选取机具前进速度、刀盘转速和耕深3个因素，进行三因素三水平正交旋转组合试验设计(表1)。分别用MATLAB方法和有限元方法进行正交旋转组合试验。

表1　正交旋转组合试验设计因素与水平

用比功来评价旋耕切土的能耗。

式中：为比功；为旋耕所需功率；为耕作幅宽；为耕深。采用Design–Expert 13对表2中的试验数据进行拟合回归，分别获得MATLAB理论计算的比功(1)和有限元仿真的比功(2)对机具前进速度(1)、刀盘转速(2)和耕深(3)的回归方程：1=10.20–4.151+0.0102–7.003，2=–1.59–14.331+0.0912–6.003。

对比发现，理论计算与仿真结果间误差较大，且随着刀盘转速增大，仿真数值的增速明显高于理论数值的，这是因为仿真中随着转速增加，应变率随之增大，土壤发生局部硬化，间接使切土阻力增大；而MATLAB计算中，土壤物理性质比较恒定，因而比功变化较小。从趋势上来看，比功与机具前进速度和耕深呈负相关，与刀盘转速呈正相关。从降低比功的角度看来，在满足农艺要求的运动参数下，应尽可能选择较大的机具前进速度和耕深及较小的刀盘转速。

表2　正交旋转组合的旋耕刀的比功

4　结论

利用ANSYS/LS–DYNA软件对刀座式旋耕刀切土过程进行仿真分析，应力云图显示，最大剪切应力与根据库伦强度理论计算结果相符。

分别对2种模型进行了正交旋转组合试验，得出了比功与机具前进速度、刀盘转速和耕深的关系，耕深与比功的相关性较低。

仿真结果依赖于边界条件设置及土壤的本构模型，在现有的用有限体积代替无限模型情况下，如何合理设置边界条件，选择合适的本构模型及提高土壤力学性能参数的测量精度是提高仿真结果精度的关键。由于耕作土壤力学性能的复杂性与差异，提高计算结果的精度仍需继续深入研究。

[1] 管春松，崔志超，高庆生，等．双轴旋耕碎土试验台设计与分层耕作试验[J]．农业工程学报，2021，37(10)：28–37．

[2] 杜兆辉，陈彦宇，张姬，等．国内外旋耕机械发展现状与展望[J]．中国农机化学报，2019，40(4)：43–47．

[3] 田耘，金亮，赵亚祥．旋耕机作业性能影响参数的试验研究[J]．中国农机化学报，2016，37(8)：1–4．

[4] 方会敏，姬长英，张庆怡，等．基于离散元法的旋耕刀受力分析[J]．农业工程学报，2016，32(21)：54–59．

[5] 郝建军，于海杰，赵建国，等．楔形减阻旋耕刀设计与试验[J]．农业工程学报，2019，35(8)：55–64．

[6] 孙志全，段洁利，杨洲．旋耕减阻降耗技术研究进展[J]．中国农机化学报，2021，42(1)：37–45．

[7] 杨玉婉，佟金，马云海，等．鼹鼠多趾结构特征仿生旋耕刀设计与试验[J]．农业工程学报，2019，35(19)：37–45．

[8] 肖茂华，汪开鑫，杨望，等．基于东方蝼蛄爪趾的仿生旋耕刀设计与试验[J]．农业机械学报，2021，52(2)：55–63．

[9] 刘春，乐天呈，施斌，等．颗粒离散元法工程应用的三大问题探讨[J]．岩石力学与工程学报，2020，39(6)：1142–1152．

[10] BOJANOWSKI C．Numerical modeling of large deforma- tions in soil structure interaction problems using FE，EFG，SPH，and MM-ALE formulations[J]. Archive of Applied Mechanics，2014，84(5)：743–755．

[11] LIU G R，LIU M B．Smoothed Particle Hydrodyna- mics：AMeshfree Particle Method[M]. Singapore：World Scientific，2003．

[12] 高建民，周鹏，张兵，等．基于光滑粒子流体动力学的土壤高速切削仿真系统开发及试验[J]．农业工程学报，2007，23(8)：20–26．

[13] LIMIDO J，ESPINOSA C，SALAÜN M，et al．SPH method applied to high speed cutting modelling[J]. International Journal of Mechanical Sciences，2007，49(7)：898–908．

[14] 韩煜杰，李云伍，赵华慧，等．基于SPH算法的立式旋耕刀土壤切削仿真模拟[J]．西南大学学报(自然科学版)，2016，38(12)：150–155．

[15] 李云伍，张国勇，张植，等．低功耗小型立轴式深耕机分段螺旋旋耕刀具的研制[J]．农业工程学报，2019，35(4)：72–80．

[16] YANG W，ZHAO W J，LIU Y D，et al．Simulation of forces acting on the cutter blade surfaces and root system of sugarcane using FEM and SPH coupled method[J]. Computers and Electronics in Agriculture，2021，180：105893．

[17] 夏俊芳，贺小伟，余水生，等．基于ANSYS/LS–DYNA的螺旋刀辊土壤切削有限元模拟[J]．农业工程学报，2013，29(10)：34–41．

[18] TAGAR A A，JI C Y，ADAMOWSKI J，et al. Finite element simulation of soil failure patterns under soil bin and field testing conditions[J]．Soil and Tillage Research，2015，145：157–170．

[19] 张克恭，刘松玉．土力学[M]．北京：中国建筑工业出版社，2001．

Numerical simulation of rotary tillage and soil cutting based on smooth particle hydrodynamics

ZHONG Jiaxiong1，REN Shuguang1,2*，WU Mingliang1,2

(1.College of Mechanical and Electrical Engineering, Hunan Agricultural University, Changsha, Hunan 410128, China; 2.Hunan Key Laboratory of Intelligent Agricultural Machinery and Equipment, Changsha, Hunan 410128, China)

In order to clarify the process of cutting heavy clay with the type IT245 rotary cultivator(knife seat type), the rotary cultivator-soil was modeled based on the smoothed particle hydrodynamics(SPH) method. Simulating by ANSYS/LS-DYNA software , the soil-cutting torque of a single rotary cutter was analyzed under the conditions of cutter head speed of 200 r/min, forward speed of 0.6 m/s, and ploughing depth of 0.12 m. The average simulated torque was 5.84 N· m. The theoretical average torque is 6.42 N·m calculated by the Ken earth pressure theory. The error between the two is 8.98%. Using the theoretical method and the finite element method respectively, the rotation-orthogonal combination test was carried out with the forward speed, the rotating speed of the cutter roller and the ploughing depth as the influencing factors, and the specific work as the evaluation index. The regression equation of the specific work of the rotary tiller was obtained dependent on the forward speed of the machine and the speed of the cutter head and the speed of the cutter head. It indicated that the larger forward speed and ploughing depth and the smaller cutter shaft speed should be selected as far as possible in the case of meeting the agronomic requirements.

rotary tiller; soil cutting; numerical simulation; smoothed particle hydrodynamics(SPH)

S222.3

1007-1032(2022)06-0744-05

钟嘉雄，任述光，吴明亮．基于光滑粒子流体动力学的旋耕切土数值模拟[J]．湖南农业大学学报(自然科学版)，2022，48(6)：744–748．

ZHONG J X，REN S G，WU M L．Numerical simulation of rotary tillage and soil cutting based on smooth particle hydrodynamics[J]．Journal of Hunan Agricultural University(Natural Sciences)，2022，48(6)：744–748．

http://xb.hunau.edu.cn