#### Abstract

The various damages of asphalt pavement are closely related to the mesomechanical gradual behavior of asphalt materials, and it is very important to study the mesoscopic response under vibration loading in order to reveal the failure mechanism of asphalt pavement. The semisinusoidal vertical load is applied to the subgrade-surface discrete element model in this paper, and we use the model to analyze the evolution behavior of microcrack generation and expansion processes, stress distribution and stress transfer, and displacement field in various structural layers of asphalt pavement. The results show that the number of cracks increases rapidly on both sides of the vibration load, the rut is generated due to repeated load on the wheel, the asphalt mixture has bulging phenomenon on both sides of the rut and formed macroscopic cracks at the ridge, the microcracks extend mainly along the weak joints of the edges of the coarse aggregate and the asphalt cement, the number of microcracks increases slowly at the initial stage of the vibration load, the microcracks increase sharply until macroscopic cracks appear with the vibration load increases, the direction of compressive stress extends parallel to the microcrack, and the direction of tensile stress extends perpendicular to the microcracks inside the asphalt pavement. The results show that the discrete element method can not only obtain the stress and displacement of each structural layer, but also reveal the microcrack gradual behavior between particle flows.

#### 1. Introduction

Asphalt pavement is generally composed of an asphalt layer, a cement stabilized layer, and a road base layer; the asphalt layer is mainly composed of asphalt, coarse aggregate, fine aggregate, and mineral powder; the cement stabilizing layer is composed of coarse aggregate, fine aggregate, and cement; the roadbed is mainly composed of earth and stone [1]. Therefore, the overall pavement structure is a heterogeneous and noncontinuous body, the stress-strain is discontinuous, and the deformation is very complicated in the pavement structure under the load of the vehicle; if continuous medium mechanics is used to analyze the internal structural stress of the pavement, the actual deformation and stress state cannot be truly reflected inside the pavement structure [2, 3]. At present, many scholars have proposed applying the micromechanical theory to solve macroscopic mechanical problems, and micromechanical research has gradually become an area of scientific focus [4].

Macroscopic mechanics often adopts the traditional continuum mechanics theory, which assumes that the structural layers inside the subgrade are the same homogeneous body, and the geometric and physical characteristics, such as aggregate shape, spatial position, aggregate texture, and aggregate gradation, are ignored; however, the internal geometry and physical properties are of great importance to the overall performance of the pavement structure layer [4–7]. The existing discrete element theory holds that each structural layer is composed of particles of different sizes and shapes; the particles are rigid bodies, and the deformation is caused by the overlapping of the particles; the particles are in contact with each other, and the force can only be transmitted by the contact of the particles; the bonding relationship can be set between the particles, and the force can be determined by Newton’s second law; the action of particles is a dynamic process, and the application of external force will cause a certain range of particles to move and rotate, it takes a long time to reach equilibrium between particles [8–12]. Some scholars and experts have made research and exploration at home and abroad; for example, Kim et al. [13, 14] took a series of deformation pictures of asphalt mixture in different time periods through CT scanning technology and established the finite element model by using mesomechanics and fracture mechanics theory to analyze the mechanism of macrocrack generation of specimens. Mahmoud et al. [3] followed the principle of constant mass ratio and specific surface area and treated less than 4.75 mm particles inside the mixture as asphalt cement material; the distribution of aggregates was obtained by CT scanning in the asphalt mixture. Eckwright et al. [15] established the AC-16 uniaxial compression model of asphalt mixture by FISH language and analyzed the regularity of particle flow stress-strain, contact force, particle displacement, contact pressure, and tensile force. At present, it can be understood from the above that most of the research on the mesomechanical behavior of asphalt pavement is concentrated on the test of small-scale asphalt layer test pieces and analyzes the mesomechanical constitutive relation and mechanical behavior of small-scale asphalt mixture specimens under static load [16–18]; however, it is still rare to analyze the mesoscopic response and microcrack grading behavior for large-scale specimens (pavement overall structure layer) and involving vibration load in the existing literature.

