新疆農(nóng)田粉土離散元仿真參數(shù)標(biāo)定

李清超，鄭 炫，劉進(jìn)寶，楊懷君，王子龍，張魯云

（1.石河子大學(xué)機(jī)械電氣工程學(xué)院，新疆石河子 832000；2.新疆農(nóng)墾科學(xué)院機(jī)械裝備研究所，新疆石河子 832000）

0 引言

【研究意義】新疆地區(qū)部分土壤屬于粉土，土壤顆粒主要以粉粒為主，比重占到80%以上[1]，由于粒徑較小，再加上新疆地區(qū)覆膜滴灌技術(shù)使用使得秋季農(nóng)田土壤含有一定水分，具有一定黏性，不利于犁地作業(yè)[2～3]。新疆地區(qū)粉土土壤農(nóng)田不適合免耕作業(yè)，耕翻成為基礎(chǔ)性措施之一。目前對(duì)于新疆地區(qū)粉土土壤的犁地作業(yè)未能有準(zhǔn)確描述，難以形成犁地翻垡作業(yè)的理論性指導(dǎo)，其主要的關(guān)鍵問(wèn)題是缺少準(zhǔn)確的參數(shù)模型進(jìn)行仿真模擬，獲取更優(yōu)的土壤翻耕離散元模型，對(duì)機(jī)械設(shè)計(jì)人員有實(shí)際意義?！厩叭搜芯窟M(jìn)展】離散元法被廣泛應(yīng)用在農(nóng)業(yè)機(jī)械觸土部件與土壤相互作用的研究上[4～8]，耕地模型及土壤與部件模型的準(zhǔn)確性是保證結(jié)果有效的前提。建立一個(gè)系統(tǒng)、準(zhǔn)確的新疆粉土土壤模型及粉土土壤與觸土部件模型，是利用離散元分析方法對(duì)耕整地機(jī)械進(jìn)行優(yōu)化設(shè)計(jì)的基礎(chǔ)。有研究通過(guò)物理試驗(yàn)與土壤仿真試驗(yàn)結(jié)合的方法，實(shí)現(xiàn)了對(duì)土壤離散元接觸參數(shù)的標(biāo)定[9～10]。Liang 等[11]以實(shí)測(cè)的土壤堆積角和土-膜的靜摩擦系數(shù)為響應(yīng)值，通過(guò)B-B試驗(yàn)，對(duì)土-土之間以及土-膜之間接觸模型參數(shù)做了標(biāo)定。戴飛等[12]用Hertz-Mindlin 無(wú)滑動(dòng)模型，通過(guò)土壤堆積角、土壤與鋼滑動(dòng)摩擦角的試驗(yàn)方法，對(duì)全膜雙莖溝覆膜土壤進(jìn)行了標(biāo)定。謝方平等[13]通過(guò)單軸密閉壓縮試驗(yàn)和無(wú)側(cè)限抗壓強(qiáng)度試驗(yàn)，標(biāo)定了基于Edinburgh Elasto-Plastic Adhesion（EEPA）模型的土壤離散元參數(shù)?！颈狙芯壳腥朦c(diǎn)】新疆粉土土壤的待耕地，在常年滴水灌溉和耕翻作用下，其土壤耕性與其他地區(qū)不同土壤類(lèi)型的耕地主要有阻力大、易與觸土部件黏附等區(qū)別，仿真模型也會(huì)存在差異，難以通用。需要技術(shù)上重新標(biāo)定，以校驗(yàn)用于耕地機(jī)械設(shè)計(jì)制造的依據(jù)。【擬解決的關(guān)鍵問(wèn)題】以新疆地區(qū)粉土土壤的待耕地為研究對(duì)象，離散元仿真試驗(yàn)與物理試驗(yàn)相結(jié)合的方法，采用標(biāo)定土壤與機(jī)具材料之間及土壤與土壤之間的離散元接觸模型參數(shù)，為優(yōu)化農(nóng)業(yè)機(jī)械觸土部件提供參考。

1 材料與方法

1.1 材料

1.1.1 土壤本構(gòu)參數(shù)

