基于EEPA接觸模型的土壤耕作特性模擬及顆粒球型影響分析

萬里鵬程 李永磊 蘇 辰 趙 虎 董向前 宋建農(nóng) 王繼承

(中國農(nóng)業(yè)大學(xué) 工學(xué)院，北京 100083)

隨著農(nóng)業(yè)現(xiàn)代化的發(fā)展，數(shù)字技術(shù)、仿真技術(shù)越來越多地應(yīng)用于田間土壤環(huán)境的模擬中。離散元法(Discrete element method, DEM)是專門用來解決不連續(xù)介質(zhì)問題的數(shù)值模擬方法。該方法適用于工作對象可轉(zhuǎn)化為離散顆粒的研究[1-3]，能夠模擬土壤多顆粒單元聚合體[4-6]及受力條件下表現(xiàn)出的動力學(xué)性質(zhì)，廣泛應(yīng)用于土壤與機具相互作用的研究領(lǐng)域[7-12]。國內(nèi)外已有研究對構(gòu)建土壤離散元模型的顆粒球型與混配方式應(yīng)用方法差異大，探究離散顆粒球型對土壤離散元模型仿真精度及計算效率的影響很有必要。

土壤耕作特性是土壤在耕作時反應(yīng)出來的性質(zhì)，包括耕作難易程度和耕后土壤表現(xiàn)出的狀態(tài)[13-15]，離散元法模擬土壤耕作特性是目前廣泛應(yīng)用的研究方式。EEPA接觸模型兼顧粘性、彈性與塑性，與實際作業(yè)區(qū)域的土壤加載力再卸載過程展現(xiàn)出宏觀土壤力學(xué)特性一致，應(yīng)用于粘性土壤模型研究，劉宏俊等[16]利用EEPA接觸模型對稻麥周年地區(qū)黏重土壤進行相關(guān)參數(shù)標(biāo)定；謝方平等[17]采用單軸密閉壓縮試驗和無側(cè)限抗壓強度試驗標(biāo)定了基于EEPA模型的土壤離散元參數(shù)。在模擬計算散體物料力學(xué)行為時，顆粒形狀是重要的影響因素[18-21]，王憲良等[22]以粒徑為10 mm的球形顆粒為基礎(chǔ)，仿照土壤顆粒形狀填充了核狀、塊體、條狀、片狀土壤顆粒單元，進行了0.95～1.05倍粒徑隨機分布的球形顆粒填充土槽，建立了土壤模型標(biāo)定方法；Hang等[23]采用粒徑12 mm的球形顆粒并填充了三球顆粒與四球顆粒為顆粒單元，分析了深松鏟齒間距對土壤擾動的作用；Sun等[24]使用粒徑8～12 mm隨機分布單球顆粒填充土槽，分析了耕作深度大于40 cm的仿生深松鏟的減阻、土壤擾動特性。國內(nèi)外已有研究應(yīng)用離散元法模擬土壤，選擇不同顆粒球型建立虛擬模型[25-30]，卻未闡明所選顆粒球型對土壤離散元模型仿真精度及計算效率的影響。

本研究擬以土槽實驗室的土槽土壤為研究對象，采用靜態(tài)堆積角試驗、單軸密閉壓縮試驗測量土壤摩擦特性、粘性、彈性和塑性參數(shù)，建立挖掘犁鏟-土壤離散元耦合仿真模型驗證參數(shù)優(yōu)化準(zhǔn)確性，進一步分析比較顆粒球型、顆粒粒徑配比對計算效率、仿真精度的影響，以期為土壤-機具互作離散元仿真的建模，EEPA接觸模型的參數(shù)選擇，顆粒球型的確定提供依據(jù)和參考。

1 EEPA接觸模型

