*Article* **Optimal Design and Analysis of Cavitating Law for Well-Cellar Cavitating Mechanism Based on MBD-DEM Bidirectional Coupling Model**

**Gaowei Xu 1,\*, Huimin Fang 2, Yumin Song <sup>1</sup> and Wensheng Du <sup>3</sup>**


**Abstract:** A variable velocity parallel four-bar cavitating mechanism for well-cellar can form the well-cellar cavitation which suits for well-cellar transplanting under a continuous operation. In order to improve the cavitating quality, this paper analyzed the structural composition and working principle of the cavitating mechanism and established the bidirectional coupling model of multi-body dynamics and the discrete element between the cavitating mechanism and soil through Recurdyn and EDEM software. Based on the model, a three-factor, five-level quadratic orthogonal rotational combination design test was conducted with the parameters of the cavitating mechanism as the experimental factors and the parameters of the cavitation as the response index to obtain the optimal parameter combination, and a virtual simulation test was conducted for the optimal parameter combination in order to study the cavitating law of the cavitating mechanism and soil. The test results showed that the depth of the cavitation was 188.6 mm, the vertical angle of the cavitation was 90.4◦, the maximum diameter of the cavitation was 76.1 mm, the minimum diameter of the cavitation was 68.5 mm, and the variance in the diameters for the cavitation was 5.42 mm2. The cavitating mechanism with optimal parameters based on the Recurdyn–EDEM bidirectional coupling mode could further improve the cavitating quality.

**Keywords:** well-cellar cavitating mechanism; MBD-DEM bidirectional coupling model; optimal design; cavitating law

#### **1. Introduction**

Well-cellar transplanting is a kind of agronomic technology which transplants the seedlings of some of the crops used for the planting of the hilly and mountainous areas, such as of tobacco [1] and pepper [2], into the well-cellar cavitation. The well-cellar cavitation is a cavity with a certain depth and height, as well as a consistent diameter. According to Darcy's law, the hydraulic gradient is generated in the well-cellar cavitation with the change in the external temperature and humidity, which forces the soil moisture in the well cellar to evaporate or condense, so as to keep the relative stability of the temperature and humidity inside the well-cellar cavitation [3]. Therefore, by taking the advantage of the relatively stable of the humidity and temperature in the well-cellar cavitation, watering, applying pesticides, and covering with a plastic film is immediately carried out after transplanting the crop seedlings into the well-cellar cavitation, which can realize the moderate early planting and deep planting inside the high ridge for the crop transplanting in hilly and mountainous areas, so as to increase the drought resistance and reduce the pests and diseases of transplanting seedlings [4]. At present, it has been partially promoted in the hilly and mountainous areas of Guizhou, Yunnan, Shandong, and Guangxi Province [5,6].

**Citation:** Xu, G.; Fang, H.; Song, Y.; Du, W. Optimal Design and Analysis of Cavitating Law for Well-Cellar Cavitating Mechanism Based on MBD-DEM Bidirectional Coupling Model. *Agriculture* **2023**, *13*, 142. https://doi.org/10.3390/ agriculture13010142

Academic Editors: Mustafa Ucgul and Chung-Liang Chang

Received: 12 December 2022 Revised: 2 January 2023 Accepted: 2 January 2023 Published: 5 January 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Nowadays, the transplanting equipment in the hilly and mountainous areas is mainly based on the duckbill transplanting machine according to the traditional planting agronomy [7–13]. It cannot form the stable well-cellar cavitation and it is difficult to adapt to the agronomic requirements of the well-cellar transplanting for the crops in the hilly and mountainous area. There are currently three types of cavitating devices for hilly and mountainous transplanters. First, the man-powered knapsack cavitating device [14] is driven by a gasoline engine and its cavitating position and depth is manually controlled, which results in a low cavitating efficiency and high labor intensity. Second, the shape of the cavitation formed by the multi-link cavitating mechanism is a "horn shape", which cannot meet the agronomic requirements of the well-cellar cavitation in the cavitating process [15,16] because it cannot match the forward velocity of the transplanter. Finally, the intermittent well-cellar cavitating device works with the intermittent powered chassis. As the powered chassis reaches the specified position, the device forms the cavitation along the direction which is perpendicular to the ridge surface [17,18]. The device avoids the matching problems between the cavitating mechanism and powered chassis, which formed the cavitation with a better quality. However, due to the intermittent movement of the chassis, the machines and tools have a low operating efficiency and reliability.

