Abstract
Dynamic penetration test by a penetrator is an important detection method for lunar exploration. The analysis of the dynamic penetration properties of the penetrator significantly helps to design the penetrator structure. In this article, the dynamic penetration mechanical model was established according to soil mechanics. Then, the discrete-element simulation model, which was established using the application program interfaces of the EDEM software, was used to simulate the penetration process. The penetration experiment demonstrated the accuracy of the dynamic penetration model and the discrete-element simulation model. This study of dynamic penetration model of a penetrator may provide technical support to the deep space exploration.
Introduction
There are many mineral resources in the Moon, and the Moon is the nearest celestial body to the Earth. Moon explorations can improve the development of science and technology in a short period of time, help to harvest abundant resources, to develop greater knowledge about the origin of the cosmos, and to establish bases on the Moon in the future.
To acquire lunar soil samples, the main sampling methods consist of surface sampling and shallow depth sampling.1,2 The main surface-sampling methods are digging sampling and hammering sampling; the main shallow-depth-sampling method is drilling. These methods have been performed by the United States and former Soviet Union in lunar exploration.3,4
These sampling methods have the limitation of a high cost and low return. In addition, the soil samples cannot maintain their original state when these sampling methods are used. However, the penetration devices have smart structures, low power consumption, and light weight, so they can be used in challenging environments. Penetration devices have the advantages of low cost and high efficiency and have become a hotspot of study.
In 1996, the Space Research Centre of the Polish Academy of Sciences began the study of automatic penetrator used for the exploration of extraterrestrial objects, and the penetrator named MUPUS was designed to explore the comet.5,6
During 2001 and 2003, the ESA and the DLR designed the penetrator used for Martian exploration based on the penetrator designed by Transmash. The outline dimension is φ27 mm × 280 mm, the mass is 890 g, the power is 3 W, and the penetration depth is 2 m. Besides, a sample collection device, which can collect soil for 0.2 cm 3 , was on the front part of the penetrator.7–9
Then, from 2003 to 2006, according to the Martian penetration exploration mission, Italian Galileo Avionica asked the DLR to research the penetrator again. And the DLR designed the IMS penetrator based on the PLUTO penetrator, which had smaller volume and less power. The outline dimension is Φ26 mm × 251 mm, the mass is 435 g, the power is 1.3 W, and the theoretical penetration depth is 3 m. Penetration speed is 0.6 m/h and the practical penetration depth is 1.25 m. 10
During 2006 and 2008, in order to carry out the penetration exploration, the SRC PAS designed the KPET penetrator based on the MUPUS penetrator.11,12
In 2008, the NASA designed the MMUM penetrator for the Martian or the lunar dynamic penetration test based on the PLUTO penetrator. The outline dimension is Φ40 mm × 600 mm, the mass is 200 g, the power is 10 W, and the theoretical penetration depth is 2 m. 13
In 2013, the State Key Laboratory of Robotics and System of the Harbin Institute of Technology (HIT) designed a kind of penetrator whose outline dimension is φ35 mm × 325 mm, the mass is 420 g, the power is 10 W, and the penetration depth is 1 m. 14
Although there have been some achievements using the penetration devices, theoretical studies regarding the use of the penetration device are lacking. No successful cases detected the planet using penetration devices. The analysis of the dynamic penetration properties of the penetrator has a significant value to help design the penetrator structure. In this article, based on the results from a previous study on two arc-shaped front noses of penetrator,15–17 the axial resistance of the penetrator with concave arc-shaped nose is smaller than the convex one. Therefore, it is necessary to study the dynamic process of the concave arc-shaped nose of the penetrator. The dynamic penetration mechanical model was established according to soil mechanics. Then, the discrete-element simulation model, which was established using the application program interfaces (APIs), was used to simulate the penetration process. The penetration experiment demonstrated the accuracy of the theory. This study may provide technical support to the penetration programs of the deep space exploration.
Mechanical model of the penetrator in the penetration process
A concave arc-shaped front nose
Because the front nose first contacts the soil, the shape and the structural parameters of the front nose have a significant influence on the penetration feature of the penetrator. To maintain a depth of penetration above 1 m, a proper shape of the front nose should be determined.
The concave arc-shaped front nose is shown in Figure 1. The penetrator is considered as a rigid body in the penetration process, which remains in its original shape in the penetration process.