In this paper, according to the actual pavement structure hierarchy, and using the combination of the random placement algorithm and porosity treatment principle, the two-dimensional subgrade-pavement discrete element model is established by using the FISH language in the PFC2D software and the mesoscopic parameters of each structural layer are obtained by comparing the standard uniaxial compression stress-strain test of materials with the uniaxial compression test of the established model. Semisinusoidal vertical load is applied to simulate vehicle loading and unloading; reveal the relationship between the microcrack propagation process, mesoscopic failure, and the macroscopic dynamic response of each structural layer; and analyze the mesoscopic response gradation behavior of asphalt pavement under vibration loading; these research results provide a reference for the mechanism of pavement damage.

#### 2. Materials and Methods

According to the literature research [19–23], a typical pavement structure model is applied and the model is divided into six layers, including upper layer, lower layer, upper base layer, lower base layer, bottom base layer, and subgrade. The roadbed mainly bears and disperses the vehicle load, the vehicle load is relatively small under the roadbed, and the depth of roadbed is 80 cm in the specification. Since the model uses large-scale test pieces, the depth of the roadbed is determined to be 100 cm in order to reduce the calculation amount of the model, the material of the pavement structure layer is composed of different particle size streams, and the particles are in contact with the particles between the layers. The pavement structure parameters are shown in Table 1.

A variety of common contact models are provided in the PFC2D software. The linear parallel bond model can reflect the internal mechanical relationship between the binding material and the coarse aggregate, and parallel bonding model calculation is relatively simple [24, 25]. The particle flow model and the mechanical relationship are shown in Figure 1.

**(a)**

**(b)**

Figures 1(a) and 1(b) illustrate contact particles; , , and illustrate the centers of particles A and B and the bonding point; is the representative of parallel bond model force; represents the parallel bond model moment. The internal forces and moments are all zero when the model is calculated to reach equilibrium; is decomposed into normal force, is decomposed into tangential force, is decomposed into normal moment, and is decomposed into tangential moment [26, 27].

We have corresponding force and bending moment increment when the particles are displaced:

In these equations, ; *I* represents moment of inertia; and represent the particle’s normal and tangential displacements; and represent the particle’s normal and tangent relative rotation angles; represents the normal contact stiffness; represents the tangential contact stiffness; represents the normal contact stiffness of particle A; represents the tangential contact stiffness of particle A; represents the normal contact stiffness of particle B; represents the tangential contact stiffness of particle B; represents the parallel key distance; represents the tensile stress between particles; represents the friction coefficient; represents the tangential contact stiffness; represents cohesion; and represents the internal friction angle [28].

Force and bending moment are updated according to the following formula where represents time:

So, we can obtain normal stress and tangential stress:

In the particle flow calculation process, the contact force and the unbalanced force are constantly balanced in the particle flow motion; microcracks are generated between particles when the tensile and shear resistance exceeds the bond strength in the model [29].

According to the theory of asphalt cement, the coarse aggregate (particle size greater than 2.36 mm) plays a skeleton role and bears most of the internal stress in the asphalt mixture, and the fine aggregate (particle size less than 2.36 mm) plays the role of filling skeleton space and binding between the coarse aggregates. In order to simplify the model and simplify the calculation data, we consider fine aggregate, water, and asphalt as binding materials, so the discrete element model only contains coarse aggregate and binding material particles [30–33].

Six spaces are generated and separated by walls using the PFC2D software, according to the actual structural material grading, the distribution area of the coarse aggregate is calculated and placed it in the wall step by step. The specific model is shown in Figure 2.

#### 3. Vibration Load Realization

The two-wheel set BZZ-100 is used as the standard load in the pavement design specification, the tire marks are actually elliptical; however, for the convenience of calculation, more circular tire marks are used in the model calculation, the tire load is regarded as the equivalent circular uniform load, and the internal pressure of the tire is used instead of the contact pressure at the bottom of the wheel. The specific values are shown in Table 2.