In the view of the above situation, the project team designed a kind of well-cellar-type cavitating mechanism for the hilly and mountainous transplanting machine [19], which is based on the noncircular gear parallel four-bar mechanism. It used the characteristics of a nearly uniform linear motion of the non-circular gear linkage to offset the forward working velocity of the mechanism. Under the continuous operation of the tool, the well-cellar cavitation which suited for transplanting the large depth and high specifications of the hilly and mountainous crops was formed, so as to improve the working efficiency and operation quality of the mechanized well-cellar cavitation and reduce the labor intensity of the artificial well-cellar cavitation.

At present, the cavitating mechanism was designed from the perspective of a mechanical design, and the interaction between the machine and soil was not studied. Accordingly, this paper used MBD-DEM bidirectional coupling technology to establish the interactive coupling model between the cavitating mechanism and soil in order to study the law of cavitation between the cavitating mechanism and the soil and optimize the parameters of the cavitating mechanism.

#### **2. Materials and Methods**

#### *2.1. Composition and Working Principle of the Cavitating Mechanism*

The structural diagram for the well-cellar cavitating mechanism is shown in Figure 1. The mechanism was mainly composed of the non-circular gear mechanism, the spur gear mechanism, the parallel four-bar mechanism, the cavitator, and the frame. One and two were the non-circular gear which constituted the non-circular gear mechanism. Three, four, and five were the straight tooth cylindrical gear which composed the spur gear mechanism. Six, seven, and eight were the double crank, frame, and linkage which composed the parallel four-bar mechanism. Nine and ten were the powershuttle mechanism and cavitator which composed the cavitating device.

The non-circular gear mechanism and the spur gear mechanism were installed on the frame, and the non-circular gear 2 was fixedly connected with the spur gear 3. The double crank 6 of the parallel four-bar mechanism was fixedly connected with the gears on both sides of the cylindrical gear mechanism 3 and 5, separately. The cavitator was fixedly connected with the connecting rod of the parallel four-bar mechanism.

**Figure 1.** Structural diagram for the well-cellar cavitating mechanism: 1. non-circular gear 1; 2. noncircular gear 2; 3. spur gear 3; 4. spur gear 4; 5. spur gear 5; 6. crank; 7. linkage; 8. frame; 9. cavitating device; 10. cavitator.

The working principle of the cavitating mechanism is shown in Figure 2. During the operation, the power of uniform rotation was transmitted to the non-circular gear 1. The non-circular gear mechanism converted the uniform rotation into the variablevelocity rotation and transmitted it to the spur gear mechanism. The spur gear mechanism transmitted the variable velocity rotation to the double crank of the parallel four-bar mechanism and drove its connecting rod to oscillate at a variable velocity, and the cavitating device oscillated with the connecting rod. In the process of invading and extracting from the soil, the cavitator in the horizontal direction formed the velocity, which was similar to the forward velocity and opposite to the direction of the mechanism, so as to offset the forward velocity of the well-cellar cavitating mechanism and form the well-cellar cavitation with a uniform diameter. The cavitating device oscillated with the connecting rod and the power of the uniform rotation was transmitted to the cavitator through the flexible shaft. It drove the cavitating device to rotate at a constant velocity so as to compact the inner wall of the cavitation during the cavitating process. The shape of the well-cellar cavitation formed by the cavitating mechanism was approximately cylindrical, the depth of the well-cellar cavitation was 180~200 mm, and the diameter was 60~80 mm. The overall diameter of the well-cellar cavitation should be uniform, and its inner wall should be solid without the collapse.

