基于離散元法的油莎豆降阻挖掘裝置設計與試驗

何曉寧 張學軍 趙 壯 尚書旗 王東偉 袁小偉

(1.新疆農業大學機電工程學院， 烏魯木齊 830052； 2.青島農業大學機電工程學院， 青島 266109)

0 引言

油莎豆是原產于非洲及地中海沿岸國家的一種草本植物[1-2]，是集油、糧、牧、飼、綠化觀賞于一體的新型經濟作物，具有適應性廣、含油量高、營養價值高等特性[3-4]。目前，我國大多數地區油莎豆收獲以正旋旋耕和挖掘鏟式作業方式為主，但是由于正旋挖掘方式存在埋果率高、壅土嚴重，挖掘鏟式存在挖掘阻力大、土壤破碎率低等問題[5-6]，無法滿足油莎豆規模化收獲要求。國內外學者近年來對旋耕機理進行深入研究，ASAF等[7]研究了旋耕刀動力預測模型，推導了旋耕刀對土壤切削時的支反作用力和扭力方程，闡明了挖掘阻力產生的原理，ROSA等[8]利用有限元方法模擬出刀具耕作時的工作情況，在以不同耕速比為變量前提下，得出刀具不同結構參數對于性能指標的影響，LEE等[9]基于旱地水稻直播機的旋耕部件，研究旋耕刀在不同結構參數和工作參數下對旋耕作業扭矩、碎土率的影響。

文獻[10-11]基于深松部件、旋耕部件作業之間的交互作用，設計了一種用于深耕的反旋深松聯合作業耕整機，相比于普通深松旋耕聯合整地機，牽引阻力降低了16.21%，方會敏等[12]建立了秸稈-土-旋耕刀離散元模型，進行EDEM仿真試驗和田間試驗對比，結果顯示秸稈位移隨著轉速增加呈現增加的趨勢，秸稈的水平運動位移總是大于同轉速下的側向位移，同時闡明了土壤的破碎機理，郭俊等[13]為研究正、反旋兩種旋耕方式在不同轉速下土壤和秸稈位移情況，運用示蹤法原理，對鋁塊和秸稈進行標記，試驗結果表明反轉旋耕破碎率優于正旋。其他相關研究主要針對土壤破碎性能設計一種旋耕裝置，未結合油莎豆種植模式和收獲模式，適應性差，無法滿足油莎豆收獲基本要求。

本文基于離散元原理與旋耕刀-油莎豆團聚體力學模型，以反旋作業方式進行運動學以及旋耕刀結構參數分析，建立油莎豆團聚體離散元模型，通過EDEM仿真試驗，分析反轉旋耕刀結構參數對試驗指標影響規律，結合試驗結果，設計一種高效破土旋耕刀，以期解決油莎豆聯合收獲過程中挖掘阻力大、土壤破碎率低、埋果率高等問題。

1 油莎豆植株特性分析

由于油莎豆植株成簇生長，收獲期時根系發達且含土量高[14]，如圖1所示。本文針對油莎豆種植模式進行試驗分析，發現河南民權地區油莎豆種植特點：行距15～20 cm，株距12～16 cm，生長深度為12～15 cm，為挖掘裝置結構參數設計提供依據。

為保證后續仿真試驗的精準性和科學性，結合河南民權地區油莎豆種植土壤特點，通過五點取樣法進行土壤性質測定，得到其不同深度的土壤含水率、密度以及土壤堅實度，如表1所示。

表1 土壤物性參數Tab.1 Soil physical parameters

2 旋耕刀反旋工作原理

2.1 旋耕刀運動分析

當反轉旋耕刀切削土壤時，旋耕刀在任何時刻的絕對速度都是由機器前進速度和刀軸轉速矢量合成，運動軌跡模型如圖2所示，設A為旋耕刀刀尖，其坐標為(X,Y)，其運動軌跡方程為

(1)

式中vm——前進速度，m/s

R——刀軸回轉半徑，mm

α——彎折角，(°)

ω——刀軸角速度，rad/s

t——刀軸運動時間，s

β——旋耕刀轉動角，(°)

λ——旋耕速比