To simulate wheel load with CLUMP unit on the surface layer above the model, five PEBBLE units are used to represent one wheel, five PEBBLE units are tightly connected (21.3 cm in diameter), and the distance is 10.65 cm between the two CLUMP units (the upper layer of Figure 2). In order to accelerate the failure behavior of various structural layers, a semisinusoidal vertical load is applied to the CLUMP unit with a peak value of 100 kN and a cycle time of 0.06 s, and the model repeatedly applies periodic loads, eventually producing macroscopic cracks, the vibration load is performed by writing a subroutine in the FISH language within the PFC. The vibration load amplitude curve is shown in Figure 3.

#### 4. Particle Flow Mesoscopic Parameters

In the literature, most experts and scholars believe that the macroscopic mechanical parameters cannot be directly applied to the discrete element model [34–36]; however, the macroscopic mechanical characteristics are closely related to the microscopic mechanics of discrete element models. At present, the stress-strain (load-displacement) curve is generally used to illustrate the relationship between macroscopic and mesoscopic mechanics, by comparing the results of compressive strength test in the laboratory with those of the discrete element model. The discrete element model mesoscopic parameters can be obtained by iterative calculation [37, 38]. Figure 4 shows the comparison between the uniaxial compression test data of the SMA-13 asphalt mixture and the model data; the model geometry is 80 mm × 40 mm. Figure 5 shows the SMA-13 asphalt mixture model and the entire process of destruction. Figure 6 shows the failure process of the uniaxial compression specimen of asphalt mixture under indoor standard conditions; the red line is the main crack extension line.

**(a)**

**(b)**

**(c)**

**(a)**

**(b)**

**(c)**

**(d)**

It is necessary to know the mesoscopic parameters (, , , ) when we calculate the mesomechanical behavior of the discrete element model [39], it can be seen from Figure 4 that the uniaxial compression stress-strain curve is basically consistent with the experimental curve, and the stress peak error is less than 5.2%. Therefore, we can get the SEM-13 asphalt mixture discrete element mesoscopic parameters; the mesoscopic parameters of the remaining pavement materials are the same as those of the SMA-13 asphalt mixture. The microscopic parameters of the linear parallel bond model are shown in Table 3.

#### 5. Results and Description

##### 5.1. Microcrack Gradual Behavior

Applying a half-wave sinusoidal vertical load to the surface of the model, and gradually generating microcracks inside the pavement, microcrack density continues to increase with the load continues; finally, the microcracks penetrate each other to form macrocrack. Figures 7–10 show the whole process of microcrack generation in the model under vibration load.

It can be seen from Figure 7 that there is no obvious microcrack in the pavement structure from the upper part to the lower part when the vibration load reaches 1 s; this stage is represented by the compaction process of the particle flow in each structural layer, and the material properties exhibit elasticity, the density of the upper layer and the lower layer increases significantly, and the contact force of the particle flow increased significantly.

It can be seen from Figure 8 that microcracks appear in the upper layer and the lower layer when the vibration load reaches 2 s, microcracks are concentrated at both ends of the wheel load; at this time, the tensile stress is obviously increased at the bottom layer of the lower layer, and the microcracks appear as bond break between the particle streams, and the surface is slightly sunken on the upper layer.

It can be seen from Figure 9 that the number of microcracks increases significantly in the pavement when the vibration load reaches 4 s (blue indicates crack), a large number of microcracks are concentrated on both ends of the wheel load, and the microcrack develops vertically downward and toward the sides. There is a slight rut in the upper layer, microcracks appear in the lower layer, and the bond key is basically destroyed at the bottom tension zone. Microcrack concentration occurs at the bottom of the upper base layer, and the microcrack width is larger than the bottom of the lower layer.

