Determination of Discrete Element Modelling Parameters of a Paddy Soil with a High Moisture Content (>40%)

Abstract

Discrete element modelling (DEM) parameters are of great importance for the accurate prediction of soil properties and disturbance. This study aimed to provide an efficient method for accurately determining the DEM parameters of a paddy soil with a high moisture content (45.66%). The DEM parameters of the paddy soil modelled using the Hertz–Mindlin with JKR contact model were determined by using the Plackett–Burman (PB), steepest ascent, and central composite tests. The accuracies of the developed DEM models were evaluated using actual slumping tests. Based on the PB test, the surface energy of soil, coefficients of soil–soil rolling friction, and coefficients of soil–steel static friction exerted larger influences on the overall relative error (δ_ZH). The optimization results showed that to achieve a minimum δ_ZH (5.96%), the surface energy of soil, the coefficients of soil–soil rolling friction, and the coefficients of soil–steel static friction should be 0.808 J m⁻², 0.11, and 0.6, respectively. The optimized DEM model had an overall relative error of 7.27% with a coefficient of variation of 1.32%, indicating that the DEM parameters of the calibrated paddy soil had good accuracy.

1. Introduction

Rice is one of the main grain crops in China, and paddy field cultivation is an important process in rice production [1,2,3]. Precision tillage in paddy fields is conducive to improving water and fertilizer utilization, irrigation uniformity, and rice yield; reducing weeds, pests, diseases and environmental pollution; and promoting sustainable agricultural development [4,5,6]. Understanding the interaction between flat shovels and paddy soil is the key to optimizing flat shovels and to improving soil disturbance in paddy fields. Previous studies on soil–tool interactions can be grouped into three types: experimental, analytical, and numerical methods [7,8,9,10,11,12]. Experimental methods (e.g., field tests and lab soil bin tests) are time-consuming and costly, while analytical methods are limited by simple geometries and soil failure assumptions [10,11]. Numerical methods are able to examine complex tool geometries and more soil dynamic attributes [12]. Discrete element modelling (DEM) is a numerical method used for simulating the mechanical behavior of discontinuous particles (e.g., soil) without limitations on the magnitudes of particle displacement and tool shapes [10,13,14,15]. In recent years, DEM has been widely used as a modelling tool in agricultural research. Huang et al. [16] developed a subsoiling tool–soil interaction model using DEM, which could predict the soil’s microscopic behaviors at different positions with relative errors less than 20%. Hang et al. [17] investigated the effects of tine spacing of a subsoiler on soil disturbance behaviors (e.g., soil particle force) using DEM and found that tine spacing of 400 mm outperformed other spacings. Ding et al. [18] studied the performance of a deep tillage tool as affected by working depth under near-quasistatic conditions (working speed: 0.1 m/s) using DEM simulations. Tanaka et al. [19] simulated the soil loosening process under the action of a vibrating subsoiling tool based on DEM. Tamas and Bernon [20] simulated sweep–soil interactions using DEM simulations and validated them using field tests. Previous DEM studies mainly focused on simulations of disturbance behaviors of soil with a lower moisture content (<25%). By contrast, the soil in paddy fields before rice transplanting has a much higher moisture content (>40%), cohesion, and fluidity, which result in quite different disturbance behaviors. However, the DEM parameters, which are a prerequisite for simulating paddy soil–tool interactions, are absent in previous literature.

