#### Abstract

Peridynamics has a great advantage over modeling the damage process of rock-like materials, which is assumed to be in a continuum interaction with each other across a finite distance. In the paper, an approach to incorporate classical elastic damage model in the nonordinary, state-based peridynamics is introduced. This method can model the dynamic damage process and stress change of rock-like materials. Then two instances about three-point bend experiment are simulated in the rock-like materials. Finally the conclusions are drawn that numerical results are close to the experimental results. So the method has a great predictable value in the geotechnical engineering.

#### 1. Introduction

The rock-like material is a quasibrittle material; it is widely applied to the geotechnical engineering, for example, tunnels and lots of other underground buildings. Especially for the mountain area, the research on the damage process of rock-like materials is essential, so the prediction for the damage process of rock-like materials becomes a focus problem gradually. But it is difficult to model because the mechanical character of rock-like materials is very complex and the precise damage prediction of rock-like materials is elusive. In the geotechnical engineering, the rock-like mass includes lots of tiny cracks, even before the load is applied. The damage process of the rock-like materials is often caused by tension, comparison, the change of the temperature, and so on, and the stress and displacement fields are influenced by the propagation and coalescence of the crack. Even through the damage process of rock-like materials has been investigated for many years, the damage mechanism and prediction of rock-like material are still not understood. In fact, the damage process of rock-like materials includes the initiation, propagation, and coalescence of the crack, which will lead to a sudden collapse due to brittle damage. How to model the damage character of the rock-like materials has brought great challenge.

Over the past decades, many methods are put forward to model the damage process of rock-like materials. In the finite element-based method, singular crack-tip elements are frequently encountered [1]. Because of the crack-tip stress singularity, an external fracture criterion must be introduced to determine propagation and bifurcation of the cracks, and the nucleation question of the crack is still not solved [2]. In order to overcome the above difficulties, the extended finite element theory [3] is proposed to simulate the propagation of cracks. Although many crack questions are solved by virtue of the extended finite element theory, external and bifurcation criterion must still be introduced when displacement is discontinuous and when interaction and bifurcation of multiple cracks are involved. Besides, a series of difficulties are encountered for the problem of the three-dimensional cracks by the XFEM. In order to solve the problems of the three-dimensional cracks, such as interactions among cracks and branching phenomenon of multiple cracks, meshless methods are developed [4]. The propagation and coalescence process of cracks can be simulated by Smooth Particle Hydrodynamics (SPH); however, the tensile instability problems are still encountered in method of SPH [5]. In order to avoid the aforementioned lacks, peridynamic theory, which is a numerical method based on the nonlocal thoughts, is introduced to model propagation and bifurcation process of cracks.

The peridynamic theory is a nonlocal meshless method; it is put forward by Silling [6], at Sandia National Laboratory. It is assumed that particles in a continuum interact with each other across a finite distance, and it formulates problems in terms of intergral equations rather than partial differential equations [7]. Therefore, the peridynamic method can be applied to model the problems of continuous or discontinuous displacements [8].

After this theory is put forward, it has been widely applied to model the damage process of different materials. Firstly, the theory is applied to model the fracture process of composite material; for example, the fracture processes in laminated composites subjected to low-velocity impact and in woven composites subject to static indentation are predicted by Askari et al. [9] and Colavito et al. [10, 11] and the fatigue crack growth analysis is done in layered heterogeneous material system using peridynamic method by Jung and Soek [12]. In addition, the notched laminated composite under biaxial loads is considered by Xu et al. [13]. Then it is applied to model the fracture process of metal material; for example, Foster et al. [14] use peridynamic viscoplastic theory to model the collapse process of metal and Wu et al. [15, 16] analyze the ductile fracture of metal materials using nonordinary state-based peridynamics. Meanwhile, Sun and Sundararaghavan [17] model the crystal material using the peridynamic plasticity theory. Immediately, it is extended to model the damage of concrete; for example, Gerstle et al. [18] model the damage process of concrete by introducing “micropolar peridynamic model.” The damage results of these materials are rather good by using the peridynamics, but for the geotechnical engineering, rock-like materials are an important material, the peridynamic theory is seldom used to model the damage process of rock-like materials, and only Ha et al. [19] use peridynamics to model the fracturing patterns of rock-like materials in compression, but the stress field description is not considered. In the paper, the application of nonordinary, state-based peridynamics in the damage process of rock-like material will be investigated. Not only the damage process of rock-like materials is described by using the method, but also the change of stress field in the damage process is depicted, so it provides a new idea to model the damage process of rock-like materials by using the peridynamics.

The paper is organized as follows. In Section 2, the state-based peridynamic theory is introduced at first, and the implementation on this method for the specific damage model is discussed. In Section 3, the state-based peridynamic numerical discretization is described. In Section 4, we represent numerical results consisting of the numerical simulation about the three-point bend test of the beam and the numerical simulation about the three-point bend test of Brazilian disk. In Section 5, conclusions are drawn.

#### 2. Model Description