#### *2.2. Establishment of Bidirectional Coupling Model for Cavitating Mechanism Based on DEM and MBD*

The bidirectional coupling of the discrete element method and multi-body dynamics could realize the interaction between complex mechanisms and discrete particles and has been widely used in the field of agricultural machinery [20,21]. The cavitating process of the cavitation under the action of the cavitating mechanism was the bidirectional interaction between the mechanism and the soil. Therefore, the method of the bidirectional coupling of DEM and MBD was used to analyze the cavitating process of the mechanism [22].

#### 2.2.1. Multibody Dynamics Modeling

The multi-body dynamics model of the cavitating mechanism was established in the Recurdyn software, and the kinematic pairs between the components are set as shown in Table 1.

**Figure 2.** Working principle of the cavitating mechanism. Notes: *v* represents the operating velocity of the mechanism.


**Table 1.** Kinematic pair configuration table of multi-body dynamics model.

In the model, the forward velocity of the mechanism was set as 0.25 m/s, the rotary velocity of the non-circular 1 gear was set as 30 r/min, and the rotary velocity of the cavitating device was set as 500 r/min.

#### 2.2.2. Discrete Element Modeling

Soil Particle Model

The basic structure of the soil particles for the plough layer mainly included block particles, core particles, flake particles, and columnar particles after tillage. In order to improve the authenticity of the model, spherical, linear, triangular, and tetrahedral particles were established in the discrete element software to approximate the massive, nuclear, flaky, and columnar soil particles, respectively [23,24].

The filling radius of the soil particles was set to 5 mm, and the model parameters of the soil are seen in Table 2 [25–27]. The soil particle parameters were obtained by referring to the relevant literature [28,29]. The soil moisture content and soil accumulation angle of the selected test site were measured. The basic parameters of the selected discrete element model are shown in Table 3.


**Table 2.** Model parameters of soil.

**Table 3.** Basic parameters of discrete element model.


Soil Particle Contact Model

The EEPA (Edinburgh Elasto-Plastic Adhesion) contact model includes the plasticity and viscosity of the particles and is suitable for simulating the farmland soils with a strong plasticity [30]. The soil particle contact model was set as the EEPA contact model, and the contact model between the cavitator and soil was set as the Hertz Mindlin with the JKR contact model with reference to the soil simulation routine of theEDEM 2020 (DEM Solutions Ltd., Edinburgh, Scotland, UK). According to the related literatures [31–33] on the discrete element parameters calibration of the physical characteristics of the soil particles, the configuration of the soil particles' simulation parameters is shown in Table 4.


**Table 4.** Contact model parameters of soil.

2.2.3. Establishment of Multi-Body Dynamics-Discrete Element Bidirectional Coupling Model of Cavitating Mechanism

The bidirectional coupling interface between Recurdyn and EDEM was set up. On the basis of the multi-body dynamics model for the cavitating mechanism and the discrete element model of the cavitating soil, the bidirectional coupling model of Recurdyn–EDEM for the cavitating mechanism was constructed; the simplified model is shown in Figure 3.

**Figure 3.** Bidirectional coupling model of Recurdyn–EDEM for the cavitating mechanism: (**a**) front view of model in Recurdyn; (**b**) axonometric drawing of bidirectional coupling model; (**c**) bidirectional coupling model in EDEM.

#### *2.3. Orthogonal Simulation and Experimental Optimization*

In this section, through the multi-body dynamics-discrete element bidirectional coupling model between the cavitating mechanism and the soil, the regression equation between the optimization variables and optimization objectives of the cavitating mechanism was obtained in combination with the experimental design optimization method. On the basis of the regression equation, the multi-objective function was optimized to obtain the optimized combination of the parameters for the cavitating mechanism.

#### 2.3.1. Determination of Test Factors and Indicators