Structural model of the concave arc-shaped front nose.
In Figure 1, R is the radius of the penetrator shell, S is the radius of the curvature of the camber, l is the length of the concave arc-shaped front nose; β is the angle between any point of the camber and the vertical direction; and β0 is the angle between the peak of the camber and the vertical direction.
The caliber-radius-head is defined according to the following equation
Model of soil failure
In the penetration process of the penetrator, the region of soil failure is shown in Figure 2. According to Figure 2, during the penetration process, the soil near concave arc-shaped front nose of penetrator consists of the compacted region I, the transition region II, and the passive region III. The soil in region I moves down with the penetrator in the penetration process. The soil in region II initially moves to region III, and the soil in region III will fail when the pressure improves to a certain extent. There are two slip lines in region II, and the equation of the slip is expressed as
where

Region of soil failure.
The mechanic state of block BFHK is shown in Figure 3 when the soil in region III begins to fail.

Mechanic state of the block BEFK.
In Figure 3, q is the soil stress on plane BK; q = γh, where γ is the bulk density of the soil.
where D is the penetration depth of the penetrator; l is the length of the concave arc-shaped front nose; pP is the normal Rankine’s passive earth pressure on the surface of FK; C is the cohesion of soil; pn is the normal stress on surface EF; pntgφ is the friction force on surface EF; pr is the joint stress of pn and pntgφ; φ is the angle between pr and pn; pF is the Rankine’s active earth pressure on the surface of BE; δ is the angle between the slip line and the horizontal plane, δ = π/4 −φ/2; and r0 is the length of BE.
From the geometric relationships in Figure 4, we easily see that

Mechanic state of the block BOE.
Based on the Mohr–Coulomb failure criteria, the normal Rankine’s passive earth pressure on the surface FK is expressed by the following equation
The surface BE can be considered to be the retaining wall because of the fact that the soil in region I moves down with the penetrator. Taking point B as the centroid, we can get the following equation according to the static equilibrium condition of block BEFK
And we can easily get that
Substituting equation (6) into equation (5), we find that the Rankine’s active earth pressure on the surface BE is as follows
where
here, Nq and Nc are the bearing capacity coefficients of the soil.
The mechanical state of the block BOE in the compacted region is to be analyzed after we get the Rankine’s active earth pressure on the surface BE pF. Its mechanical state is shown in Figure 4.
In Figure 4, τF is the friction stress on the surface BE, τF = c + pFtgφ; P0 is the normal stress on the surface OE; τ0 is the friction stress on the surface OE, τ0 = c + p0tgφ; p is the normal stress on the surface OB; and f is the friction stress on the surface OB.
The force between the concave arc-shaped front nose of the penetrator and the block BOE is expressed by the following equation
According to the static equilibrium condition of block BOE, we find
Substituting equations (7) and (8) into equation (9), we obtain the following radial stress
According to the above analysis, the normal stress on the surface OB in the penetration is obtained. We can use its value to help analyze and calculate the axial resistance of the penetrator.
Axial resistance model of the penetrator
Assume that the penetrator will always travel along the vertical direction during penetration. The penetration process can be divided into three stages before the penetrator completely penetrates into the soil.
1. D < l: in this stage, resistance of the penetrator focuses on concave arc-shaped front nose, and the mechanical analysis is shown is Figure 5.

Mechanical analysis of the concave arc-shaped nose.
According to Figure 5, the infinitesimal axial resistance of the penetrator is expressed by the following equation
Then, the total axial resistance is derived after the integral operation of equation (11)
where
2. l < D < (L + l): L is the noumenon length of the penetrator. In this state, the concave arc-shaped front nose of the penetrator has completely penetrated into the soil. The total axial resistance can be obtained by substituting β = β0 into equation (15). Besides, in this state, the soil also contacts the lateral surface of the penetrator. The mechanical model of the lateral surface is shown in Figure 6.