土壤樣品來(lái)源于2021年秋季作物收獲后，采自新疆石河子市小麥屯村農(nóng)場(chǎng)（44°21′31″N，85°55′59″E），前茬作物為棉花，土壤質(zhì)地為粉土。采用5點(diǎn)取樣法取耕翻作業(yè)0 ～30 cm深度的土壤，利用環(huán)刀（100 cm2）、YH-M50002電子天平（浙江東陽(yáng)英衡智能設(shè)備廠，精度0.01 g）、DHG-9070A干燥機(jī)（上海一恒科技有限公司）、ZJ型應(yīng)變控制式直剪儀（南京土壤儀器廠）、BT-2001 型激光粒度分布儀（濕法），對(duì)土壤的密度、含水率、泊松比、土壤粒徑進(jìn)行測(cè)量。密度、含水率、粒徑的測(cè)量按照土壤試驗(yàn)標(biāo)準(zhǔn)進(jìn)行；土壤泊松比的測(cè)量方法參照文獻(xiàn)[14]。土壤本構(gòu)參數(shù)：密度1 480 kg/m3，含水率14.5%，泊松比0.37，粒徑分布＜0.001 mm、0.001 ～0.05mm和＞0.05mm分別占1.82%，85.55%和12.63%。計(jì)算泊松比。

表1 土壤本構(gòu)參數(shù)Table 1 Soil constitutive parameters

1.1.2 土壤與65 Mn靜摩擦系數(shù)

土壤與65 Mn之間的靜摩擦系數(shù)是兩者之間的特有屬性，不會(huì)因形狀改變而改變。土壤與65 Mn 之間的靜摩擦系數(shù)由靜摩擦試驗(yàn)臺(tái)測(cè)得。將土壤用方刀制成1 cm×1 cm×1 cm的土塊，放于65 Mn 斜面之上，旋轉(zhuǎn)推桿來(lái)逐漸提高斜面角度，通過(guò)斜面儀來(lái)記錄土塊下滑瞬間的角度，由公式（2）計(jì)算得到靜摩擦系數(shù)，重復(fù)10次，取平均值作為最終結(jié)果，得到土壤與65 Mn 之間的靜摩擦系數(shù)為0.56。圖1

圖1 實(shí)際物理試驗(yàn)Fig.1 Actual physical test

式中：為靜摩擦系數(shù)，為土塊下滑瞬間角度。

1.1.3 土壤滾動(dòng)距離

用65Mn 板搭建的斜坡試驗(yàn)臺(tái)，標(biāo)定土壤與65 Mn 之間的接觸模型參數(shù)。土球顆粒由3D 打印的球形模具制得，半徑為5 mm。在預(yù)試驗(yàn)的基礎(chǔ)上，取斜面的高度為30 mm，斜面長(zhǎng)度為100 mm，水平面長(zhǎng)度為400 mm。將土球放于斜面上端中間位置處，使其沿斜面自然滾下，至土球停止運(yùn)動(dòng)。用直尺測(cè)量土球沿水平方向的滾動(dòng)距離，重復(fù)10 次，平均值為227.29 mm，作為最終的目標(biāo)值。

1.1.4 土壤堆積角

土壤堆積角試驗(yàn)裝置，主要有萬(wàn)能試驗(yàn)機(jī)（AI-7000LAU30，高鐵檢測(cè)儀器有限公司）、65 Mn圓筒、65 Mn板組成。圓筒[15]內(nèi)徑為90 mm，高為270 mm，下方與65 Mn 板（500 mm × 500 mm）接觸，內(nèi)部裝填土壤試樣1 kg。采用萬(wàn)能試驗(yàn)機(jī)，以0.05 m/s 的速度向上垂直升起圓筒，土壤向下滑落在65 Mn板上，待土壤靜置穩(wěn)定后，用相機(jī)在垂直的方向上拍照，得到土樣堆積角的正視圖。后續(xù)使用matlab圖像識(shí)別功能，提取邊界，通過(guò)線性回歸來(lái)獲得土樣堆積角的大小。重復(fù)10 次取平均值，得到土壤的堆積角為37.64°。圖1

1.1.5 犁耕模型

選取顆粒形狀為EDEM 軟件自帶的球狀顆粒，半徑為5 mm。采用EDEM 的Geometries 功能，生成2 500 mm×3 600 mm×400 mm 的土槽，顆?？倲?shù)為3617032個(gè)，將用solidworks繪制好的犁體曲面模型保存成x_t 格式并導(dǎo)入到EDEM中。前犁體和后犁體設(shè)置耕寬為450 mm，耕深為300 mm，水平前進(jìn)速度分別設(shè)置為8、9 和10 km/h，總仿真時(shí)間為2.3 s，其中前0.3 s 生成土壤，后2s犁體模型完成在土壤中的運(yùn)動(dòng)。圖2

圖2 犁體仿真耕作過(guò)程Fig.2 Plow simulation tillage process.

1.2 方法