The DEM parameters of soil particles have been calibrated using angle of repose (AOR) tests in many previous studies. Wang et al. [21] determined the DEM parameters of soils with various moisture contents (from 0.27% to 22%) using both funneling angle of repose tests and slumping angle of repose tests. Qi et al. [22] obtained the AOR of a sandy loam soil using the lifting cylinder method and calibrated the soil DEM parameters by matching simulated and measured AORs. Barr et al. [23] calibrated the soil rolling friction coefficient and cohesive energy density when the AOR in DEM reached a good match with that measured. Ucgul et al. [24] determined the restitution coefficient between soil particles and yield strength based on the smallest relative error between forces from DEM simulations and measurements. Yang et al. [25] calibrated soil properties by listing a series of values of various soil properties and determined that the relative error between simulated and measured AORs was minimal. The angle of repose test is generally used to calibrate the DEM parameters of granular dryland soil with low moisture content (≤22%) [21,22,23]. However, the moisture content of paddy field soil is generally larger than 40% according to actual drying experiments. The moving behavior of paddy soil is very similar to that of concrete, and soil cohesion plays a key role in soil movement or flow. To characterize the moving behaviors of highly mobile paddy soil (moisture content: 45.66%) and determine the DEM parameters, a slumping test was proposed using a slumping cylinder with performance indices of the slumping value and slumping expansion. In the slumping test, both the horizontal and vertical soil disturbance behaviors could be quantified.

In this study, the DEM parameters’ interval and slumping values and the slumping expansion of paddy soil were measured with physical measurements. On this basis, the Plackett–Burman test, steepest ascent test, and Box–Behnken test were used to calibrate the specific values of the DEM parameters. The accuracy of the DEM parameters after calibration was validated by comparing the simulated and measured slumping values and slumping expansion values. It is expected that this result will provide a basis for the selection of the DEM parameters of paddy soil.

2. Materials and Methods

2.1. Laboratory Test

To characterize the moving behaviors of highly mobile paddy soil (moisture content: 45.66%) and determine the DEM parameters, a slumping test was carried out using a slumping cylinder in the laboratory, as shown in Figure 1.

The measuring method was based on the Chinese standard “GBT50080-2016: Standard test method for performance of ordinary concrete mixtures”. The dimensions of the funnel are shown in Figure 1b. Initially, the paddy field soil was poured into a funnel, and the funnel was then raised vertically at a speed of 0.06 m/s (manually controlled) (Figure 1a). Eventually, soil flowed out of the funnel under the action of gravity. The slumping value (H) and slumping expansion (W) were measured when there was no obvious movement of soil particles. The W was measured according to the maximum diameter of the soil bottom, while the H was calculated by the difference between the funnel height (H_f) and final soil height (h).

2.2. DEM Simulations

2.2.1. DEM Contact Model

The DEM simulations were performed using EDEM2.7 software (DEM Solutions Ltd. Edinburgh, Scotland, UK). The soil tested was a sandy clay soil that was highly cohesive. The Hertz–Mindlin with JKR (Johnson–Kendall–Roberts) contact model is a cohesive contact model that incorporates the influence of van der Waals forces in the contact area into the calculation. The JKR model uses the same calculation method of the Hertz–Mindlin model to calculate tangential elastic force, normal dissipative force, and tangential dissipative force. The realization of normal elastic contact force is based on the Johnson–Kendall–Roberts theory. Therefore, the tested soil was modelled using the Hertz–Mindlin with JKR model.

2.2.2. Model Parameter Calibration

Soil parameters in DEM mainly consist of two types: material and interaction properties. Material properties are the soil’s intrinsic characteristics, e.g., soil bulk density, shear modulus, shape and size of particles, and Poisson’s ratio. Interaction properties are characteristics exhibited by the particle in contact with boundaries, surfaces, and other (or similar) particles, e.g., surface energy, coefficients of restitution, and static and rolling friction [21]. The nominal soil particle radius was 2 mm, and the size range was 0.8–1.2 times the nominal particle size, as determined by a combination of preliminary tests and the literature [7,10].

Soil bulk density was determined by measurement. The shear modulus, Poisson’s ratio, and density of steel were determined according to previous studies [7,16]. The ranges of shear modulus, surface energy, and Poisson’s ratio of soil, as well as the coefficients of restitution, static friction, and rolling friction of soil–soil and soil–steel, were selected from the literature (Table 1) [26,27,28,29,30,31]. The parameters that significantly affect slumping value and slumping expansion of the paddy soil were screened out using the Plackett–Burman (PB) test.