EEPA接觸模型以7個模型參數(shù)表征法向粘結(jié)力加載、卸載、再加載、再卸載過程，其加載/卸載路徑模型[30]見圖1：最初加載時，顆粒間法向粘結(jié)力沿初始加載路徑①，卸載后，切換到卸載/再加載路徑②并在特定重疊處達到零，描述為塑性重疊區(qū)長度δp。若發(fā)生再加載，顆粒間法向粘結(jié)力加載方式將保持在路徑②上，但若超過歷史最大加載力，則切換到初始加載路徑①。卸載超過塑料重疊區(qū)長度δp，將導(dǎo)致粘合力的變化，直到最大粘合力fmin達到f0。在這一點上進一步卸荷，將使法向重疊量和吸引力都減小，直至在重疊區(qū)長度δ=0處發(fā)生分離。如果再次發(fā)生加載，接觸力將沿著一條與初始卸載/再加載路徑平行的路徑②，直至到達最初加載路徑①重合點后，在初始加載路徑①上進一步加載。

、 法向粘結(jié)力加載/卸載方向， 法向粘結(jié)力加載/卸載路徑。 fn為法向粘結(jié)力，N；f0為恒定拖拉力， N；fmin為顆粒間最大粘結(jié)力，N；δ為重疊區(qū)長度；δp為顆粒塑性變形量；δmin為最大粘結(jié)力的顆粒重疊量；δmax為重疊區(qū)長度最大值；K1為原始加載剛度；K2為再加載剛度；Kadh為衰減剛度；x為黏附分支曲線冪指數(shù)；n為非線性曲線冪指數(shù)。 and are loading/unloading direction of normal bond force, is the loading/unloading path of normal adhesive force. fn is the normal bond force, N; f0 is constant drag force, N. fmin is the maximum bonding force between particles, N. δ is the length of overlap. δp is the plastic deformation of particles. δmin is the particle overlap of the maximum bonding force. δmax is the maximum length of overlap region. K1 is the original loading stiffness. K2 is reloading stiffness. Kadh is attenuation stiffness. x is the power exponent of adhesion branch curve. n is the power exponent of nonlinear curve.圖1 EEPA-加載/卸載路徑模型Fig.1 EEPA-Load/Unload path model

加載路徑下的原始加載剛度為K1，表達式為：

(1)

式中：E*為楊氏模量；R*為顆粒間接觸距離；n1為加載分支指數(shù)；塑性變形比λp由初始加載剛度、再加載剛度決定，表達式為：

(2)

式中：K2為再加載剛度。塑性重疊區(qū)變形量δp與λp、加載分支指數(shù)n1密切相關(guān)，表達式為：

(3)

式中：δ為重疊區(qū)長度。顆粒分離前抵抗拉脫的最大粘結(jié)力fmin的表達式為：

fmin=πΔγψr

(4)

式中：Δγ為粘附能；ψ為粘著常數(shù)；r為顆粒接觸半徑。

綜上，顆粒間接觸力主要由EEPA模型參數(shù)、顆粒間重疊區(qū)長度計算得出。耕作土壤的宏觀力為定值，應(yīng)用不同粒徑的顆粒單元建立對應(yīng)土壤模型時，模型顆粒間重疊區(qū)長度差異大，導(dǎo)致模型參數(shù)變化大，模型宏觀表現(xiàn)差異顯著；而等效粒徑下的土壤顆粒單元將減弱顆粒間重疊區(qū)長度變化差異，該條件下的模型優(yōu)化參數(shù)組合在表征土壤宏觀力學(xué)性質(zhì)具有通用性。

2 土壤力學(xué)特性參數(shù)測定

室內(nèi)土槽試驗?zāi)軌蛴行M田間工況和較為準(zhǔn)確便捷的獲得試驗結(jié)果，且具有較強的可重復(fù)性和可對比性，是研究土壤-機具-植物互相作用關(guān)系的重要方法。為便于獲取土壤耕作特性參數(shù)，選擇中國農(nóng)業(yè)大學(xué)工學(xué)院土壤—機器—植物系統(tǒng)技術(shù)實驗室土槽土壤為對象開展相關(guān)研究。

2.1 土壤物理參數(shù)