針對(duì)常用的耕翻機(jī)具材料65 Mn，采用靜摩擦試驗(yàn)測(cè)定土壤與65 Mn 材料的靜摩擦系數(shù)，運(yùn)用斜坡物理試驗(yàn)與斜坡EDEM仿真試驗(yàn)相結(jié)合的方法來(lái)標(biāo)定土壤與65 Mn鋼之間的離散元接觸模型參數(shù)。運(yùn)用堆積角物理試驗(yàn)與堆積角EDEM仿真試驗(yàn)相結(jié)合的方法建立粉土土壤待耕地的離散元仿真模型。應(yīng)用Design-Expert 軟件，根據(jù)Central Composite 試驗(yàn)原理建立粉土土壤的斜坡滾動(dòng)距離和堆積角與相關(guān)參數(shù)的回歸模型，并找到最優(yōu)參數(shù)組合。在最優(yōu)參數(shù)組合下，比較堆積角的仿真值與實(shí)際值，以及耕翻仿真試驗(yàn)與耕翻田間試驗(yàn)的對(duì)照，驗(yàn)證模型的準(zhǔn)確性。

1.2.1 接觸模型選取

秋季待耕期間，由于膜下滴灌的影響，粉土土壤待耕地含有一定水分，顆粒間存在明顯的黏附現(xiàn)象，且在耕翻時(shí)，土壤較容易黏附在觸土部件上。Hertz-Mindlin with JKR 接觸模型克服了經(jīng)典的Hertz接觸模型不涉及黏結(jié)力的缺點(diǎn)，通過(guò)表面能來(lái)表現(xiàn)顆粒間的粘結(jié)力，適用于有明顯黏結(jié)和團(tuán)聚的濕黏物料[16～19]。該模型的JKR 法向彈性接觸力計(jì)算公式為：

式中FJKR：JKR法向彈性接觸力，N；

E*：等效彈性模量（Pa）；

R*：等效接觸半徑（m）；

α：兩顆粒切向重疊量（m）；

γ：表面能（J/m3）；

δ：兩顆粒法向重疊量（m）。

其中，等效彈性模量和等效接觸半徑由式（5）和（6）定義：

式中E1、E2：兩顆粒的彈性模量（Pa）；

V1、V2：兩顆粒的泊松比；

R1、R2：兩顆粒的接觸半徑（m）。取表面能γ為0時(shí)，此時(shí)JKR模型與Hertz接觸模型一致，即：

即使兩顆粒沒(méi)有直接接觸，該模型也會(huì)在一定距離范圍內(nèi)，存在黏結(jié)力，顆粒間存在黏結(jié)力的最大間隙為：

式中，δC和αC分別為顆粒間存在粘結(jié)力的法向和切向最大間隙。

當(dāng)顆粒間隙大于0，小于δC時(shí)，兩者間黏結(jié)力達(dá)到最大值，計(jì)算公式為：

JKR接觸模型中的表面能，能較接近地模擬土壤顆粒之間以及土壤與觸土部件之間地黏結(jié)力，選用Hertz-Mindlin with JKR 接觸模型來(lái)建立新疆粉土土壤類(lèi)型待耕地的離散元仿真模型。

1.2.2 EDEM仿真試驗(yàn)

應(yīng)用離散元軟件EDEM 2018 版本建立仿真模型，所有仿真試驗(yàn)的土壤顆粒形狀選為球狀，半徑設(shè)置為5 mm，生成方式為隨機(jī)分布，三維結(jié)構(gòu)模型的材料均為65 Mn。離散元仿真試驗(yàn)中的土壤和65 Mn 的本構(gòu)參數(shù)，通過(guò)測(cè)量或查閱參考文獻(xiàn)[20]獲得。表2

表2 土壤與65Mn的本構(gòu)參數(shù)Table 2 Constitutive parameters of soil and 65Mn

斜坡仿真試驗(yàn)?zāi)Ｐ团c物理試驗(yàn)保持一致，將solidworks 繪制的斜坡模型保存成x_t 格式，導(dǎo)入到EDEM軟件中，坐標(biāo)原點(diǎn)位置設(shè)置在斜坡頂端，中間位置處，并在此位置設(shè)置顆粒工廠，顆粒生成方式為dynamic，顆粒數(shù)量為1，數(shù)據(jù)保存間隙為0.01 s，重力加速度為9.81 m/s2。顆粒球順著斜坡滑下，到靜止時(shí)，將水平坐標(biāo)導(dǎo)出，再減去95.39 mm（斜坡水平長(zhǎng)度），即為顆粒的水平滾動(dòng)距離。表3