The steepest ascent (SA) test was used to further narrow the value ranges of significant parameters [21]. A second-order regression model between indices (i.e., slumping value, slumping expansion, and overall relative error) and significant parameters was then established using a central composite test (CCT). The overall relative error (δ_ZH) was used to evaluate the overall error between simulated and experimental slumping values and slumping expansion. It was calculated as follows:

(1)${\delta }_1 = \frac{H_{\mathrm{s}}{ - H}_{\mathrm{e}}}{H_{\mathrm{s}}} \times 100\%$

(2)${\delta }_2 = \frac{W_s{ - W}_{\mathrm{e}}}{W_{\mathrm{s}}} \times 100\%$

(3)${\delta }_{ZH} = \frac{\left|{\delta }_1\right| + \left|{\delta }_2\right|}{2}$

where H_s and H_e represent the simulated and experimental slumping values in units of mm, and W_s and W_e represent simulated and experimental slumping expansion in units of mm.

Finally, the optimal combination of significance parameters concerning the indices’ optimization was calculated according to the second-order regression model.

3. Results and Discussion

3.1. Laboratory Test Results

The slumping test in the laboratory was repeated five times. As shown in Table 2, the slumping values ranged from 230.17 mm to 238.67 mm with a mean value of 233.5 mm. By contrast, slumping expansions ranged from 344.17 mm to 360.00 mm with a mean value of 350.8 mm. The coefficients of variation of slumping values and slumping expansion from different repetitions were 1.24% and 1.77%, respectively, indicating that the experimental results were reliable.

3.2. Plackett–Burman (PB) Test

The results of the Plackett–Burman (PB) test are shown in Table 3. According to the ANOVA outputs (Table 4), the contributions of X₂, X₆, and X₈ to overall relative error (δ_ZH) were relatively higher than those of the other parameters, so these three parameters were calibrated using the following steepest ascent (SA) and central composite tests.

3.3. Steepest Ascent (SA) Test

To further investigate the feasible range of soil parameters (i.e., X₂, X₆, and X₈), the SA test was designed and carried out in accordance with their effects on overall relative error (δ_ZH) (see Table 4). As shown in Table 5, the lowest relative error (9.91%) was at the fourth simulation, implying that more reasonable ranges were 0.60–1.10 for X₂, 0.05–0.09 for X₆, and 0.60–1.10 for X₈. The obtained feasible ranges of soil parameters were then utilized in the following central composite test.

3.4. Central Composite Test (CCT)

The specific effects of the surface energy of soil (X₂), coefficients of soil–soil rolling friction (X₆), and coefficients of soil–steel static friction (X₈) on the slumping value (H) and slumping expansion (W), as well as the relative errors between the simulated and measured values of both H and W (i.e., δ₁ and δ₂), were investigated using a central composite test. The test design and results are shown in Table 6. The regression equations between the indicators (i.e., H, W, δ₁, and δ₂) and X₂, X₆, and X₈ were obtained as follows:

H = 241.28 − 60.26X₂ − 242.37X₆ − 8.87X₈ + 37.77X₂² + 4002.61X₆² + 48.88X₈² + 200.69X₂X₆ − 52.7X₂X₈ − 574.31X₆X₈(4)

W = 556.59 − 263.69X₂ − 1342.58X₆ + 21.62X₈ + 94.23X₂² + 2397.09 X₆² − 2.94X₈² + 667.01X₂ X₆ − 39X₂X₈ + 208.68 X₆X₈(5)

δ₁ = −3.26 + 25.77X₂ + 103.29X₆ + 3.74X₈ − 85.76X₂X₆ + 22.57X₂X₈ + 245.49X₆X₈ − 16.16X₂² − 1710.1X₆² − 20.88X₈²(6)

δ₂ = 118.12 − 170.20X₂ − 666.68X₆ − 21.31X₈ + 289.58X₂X₆ + 0.69X₂X₈ + 154.86 X₆X₈ + 76.65X₂² + 1437.9X₆² + 4.98X₈²(7)