It can be seen from Figure 10 that there are more and more microcracks in the pavement structure when the vibration load reaches 12 s, and the microcracks density is getting larger and larger. The depth of the microcrack gradually extends to the subgrade, and the overall width of the microcrack and the surface subsidence of the pavement continues to increase. The functional diseases have occurred in the overall structure of the pavement when the vibration load reaches the end, and the microcracks of each structural layer are merged and penetrated to form macrocracks in addition to the subgrade. The macrocracks are dominated by the vertical direction. The width of the microcracks is increased by more than 2 times at the bottom of the lower layer, and these microcracks are the initial cause of the pavement structure damage.

By studying the whole process of microcrack propagation, it can be known that most of the microcracks expand along the shear stress direction under the double-wheel load, the number of microcracks is the highest at both ends of the load, the shape of microcracks is symmetrically distributed overall, and the microcrack density in the middle strong stress region of the load is much larger than the microcrack density at both ends of the load.

Figure 11 shows partial enlargement of the calculation time 4 s. It can be seen from Figure 11 that the wheel load portion has obvious rutting (middle part of the figure) with the calculation time increases, and the rut depth is about 0.043 mm. There is an asphalt mixture bulging on the outside of the wheel, and the phenomenon of bulging is more obvious on both sides with the load increases. There are two obvious macrocracks on both sides of the rut, the microcrack on the left side maintains 40° curve extending from the horizontal direction among them and penetrates the upper base layer and gradually returns to the lower layer; the microcrack on the right side also maintains a 40° curve extending from the horizontal direction and passes through the surface of the lower layer; the model produces microcracks in different directions during the macrocrack propagation process.

Figures 12 and 13 are enlarged views of the internal macrocracks of the upper layer and the lower layer. It can be seen from the figure that macrocracks will be generated inside the upper layers and lower layers with the load increases, the initial macrocrack expands along the rough aggregate surface and the weakest edge of the asphalt cement, and the coarse aggregate particles are almost not broken and the macrocrack gradually penetrates the asphalt cement, eventually leading to macrocracks. The second crack will be generated at the crack extension stress concentration with the load continuing to increase, and the remaining cracks will develop like this; eventually, the microcracks are distributed throughout the entire area, and microcracks are merged in different directions.

From the trend of microcrack generation and expansion, it is known that the combination is the main weak part of asphalt mixture between coarse aggregate and asphalt cement, which is easy to produce early microcracks; microcracks eventually become macroscopic cracks under continuous vehicle loading. Functional conditions have already appeared in the pavement structure in this state, and we should protect it.

##### 5.2. Microcracks Number Gradual Behavior

Figure 14 shows the relationship between the number of microcracks and the load time. It can be seen from the figure that the number of microcracks is continuously increasing in the model, the growth rate of microcracks is relatively large in the 0–8 s stage, which shows a linear relationship, and the growth rate of microcracks is relatively small in the 8–12 s stage, which presents a nonlinear relationship.

It can be seen from the above figure that there is almost no microcrack inside the 0-1 s phase discrete element model, stress is continuously transmitted between the particles, and the stress continues to accumulate and self-adjust, eventually reaching dynamic equilibrium. A small number of microcracks appear inside the 1 s-2 s stage discrete element model, and the number of microcracks is very small, most of the macroscopic physical and mechanical performance indicators are not damaged in this state. The number of microcracks increases rapidly inside the 2 s–8 s stage, and the length of microcracks continued to increase. The number of microcracks is linear with the load time. The dynamic balance of the model is broken at this stage. Microcracks do not appear randomly, but spread rapidly downwards and outwards on both sides of the tire. Most of the microcracks cross each other and merge in this state, and the surface layer shows macroscopic cracks. The ability to resist external loads is much reduced in the model, and the rutting is gradually increasing. The number of microcracks is in a stationary phase inside 8 s–10 s stage, and the number of microcracks does not increase significantly. The number of microcracks has a certain increasing trend during the 10 s–12 s stage, but the rate of growth is small.