表3 土-65Mn仿真試驗(yàn)因素水平編碼Table 3 Factor level coding table of soil-65Mn simulation test

參考文獻(xiàn)[14～15，20～22]的接觸參數(shù)選取范圍，采用二次回歸正交旋轉(zhuǎn)中心組合試驗(yàn)的方法，確定土壤與65 Mn接觸參數(shù)的范圍。表4

表4 土-土仿真試驗(yàn)因素水平編碼表因素水平Table 4 Factor level coding table of soil-soil simulation test

將solidworks 繪制的堆積角試驗(yàn)?zāi)Ｐ捅４娉蓌_t 格式，導(dǎo)入到EDEM 軟件中，圓筒和65 Mn 板各尺寸與物理試驗(yàn)保持一致。在圓筒上方創(chuàng)建顆粒工廠，生成顆粒方式為dynamic，數(shù)量為1 kg，待落在圓筒中的顆粒穩(wěn)定后，圓筒以0.05 m/s 的速度垂直向上運(yùn)動(dòng)，顆粒在重力作用下向下滑落，最終顆粒停止運(yùn)動(dòng)形成堆積角，測(cè)量。圖3

圖3 仿真試驗(yàn)Fig.3 simulation test.

參考文獻(xiàn)[14～15，20～22]的接觸參數(shù)選取范圍，采用二次回歸正交旋轉(zhuǎn)中心組合試驗(yàn)的方法，確定土壤接觸參數(shù)的范圍。表5

表5 斜坡仿真試驗(yàn)設(shè)計(jì)及結(jié)果Table 5 Slope simulation test design and results

1.2.3 阻力對(duì)比

田間試驗(yàn)在新疆石河子市小麥屯村農(nóng)場(chǎng)（44°21′31″N，85°55′59″E）進(jìn) 行，土 壤 含 水 率為15.51%，土壤堅(jiān)實(shí)度為2152.8 kPa。使用的儀器有：雷沃歐M904 拖拉機(jī)、翻轉(zhuǎn)犁、遙測(cè)儀（Autobona，黑龍江省農(nóng)業(yè)機(jī)械工程科學(xué)研究院制），直尺、卷尺、秒表等。量取地塊中間100 m作為測(cè)試區(qū)間，在保證相應(yīng)耕作速度的前提下，通過(guò)遙測(cè)儀測(cè)得耕作阻力，并對(duì)區(qū)間內(nèi)的試驗(yàn)結(jié)果取平均值。

2 結(jié)果與分析

2.1 土壤與65Mn接觸參數(shù)標(biāo)定

研究表明，二次回歸方程如下：

該回歸模型P＜0.000 1，失擬項(xiàng)P=0.132 0 ＞0.05，決定系數(shù)R2=0.983 5，校正系數(shù)Radj2=0.968 1，均＞0.9，模型具有良好的準(zhǔn)確度，該模型可較好地表達(dá)滾動(dòng)距離與土-65 Mn 恢復(fù)系數(shù)、土-65 Mn 靜摩擦系數(shù)、土-65 Mn 動(dòng)摩擦系數(shù)、土-65 MnJKR表面能之間的關(guān)系。4個(gè)因素對(duì)滾動(dòng)距離影響的主次順序?yàn)橥?65 Mn 動(dòng)摩擦系數(shù)、土-65 MnJKR 表面能、土-65Mn 恢復(fù)系數(shù)、土-65 Mn 靜摩擦系數(shù)。其中，C、D對(duì)滾動(dòng)距離影響極顯著，A、B對(duì)滾動(dòng)距離影響不顯著。一次交互作用均對(duì)滾動(dòng)距離影響不顯著。二次交互作用中，A2、B2、D2對(duì)滾動(dòng)距離影響顯著。C2對(duì)滾動(dòng)距離影響極顯著。選定靜摩擦系數(shù)為0.56 的一組最優(yōu)解。選擇土球在65 Mn 板上滾動(dòng)距離相近的優(yōu)化解。表6，表7

表6 斜坡仿真試驗(yàn)方差Table 6 Variance analysis of slope simulation test

表7 土-65Mn接觸模型參數(shù)Table 7 Parameters of soil-65mn contact model

2.2 土壤與土壤接觸參數(shù)標(biāo)定

研究表明，建立堆積角Y2與土-土恢復(fù)系數(shù)X1、土-土靜摩擦系數(shù)X2、土-土動(dòng)摩擦系數(shù)X3、土-土JKR表面能X4之間的二階回歸模型。表8

表8 土-土參數(shù)標(biāo)定試驗(yàn)設(shè)計(jì)與結(jié)果Table 8 Design and results of soil-soil parameter calibration test