Mechanical model of the lateral surface.
The equation of the total axial resistance on the concave arc-shaped front nose in this state is like equation (12), and it is as follows
where
Besides, the lateral surface is considered as the retaining wall, and then the passive earth pressure is expressed by the following equation according to Rankine’s earth pressure theory
where KP is the coefficient of the passive earth pressure, and
So, the horizontal positive pressure on the infinitesimal surface is
The friction force between the soil and the lateral surface is
So, the axial resistance force of the penetrator is expressed by the following equation
3. D ≥ (L + l): in this state, the axial mechanical analysis of the penetrator is shown in Figure 7.

Axial mechanical analysis of the penetrator.
In this state, the total axial force of the concave arc-shaped front nose F″ is
The friction force between the soil and the lateral surface is given by
So, the axial force that resists the motion of the penetrator can be obtained as follows
According to the above analysis, the axial resistance model of the penetrator in the penetration is obtained. And then, we can use its value to establish the dynamical model of the penetrator in the penetration process.
Dynamic model in the penetration process
After the impact hammer impacts the cone, they will move down together. In addition, the concave arc shape of the penetrator will rapidly penetrate into the soil. Then, the downward movement speed will descend with the increase in penetration depth; so, the study of penetration dynamics focuses on the second and third states of penetration.
When the impact hammer and cone move down together, the following equation is obtained according to Newton’s Second Law of Motion
where M is the total mass of the impact hammer and cone (kg) and v is the downward motion speed (mm/s).
The impact hammer and the cone move down at an identical initial speed V after they impact at any time, and every penetration process finishes when their speed is zero. In this article, we assume that the penetration depth DN−1 increases to DN after the Nth penetration.
The following equation is obtained by integrating equation (21)
In the second penetration state of the penetrator, l < D < l + L; we substitute equation (17) into equation (22); then, DN and DN−1 are expressed as follows
where
Assuming that the initial state in the second penetration state of the penetrator is that the concave arc-shaped front nose is correctly buried by the soil, the initial depth is
where
In the third penetration state of the penetrator, D > l + L. Substituting equation (20) into equation (23), DN and DN−1 can be expressed as follows
where
Assuming that the initial state in the second penetration state of the penetrator is that the concave arc-shaped front nose is correctly buried by the soil, the initial depth is
where
Theoretical penetration depth according to the penetration dynamics
According to the established theoretical equation, the program chart to calculate the penetration depth is shown in Figure 8.

Program chart to calculate the penetration depth.
In this article, the impact work in every impact state is set to be W = MV2/2 = 1 J, and the penetration depths of the penetrator with different structural parameters can be obtained by changing the structural parameter ψ. Figure 9 shows the change in penetration depth with the number of impact times.

Change in penetration depth with the impact times.
According to Figure 9, in the simulation process, the penetration depth increases with the increase in the number of impact times. The penetration depth of the cone is 160 mm when the number of impact times is 100.
Discrete-element simulation model
To demonstrate the accuracy of the provided theory, the discrete-element simulation model was established to simulate the interaction process between the soil and the penetrator.
Discrete-element simulation model of the lunar soil
Lunar soil is composed of numerous micro particles, and the shapes of the particles are complicated and fickle. In the simulation process, the irregular lunar soil particles are alternated by spherical particles after the stacking angle test. In this article, to decrease the simulation time based on the simulation accuracy, the discrete-element simulation model of the soil has three types of particles with different diameters: 0.5, 1, and 1.5 mm. The simulation model of the lunar soil is shown in Figure 10.