According to the kinematic equation, the crank's length *l*<sup>2</sup> of the parallel four-bar, the eccentricity *e* of the non-circular gear, and the deformation coefficient *m*<sup>1</sup> of the non-circular gear in the kinematic equation were the significant factors which affected the quality of the cavitation. The standard of the cavitating quality was based on the results of the agronomic analysis for the cavitation, the depth of cavitation *y*1, the vertical angle between the line connecting the top midpoint to the bottom midpoint of cavitation and the line of the ridge surface *y*2, the maximum diameter of cavitation *y*3, the minimum diameter of the cavitation *y*4, the difference between the maximum and minimum value diameter of the cavitation *y*5, and the variance in the diameters of the cavitation *y*<sup>6</sup> were used as the evaluation indicators [19]. In the above indicators, the perpendicularity for the cavitation was measured by the vertical angle of cavitation, the maximum diameter of cavitation, the minimum diameter of cavitation, the difference between the maximum and minimum diameter of cavitation, the variance in the diameters which measured the overall range of the diameter, and the consistency of the contour diameter for cavitation.

#### 2.3.2. Test Contents and Methods

According to the single factor pre-test combined with the test index, the length range of the crank for the cavitating mechanism was 175~195 mm, the eccentricity of non-circular gear was 0.35~0.55, the deformation coefficient of the non-circular gear was 1.2~1.4, and the level of the factors was determined as five levels. On this basis, the optimized parameter combination of the cavitating mechanism was determined by using the bidirectional coupling model of Recurdyn–EDEM for the cavitating mechanism which was combined with the combination test method of the quadratic orthogonal rotation center. The diameter data of the 11 groups for the cavitation were measured through the bidirectional coupling model and the data were calculated to obtain the evaluation index. The test factor codes are shown in Table 5.

**Table 5.** Experimental factors and codes.


#### 2.3.3. Test Results and Analysis

The test scheme and results are shown in Table 6, the analysis of variance for each response index is shown in Table 7, and *x*1, *x*2, and *x*<sup>3</sup> were the factor coding values.




**Table 7.** Variance analysis of regression model.

Notes, \*\* indicates that the difference is extremely significant (*p* < 0.01),\* indicates that the difference is significant (0.01 ≤ *p* ≤ 0.05).

Regression analysis and factor variance analysis were carried out on the experimental data by Design-Expert 8.0 software (Stat-Ease Ltd., Minneapolis, MN, USA). After eliminating the insignificant items in the interactive and quadratic term, the regression equation between the index and the factor coding value was obtained:

$$y\_1 = 183.60 + 5.33 \mathbf{x}\_1 + 0.36 \mathbf{x}\_2 - 0.25 \mathbf{x}\_3 + 1.76 \mathbf{x}\_1 \mathbf{x}\_3 + 1.21 \mathbf{x}\_3^2$$

$$y\_2 = 92.15 + 0.16 \mathbf{x}\_1 + 1.51 \mathbf{x}\_2 + 2.17 \mathbf{x}\_3 + 1.46 \mathbf{x}\_2 \mathbf{x}\_3$$

$$y\_3 = 73.57 + 0.38 \mathbf{x}\_1 + 1.53 \mathbf{x}\_2 - 0.81 \mathbf{x}\_3 - 6.95 \mathbf{x}\_1 \mathbf{x}\_2 + 5.33 \mathbf{x}\_1^2 + 11.75 \mathbf{x}\_2^2$$

$$y\_4 = 64.84 - 3.32 \mathbf{x}\_1 + 4.13 \mathbf{x}\_2 + 0.19 \mathbf{x}\_3 - 5.29 \mathbf{x}\_1 \mathbf{x}\_2 + 5.54 \mathbf{x}\_2^2 - 1.44 \mathbf{x}\_3^2$$

$$y\_5 = 9.22 + 3.84 \mathbf{x}\_1 - 2.44 \mathbf{x}\_2 - 1.21 \mathbf{x}\_3 - 1.79 \mathbf{x}\_1 \mathbf{x}\_2 + 4.67 \mathbf{x}\_1^2 + 6.03 \mathbf{x}\_2^2 + 1.06 \mathbf{x}\_3^2$$