模型P值小于0.000 1，失擬項(xiàng)P值大于0.05，決定系數(shù)R2均大于0.9，具有良好的準(zhǔn)確度，此模型可較好地表達(dá)堆積角與土-土恢復(fù)系數(shù)、土-土靜摩擦系數(shù)、土-土動(dòng)摩擦系數(shù)、土-土JKR 表面能之間的關(guān)系。其中，有4個(gè)回歸項(xiàng)影響極顯著（P＜0.01），分別為X4、X1X4、X3X4、X22，2個(gè)回歸項(xiàng)對(duì)模型影響顯著（P＜0.05），分別為X2、X42。表9

表9 土-土參數(shù)標(biāo)定試驗(yàn)方差Table 9 Analysis of variance of soil-soil parameter calibration test

取堆積角37.64°為目標(biāo)，選擇與目標(biāo)值相近的優(yōu)化解。選擇軟件給出的多組最優(yōu)解中的第一個(gè)組解作為土-土接觸模型的參數(shù)，最終確定土-土接觸模型參數(shù)。表10

表10 土-土接觸模型參數(shù)Table 10 Parameters of soil-soil contact model

2.3 物理試驗(yàn)驗(yàn)證

研究表明，土壤參數(shù)標(biāo)定后的滾動(dòng)距離、堆積角的相對(duì)誤差分別為6.05%、1.28%。標(biāo)定后堆積角的形狀與實(shí)際形狀更為接近。表11，圖4

圖4 土壤堆積角仿真試驗(yàn)與物理試驗(yàn)對(duì)比Fig.4 Comparison between soil accumulation angle simulation test and physical test.

表11 試驗(yàn)結(jié)果對(duì)比Table 11 comparison of test results

2.4 犁體阻力試驗(yàn)驗(yàn)證

研究表明，仿真的滾動(dòng)距離比實(shí)際滾動(dòng)距離大了6.05%，堆積角小了1.28%，在耕速10 km/h情況下，犁體耕作阻力小了6.49%；土壤含水率在田間試驗(yàn)時(shí)，與取土階段相比發(fā)生了輕微的變化，也是造成誤差的原因之一。所建立的新疆粉土土壤待耕地離散元仿真模型，在誤差范圍內(nèi)是可以替代土壤的實(shí)際參數(shù)的，適合秋季進(jìn)行耕翻的干旱地區(qū)。表12，圖5

表12 阻力對(duì)比Table 12 Resistance comparison

圖5 田間試驗(yàn)Fig.5 field experiment

3 討論

研究標(biāo)定的是新疆地區(qū)易于粘附的粉土土壤離散元模型，雖然與東北地區(qū)黑土的離散元模型選用的接觸模型是一致的[20]，但具體接觸參數(shù)標(biāo)定結(jié)果不同，新疆農(nóng)田粉土的JKR 表面能要明顯小于東北黑土；北方地區(qū)的葡萄藤防寒土所選用的基礎(chǔ)模型為整合延遲彈性模型和線性粘附模型，與研究的接觸模型有一定區(qū)別[21]；與新疆棉田耕后土壤離散元模型相比，研究所選的接觸模型更適合有一定粘附特性的土壤[22]。研究與黑土區(qū)玉米秸稈-土壤混料離散元模型相比，只考慮土-土之間以及土-65Mn 之間的接觸參數(shù)，并未考慮殘茬的相關(guān)作用[15]。

4 結(jié)論

4.1 利用二次回歸擬合方法獲得土-65 Mn之間最優(yōu)接觸參數(shù)組合：恢復(fù)系數(shù)0.51、靜摩擦系數(shù)0.56、動(dòng)摩擦系數(shù)0.08、JKR表面能4.12。

4.2 土-土之間最優(yōu)接觸參數(shù)組合為：恢復(fù)系數(shù)0.57、靜摩擦系數(shù)0.65、動(dòng)摩擦系數(shù)0.23、JKR表面能4.49。

4.3 對(duì)比仿真試驗(yàn)和實(shí)際試驗(yàn)的滾動(dòng)距離被堆積角，誤差分別為：6.05%、1.25%。

4.4 犁體8、9、10 km/h 作業(yè)速度下，仿真試驗(yàn)和實(shí)際試驗(yàn)的阻力相對(duì)誤差為：5.69%、5.95%、2.74%，誤差在可接受范圍內(nèi)，驗(yàn)證了標(biāo)定的新疆粉土土壤離散元接觸參數(shù)的準(zhǔn)確性。