對式(1)求導可得刀尖A(X，Y)處在x和y方向的分速度為

(2)

刀尖A處絕對速度為

(3)

聯立式(2)、(3)可得

(4)

其中

λ直接影響旋耕刀的運動軌跡、切土節距以及機器功率損等工作性能。

切土節距計算式為

(5)

式中S——切土節距，mm

Z——刀軸上同一徑向旋耕刀數量

n——旋耕刀轉速，r/min

與正旋挖掘相比，機器以反旋的工作方式作業時，可以在低速下進行挖掘作業，保證機器在較低功耗下實現較高的土壤破碎率，提高機具的作業性能，同時土垡切削量小且穩定，作業后土壤顆粒層次分布，通透性好，利于后續耕整地作業[15]。

2.2 旋耕刀結構參數分析

在油莎豆反轉挖掘收獲過程中，旋耕刀正切刃和刀端彎曲部位首先開始切削土壤，然后正切刃由下至上切削土壤，正切面的主要功能為切土、碎土、拋土等，是反轉旋耕刀重要組成結構，其結構形狀和參數對于挖掘性能有著重要影響，其中影響正切面的結構參數主要有回轉半徑R、彎折角α以及單刀工作幅寬L[16]，如圖3所示。本文結合油莎豆種植模式研究，將回轉半徑設為150 mm，通過改變彎折角α以及單刀工作幅寬L來改變正切面參數，進一步提高正切面的拋土性能，通過EDEM仿真探究其最佳結構參數組合。

3 離散元模型建立

3.1 土壤模型

結合前面對于土壤的特性分析，本文以參數平均值建立土壤顆粒床模型，在確保試驗精度的前提下，將土壤顆粒進行適當放大，將土壤顆粒半徑設為5 mm，同時為了確保土壤環境的真實性，采用顆粒組合的方式將土壤顆粒模型設置為以下4種模型[17-19]，即粒狀顆粒、團粒狀顆粒、塊狀顆粒、桿狀顆粒，如圖4所示。

3.2 油莎豆塊莖建模

在實際油莎豆收獲過程中，油莎豆塊莖是挖掘過程中油莎豆團聚體主要組成部分，很多學者在模擬試驗時，基本都忽略油莎豆塊莖的影響[20]。為精準模擬在收獲過程中油莎豆塊莖與根系連接力破壞的影響規律，本文利用三維激光掃描技術確定油莎豆塊莖的形狀特征及三維尺寸，對油莎豆塊莖的形狀結構進行掃描、曲面重構，得出塊莖的幾何體模型，利用EDEM軟件中的填充功能完成塊莖離散元模型的建立，對油莎豆塊莖三維尺寸進行試驗分析，油莎豆塊莖長為12 mm，寬為10 mm，高為8 mm，如圖5所示。

3.3 油莎豆根系模型

為精準模擬油莎豆收獲機在工作過程中旋耕刀與油莎豆團聚體的接觸規律，建立油莎豆根系離散元模型。由于油莎豆根系具有柔軟、細長以及方向不一致等特點[21]，本文對油莎豆根系建模進行簡化處理，建立根系主要須莖，通過試驗統計與分析，油莎豆根系長度范圍為90～120 mm，半徑為0.05～0.2 mm，結合對油莎豆根系長度的試驗分析，通過顆粒組合的形式將油莎豆根系長度設為120 mm，如圖6所示。

3.4 接觸模型及參數設定

3.4.1接觸模型

為精準模擬河南民權地區粘性土壤性質，土壤顆粒之間的接觸模型選擇Hertz-Mindlin with Bonding模型，如圖7所示，該模型可以通過結合力將相鄰2個土壤顆粒粘結在一起，結合力可以承受切向和法向位移[22-24]。在實際田間作業過程中，土壤與土壤之間存在粘附力，挖掘鏟與土壤存在受力關系，該模型可以模擬土壤顆粒之間的粘結作用和土壤顆粒破碎的現象。

土壤顆粒間Hertz-Mindlin with Bonding接觸模型計算式為

(6)

其中

(7)

式中Fn——法向接觸力，N

Ft——切向接觸力，N

Tn——法向力矩，N·m

Tt——切向力矩，N·m