For completeness, nonordinary state-based peridynamics is reviewed briefly; it includes a summary of basic peridynamic equation, the idea of constitutive relation, and its property.

##### 2.1. Basic Theory

Peridynamic theory is put forward by Silling [6] at Sandia laboratory in the United States; its basic equation of motion is shown as follows:where is the density of material point , is a neighborhood (Figure 1), and is a prescribed body force density field of material point at the instance , which represents the external force per unit reference volume square. and are the force density vector of material points and , respectively.

In the neighborhood of point , the relative position of any point, and , is , the deformation state of the bond is defined as , and they are expressed, respectively, as follows:

The nonlocal deformation gradient of material point is shown as the following expression:where is the influence function of the bond. It is the function of the relative origin position between points and , and is a nonlocal shape tensor defined by

Alternatively, this integral can be defined without recourse to peridynamic states as follows:

The Green-Lagrange strain sensor can be expressed as

The elastic strain energy is written aswhere is the fourth-order elastic tensor and is assumed as anisotropic with cubic symmetry.

From (7), the second Piola-Kirchhoff stress is written as

And Cauchy stress is found:

The first Piola-Kirchhoff stress can be written from the Cauchy stress as

Equation (9) is substituted into (10); the expression of the first Piola-Kirchhoff stress is rewritten as

Peridynamic force vector state can be expressed by the traditional stress state as

Substituting (12) into (13), (14) can be obtained as follows:

From equations (1) and (14), the basic equation in the nonordinary, state-based peridynamics can be written as

##### 2.2. Damage Correspondence

To describe the damage correspondence within the peridynamic framework, the influence function is introduced to represent the damage; it can be expressed as follows:where is the horizontal radius of material point , is the relative displacement of any point and point , is a scalar function, and its expression is given:where is the stretch of the bond; it is defined as is critical stretch; in two dimensions it is defined as [20]where is volume modulus, is the shear modulus, is the critical energy dissipating ratio, and it is related to the fracture toughness .

To model the problem of damage, the concept of local damage value is introduced; it is defined aswhere is local damage value; its limitation is , where 0 represents original materials and 1 represents complete disconnection of a point from all of the points with which it initially interact.

#### 3. The Discretization

The region is discretized into nodes, each with a known volume in the reference configuration; taken together, the nodes form a grid. So the deformation gradient is discretized into the following form:where a nonlocal shape sensor is discretized aswhere is the volume of nodes ; is the magnitude of node in the range of horizontal radius of the node .

Likewise, the discretization form of basic equation in the nonordinary, state-based peridynamics is shown as follows:where and are the first Piola-Kirchhoff stress of nodes and , respectively, and are the shape sensor of the nodes and .

The acceleration can be shown aswhere the superscript represents the time step and the subscript represents the number of nodes, so represents the acceleration of the node at time step .

#### 4. The Numerical Example

To evaluate the functionality of the proposed state-based peridynamic formulation, the discretized equations were implemented in a Fortran computer code for two-dimensional simulation. The following two examples are investigated, respectively.

##### 4.1. The Three-Point Bend Test of the Beam

Three-point bend test [21] has been widely used for characterizing the dynamic response of materials as well as for validating numerical methods in the case of the state-based peridynamic formulation with damage. This test consists of a rectangular concrete beam with a preexisting vertical crack; its geometric configuration is shown in Figure 2. Its length is 0.32 m, the height is 0.07 m, the length of the crack in the beam is 0.0233 m, and the mechanic parameters of the concrete beam are listed as follows: the elastic modulus GPa, Poisson’s ratio , and center point of coordinate axis is positioned in the bottom center of concrete beam. -axis is vertical to the pressure, tension is positive, and compression is negative. The concentrated load is positioned at the upper center of the beam. the magnitude of load is 0.075 mm/min, the concrete beam model is discretized into particles, the distance between the adjacent two particles is m, time step s, the critical stretch , the density kg/m^{3}, and is adopted as . When it arrives at 25000 time steps, the concrete beam is fractured fully. The computed damage contour of concrete beam is shown in Figure 3, the maximum principal stress contour in the damage process of concrete beam is shown in Figure 4, the damage result in specific test is given in Figure 5, and the comparison curve between test and numerical solution is shown in Figure 6.

**(a)**s

**(b)**s

**(c)**s

**(d)**s

**(a)**s

**(b)**s

**(c)**s