$$y\_6 = 8.73 + 23.07 \mathbf{x}\_1 - 9.94 \mathbf{x}\_2 - 5.42 \mathbf{x}\_3 - 6.16 \mathbf{x}\_1 \mathbf{x}\_2 + 19.63 \mathbf{x}\_1^2 + 27.48 \mathbf{x}\_2^2$$

In order to analyze the influence of each factor on each response index, the response surface was obtained by using Design-Expert 8.0 software, as shown in Figure 4.

Figure 4a was the influence of the crank's length and non-circular gear deformation coefficient on the depth of the cavitation at the zero level of other factors. The deformation coefficient for the non-circular gear was below the zero level; the depth of cavitation showed an upward trend with the increase in the crank's length. The deformation coefficient *m*<sup>1</sup> for the non-circular gear was above the zero level, the depth of cavitation showed an upward trend with the increase in the crank's length, and the upward trend was significantly larger than that of the former. The reason was that the longer the crank's length, the lower the position of the cavitator in the vertical direction, and the greater the depth of the cavitation. The deformation coefficient for the non-circular gear changed the velocity of the cavitator in the horizontal direction, resulted in the change in the angle for the cavitation, and then affected the vertical depth of the cavitation. The larger the deformation coefficient for the non-circular gear, the smaller the influence on the verticality of the cavitation, and the greater the vertical depth of the cavitation. The crank's length was below the zero level; the depth of the cavitation decreased with the rise of the non-circular gear deformation coefficient. The reason was that the crank's length changed below the zero level, the interaction between the non-circular gear deformation coefficient and the crank's length reduced the velocity of the cavitator in the horizontal direction, and then increased the inclined angle for the cavitation, which resulted in the decrease in the vertical depth for the cavitation. The crank's length was above the zero level; the length of the crank increased with the rise of the non-circular gear deformation coefficient. The reason was that the crank's length changed above the zero level, the interaction with the deformation coefficient of the non-circular gear continuously increased the velocity of the cavitator in the horizontal direction and gradually approached the forward velocity of the mechanism so that the inclined angle of the cavitation was continuously reduced and the vertical depth of the cavitation was increased.

**Figure 4.** Response surface of factor to test index: (**a**) influence of crank's length and non-circular gear deformation coefficient on the depth of cavitation at the zero level of other factors; (**b**) influence of the deformation coefficient and the eccentricity for the non-circular gear on the vertical angle of cavitation at the zero level of other factors; (**c**) influence of crank's length and non-circular gear eccentricity on the maximum diameter of cavitation at the zero level of other factors; (**d**) influence of crank's length and non-circular gear eccentricity on the minimum diameter of cavitation at the zero level of other factors; (**e**) influence of crank's length and non-circular gear eccentricity on the difference between the maximum and minimum diameter of cavitation at the zero level of the other factors; (**f**) influence of the length for the crank and the eccentricity for the non-circular gear on the variance in diameters for cavitation at the zero level of other factors.

Figure 4b was the influence of the deformation coefficient and the eccentricity for the non-circular gear on the vertical angle between the line connecting the top midpoints to the bottom midpoint of cavitation and the line of the ridge surface at the zero level of the other factors. The deformation coefficient for the non-circular gear was below the zero level and the vertical angle of the cavitation and the line of the ridge surface showed a significant upward trend with the increase in the eccentricity for the non-circular gear. The deformation coefficient for the non-circular gear was above the zero level and the vertical angle of cavitation and the line of the ridge surface showed an upward trend significantly with the increase in the eccentricity of the crank's non-circular gear. The eccentricity of the non-circular gear was below the zero level and the vertical angle of cavitation and the line of the ridge surface showed an upward trend with the increase in the deformation coefficient for the non-circular gear. The eccentricity of the non-circular gear was above the zero level and the vertical angle of cavitation showed the upward trend with the increase in the deformation coefficient of the non-circular gear. The reason was that the increase in the deformation coefficient and the eccentricity for the non-circular gear constantly changed the velocity of the cavitator in the horizontal direction, which resulted in the continuous increase in the vertical angle of cavitation.