vn——法向速度，m/s

vt——切向速度，m/s

kn——法向剛度，N/m

kt——切向剛度，N/m

ωn——法向角速度，rad/s

ωt——切向角速度，rad/s

A——單位接觸面積，mm2

J——慣性矩，mm4

δt——時間步長，s

Ra——顆粒粘結半徑，mm

土壤顆粒間法向和切向應力達到一定極值時粘結鍵將會斷裂，設斷裂時法向應力和切向應力分別為σmax、τmax，當粘結鍵斷裂時，有

(8)

3.4.2離散元參數設定

為了確保仿真試驗的準確性，本文通過參數試驗以及查閱文獻[25]，確定仿真試驗所涉及的接觸參數和基本物理參數，如表2所示。

3.5 虛擬仿真試驗過程

為模擬油莎豆挖掘裝置實際作業狀態，根據河南民權地區油莎豆種植模式，運用EDEM建立長1 200 mm、高400 mm虛擬土槽，其中油莎豆株距120 mm，行距150 mm，如圖8所示。

在虛擬仿真試驗過程中，將旋耕刀前進速度設置為0.5 m/s，刀軸轉速30 r/min，旋耕刀彎折角110°，旋耕刀工作幅寬30 mm；為保證仿真過程精確性和連續性，設置仿真步長為3.14×10-6[26-28]，網格單元尺寸9 mm，為最小顆粒半徑的3倍，顆粒粘結半徑為5.5 mm，顆粒生成時間0.4 s，仿真時長2 s，通過仿真試驗可得到旋耕刀對土壤的破碎情況，從EDEM后處理工具中導出旋耕刀所受阻力以及土壤破碎率的實時數據。

表2 離散元參數Tab.2 Discrete element parameter values

4 仿真試驗

4.1 挖掘過程動態分析

工作阻力是衡量聯合收獲機綜合性能的主要技術參數，在旋耕刀工作過程中，通過EDEM后處理分析不同時刻下挖掘鏟與油莎豆團聚體破壞動態過程，如圖9、10所示。

油莎豆反轉旋耕作業過程中(圖9、10)，旋耕刀正切刃和刀端彎曲部位首先開始切削底層土壤，土壤垂直于旋耕刀尖上向上方及兩側運動，深層土壤被打破，隨著旋耕刀不斷對中部土層進行剪切、撕裂、翻轉、拋擲等作用，土壤失效區域范圍變大，呈松散狀態，油莎豆塊莖和油莎豆根系受到擾動，如圖10所示，土壤呈紅色狀態，逐漸沿旋耕刀尖切線方向運動，在旋耕刀往復作用下油莎豆團聚體在被破壞的同時被拋起，完成反轉挖掘過程。

4.2 單因素試驗設計

為探究影響挖掘作業質量的主要因素取值范圍，以彎折角α以及單刀工作幅寬L為試驗因素進行單因素試驗，以工作阻力和刀軸扭矩為試驗指標，由上述旋耕刀結構參數分析可知，油莎豆反轉作業質量與彎折角α以及單刀工作幅寬L因素密切相關，根據河南油莎豆種植農藝要求，選取彎折角為90°～170°，工作幅寬為25～65 mm研究其對工作阻力和刀軸扭矩的影響規律，為后續多因素試驗提供理論基礎。

4.2.1單因素試驗

根據EDEM仿真數據，導出工作阻力和刀軸扭矩在不同仿真時間段的變化規律，如圖11、12所示。在旋耕刀工作過程中，當旋耕刀未接觸土壤時工作阻力和刀軸扭矩為0，隨著刀軸的轉動，旋耕刀開始切削土壤，所受工作阻力和刀軸扭矩開始逐漸增大，隨著切削運動的進行，旋耕刀切削土壤的體積不斷增加，當旋耕刀開始接觸油莎豆團聚體，對油莎豆團聚體進行破壞，旋耕刀所受阻力瞬間增大，在旋耕刀與油莎豆團聚體破壞過程中耕深達到最大值時，所受阻力達到最大值；當旋耕刀越過最大耕深時，切土體積逐漸減少，所受阻力相應減少；由于旋耕刀在切削土壤時每轉過一定的角度都有相同的4把旋耕刀進入土壤，因此刀軸在工作時受到的總阻力呈現周期性變化。