Discrete-element simulation model of the lunar soil.
The physical and mechanical properties of the lunar soil in the discrete-element simulation model are determined by parameter matching. The concrete method is as follows: first, the mechanical properties of the soil are determined by three-triaxial compression experiment, and then, the parameters of the discrete-element simulation model are changed until the results of the three-triaxial compression simulation in the EDEM software are similar to the result of the three-triaxial compression experiment. According to this parameter matching method, the final material parameters and contact parameters of the soil and the penetrator are shown in Tables 1 and 2, respectively.
Material parameters.
Contact parameters.
Discrete-element simulation model of the interaction process between the lunar soil and the penetrator
Because of the motion limitation of the geometry in the EDEM software, the traditional simulation model cannot simulate the dynamic penetration process of the penetrator. Thus, we wrote simulation programs and imbedded them into the EDEM software through the APIs of EDEM. The simulation program flow chart is shown in Figure 11.

Simulation program flow chart.
In the simulation, the impact work of every time is 1 J. The dynamic interaction process between the concave arc-shaped nose and the lunar soil in one time of the penetration as shown in Figure 12.

Dynamic penetration process in one time.
Figure 13 shows the soil failure condition when the number of impact times between the hammer element and the concave arc-shaped nose is 100.

Soil failure condition.
With the increase in number of impact times, the failure of the lunar soil becomes increasingly more serious according to Figure 13. In addition, the red and blue lunar soil particles mix together, particularly in the area of the concave arc-shaped nose. Figure 14 shows the change in penetration depth with the number of impact times when the structural parameter ψ is 2.

Change in penetration depth with the number of impact times.
According to Figure 14, in the simulation process, the penetration depth first rapidly increases with the increase in the number of impact times, but the increase in the penetration depth with the increase in the number of impact times subsequently becomes slow. The penetration depth of the cone is 102 mm when the number of impact times is 100.
Dynamical penetration experiment
Except for the penetration simulation, we also carried out the penetration experiment to prove the accuracy of our theory. The experiment and facilities are as follows.
Experiment table
To perform the penetration experiment, the experiment facilities were produced as shown in Figure 15. The penetration experiment table consists of the drive rod, guide rod, impact plate A, impact plate B, sleeve, and trestle table. The parameters of the experiment facilities are shown in Table 3.

Penetration experiment facilities: (a) lunar soil container, (b) hammer element, (c) housing element, and (d) impact test device.
Parameters of the penetration experiment.
Penetration experiment
When we performed the penetration experiment, the experiment parameters were set according to Table 4, and the penetration experiment is shown in Figure 16. The experiment result is shown in Figure 17.
Experiment parameters.

Penetration experiment and the concave arc-shaped nose.

Change in penetration depth with the impact times.
According to Figure 17, in the simulate regolith, the penetration depth is approximately 90 mm when the number of impact times is 100.
To demonstrate the accuracy of our impact dynamic theory, the theoretical, simulation, and experiment results are listed in Figure 18.

Different results.
According to Figure 18, the simulation result curve is notably close to the experiment result curve, which demonstrates the accuracy of the discrete-element simulation model. Because of some of the simplified treatments and many empirical parameters provided by soil mechanics, there are some errors between the theory result and the other results, but the entire variation trend of the theory result curve is close to the others. The established dynamic penetration model in this article can be used to analyze the dynamic penetration properties of the penetrator. The parameters provided by the soil mechanics will be calibrated to improve the accuracy trend of the theory result.
Conclusion and outlook
In this article, the dynamic penetration model of the penetrator was established. To demonstrate the accuracy of the theoretical model, the discrete-element simulation model was established to simulate the penetration process. Finally, the penetration experiment was also performed. The theoretical, simulation, and experiment results show a good agreement, and the study method proposed in this article might be a kind of method that can be used for planetary penetration exploration. In the future, the parameters given by soil mechanics will be calibrated to improve the accuracy trend of the theory result.
Footnotes
Academic Editor: ZW Zhong
Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This project was supported by the National Nature Science Foundation of China (nos 51575123 and 51575122), Self-Planned Task of State Key Laboratory of Robotics and System, Harbin Institute of Technology (HIT; no. SKLRS201616B), the Fundamental Research Funds for the Central Universities (HIT.NSRIF.2017028), and the International Science and Technology Cooperation Program of China (no. 2014DFR50250).