**(d)**sIt can be seen from Figure 3 that, when time arrives at s, the crack begins to initiate along the upper tip; when time arrives at s, the crack begins to propagate along the position where the vertical load is applied; when time arrives at s, the crack continues to propagate and elongate along the position where the vertical load is applied, and two support positions in the bottom of concrete beam begin to be damaged; when time arrives at s, the crack propagates to the position where the concentrated load is applied, and the whole concrete beam broke down fully. It can be seen from Figure 4 that, when time arrives at s, the stress concentration begins to occur in the tip of crack, the stress concentration phenomenon is found at the position where the concentrated load is applied, and the magnitude of maximum principal stress in the tip of the crack reaches Pa; when time arrives at s, the stress concentration in the tip of concrete beam continues to propagate to the position where the concentrated load is applied, and the magnitude of maximum principal stress reaches Pa; when time arrives at s, the magnitude of the stress in the tip of concrete beam reaches the maximum, whole concrete beam is destroyed fully, the magnitude of the maximum principal stress reaches Pa, the numerical result in Figure 3 is close to the experimental observation in Figure 5, and from Figure 6, it can be found that the numerical process in the curve of load versus displacement is consistent with the trend of experimental results. These conclusions demonstrate that nonordinary, state-based peridynamics can provide good prediction for the damage process of rock-like materials. In the example, the difference of the damage result between the numerical prediction and experimental observation exists, and the constitutive relation of material has important influence on the difference, so the improvement for the constitutive relation of rock-like materials is necessary in my future research. It can not only improve the accuracy of numerical prediction, but also model the large deformation problem by altering the constitutive relation of materials, such as impact damage problem and blast. A critical extension is furnished to solve the problems involving severe deformation and damage by using the method.

##### 4.2. Failure Simulation of a Brazilian Disk under Compressive Loading

This test consists of a circular rock sample with a preexisting vertical crack, its geometric configuration is shown in Figure 7, its diameter is 100 mm, there is a crack in the middle of sample, and its inclination angle is . The length of crack is 2 = 30 mm, the mechanic parameters in this rock sample are listed as follows: the elastic modulus GPa, Poisson’s ratio , and center point of coordinate axis is positioned in the center of circular disk. -axis is vertical to the pressure, tension is positive, and compression is negative. The compressive load is located at the upper and lower tips of the circular disk and their magnitude is 0.05 m/s; rock sample is discretized into particles, the distance between the adjacent two particles is m, time step s, the critical stretch , the density kg/m^{3}, and is adopted as . When it arrives at 700 time steps, the concrete beam is fractured fully. The numerical damage contour of rock sample is shown in Figure 8, the maximum horizontal principal stress contour is shown in Figure 9, and the damage result in specific test is given in Figure 10.

**(a)**s

**(b)**s

**(c)**s

**(d)**s

**(a)**s

**(b)**s

**(c)**s

**(d)**sIt can be seen from Figure 8 that, when time arrives at s, the upper and lower tips of the crack begin to initiate; when time arrives at s, the crack along its tips continues to propagate; when time arrives at s, the crack begins to penetrate along the direction where the concentrated load is applied. The whole rock sample fractures fully. Likewise, it can be seen from Figure 9 that, when time arrives at s, the stress concentration in the tips of crack has taken place before the tips of crack begin to initiate, and the magnitude of the maximum principal stress in the tips of crack reaches Pa; when time arrives at s, the crack begins to initiate, and the magnitude of the maximum principal stress reaches Pa; when time arrives at s, the horizontal stress in the tips of crack continues to increase, and its magnitude reaches Pa; when time arrives at s, the magnitude of the maximum horizontal stress reaches Pa; rock sample is fractured fully. The conclusions can be drawn from the comparison between Figures 8 and 10, the numerical results are similar to the experimental ones, and good effect can be obtained for the numerical simulation of the damage process in the rock material by using nonordinary, state-based peridynamic method, so it provides great instructional significance for the prediction of damage process in the future.

#### 5. Conclusions

An extension of the state-based peridynamic constitutive correspondence framework to incorporate into elastic damage has been proposed. Then three-point bending simulation is performed with damage process. In two examples, the numerical results are compared with the experimental observation. In conclusion, the suggested method provides a new thought for the damage process prediction of the rock-like materials.

#### Nomenclature

: | The density of material point |

: | The applied body force density of material point at the time |

: | The displacement of material point |

: | The relative origin position between material points and |

: | The force density vector of material point |

: | The force density vector of material point |

: | The relative position of points and |

: | The deformation state of the bond between material points and |

: | The deformation gradient of material point |

: | The influence function of the bond |

: | The volume of material point |

: | The elastic strain energy |

: | The second Piola-Kirchhoff stress |

: | The horizontal radius of material point |

: | A scalar function of the bond between material points and at time |

: | The stretch of the bond |

: | The volume modulus |

: | The critical energy dissipating ratio |

: | The local damage value |

: | The neighborhood of material points |

: | Poisson’s ratio |

: | The Green-Lagrange strain sensor |

: | The elastic tensor |

: | The first Piola-Kirchhoff stress |

: | The critical stretch of the bond |

: | The shear modulus |

: | The fracture toughness |

: | The acceleration of the node at time step . |

#### Competing Interests

The authors declare that they have no competing interests.

#### Acknowledgments

This work is supported by the first batch of Talent Introduction Projection in 2016 in Sichuan University of Science & Engineering (no. 2016RCL19), the projection in the Sichuan Provincial Department of Science (no. 2016GFW0137), and financial contributions from the Sichuan Province Science and Technology Plan Project (2013JY0119).