4.2.2彎折角α對試驗指標的影響

固定旋耕刀前進速度0.5 m/s，刀軸轉速30 r/min，旋耕刀工作幅寬30 mm，彎折角α分別取90°、110°、130°、150°和170°，進行5組仿真試驗，以工作阻力和刀軸扭矩作為試驗指標，觀察旋耕刀所受三向阻力最大值的變化規律，如圖13、14所示。在彎折角90°～150°范圍內，彎折角的變化對于工作阻力影響較為顯著，旋耕刀所受三向阻力逐漸增加，在150°～170°之間，彎折角的變化對于工作阻力影響較低，旋耕刀所受三向阻力變化較為平穩；當彎折角在90°～170°范圍內，對于x方向刀軸扭矩先增加至51 N，然后逐漸降低趨于平穩，當彎折角在90°～150°范圍內，對于y、z方向刀軸扭矩變化影響顯著，在150°～170°區間范圍內其刀軸扭矩變化趨于平穩。

4.2.3工作幅寬L對工作阻力的影響

固定旋耕刀前進速度0.5 m/s，刀軸轉速30 r/min，彎折角90°，旋耕刀工作幅寬L分別取25、35、45、55、65 mm，進行5組仿真試驗，旋耕刀所受三向阻力最大值隨轉速的變化情況如圖15、16所示。當工作幅寬在25～55 mm區間范圍內，x、z方向工作阻力逐漸增加，工作幅寬的變化對于x方向工作阻力影響較為顯著，在55～65 mm之間，x方向工作阻力變化較為平穩；工作幅寬在25～55 mm區間范圍內，對于y方向工作阻力變化影響顯著，在55～65 mm范圍其工作阻力有略微下降趨勢，但阻力波動范圍較小；當工作幅寬在25～55 mm區間范圍內，對于x方向刀軸扭矩變化影響顯著，但工作幅寬在55～65 mm內其刀軸扭矩略微下降，當工作幅寬在25～55 mm區間范圍內，y、z方向刀軸扭矩逐漸增加，在55～65 mm之間，刀軸扭矩變化較為平穩。

根據以上分析，為得到最佳旋耕刀結構參數，選取彎折角110°～150°，工作幅寬35～55 mm為試驗因素取值范圍。

4.3 多因素仿真試驗

在前期油莎豆團聚體離散元模型建立及旋耕刀動力學分析的基礎上，為得到理想的旋耕刀結構參數組合，進行EDEM正交旋轉虛擬仿真試驗研究，結合前期理論分析和單因素試驗結果，設定試驗因素編碼如表3所示。

表3 試驗因素編碼Tab.3 Test factors coding

以彎折角α、工作幅寬L為試驗因素，采用二因素五水平二次正交旋轉組合設計進行旋耕刀虛擬仿真試驗，根據所得數據進行結果分析，并對影響指標的主要因素進行顯著性分析，試驗結果如表4所示，表中A、B為編碼值。

表4 仿真試驗結果Tab.4 Simulation test results

對表4數據進行二次多元回歸擬合，得到工作阻力和刀軸扭矩對編碼自變量的二次多元回歸方程

y1=703.22+103.48A+33.72B+5.42AB+
222.50A2+47.92B2

(9)

y2=31.06+2.09A+3.19B+0.15AB-
1.37A2+1.96B2

(10)

同時經過Design-Expert軟件處理后得到關于工作阻力和刀軸扭矩的方差分析結果如表5、6所示。

表5 工作阻力方差分析Tab.5 Working resistance analysis of variance

表6 刀軸扭矩方差分析Tab.6 Analysis of variance of cutter shaft torque

由表5可知，模擬的一次項彎折角A對工作阻力的影響顯著，二次項A2對工作阻力影響極其顯著。從單因素水平分析，各項因素對工作阻力的影響由大到小為彎折角A、工作幅寬B。響應面分析如圖17所示。

由表6可知，模擬的一次項工作幅寬B對刀軸扭矩的影響極其顯著，一次項彎折角A和二次項B2對刀軸扭矩影響顯著。從單因素水平分析，各項因素對工作阻力的影響由大到小為工作幅寬B、彎折角A，響應面分析如圖18所示。