室內(nèi)土槽土壤來源為北京地區(qū)棕壤土，質(zhì)地比較黏重，為獲得實際土壤的物理參數(shù)，取試驗土樣應(yīng)用烘干法測得土壤實際含水率為16.09%，采用環(huán)刀法測得土壤實際密度為2 130 kg/m3，三軸剪切試驗測得土壤剪切模量為0.96 MPa，取土壤泊松比為0.36。

2.2 土壤靜態(tài)堆積角

靜態(tài)堆積角試驗[4]常用來表征土壤的流動性、粘性與摩擦特性。試驗采用內(nèi)徑為200 mm，高度400 mm的有機玻璃材料圓筒，將盛滿土壤的圓筒開口端向下置于底板上，沿豎直方向0.1 m/s提起圓筒得到土壤堆積角。

采用MATLAB 2017R軟件對靜態(tài)堆積角試驗的原始圖像進行二值化處理得到二值圖像，借助Canny階梯型邊緣檢測算子進行邊緣檢測，最小二乘法將邊緣輪廓擬合成直線，添加水平線后通過霍夫變換函數(shù)提取二值圖像中直線夾角(靜態(tài)堆積角β)數(shù)值為33.52°。具體測量過程見圖2。

圖2 土壤靜態(tài)堆積角(β)的測量過程Fig.2 Measurement process of soil static accumulation angle (β)

2.3 土壤單軸密閉壓縮試驗

土壤在實際作業(yè)過程中表現(xiàn)出松散性、可壓實性，在受壓條件下土壤顆粒將產(chǎn)生塑性變形并對施壓部件產(chǎn)生反作用力(包含彈力)。單軸密閉壓縮試驗?zāi)軌驅(qū)崿F(xiàn)土壤顆粒群的單向壓縮，測得土壤軸向壓力與軸向應(yīng)變的變化關(guān)系曲線，表征出土壤彈性與塑性。試驗用圓筒內(nèi)徑為95 mm，注入土壤高度150 mm，使用美國FTC公司的TMS-Touch物性分析質(zhì)構(gòu)儀，設(shè)置圓盤壓板加載速度為60 mm/min，向下壓縮圓筒內(nèi)土壤12 mm，測得加載壓力與土壤的回彈高度。重復(fù)10次試驗，計算得軸向壓力均值為653.43 N(圖3)。

圖3 單軸密閉壓縮試驗結(jié)果Fig.3 Result of uniaxial closed compression test

2.4 土壤耕作阻力與耕作堆積角試驗

土壤耕作阻力與耕作堆積角試驗在中國農(nóng)業(yè)大學(xué)土槽實驗室進行。物理土槽(標(biāo)記為TC0)中土槽土壤制備過程見圖4：采用TCC電力四驅(qū)土槽試驗臺車的灑水功能對試驗所需的土槽土壤進行均勻的灑水作業(yè)，重復(fù)3遍灑水作業(yè)，待土壤存水濕透，靜置3天；檢查土壤含水率，確保土壤含水濕透，保證土壤含水率在17%左右，并使用1GQN-125型旋耕機將土壤旋耕，將表面土壤全面松碎并打勻，保證旋耕深度穩(wěn)定在20 cm；采用鎮(zhèn)壓輥保持額定前進速度對旋耕后的土壤進行壓實，鎮(zhèn)壓3遍，使用堅實度儀測量0～200 mm平均土層堅實度為247 kPa，烘干法測得土壤的含水率為16.09%。

耕作阻力測定試驗綜合考慮了挖掘犁鏟工作參數(shù)許用范圍與土槽工作環(huán)境，設(shè)置挖掘鏟耕作深度為150 mm，挖掘犁鏟前進速度為0.2 m/s，單次試驗前進距離為2 m。共進行6組試驗，耕作阻力值取穩(wěn)定區(qū)間數(shù)據(jù)點均值。耕作試驗過程見圖5。