Figure 4c,d were the influence of the crank's length and non-circular gear eccentricity on the maximum and minimum diameter of cavitation at the zero level of the other factors. The eccentricity for the non-circular gear was above the zero level, the maximum and minimum diameters of cavitation both showed the downward trend with the increase in the crank's length. The non-circular gear coefficient was below the zero level and the maximum and minimum diameters of cavitation showed a rising trend with the increase in the crank's length. The reason was that the increase in the crank's length led to the increase in the difference between the velocity of the crank in the horizontal direction and the forward velocity of the mechanism, which made the maximum and minimum value of the cavitation diameter show a rising trend. With the increasing eccentricity for the non-circular gear (above the zero level), the influence on the velocity of the crank's length in the horizontal direction was increasing, which was reducing the velocity of the crank's length in the horizontal direction constantly, so that the maximum and minimum value of the cavitation diameter showed a rising trend. The length of the crank was above the zero level and the maximum and minimum diameter of cavitation showed a downward trend first and then an upward trend with the increase in the coefficient for the non-circular gear. The crank's length was below the zero level and the maximum diameter of the cavitation diameter showed a rising trend with the increase in the crank's length. The reason was that the increase in the non-circular gear eccentricity led to the increase in the difference between the horizontal velocity of the crank and the forward velocity of the mechanism, which made the maximum and minimum diameter of the cavitation diameter show a rising trend. With the increasing length of the crank (above the zero level), the increase in the eccentricity for the non-circular gear led to the difference between the velocity of the crank in the horizontal direction and the forward velocity of the mechanism, which decreased first and then increased, so that the maximum and minimum diameter of cavitation showed a downward trend first and then an upward trend.

The eccentricity of the non-circular gear was below the zero level, which resulted in there being less influence on the velocity of the crank's length in the horizontal direction, and the diameter of cavitation showed an upward trend. As the eccentricity of the noncircular gear increased continuously (above the zero level), the influence on the velocity of the crank's length in the horizontal direction also increased continuously, which reduced the velocity of the crank's length in the horizontal direction, and the diameter of cavitation showed a downward trend. The eccentricity of the non-circular gear was above the zero level and the maximum and minimum diameter of cavitation showed a downward trend with the increase in the crank's length. The non-circular gear coefficient was below the zero level and the maximum and minimum diameter of cavitation showed an upward trend with the increase in the crank's length. This was because the velocity of the crank's length *l*<sup>2</sup> in the horizontal direction increased and the difference between the velocity of the crank's length and the velocity of the forward direction of the mechanism also increased, which resulted in the increase in the diameter for the formed cavitation. The eccentricity of the

non-circular gear was below the zero level, the effect on the velocity of the crank's length in the horizontal direction was less, and the diameter of cavitation showed an upward trend. As the eccentricity of the non-circular gear increased continuously (above the zero level), the influence on the velocity of the crank's length in the horizontal direction increased continuously, which reduced the velocity of the crank in the horizontal direction, and the diameter of cavitation showed a downward trend.

Figure 4e was the influence of the crank's length and non-circular gear eccentricity on the difference between the maximum and minimum diameter of cavitation at the zero level of the other factors. The eccentricity for the non-circular gear was constant and the difference between the maximum and minimum diameter of cavitation showed a downward trend first and then an upward trend with the increase in the length of the crank. The length of the crank was constant and the difference between the maximum and minimum diameter of cavitation showed a downward trend first and then an upward trend with the increase in the eccentricity of the non-circular gear. It meant that the crank's length and the eccentricity of the non-circular gear were near the zero level and the velocity of the cavitating mechanism in the horizontal direction was close to the forward velocity of the mechanism, so as to minimize the difference between the maximum and minimum diameter of cavitation. The eccentricity of the non-circular gear was below the zero level and the maximum and minimum diameter of cavitation showed an upward trend with the increase in the crank's length. The crank's length and the eccentricity of the non-circular gear was near the zero level, the velocity of the cavitating mechanism in the horizontal direction was closed to the forward velocity of the mechanism, and the difference between the maximum and minimum value of the hole diameter was the least.