The coefficients of correlation were 0.9778, 0.9904, 0.9778, and 0.9873 for Equations (4)–(7), respectively. Moreover, the residual standard deviations were 1.4993, 2.3696, 0.6421, and 0.5992 for Equations (4)–(7), respectively. The ANOVA outputs (Table 7) showed that Equations (4)–(7) can be used to describe the relationships between indicators (H, W, δ₁ and δ₂) and X₂, X₆, and X₈.

3.5. Parameter Optimization and Validation

To obtain the optimal combination of soil DEM parameters, the minimum overall relative error (δ_ZH) was used as the target; the constraints were the feasible ranges of soil parameters determined from the SA test. The final optimization results showed that to achieve a minimum δ_ZH (5.96%), X₂, X₆, and X₈ should be 0.808 J m⁻², 0.11, and 0.6, respectively. To validate the optimization results, five slumping tests were carried out in DEM simulations with the optimized combination of X₂, X₆, and X₈ (Figure 2). The results showed that the overall relative error (δ_ZH) was 7.27% with a coefficient of variation of 1.32% (Table 8), indicating that the calibrated soil DEM parameters had good accuracy.

4. Conclusions

In this paper, the DEM parameters of a paddy field soil modelled using JKR were determined by a combination of the Plackett–Burman (PB), steepest ascent (SA), and central composite tests. The accuracies of DEM models were evaluated using slumping tests. The following conclusions were drawn:

(1) For a given paddy field soil, the surface energy of soil (X₂), coefficients of soil–soil rolling friction between (X₆), and coefficients of soil–steel static friction (X₈) have a greater influence on the overall relative error.

(2) An optimal combination of X₂ of 0.808 J m⁻², X₆ of 0.11, and X₈ of 0.6 was obtained given overall relative error minimization (5.96%).

(3) A low relative error (7.27%) was found for simulated overall relative error using a developed DEM model of paddy field soil, indicating that the calibrated DEM parameters were reliable. Paddy soil–flat shovel interactions using DEM simulations will be investigated in our future research.

Author Contributions

Methodology, H.Z., T.Z., X.W., L.H. and S.W.; software, H.Z., T.Z. and X.W.; validation, H.Z., T.Z. and X.W.; formal analysis, X.W., H.Z., X.L. and J.J.; investigation, H.Z., T.Z., X.W. and S.W.; resources, H.Z. and J.J.; data curation, H.Z. and T.Z.; writing—original draft preparation, H.Z. and T.Z.; writing—review and editing, H.Z., X.W., L.H. and X.L.; visualization, H.Z. and T.Z.; supervision, L.H. and X.L.; project administration, H.Z. and J.J.; funding acquisition, H.Z., X.W. and J.J. All authors have read and agreed to the published version of the manuscript.

Institutional Review Board Statement

Not applicable.

Data Availability Statement

The data reported in this study are contained within the article.

Conflicts of Interest

The authors of this manuscript declare that they have no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. The slumping test and collection of slumping value (H) and slumping expansion (W) in experiments (a) and DEM simulations (b). (Hf, h, D1, and D2 stand for funnel height, final soil height, lower opening diameter, and upper opening diameter of the funnel, respectively; vf stands for travel speed of the funnel).

Figure 2. (a) Simulated slumping value (200-99.80 mm) and slumping expansion (351.79 mm) in the DEM model with optimized parameters; (b) measuring slumping value and slumping expansion in the lab.

Table 1

Parameters and corresponding value ranges in the Plackett–Burman test.