耕作阻力表征土壤耕作特性的耕作難易程度，當(dāng)挖掘鏟前進速度為0.2 m/s，耕深150 mm時，耕作阻力數(shù)值為1275.05 N；耕作堆積角為動態(tài)堆積角，表征耕作中土壤狀態(tài)，在上述工況下測得耕作堆積角數(shù)值為32.29°。

圖4 土槽土壤制備過程Fig.4 Preparation process of soil

圖5 耕作試驗過程Fig.5 Tillage test process

3 土壤離散元模型構(gòu)建

3.1 土壤離散元模型參數(shù)確定

根據(jù)項目團隊前期研究結(jié)果并參考文獻[16]-[24]，以粒徑為13 mm的單球模型填充虛擬土槽(標(biāo)記為TC1)，土槽土層外形尺寸為1 500 mm×1 000 mm×220 mm，顆粒數(shù)量約20萬個。土壤本征參數(shù)：密度、剪切模量和泊松比。土壤顆粒接觸參數(shù)：碰撞恢復(fù)系數(shù)е取值0.2～0.6，顆粒間靜摩擦因數(shù)μs取值0.3～0.9，顆粒間滾動摩擦因數(shù)μr取值0.2～0.8。土壤顆粒EEPA接觸模型參數(shù)：粘附能Δγ取值12～36 J/m2，塑性變形比λp取值0.2～0.6，黏結(jié)分支指數(shù)Χ取值2.5～5.5，切向剛度因子Ktm取值0.3～0.7，參考文獻[17]確定恒定拖拉力f0為0，加載分支指數(shù)n1為1.5。以土壤靜態(tài)堆積角、土壤軸向壓力響應(yīng)面分析確定土壤離散元模型仿真參數(shù)優(yōu)化組合。

3.2 仿真參數(shù)顯著性分析

以土壤靜態(tài)堆積角和軸向壓力為指標(biāo)，開展Plackett-Burman試驗，分析碰撞恢復(fù)系數(shù)е、靜摩擦因數(shù)μs、滾動摩擦因數(shù)μr、粘附能Δγ、塑性變形比λp、 黏結(jié)分支指數(shù)Χ、切向剛度因子Ktm等7個主要仿真參數(shù)的顯著性。Plackett-Burman試驗結(jié)果見表1。

以α=0.05為顯著性水平，靜態(tài)堆積角方差分析見表2，土壤靜態(tài)堆積角模擬顯著性影響因素為塑性變形比λp、粘附能Δγ、滾動摩擦因數(shù)μr、切向剛度因子Ktm的P值為0.035 8，顯著性較低；軸向壓力方差分析見表3，軸向壓力模擬顯著性影響因素為碰撞恢復(fù)系數(shù)е、切向剛度因子Ktm和靜摩擦因數(shù)μs。

3.3 土壤靜態(tài)堆積角響應(yīng)面分析

采用Box-Behnken試驗設(shè)計方法，以塑性變形比λp、粘附能Δγ、滾動摩擦因數(shù)μr為試驗因素開展圓筒提升土壤靜態(tài)堆積角模擬試驗。EEPA接觸模型仿真參數(shù)中切向剛度因子Ktm表征模型彈性的參數(shù)，為簡化響應(yīng)面分析試驗方案，確定切向剛度因子Ktm為區(qū)間中值0.5。碰撞恢復(fù)系數(shù)е、靜摩擦因數(shù)μs、黏結(jié)分支指數(shù)Χ等參數(shù)分別取區(qū)間中值為0.4、0.6和4.0。圓筒提升模擬試驗采用內(nèi)徑200 mm、高度400 mm圓筒模型，粒徑13 mm單球顆粒自由填充形成土層。圓筒提升速度0.1 m/s，仿真時間設(shè)置為30 s。土壤靜態(tài)堆積角模擬結(jié)果見表4。

表1 Plackett-Burman試驗結(jié)果Table 1 Results Table of Plackett-Burman test

表2 靜態(tài)堆積角試驗方差分析Table 2 Variance analysis of static packing angle test

表3 軸向加載試驗方差分析Table 3 Variance analysis of axial loading test