Figure 4f was the influence of the length for the crank and the eccentricity for the non-circular gear on the variance in the diameters for cavitation at the zero level of the other factors. The length of the crank was constant and the variance in the diameters for cavitation showed a downward trend first and then an upward trend with the increase in the eccentricity *e* of the non-circular gear. The eccentricity of the non-circular gear was constant and the variance in the diameters for cavitation showed a downward trend first and then an upward trend with the increase in the crank's length. It meant that the change in the crank's length and non-circular gear eccentricity led to the fluctuation of the difference between the velocity of the cavitator in the horizontal direction and the forward velocity of the mechanism, which resulted in a large change in the diameter of each part of the cavitation and increased the variance in the diameters for cavitation.

#### 2.3.4. Parameters Optimization

In order to obtain the best level combination of the experimental factors, the experimental factors of the optimal design were carried out, and the mathematical model of the parameter's optimization was established. The optimal target of the evaluation index was determined according to the agronomic requirements for the well-cellar cavitation with the depth of 180~200 mm, diameters of 60~80 mm, uniform overall diameter, and a high verticality. For the combination with the boundary conditions of the experimental factors, the regression equation of the response index for the cavitating mechanism was analyzed to obtain the optimization model of the nonlinear programming as follows:

$$\begin{cases} \begin{array}{l} 180 \text{ mm} \le y\_1 \le 185 \text{ mm} \\ 85^\circ \le y\_2 \le 95^\circ \\ 60 \text{ mm} \le y\_3 \le 80 \text{ mm} \\ 60 \text{ mm} \le y\_4 \le 80 \text{ mm} \\ \text{min} y\_5 \\ \text{min} y\_6 \\ \text{s.t.} \begin{cases} 175 \text{ mm} \le l\_1 \le 195 \text{ mm} \\ 0.35 \le e \le 0.55 \\ 1.2 \le m\_1 \le 1.4 \end{cases} \end{array} \end{cases}$$

This section used the Design-Expert 8.0 software optimization module to optimize the multi-objective parameters of the regression equation. The optimized parameters of the cavitating mechanism were obtained, the length of the crank *l*<sup>2</sup> was 185.71, the eccentricity *e* of the non-circular gear was 0.456, and the deformation coefficient *m*<sup>1</sup> was 1.2857.

#### **3. Results and Analysis**

The bidirectional coupling model of Recurdyn and EDEM was used to simulate the cavitating mechanism under the optimal parameters. The cavitating process of the cavitating mechanism is shown in Figure 5. The cavitator was driven by the cavitating mechanism to invade the soil and the cavitator rotated under the drive of the driving device to interact with the soil simultaneously.

**Figure 5.** Cavitating process of the cavitating mechanism: (**a**) process of the cavitator starting to interact with the soil driven by the cavitating mechanism until reaching the maximum depth of the cavitation; (**b**) process of the cavitator extracting from the cavitation driven by the cavitating mechanism after the cavitator reached the maximum cavitating depth.

Figure 5a showed the process of the cavitator starting to interact with the soil until it reached the maximum depth of cavitation. As there was lack of the compaction and restriction of the other soil layers on the top soil layer, the gap between the soil particles was large, but the force between the soil particles was opposite. After the cavitator invaded the soil, the soil of the top layer near the cavitator rose, and the rotation of the cavitator would rotate and squeeze the soil particles, which made the soil particles generate the centrifugal force. Under the action of the centrifugal force, the soil particles on the top layer were thrown out and compacted to leave the cavitator, which resulted in a larger cavitating diameter for the soil of the top layer. With the increase in the cavitating depth, the force between the soil particles gradually increased, and the diameters of cavitation gradually decreased until it was the same as the diameter of the cavitator.