Symbol	Parameter	Unit	Low Level(−1)	High Level(+1)
X 1	Poisson’s ratio of soil	-	0.32	0.42
X 2	Surface energy of soil	J/m2	0.1	0.9
X 3	Shear modulus of soil	Mpa	0.1	1
X 4	Coefficients of soil–soil restitution	-	0.01	0.15
X5	Coefficients of soil–soil static friction	-	0.1	0.9
X 6	Coefficients of soil–soil rolling friction	-	0.01	0.1
X 7	Coefficients of soil–steel restitution	-	0.01	0.15
X 8	Coefficients of soil–steel static friction	-	0.1	0.9
X 9	Coefficients of soil–steel rolling friction	-	0.01	0.1
X10, X11	Virtual parameters	-	-	-

Table 2

Results of the slumping test in the laboratory.

Test No.	Slumping Value (mm)	Slumping Expansion (mm)
1	230.17	360.00
2	238.67	356.33
3	234.00	347.50
4	233.00	346.00
5	231.67	344.17
Mean value	233.5	350.8
CV (%)	1.24	1.77

“CV” in the table represents coefficient of variation.

Table 3

Plackett–Burman test scheme and results.

No.	X 1	X 2	X 3	X 4	X 5	X 6	X 7	X 8	X 9	X 10	X 11	H (mm)	W (mm)	δZH (%)
1	0.32	0.9	0.3	0.15	0.9	0.01	0.15	0.9	0.1	−1	−1	220.45	469.56	19.72
2	0.32	0.9	1	0.15	0.1	0.01	0.01	0.9	0.01	1	1	208.75	413.83	14.28
3	0.42	0.1	0.3	0.01	0.9	0.01	0.15	0.9	0.01	1	1	231.94	569.14	31.45
4	0.42	0.1	1	0.15	0.1	0.1	0.15	0.9	0.01	−1	−1	227.16	532.57	27.27
5	0.42	0.9	0.3	0.15	0.9	0.1	0.01	0.1	0.01	1	−1	259.10	752.85	62.79
6	0.42	0.9	1	0.01	0.1	0.01	0.15	0.1	0.1	1	−1	263.53	719.80	59.02
7	0.42	0.1	1	0.15	0.9	0.01	0.01	0.1	0.1	−1	1	272.44	1244.42	135.71
8	0.42	0.9	0.3	0.01	0.1	0.1	0.01	0.9	0.1	−1	1	214.03	373.69	7.43
9	0.32	0.1	0.3	0.01	0.1	0.01	0.01	0.1	0.01	−1	−1	274.64	1243.41	136.03
10	0.32	0.1	1	0.01	0.9	0.1	0.01	0.9	0.1	1	−1	209.40	473.07	22.59
11	0.32	0.1	0.3	0.15	0.1	0.1	0.15	0.1	0.1	1	1	267.12	294.14	109.27
12	0.32	0.9	1	0.01	0.9	0.1	0.15	0.1	0.01	−1	1	262.11	107.87	60.06

Table 4

Significance analysis of parameters in PB test.

Parameter	Effect	Sum of Squares	Contribution/%	Significance Rank
X 1	−127.60	22,905.52	0.52	7
X 2	−796.73	122.11	20.16	2
X 3	−159.20	4760.88	0.81	6
X 4	174.87	190.08	0.97	5
X 5	−69.93	229.34	0.16	9
X 6	−356.00	36.68	4.03	3
X 7	−240.13	950.52	1.83	4
X 8	−1467.13	432.48	68.37	1
X 9	72.87	16,143.60	0.17	8

Table 5

Scheme and results of the steepest ascent test.

No.	X2 (J m−2)	X 6	X 8	δZH (%)
1	0.10	0.01	0.10	132.70
2	0.35	0.03	0.35	21.42
3	0.60	0.05	0.60	14.31
4	0.85	0.07	0.85	9.91
5	1.10	0.09	1.10	12.28

Table 6

Scheme and results of the test.