表4 土壤靜態(tài)堆積角模擬結(jié)果Table 4 Simulation results of soil static accumulation angle

表5 靜態(tài)堆積角響應(yīng)面方差分析Table 5 Variance analysis of static packing angle response surface

3.4 土壤軸向壓力響應(yīng)面分析

采用Box-Behnken試驗設(shè)計方法，以碰撞恢復(fù)系數(shù)е、切向剛度因子Ktm、靜摩擦因數(shù)μs為試驗因素開展土壤單軸密閉壓縮仿真試驗。EEPA接觸模型仿真參數(shù)中塑性變形比λp、粘附能Δγ、滾動摩擦因數(shù)μr分別取優(yōu)化值：0.31、26.30和0.52，黏結(jié)分支指數(shù)Χ等參數(shù)分別取值為4。

單軸密閉壓縮模擬試驗采用內(nèi)徑為95 mm、高度150 mm圓筒模型，粒徑13 mm單球顆粒自由填充形成土層。下壓盤實體下降速度為60 mm/min，設(shè)置仿真時間12 s。土壤軸向壓縮仿真試驗結(jié)果如表6所示。

軸向加載力響應(yīng)面方差分析結(jié)果見表7：試驗因素Ktm、μs達到極顯著水平，因素е達到顯著水平，其他因素均不顯著。以實測軸向壓力653.43 N為目標(biāo)值，采用Design-Expert 12軟件分析，優(yōu)選顯著影響彈塑性的因素參數(shù)組合為：碰撞恢復(fù)系數(shù)е為0.58、切向剛度因子Ktm為0.39、靜摩擦因數(shù)μs為0.67。

表6 土壤軸向壓縮仿真試驗結(jié)果Table 6 Results of axial compression simulation test

表7 軸向壓力響應(yīng)面方差分析Table 7 Variance analysis of axial pressure response surface

3.5 土槽土壤-挖掘犁鏟離散元仿真試驗

采用優(yōu)化后的仿真參數(shù)構(gòu)建土槽土壤-挖掘犁鏟離散元仿真模型，挖掘犁鏟材料為65 Mn根據(jù)參考文獻[4]確定本征參數(shù)：泊松比ν1為0.3，密度ρ1為7 865 kg/m3，剪切模量G為79 GPa；土壤與挖掘犁鏟接觸參數(shù)：碰撞恢復(fù)系數(shù)е1為0.5、靜摩擦因數(shù)μs1為0.3、動摩擦因數(shù)μr1為0.1。設(shè)置挖掘犁鏟作業(yè)速度0.2 m/s、挖掘深度150 mm，仿真時間5 s，采樣間隔0.05 s。

挖掘犁鏟模擬耕作阻力與實測耕作阻力對比見圖6。在某一穩(wěn)定作業(yè)區(qū)間內(nèi)模擬耕作阻力均值為1 253.77 N，實測耕作阻力均值為1 275.05 N，誤差約為1.70%，此離散元模型較好的模擬了挖掘犁鏟的耕作阻力。

圖6 實測耕作阻力-虛擬耕作阻力對比Fig.6 Comparison between measured and virtual tillage resistance

4 顆粒球型對仿真過程的影響分析

4.1 顆粒球型結(jié)構(gòu)尺寸確定

在前述研究基礎(chǔ)上，以13 mm為等效粒徑確定雙球、三球、四球顆粒結(jié)構(gòu)尺寸，簡圖見圖7。雙球模型顆粒半徑為4.5 mm，(中心點坐標(biāo)長度單位為mm，下同)中心點坐標(biāo)分別為(2，0，0)、(-2，0，0)；三球模型顆粒半徑為4 mm，中心點坐標(biāo)分別為(0，0，2.5)、(-2.165，0，-1.25)、(2.165，0，-1.25)；四球模型的顆粒半徑為3.5 mm，中心點坐標(biāo)分別為(2.121，0，2.121)、(-2.121，0，2.121)、(2.121，0，-2.121)、(-2.121，0，-2.121)。