Due to external wheel vibration load, increasing microcracks have seriously changed the ability to withstand external loads; if the model is subjected to wheel loads for a long time, the particle flow cannot be completely restored to its original state, and the asphalt pavement is already in an initial destruction state.

##### 5.3. Structural Layer Stress Response

Since the discrete element model is a point contact between particles, the particle flow will continue to shift and rotate under the continuous vibration load of the wheel, stress is transmitted between the particles, and the change of stress reflects the macroscopic structural mechanical properties in the discrete element model. The stress is closely related to the microcracks generation, microcracks expansion, and microcracks penetration in particles.

Figure 15 shows the overall force chain diagram for the asphalt pavement model. The blue color represents the compressive stress between the particles in the figure, and the red color represents the tensile stress between the particles in the figure; the color thickness represents the magnitude of the stress. It can be seen that most of the pavement structure exhibits compressive stress from Figure 15, and the compressive stress is arranged in a layer pattern (blue layering) in the overall model, and the pressure layer style is “U” type. The direction of the compressive stress is vertically downward in the middle, and the direction of the compressive stress is outwardly distributed on both sides. The red color part (tensile stress) exists in the form of dots, and the tensile stress of individual red parts is relatively dense, the red part accounts for about 43% in the upper layer load area, and the interaction force between particles shows both tensile stress and compressive stress in the upper layer. The red color is relatively heavy and forms a “U” shape at the bottom of the lower layer, and the red part accounts for about 60–70%. The tensile stress is distributed outward and downward at 40 degrees from the horizontal direction; this is closely related to the number of microcracks generated at the bottom of the lower layer. The tensile stresses are distributed in a point in the upper base layer and the lower base layer, and the tensile stress accumulation is not obvious. There is tensile stress concentration on the contact surface between the bottom base layer and the subgrade, which is closely related to the density of the microcracks on the contact surface. The tensile stress is a uniformly dispersed arrangement inside the subgrade, and the tensile stress concentration is rare, and the particles mainly show compressive stress.

Figure 16 shows the force chain diagram at the local crack in the upper layer. It can be seen from Figure 16 that the force chain is completely broken on both sides of the crack, and there is no relationship between the particles. The direction of the compressive stress extends parallel to the crack at the upper part of the crack, the starting load section (upper left corner) shows the tensile stress, and the tensile stress extends to the lower right; the rest of the area exhibits compressive stress and compressive stress concentration at some locations (blue line concentration). The particles are mainly subjected to compressive stress under the crack, and the direction of the compressive stress extends parallel to the crack, the tensile stress is dispersed and produces tensile stress concentration phenomenon in some places of the model (red colour concentration point).

Figure 17 shows the force chain diagram at the local crack in the lower layer. It can be seen from Figure 17 that the force chain is almost completely broken on both sides of the crack, and only a small part of the crack is connected together (2 positions). The compressive stress and the tensile stress exist simultaneously between the particles at the upper part of the crack, and tensile stress concentration occurs at the upper right layer near the upper layer (1 positions); the direction of the compressive stress is collected to the upper right at the lower part of 1 position, and finally concentrated at the crack junction (2 positions), the tensile stress concentration (red position) appears in the lower right corner of the joint. The direction of the tensile stress is concentrated to the lower right corner under the crack, and the direction of the compressive stress is parallel to the direction of the crack, and the direction of compressive stress and direction tensile stress are perpendicular to each other.

##### 5.4. Structural Layer Displacement Response

The particle flow will continue to shift and rotate under the continuous vibration load of the wheel, and the appearance of the model shows that the particles have a tendency to move around.

Figure 18 shows the overall displacement nephogram of the discrete element model. It can be seen from Figure 18 that the overall displacement nephogram shape is “U” type under vibration load, the displacement gradually decreases with the depth of the structural layer increases, the upper layer produces the largest displacement among them, and the displacement of the lower layer is smaller than the upper layer and the subgrade has almost no displacement. The surface of the upper layer continues to produce approximately 5.5 mm ruts, and particle bulging occurs on both sides of the rut; these phenomena are closely related to the formation of macroscopic cracks and force accumulation on both sides of the rut.