4.4 試驗優化

為得到試驗因素最佳旋耕刀結構參數組合，對試驗進行優化設計。建立參數化數學模型，結合試驗因素的邊界條件，以降低工作阻力和刀軸扭矩為目標，對旋耕刀工作阻力和刀軸扭矩進行回歸分析，建立非線性規劃參數模型為

(11)

利用Design-Expert對數學模型進行分析求解，彎折角130°、工作幅寬45 mm為最優參數組合，此時工作阻力為703.22 N，刀軸扭矩32.06 N·m，將最優參數進行EDEM虛擬仿真試驗驗證，其工作阻力為928.63 N，刀軸扭矩34.56 N·m，與優化結果基本一致。在相同參數設置下與目前正旋旋耕挖掘裝置進行試驗對比，結果表明，正旋挖掘裝置工作阻力為1 015.52 N，刀軸扭矩為36.72 N·m，工作阻力降低了8.6%，刀軸扭矩降低了5.9%，驗證了反旋式旋耕刀片具有一定的降阻作用。

5 田間試驗

5.1 田間試驗條件與方案

旋耕刀挖掘阻力和刀軸扭矩并不是衡量作業效果主要標準，由于油莎豆目前無相關收獲作業標準，本文參考NY/T 502—2016中的收獲標準以及GB/T 5668—2017中的土壤破碎率計算方法，將埋果率和土壤破碎率作為田間試驗的指標，目的是檢驗優化設計的旋耕刀作業性能，結合理論分析與仿真模擬結果進行田間驗證，試驗地點選擇河南省民權縣油莎豆種植基地，該地區油莎豆種植模式為平作，行距為140 mm，株距為120 mm，油莎豆生長深度160 mm左右，如圖19所示。

試驗劃為5個區域，每個試驗區域長50 m，將新型旋耕刀和IT245旋耕刀分別安裝在油莎豆聯合收獲機上，進行油莎豆收獲性能對比分析。機器作業速度設置為0.5 m/s，將埋果率和土壤破碎率作為試驗指標，進行試驗分析，如圖20所示。

5.2 田間試驗結果與分析

試驗結果如表7所示。

表7 田間性能試驗對比Tab.7 Field performance test comparison %

由表7可知，新型旋耕刀作業時埋果率最高值為1.23%，最低值為0.85%，平均值為1.07%，土壤破碎率最高值為94.61%，最低值為92.67%，平均值為93.48%，標準旋耕刀埋果率最高為2.85%，最低為1.75%，平均值為2.27%，土壤破碎率最高為92.87%，最低為91.54%，平均值為92.18%，與標準旋耕刀相比，新型旋耕刀油莎豆塊莖埋果率減低了1.20個百分點，土壤破碎率提高了1.3個百分點，參考NY/T 502—2016的規定，進一步提高了油莎豆的收獲效率。

6 結論

(1)通過對旋耕刀進行運動學分析和結構分析，結合虛擬仿真試驗，得到旋耕刀結構參數彎折角取值范圍為110°～150°，工作幅寬取值范圍為35～55 mm。

(2)建立了旋耕刀-土壤-油莎豆團聚體離散元模型，進行了油莎豆虛擬仿真挖掘試驗，從微觀角度分析了土壤破碎的動態過程和油莎豆團聚體的速度關系，根據Design-Expert試驗分析，確定了旋耕刀的最佳結構參數：彎折角130°，工作幅寬45 mm，在相同參數設置下與正旋旋耕方式進行對比試驗，工作阻力降低了8.6%，刀軸扭矩降低了5.9%，證明反旋挖掘具有降阻作用。

(3)為驗證反旋刀片設計的合理性，以埋果率和土壤破碎率進行油莎豆收獲性能試驗，同時與標準旋耕刀進行試驗對比，結果表明：反旋作業方式埋果率1.07%，土壤破碎率93.48%，與標準旋耕刀相比，新型旋耕刀油莎豆塊莖埋果率減低了1.2個百分點，土壤破碎率提高了1.3個百分點。