4.2 不同球型顆粒填充條件下虛擬土槽仿真試驗

4.2.1虛擬土槽建立與仿真平臺配置

此外，為了使新的生產(chǎn)系統(tǒng)能夠靈活、可靠地工作，并且能夠順利運轉(zhuǎn)物流過程，Leadec既使用現(xiàn)有的技術(shù)，同時也研發(fā)了自己的數(shù)字化解決方案，從而優(yōu)化自身以及客戶的內(nèi)部流程。例如，Leadec的專業(yè)人員正在為生產(chǎn)、物流和原材料的控制開發(fā)個性化的軟件解決方案。因為無論是在產(chǎn)品中還是在工廠基礎(chǔ)設(shè)施中的傳感器和執(zhí)行器相互聯(lián)系起來，都可以大大改變維護保養(yǎng)的措施。通過個性化的軟件解決方案，Leadec可以更好地預(yù)見、規(guī)劃這些改善措施，最終為客戶創(chuàng)造更好的經(jīng)濟效益。

采用單球、單球變粒徑、雙球、三球、四球、多球混合等填充方式建立6種虛擬土槽，標(biāo)記為TCi，各虛擬土槽的填充方式表示為：i=1，單球；i=2，單球變粒徑；i=3，雙球；i=4，三球；i=5，四球；i=6，多球混合。單球變粒徑填充時以粒徑13 mm為基準(zhǔn)，尺寸縮放比為0.5∶1.0∶1.5，填充數(shù)量比為3∶4∶3；多球混合填充時單球、雙球、三球、四球四種顆粒球型的尺寸縮放比和填充數(shù)量比為1∶1∶1∶1。

仿真平臺硬件配置：處理器為Intel(R) Xeon(R) Gold 6226R CPU @ 2.90 GHz(雙CPU，64核)，顯卡(GPU)為GeForce RTX 3080。在EDEM2020軟件中設(shè)置固定時間步1.2×10-5s、仿真時間5 s、寫出時間間隔0.01 s、Grid Size(網(wǎng)格大小)為3Rmin。為有效評估顆粒球型對計算效率的影響，仿真過程不使用動態(tài)計算域。

圖7 顆粒球型結(jié)構(gòu)尺寸簡圖Fig.7 Size diagram of granular spherical structure

4.2.2評價指標(biāo)

為合理評估離散元仿真模型仿真精度與計算效率，采用耕作堆積角相對誤差φi、耕作阻力相對誤差σi評價模型仿真精度；采用計算時間增長率τi評價模型仿真計算效率。

耕作堆積角相對誤差φi(i=1，2，…，6)由式(5)計算：

(5)

式中：αi為虛擬土槽TCi中耕作堆積角的模擬值，(°)；α0為土槽土壤耕作堆積角實測值，(°)。耕作阻力相對誤差σi由式(6)計算：

(6)

式中：Fi為虛擬土槽TCi中挖掘犁鏟耕作阻力模擬值，N；F0為挖掘犁鏟耕作阻力實測值。計算時間增長率τi由式(7)計算：

(7)

式中：Ti為虛擬土槽TCi仿真計算時間，s；T1為單球填充虛擬土槽TC1仿真計算時間，s。

4.3 試驗結(jié)果分析

4.3.1模型仿真精度

土壤耕作阻力與耕作堆積角試驗結(jié)果表明，物理土槽(TC0)土壤的耕作堆積角為32.29°、耕作阻力1 275.05 N，各虛擬土槽類型TCi中耕作堆積角模擬值αi測定見圖8。

圖8 不同土槽類型虛擬耕作堆積角模擬值的測定Fig.8 Determination of simulation values of virtual tillage accumulation angle for different soil groove types

各虛擬土槽條件下耕作堆積角模擬值與實測值的相對誤差值由大到小依次為：TC1>TC3>TC4>TC2>TC6>TC5(表8)。虛擬土槽TC2、TC6、TC5耕作堆積角模擬值相對誤差小于10%，能夠較好模擬土槽土壤的耕作堆積角。單球、雙球、三球顆粒模型結(jié)構(gòu)相對簡單，顆粒間互鎖作用較弱，虛擬土槽TC1、TC3及TC4的耕作堆積角模擬值相對誤差分別為20.32%、16.44%、13.04%，耕作堆積角仿真精度較低。