No.	X2 (J m−2)	X 6	X 8	H (mm)	δ1 (%)	W (mm)	δ2 (%)
1	1.08	0.1	1.08	182.50	21.84	338.68	3.45
2	1.08	0.1	0.72	193.07	17.32	343.78	2.00
3	1.08	0.06	1.08	180.04	22.90	343.61	2.05
4	1.08	0.06	0.72	188.14	19.43	343.89	1.97
5	0.72	0.1	1.08	189.80	18.71	366.29	4.42
6	0.72	0.1	0.72	199.34	14.63	358.51	−2.20
7	0.72	0.06	1.08	196.03	16.05	373.00	6.33
8	0.72	0.06	0.72	191.50	17.99	376.05	7.20
9	0.6	0.08	0.9	197.07	15.60	385.59	9.92
10	1.2	0.08	0.9	182.20	21.97	336.92	−3.96
11	0.9	0.05	0.9	187.77	19.59	359.58	2.50
12	0.9	0.11	0.9	192.07	17.74	350.62	−0.05
13	0.9	0.08	0.6	193.97	16.93	353.84	0.87
14	0.9	0.08	1.2	187.30	19.79	351.18	0.11
15	0.9	0.08	0.9	185.57	20.53	352.65	0.53
16	0.9	0.08	0.9	186.30	20.21	351.57	0.22
17	0.9	0.08	0.9	185.47	20.57	350.80	0.00
18	0.9	0.08	0.9	185.83	20.41	354.11	0.94
19	0.9	0.08	0.9	185.60	20.51	350.65	−0.04
20	0.9	0.08	0.9	186.23	20.24	351.91	0.32

Table 7

Variance analysis of the regression equations.

Indicators	Source of variation	SS	f	MS	F	Significance
H	Regression	491.64	9	54.63	25.73	Significant
	Residual	21.23	10	2.12
	Sum	512.87	19
δ 1	Regression	90.16	9	10.02	25.73	Significant
	Residual	3.89	10	0.39
	Sum	94.05	19
Y 2	Regression	2878.03	9	319.78	58.26	Significant
	Residual	54.87	10	5.49
	Sum	2932.91	19
δ 2	Regression	138.16	9	15.35	37.50	Significant
	Residual	4.09	10	0.41
	Sum	142.26	19

Table 8

Validation test.

No.	X 2	X 6	X 8	H (mm)	W (mm)	δZH (%)
1	0.808	0.11	0.6	199.17	351.09	7.39
2				200.37	350.85	7.1
3				200.62	352.78	7.32
4				200.50	352.17	7.26
5				200.32	352.06	7.29
Average				200.20	351.79	7.27

References