Figure 19 shows the displacement of the horizontal displacement of the asphalt pavement. The red color indicates a positive direction displacement (positive in the right direction), and the blue color indicates a negative direction displacement (negative in the left direction). It can be seen from Figure 18 that the structural layers of the pavement have a certain horizontal displacement under vibration load; there is a clear separation between the layers of the pavement structure. The area of the layer contact surface is getting smaller and smaller with the depth of the structural layer increases, and the horizontal displacement is almost symmetrically arranged. The displacement is the largest (red color area) in the horizontal direction of the middle compression zone, and the displacement gradually decreases with the structural depth increases, The horizontal displacement of the particles moves oppositely on the right side of the rutting area (the lower part is red and the upper part is blue), and the horizontal displacement of the particles moves in the same direction on the left side of the rut area (blue color). The horizontal displacement values are almost the same at the upper base layer and the lower base layer, and the particles move in opposite directions (red color and blue color). The horizontal displacement is the largest at the contact layer position, which is consistent with the maximum number of cracks and the maximum tensile stress at the contact layer, and the subgrade has almost no horizontal displacement.

Figure 20 shows the vertical displacement nephogram of the discrete element model. The red color indicates that the displacement direction is upward. The red color indicates an upward displacement. The blue color indicates a downward displacement. It can be seen from Figure 20 that the particle flow has an upward movement tendency on both sides of the rutting area under vibration load. The red color area mainly occurs in the upper layer and the lower layer, and red color area symmetrically arranged. There is a clear dividing line between the red color and the green color, and the dividing line is a macrocrack. Most of the particle flow has a downward trend at the bottom of the rutting area, the displacement of the particle flow gradually decreases with road depth increases, and the particle flow also has a tendency to move to both sides under the action of wheel load.

#### 6. Conclusion

By comparing the indoor standard uniaxial compression test with the discrete element model, we can obtain the mesoscopic parameters of each structural layer and establish a discrete element model of the pavement structure and the mesoscopic response gradual behavior of the pavement under the given vibration load. The results show that the initial crack growth is slow in the particle flow under the action of vibration load and the crack growth is faster in the later stage, which eventually leads to rutting on the surface of the model. Most of the cracks extend along the shearing direction in the particle streams, and there are many cracks at the bottom of the wheel, and the cracks are almost symmetrically arranged. Most of the fine cracks extend along the edge of the coarse aggregate and the weakness of the asphalt cement in the particles under the action of vibration load; eventually, microcracks merge into macroscopic cracks in different directions. The force chain pattern is “U” type in the particles under the action of vibration load, which has both compressive and tensile stresses. The upper layer and the lower layer mainly exhibit tensile stress, and the other regions exhibit compressive stress; the particles are completely separated on both sides of the crack. The upper layer particles have the largest displacement under the action of vibration load, the displacement gradually decreases with road depth increases, and the subgrade has almost no displacement. The horizontal displacement layer has obvious boundary and has a tendency toward both sides, and the horizontal displacement is symmetrically arranged; the particle flow has a bulging phenomenon on both sides of the rut area. The discrete element model can not only obtain the stress and displacement under the action of vibration load, but also obtain the variation trend between the particles; the discrete element model can express both macroscopic mechanical behavior and mechanical mesoscopic behavior.

#### Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare no conflicts of interest.

#### Authors’ Contributions

Z.Y. and E.C. conceived the algorithm and designed experiments; Z.W. implemented the experiments and processed test data; and C.S. analyzed the results. All authors read and revised the manuscript.

#### Acknowledgments

This research was funded by the National Natural Science Foundation of China (grant 11172183) and Central Guided Local Science and Technology Development Project of China (18242219G).