表8 耕作堆積角模擬值與實測值誤差分析Table 8 Error analysis between simulated and measured values of tillage accumulation angle

各虛擬土槽條件下耕作阻力模擬值與實測值的相對誤差值由大到小依次為：TC3>TC4>TC2>TC5>TC6>TC1(表9)，各虛擬土槽挖掘犁鏟耕作阻力模擬值的相對誤差值均小于8%，能較好模擬土槽土壤的耕作阻力。

表9 耕作阻力模擬值與實測值誤差分析Table 9 Error analysis between simulated and measured values of tillage resistance

4.3.2模型仿真計算效率

采用計算時間增長率τi評價模型仿真計算效率，不同虛擬土槽類型(TCi)仿真計算效率分析見表10。可知：同等尺寸虛擬土槽顆粒填充數(shù)量由大到小依次為：TC3>TC4>TC6>TC1>TC5>TC2。隨著顆粒表面凸起數(shù)量增加，接觸頻次增多，計算時間增長率由大到小依次為：TC5>TC6>TC4>TC2>TC3>TC1。雙球填充方式仿真時間增加率為140.01%、顆粒填充數(shù)量約為31萬個，四球填充方式仿真時間增長率為328.68%、顆粒填充數(shù)量增

表10 不同虛擬土槽類型仿真計算效率分析Table 10 Simulation calculation efficiency analysis of different virtual soil bin types

長率為4.23%。仿真效率受土壤顆粒填充數(shù)量、接觸關(guān)系及顆粒球型等因素影響，大計算規(guī)模時顆粒球型對仿真效率影響更顯著。

綜合考慮仿真精度及計算效率，土壤耕作阻力模擬宜選用單球顆粒填充方式、耕作堆積角模擬宜選用四球或多球混合填充方式建立虛擬土槽。

5 結(jié) 論

本研究以室內(nèi)土槽土壤為研究對象，采用EEPA接觸模型模擬土壤耕作特性，建立了6種虛擬土槽，研究顆粒球型對土壤離散元模型仿真精度及計算效率的影響，得到結(jié)論如下：

1)通過試驗建立起耕層土壤離散元模型，確定并優(yōu)化了基于EEPA接觸模型的單球顆粒土槽土壤離散元模型仿真參數(shù)：塑性變形比λp為0.31、粘附能Δγ為26.30 J/m2、滾動摩擦因數(shù)μr為0.52、碰撞恢復(fù)系數(shù)е為0.58、切向剛度因子Ktm為0.39、靜摩擦因數(shù)μs為0.67、黏結(jié)分支指數(shù)Χ為4。

2)依次進行了6種虛擬土槽土壤-挖掘犁鏟仿真試驗，耕作堆積角模擬值與實測值的相對誤差值由大到小為：TC1>TC3>TC4>TC2>TC6>TC5，其中TC2、TC6、TC5耕作堆積角模擬值相對誤差小于10%，適用于土壤耕作堆積角等耕作土壤表現(xiàn)狀態(tài)的研究；耕作阻力模擬值與實測值的相對誤差值由大到小為：TC3>TC4>TC2>TC5>TC6>TC1，6種虛擬土槽的耕作阻力模擬值與實測值相對誤差均小于8%，均可用于土壤耕作阻力的研究。

3)填充同等尺寸虛擬土槽所需顆粒填充數(shù)量由大到小依次為：TC3>TC4>TC6>TC1>TC5>TC2；土槽土壤-挖掘犁鏟仿真計算時間增長率由大到小依次為：TC5>TC6>TC4>TC2>TC3>TC1。顆粒填充方式顯著影響仿真效率，其中單球土槽仿真計算時間短，具有較優(yōu)的計算效率。