1. Aryal, J.P.; Mehrotra, M.B.; Jat, M.L.; Sidhu, S.H. Impacts of laser land leveling in rice–wheat systems of the north–western indo-gangetic plains of India. Food Secur.; 2015; 7, pp. 725-738. [DOI: https://dx.doi.org/10.1007/s12571-015-0460-y]

2. Jat, M.L.; Gathala, M.K.; Ladha, J.K.; Saharawat, Y.S.; Jat, A.S.; Kumar, V.; Sharma, S.K.; Kumar, V.; Gupta, R. Evaluation of precision land leveling and double zero-till systems in the rice-wheat rotation: Water use, productivity, profitability and soil physical properties. Soil Tillage Res.; 2009; 105, pp. 112-121. [DOI: https://dx.doi.org/10.1016/j.still.2009.06.003]

3. Agarwal, M.C.; Goel, A.C. Effect of field leveling quality on irrigation efficiency and crop yield. Agric. Water Manage.; 1981; 4, pp. 89-97. [DOI: https://dx.doi.org/10.1016/0378-3774(81)90033-0]

4. Towa, J.J.; Guo, X.P. Effects of irrigation and weed-control methods on growth of weed and rice. Int. J. Agric. Biol. Eng.; 2014; 7, pp. 22-33. [DOI: https://dx.doi.org/10.3965/j.ijabe.20140705.003]

5. Finney, C. The benefits of land leveling on irrigation schemes in Turkey and Sindh Province, Pakistan. ICID J.; 1996; 45, pp. 1523-1539.

6. Onodera, T.; Tsuruta, M.; Osari, H. Leveling rice paddies using a laser land leveler when tilling. Trans. Jpn. Soc. Irrig. Drain. Reclam. Eng. Jpn.; 2002; pp. 1-8.

7. Kan, Z.; Jin, H.; Hongwen, L.; Peisong, D.; Qingjie, W.; Hongbo, Z. Research on polyline soil-breaking blade subsoiler based on subsoiling soil model using discrete element method. Trans. Chin. Soc. Agric. Mach.; 2016; 47, pp. 62-72.

8. Zhang, L.; Zhai, Y.; Chen, J.; Zhang, Z.; Huang, S. Optimization design and performance study of a subsoiler underlying the tea garden subsoiling mechanism based on bionics and EDEM. Soil Tillage Res.; 2022; 220, 105375. [DOI: https://dx.doi.org/10.1016/j.still.2022.105375]

9. Wang, Y.; Li, N.; Ma, Y.; Tong, J.; Pfleging, W.; Sun, J. Field experiments evaluating a biomimetic shark-inspired (BioS) subsoiler for tillage resistance reduction. Soil Tillage Res.; 2020; 196, 104432. [DOI: https://dx.doi.org/10.1016/j.still.2019.104432]

10. Chen, Y.; Munkholm, L.J.; Nyord, T. A discrete element model for soil-sweep interaction in three different soils. Soil Tillage Res.; 2013; 126, pp. 34-41. [DOI: https://dx.doi.org/10.1016/j.still.2012.08.008]

11. Zeng, Z.W.; Chen, Y.; Zhang, X.R. Modelling the interaction of a deep tillage tool with heterogeneous soil. Comput. Electron. Agric.; 2017; 143, pp. 130-138. [DOI: https://dx.doi.org/10.1016/j.compag.2017.10.005]

12. Wang, X.; Li, P.; He, J.; Wei, W.; Huang, Y. Discrete element simulations and experiments of soil-winged subsoiler interaction. Int. J. Agric. Biol. Eng.; 2021; 14, pp. 50-62. [DOI: https://dx.doi.org/10.25165/j.ijabe.20211401.5447]

13. Sadek, M.A.; Chen, Y.; Zeng, Z.W. Draft force prediction for a high-speed disc implement using discrete element modelling. Biosyst. Eng.; 2021; 202, pp. 133-141. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2020.12.009]

14. Song, W.; Jiang, X.H.; Li, L.K.; Ren, L.L.; Tong, J. Increasing the width of disturbance of plough pan with bionic inspired subsoilers. Soil Tillage Res.; 2022; 220, 105356. [DOI: https://dx.doi.org/10.1016/j.still.2022.105356]

15. Wang, Y.; Zhang, D.; Yang, L.; Cui, T.; Jing, H.; Zhong, X. Modeling the interaction of soil and a vibrating subsoiler using the discrete element method. Comput. Electron. Agric.; 2020; 174, 105518. [DOI: https://dx.doi.org/10.1016/j.compag.2020.105518]

16. Huang, Y.; Hang, C.; Yuan, M.; Wang, B.; Zhu, R. Discrete element simulation and experiment on disturbance behavior of subsoiling. Trans. Chin. Soc. Agric. Mach.; 2016; 47, pp. 80-88. [DOI: https://dx.doi.org/10.6041/j.issn.1000-1298.2016.07.012]

17. Hang, C.; Gao, X.; Yuan, M.; Huang, Y.; Zhu, R. Discrete element simulations and experiments of soil disturbance as affected by the tine spacing of subsoiler. Biosyst. Eng.; 2018; 168, pp. 73-82. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2017.03.008]

18. Ding, Q.; Ren, J.; Adam, B.; Zhao, J.; Ge, S.; Li, Y. DEM analysis of subsoiling process in wet clayey paddy soil. Trans. Chin. Soc. Agric. Mach.; 2017; 48, pp. 38-48.

19. Tanaka, H.; Oida, A.; Daikoku, M.; Lnooku, K.; Sumikawa, O.; Nagasaki, Y.; Miyazaki, M. DEM simulation of soil loosening process caused by a vibrating subsoiler. Agric. Eng. Int. CIGR J.; 2007; 9, pp. 1-18.

20. Tamas, K.; Bernon, L. Role of particle shape and plant roots in the discrete element model of soilesweep interaction. Biosyst. Eng.; 2021; 211, pp. 71-96. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2021.09.001]

21. Wang, X.; Zhang, Q.; Huang, Y.; Ji, J. An efficient method for determining DEM parameters of a loose cohesive soil modelled using hysteretic spring and linear cohesion contact models. Biosyst. Eng.; 2022; 215, pp. 283-294. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2022.01.015]

22. Qi, L.; Chen, Y.; Sadek, M. Simulations of soil flow properties using the discrete element method (DEM). Comput. Electron. Agric.; 2019; 157, pp. 254-260. [DOI: https://dx.doi.org/10.1016/j.compag.2018.12.052]

23. Barr, J.B.; Ucgul, M.; Desbiolles, J.M.A.; Fielke, J.M. Simulating the effect of rake angle on narrow opener performance with the discrete element method. Biosyst. Eng.; 2018; 171, pp. 1-15. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2018.04.013]

24. Ucgul, M.; Fielke, J.M.; Saunders, C. Three-dimensional discrete element modelling (DEM) of tillage: Accounting for soil cohesion and adhesion. Biosyst. Eng.; 2015; 129, pp. 298-306. [DOI: https://dx.doi.org/10.1016/j.biosystemseng.2014.11.006]

25. Yang, Y.; Li, M.; Tong, J.; Ma, Y. Design and experiment of biomimetic rotary tillage blade based on multiple claws characteristics of mole rats. Trans. Chin. Soc. Agric. Eng.; 2019; 35, pp. 37-45. [DOI: https://dx.doi.org/10.11975/j.issn.1002-6819.2019.19.005]

26. Shi, L.; Zhao, W.; Sun, W. Parameter calibration of soil particles contact model of farmland soil in northwest arid region based on discrete element method. Trans. CSAE; 2017; 33, pp. 181-187. [DOI: https://dx.doi.org/10.11975/j.issn.1002-6819.2017.21.021]

27. Wu, T.; Huang, W.; Chen, X.; Ma, X.; Han, Z.; Fan, T. Calibration of discrete element model parameters for cohesive soil considering the cohesion between particles. J. South China Agric. Univ.; 2017; 38, pp. 93-98. [DOI: https://dx.doi.org/10.7671/j.issn.1001-411X.2017.03.015]

28. Wang, X.; Hu, H.; Wang, Q.; Li, H.; He, J.; Chen, W. Calibration method of soil contact characteristic parameters based on DEM theory. Trans. CSAM; 2017; 48, pp. 78-85. [DOI: https://dx.doi.org/10.6041/j.issn.1000-1298.2017.12.009]

29. Xiang, W.; Wu, M.; Lü, J.; Quan, W.; Ma, L.; Liu, J. Calibration of simulation physical parameters of clay loam based on soil accumulation test. Trans. CSAE; 2019; 35, pp. 116-123. [DOI: https://dx.doi.org/10.11975/j.issn.1002-6819.2019.12.014]

30. He, Y.; Xiang, W.; Wu, M.; Quan, W.; Chen, C. Parameters calibration of loam soil for discrete element simulation based on the repose angle of particle heap. J. Hunan Agric. Univ.; 2018; 44, pp. 216-220. [DOI: https://dx.doi.org/10.13331/j.cnki.jhau.2018.02.019]

31. Xing, J.; Zhang, R.; Wu, P.; Zhang, X.; Dong, X.; Chen, Y.; Ru, S. Parameter calibration of discrete element simulation model for latosol particles in hot areas of Hainan Province. Trans. CSAE; 2020; 36, pp. 158-166. [DOI: https://dx.doi.org/10.11975/j.issn.1002-6819.2020.05.018]