Figure 5b shows the process of the cavitator extracting from the cavitation driven by the cavitating mechanism after the cavitator reached the maximum cavitating depth. As the cavitator extracted from the cavitation, the velocity of the endpoint for the cavitator in the horizontal direction could not offset the forward velocity of the cavitating mechanism completely, which resulted in the small forward movement of the cavitator. At this time, the rotating cavitator acted on the front side of the cavitation, which changed the shape of the hole for a well-cellar cavitation from the circular to the approximate elliptical, as shown in Figure 6a. As the cavitator only acted on the front side of the well-cellar cavitation, the soil of the top layer on the right side for the well-cellar cavitation was higher than that of the left side. The sectional shape of the well-cellar cavitation is shown in Figure 6b.

**Figure 6.** Three-dimensional shape of the well-cellar cavitation: (**a**) shape of the hole for the well-cellar cavitation; (**b**) sectional shape of the well-cellar cavitation.

According to the optimal parameters, the discrete element model of cavitation was obtained through the bidirectional coupling model of the cavitating mechanism and the soil (see Figure 6). The measuring tools in the EDEM software were used to measure the parameters of cavitation, as shown in Table 8.

**Table 8.** Parameter comparison before and after optimization of mechanism parameters.


The depth of the cavitation *y*<sup>1</sup> after the optimization increased by 4.3% compared with that obtained before the optimization; the vertical angle of the cavitation *y*<sup>2</sup> before and after the optimization was basically unchanged; the maximum diameter of the cavitation *y*<sup>3</sup> after the optimization increased by 5.9% compared with that obtained before the optimization; the minimum diameter of the cavitation *y*<sup>4</sup> after the optimization increased by 6.3% compared with that obtained before the optimization; the difference between the maximum and minimum diameter of the cavitation *y*<sup>5</sup> before and after the optimization was basically unchanged; and the variance in the diameter for the cavitation *y*<sup>6</sup> decreased by 10.4% compared with that obtained before the optimization. The cavitating mechanism with optimal parameters based on the Recurdyn–EDEM bidirectional coupling mode could further improve the quality of the cavitation.

#### **4. Conclusions**

The work reported in this article established the multi-body dynamics and discrete element bidirectional coupling model of the cavitating mechanism and soil through the software of Recurdyn and EDEM. With the aid of the model for the cavitating mechanism and soil, combined with the quadratic orthogonal rotation center combination test method, the regression equation between the parameters of the cavitating mechanism and cavitation was established through the Design Expert 8.0 software, and the influence trend and interaction relationship between the parameters of the cavitating mechanism and the cavitation were obtained through the response surface. On the basis of the regression equation, the optimal parameters combination of the cavitating mechanism was obtained by using multi-objective function optimization. According to the optimal parameters combination of the cavitating mechanism, the multi-body dynamics and discrete element bidirectional coupling model of the cavitating mechanism and soil was used for virtual simulation tests, and the cavitating law of the mechanism was analyzed. The test results showed that the cavitating mechanism with optimal parameters based on the Recurdyn– EDEM bidirectional coupling mode could further improve the cavitating quality.

**Author Contributions:** Conceptualization, G.X. and H.F.; data curation, G.X. and Y.S.; formal analysis, W.D.; investigation, G.X. and Y.S.; methodology, H.F. and Y.S.; project administration, H.F.; resources, G.X. and H.F.; software, G.X. and W.D.; supervision, G.X. and Y.S.; validation, G.X., H.F. and W.D.; visualization, G.X. and W.D.; writing—original draft, G.X. and H.F.; writing—review and editing, Y.S. and W.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was supported by the Optimal design and experiment of parallel Combinationplanting mechanism (202211510017), which is the National College Students Innovation and Entrepreneurship Training.

**Institutional Review Board Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available on demand from the first author at (202107@sdjtu.edu.cn).

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
