Numerical study on failure characteristics for polymeric plates under oblique impacting with an improved peridynamic model

Abstract

To investigate the ricochet phenomenon and failure characteristics for PMMA (Polymethyl-methacrylate) plates subjected to oblique penetration, an improved peridynamic contact model is proposed in this study. The proposed model incorporates the dynamic friction coefficient and shank friction force, complemented by a hybrid contact search strategy that combines grid search and KD-tree methods to enhance computational efficiency. To validate the model and approach, numerical simulations for the oblique impact tests on PMMA are conducted, where the peridynamic results have a satisfactory agreement with the corresponding experimental data. The ricochet phenomenon is accurately reproduced in the peridynamic simulations, underscoring that the proposed model and approach have the capability to effectively capture the intricate failure characteristics in PMMA under oblique penetration. Additionally, a systematic analysis is performed to investigate the influence of critical parameters, including the incidence angle, initial velocity, and attack angle, on the overall oblique penetration performance. These findings highlight the robustness and practical applicability of the proposed model and approach in predicting and analyzing intricate failure responses in PMMA under dynamic loading conditions.

Keywords

peridynamics contact model oblique penetration ricochet phenomenon impact failure

Introduction

Polymethyl-methacrylate (PMMA) is extensively utilized in protection applications due to its commendable resistance to penetration and perforation (Esfahlani, 2021; Faye et al., 2016). Consequently, considerable attention has been devoted to understanding the failure mechanisms of PMMA under impact loading conditions. In practical scenarios, achieving a perfectly perpendicular impact during projectile-target interactions is often challenging, resulting in projectiles frequently penetrating targets at oblique angles (Rosenberg and Dekel, 2012; Stronge, 2018). Thus, investigating the oblique penetration characteristics of PMMA is of critical importance.

Traditional solution methods for analyzing oblique penetration typically encompass empirical and analytical approaches (Forquin et al., 2012; Kazarinov et al., 2020; Wu et al., 2004). The empirical method integrates theoretical analysis with experimental data to formulate algebraic equations that reflect the principles of momentum and energy conservation. By solving these algebraic equations, key penetration parameters such as penetration depth, residual velocity, and ballistic limits can be derived. Nevertheless, these methods often provide a limited and oversimplified representation of the constitutive relations and penetration mechanisms, failing to fully capture the complexities inherent in the actual process. As a result, comprehensively characterizing the entire physical process of oblique penetration remains a significant challenge.

In recent years, numerical simulation has emerged as an indispensable tool for the analysis of oblique penetration. Through the numerical analysis, the oblique penetration process can be effectively visualized, and various detailed information can be captured, which is always difficult to observe in experimental tests. For instance, Liu et al. (2009) conducted three-dimensional numerical simulations of oblique penetration into concrete, and the final penetration depths under different oblique angles were analyzed. Børvik et al. (2011) proposed a modified Johnson-Cook model and Cockcroft–Latham failure criterion. Numerical simulations on the oblique impact of aluminum plates were then conducted, and the numerical results agreed well with experimental data. Additionally, numerical analyses on the dynamic behavior of functionally graded plates under oblique impact loadings were performed in Hakan et al. (2023), in terms of damage, deformation, and failure mechanisms. Despite these advances, the numerical methods currently employed in studies on oblique penetration analysis primarily adhere to the traditional continuum mechanics theory, where the continuity hypothesis inherently contradicts the presence of discontinuities. This discrepancy poses challenges in effectively characterizing the failure mechanism associated with oblique penetration, thereby affecting both numerical efficiency and accuracy.

Peridynamics, as a non-local theory, has demonstrated remarkable potential in addressing discontinuous problems (Silling, 2000; Silling and Lehoucq, 2008). Unlike the continuum mechanics theory, peridynamics adopts spatial integral equations rather than partial spatial derivatives to describe the motion equations. Thus, the above contradiction between the continuity hypothesis and intrinsic discontinuity is totally addressed in peridynamics. Moreover, the interaction between material points in peridynamic theory is conceptualized as a bond. These bonds represent the forces exchanged between pairs of material points within a defined horizon, capturing the nonlocal nature of the interactions. Damage and failure are elegantly described through the progressive breakage of these bonds (Silling and Askari, 2005). This bond-breaking mechanism provides a robust and intuitive framework for modeling the spontaneous initiation and propagation of cracks, even in complex scenarios involving branching, merging, or other discontinuous behaviors.

The initial formulation of peridynamics is known as the bond-based peridynamics. It owns the advantage of being simple and easy to use, while still facing the restriction of a fixed Poisson's ratio. To overcome this limitation, state-based peridynamics is developed as a more general framework (Silling et al., 2007), further divided into ordinary and non-ordinary types. Nowadays, promising applications and achievements of peridynamics have been gained in various areas, including quasi-static or dynamic fracture of diverse materials (He et al., 2021; Li et al., 2021; Ma et al., 2023; Zhang and Dai, 2025), impact and explosion failure (Wu et al., 2022, 2023; Zhu and Zhao, 2021), multi-scale modeling (Bie et al., 2024c; Dai et al., 2021; Gu et al., 2017; Guo et al., 2020), multi-field coupling study (Amani et al., 2016; Bie et al., 2024a, 2024b; Hu et al., 2018; Zhang and Dai, 2024), and the coupling between peridynamics and other numerical methods (Chen and Yu, 2022; Dai et al., 2023; Han et al., 2016; Zhong et al., 2024).

In the context of impact failure, peridynamics has yielded substantial achievements over the past few decades. At relatively low impact velocities, the predominant failure mode is impact-induced fracture. Research in this domain has focused on brittle (Wu and Huang, 2022), ductile (Wang et al., 2019), and composite materials (Wu et al., 2020). Peridynamic simulations can satisfactorily capture the crack initiation and propagation, including the crack branching and multiple crack growth. Conversely, at higher impact velocities, penetration failure becomes the primary concern, with analyses performed for concrete (Liu et al., 2025), glass (Bobaru et al., 2012), steel (Pathrikar et al., 2019), PMMA (Wu et al., 2021), and soft materials (Huang et al., 2023), respectively. Peridynamic simulations provide detailed insights into the penetration process, enabling the characterization of intricate failure mechanisms.

In previous studies, research emphases mainly focused on the fracture and failure analysis related to perpendicular impact issues, with limited attention given to oblique penetration. Different from perpendicular penetration, the oblique penetration process is inherently more complex, as projectiles may deflect due to asymmetric forces, ultimately penetrating the target at a negative angle. This ricochet phenomenon has been observed and investigated in several studies, yet it remains an area that requires further investigation. In this work, we propose an improved peridynamic contact model that incorporates the dynamic friction coefficient and shank friction force, along with a hybrid contact search strategy integrating grid search and KD-tree methods. The failure characteristics of PMMA under oblique penetration are comprehensively analyzed, with particular attention given to the ricochet phenomenon. Additionally, the effects of incidence angle, initial velocity, and attack angle on the penetration process are systematically investigated, focusing on failure patterns, ballistic trajectories, and projectile deflection angles. The proposed model and approach effectively capture the intricate failure mechanisms of PMMA under oblique penetration, offering valuable insights and advancing the understanding of its dynamic response.

Methodology

Basic brief to the ordinary state-based peridynamics

The motion equation in state-based peridynamics is formulated to describe the dynamic behavior of a material point by considering the nonlocal interactions with its surrounding points within a defined horizon. It is expressed as:

ρ \ddot{u} (x, t) = \int_{H_{x}} {\underline{T} [x, t] ⟨ x^{'} - x ⟩ - \underline{T} [x^{'}, t] ⟨ x - x^{'} ⟩} d V_{x^{'}} + b (x, t),

(1)

where ρ represents the material density, b denotes the external force density, and

u

is the displacement vector of a material point, respectively. The term

\underline{T}

refers to the force vector state, which quantifies the interaction between material points within the defined domain. The integral in the equation is evaluated over the volume element

d V_{x^{'}}

, corresponding to the surrounding material points within the horizon. Additionally, Figure 1 illustrates the interaction between material points. The domain of interaction for any material point is defined as a spherical region, denoted as

H_{x}

H_{x} = {0 < | x - x^{'} | < δ}

, and

δ

is referred to as a peridynamic horizon.

In ordinary state-based peridynamics (Silling et al., 2007), the force vector state $\underline{T} [x, t]$ is aligned with the deformation vector state,

\underline{T} = \underline{t} \underline{M} .

(2)

where

\underline{t}

is the force scalar state, and

\underline{M}

is a unit direction vector along the deformation vector state,

\underline{M} ⟨ ξ ⟩ = {\begin{matrix} \frac{ξ + η}{| ξ + η |} & | ξ + η | \neq 0 \\ 0 & or else \end{matrix}

(3)

The strain energy density of ordinary state-based peridynamics is given by:

W (θ, {\underline{e}}^{d}) = \frac{k θ^{2}}{2} + \frac{α}{2} (\underline{ω} {\underline{e}}^{d}) ∙ {\underline{e}}^{d}

(4)

\underline{e} = \underline{y} - \underline{x}, \underline{y} = | \underline{Y} |, \underline{x} = | \underline{ξ} |

(5)

where

\underline{e}

is the extension scalar state, and

{\underline{e}}^{d}

represents the deviatoric part of

\underline{e}

\underline{ω}

is the influence function.

$θ$ and m denote the scalar value of dilatation and weighted volume, respectively .

θ = \frac{3 (\underline{ω} \underline{x}) ∙ \underline{e}}{m} = \frac{3}{m} \int_{H_{x}} (\underline{ω} \underline{x}) \cdot \underline{e} d V_{x^{'}},

(6)

m = (\underline{ω} \underline{x}) ∙ \underline{x} = \int_{H_{x}} (\underline{ω} \underline{x}) \cdot \underline{x} d V_{x^{'}} .

(7)

Using the chain rule, $\underline{t}$ in equation (2) can be written as:

\begin{aligned} \underline{t} = t^{i} + t^{d} \\ = \frac{\partial W}{\partial e^{i}} + \frac{\partial W}{\partial e^{d}} \\ = \frac{3 k θ}{m} \underline{ω} \underline{x} + α \underline{ω} {\underline{e}}^{d} \end{aligned}

(8)

Visco-elastic model in peridynamics

The dynamic mechanical behavior of PMMA is represented by the dynamic viscoelasticity model (Wu et al., 2021), as shown in Figure 2, which can describe the viscoelastic behaviors under low and high stretch rates, respectively.

The total deviatoric extension state can be decomposed and represented by:

e^{d} = e^{d e} + e^{d b_{1}} + e^{d b_{2}}

(9)

Given equation (15), equation (4) is rewritten as:

\begin{aligned} W (θ, e^{d}, e^{d b}) & = \frac{k θ^{2}}{2} + \frac{α_{1}}{2} e^{d} ∙ ω e^{d} + \frac{α_{2}}{2} (e^{d} - e^{d b_{1}}) ∙ ω (e^{d} - e^{d b_{1}}) \\ + \frac{α_{3}}{2} (e^{d} - e^{d b_{2}}) ∙ ω (e^{d} - e^{d b_{2}}) . \end{aligned}

(10)

Finally, we can derive the force scalar state,

\begin{aligned} t & = t^{i} + t^{d} \\ = - \frac{3 p}{m} ω \underline{x} + α_{1} e^{d} + α_{2} (e^{d} - e^{d b_{1}}) + α_{3} (e^{d} - e^{d b_{2}}) \\ = - \frac{3 p}{m} ω \underline{x} + α_{1} e^{d} + α_{2} \int_{0}^{t} {\dot{e}}^{d} \exp [- \frac{t - τ}{θ_{b 1}}] d τ + α_{3} \int_{0}^{t} {\dot{e}}^{d} \exp [- \frac{t - τ}{θ_{b 2}}] d τ, \end{aligned}

(11)

Where $α_{1}$ , $α_{2}$ , and $α_{3}$ are elastic parameters, following $α_{1} + α_{2} + α_{3} = \frac{15 μ}{m}$ . To prevent unbounded creep of the model, it is assumed that $α_{1} > 0$ , $α_{2} > 0$ , and $α_{3} > 0$ . $θ_{b 1}$ and $θ_{b 2}$ are relevant material parameters associated with material response under low and high extension rates in the viscoelastic model (Wu et al., 2021).

Numerical implementation

Numerical discretization

For numerical analysis, the governing equations of motion equation (1) after the discretization are described as:

ρ {\ddot{u}}_{i}^{n} = \sum_{p} {\underline{T} [x_{i}^{n}, t^{n}] ⟨ x_{p}^{n} - x_{i}^{n} ⟩ - \underline{T} [x_{p}^{n}, t^{n}] ⟨ x_{i}^{n} - x_{p}^{n} ⟩} V_{p} + b (x_{i}^{n}, t^{n}),

(12)

where n is the number of time steps,

V_{p} = | Δ x |^{3}

stands for the involved volume of the material point

x_{p}

The displacement of point $x_{i}$ at time $t^{n + 1}$ can be obtained by solving the following equations based on an explicit Verlet-Velocity algorithm (Zhang and Dai, 2024).

\begin{aligned} {\dot{u}}_{i}^{n + 0.5} = {\dot{u}}_{i}^{n} + \frac{Δ t}{2} {\ddot{u}}_{i}^{n} \\ {\dot{u}}_{i}^{n + 1} = {\dot{u}}_{i}^{n + 0.5} + \frac{Δ t}{2} {\ddot{u}}_{i}^{n + 1} \\ u_{i}^{n + 1} = u_{i}^{n} + Δ t {\dot{u}}_{i}^{n + 0.5}, \end{aligned}

(13)

Where u, $\dot{u}$ , and $\ddot{u}$ denote displacement, velocity, and acceleration, respectively.

Damage and failure criterion

In peridynamics, one of its remarkable advantages is the description of the natural initiation and propagation of cracks without extra complicated failure criteria. The most widely used criterion is the bond stretch criterion. Once any bond stretch reaches a given critical value $s_{0}$ , the bond breaks and never recovers. The bond stretch is formulated as:

s (ξ, t) = \frac{‖ ξ + η ‖ - ‖ ξ ‖}{‖ ξ ‖} .

(14)

The current state of each bond is represented by a Boolean function, that $μ = 1$ means the bond is active and intact, otherwise, the bond is broken as $μ = 0$ ,

μ (x, ξ, t) = {\begin{matrix} 1 & s < s_{0} \\ 0 & e l s e \end{matrix}

(15)

The local damage of each point can be finally obtained by the ratio:

D (x, t) = 1 - \frac{\int_{H_{x}} μ (x, ξ, t) d V_{x^{'}}}{\int_{H_{x}} d V_{x^{'}}} .

(16)

Contact model

In the peridynamic original configuration, material points in separated bodies are much greater than the grid spacing, meaning there is no interaction among them. When external loadings are so large as to make such separated bodies come in contact, or a moving body finally impacts on another one at a certain initial velocity, it necessitates the contact model to maintain the reasonable interaction between material points rather than suffering from the issue of interpenetration.

Aiming at the contact problems in peridynamics, considerable efforts have focused on this field, and various contact models in peridynamics have been proposed and utilized in dynamic contact, impact and penetration, sticking and frictional contact, and fracture analysis associated with contact (Guan et al., 2023; Kamensky et al., 2019; Lu et al., 2021; Tian et al., 2024; Wang et al., 2024; Zhang et al., 2022). In this study, the oblique penetration is modeled and analyzed. It is a fairly complex process, and the contact domain needs to be updated in each calculation step. Therefore, this work is based on the most commonly used short-range force contact model (Macek and Silling, 2007), as shown in Figure 3. Meanwhile, the dynamic friction coefficient and shank friction force are considered in the contact model to have a better description of the failure characteristic of oblique penetration.

$f_{c o n}$ denotes the contact force here and decomposes the normal contact force $f_{n}$ and friction force $f_{t}$ . The normal contact force $f_{n}$ can be presented as :

f_{n} (y_{i}, y_{j}) = min {0, \frac{c_{s h}}{δ} (‖ y_{i} - y_{j} ‖ - d_{s h})} e_{n}

(17)

where

y_{i}

and

y_{j}

are the position of material points

x_{i}

and

x_{j}

, respectively.

c_{s h}

denotes the contact stiffness related to the bulk moduli of the bodies, and

δ

is the peridynamic horizon.

d_{s h}

is the contact radius to determine whether one point touches another point or not, which is given by:

d_{s h} = min {0.9 ‖ x^{'} - x ‖, 1.35 | Δ x |} .

(18)

And $e_{n}$ is the unit vector pointing at $x_{i}$ from $x_{j}$ ,

e_{n} = \frac{y_{i} - y_{j}}{‖ y_{i} - y_{j} ‖} .

(19)

Moreover, it is notable that the short-range force is validated if and only if the distance between material points in separated bodies is greater than $d_{s h}$ , it equals 0.

According to Coulomb's law, the friction force is formulated as:

f_{t} = - μ | f_{n} | e_{t},

(20)

where

μ

is the friction coefficient, and

e_{t}

is the direction of the friction force, which is given as:

e_{t} = [I - e_{n} \otimes e_{n}] \frac{{v^{'}}_{i} - v_{j}}{‖ {v^{'}}_{i} - v_{j} ‖},

(21)

where

I

is an identity tensor, and

\otimes

is the tensor product.

\frac{{v^{'}}_{i} - v_{j}}{‖ {v^{'}}_{i} - v_{j} ‖}

represents the unit vector along the relative velocity of material points

x_{i}

and

x_{j}

Through numerous experimental observations, it is known that there is a nonlinear relationship between the friction coefficient and projectile velocity (Ben-Dor et al., 2007; Montgomery, 1976). The dynamic friction coefficient is not constant during the penetration process, especially in the oblique impact, where the friction force has a crucial effect on the ballistic trajectory. In this work, a velocity-dependent friction law (Chen, 1989) is adopted, and the dynamic friction coefficient is formulated as:

μ = μ_{\infty} + (μ_{0} - μ_{\infty}) e^{- γ v},

(22)

where

μ_{0}

and

μ_{\infty}

denote the friction coefficients under the low and high velocities, respectively.

γ

is the decay coefficient. Note that as the velocity v approaches 0,

μ

is nearly equal to

μ_{0}

In addition, many researchers ignore the shank friction in the process of projectile penetration, however, there are experimental tests that indicate that the shank friction should be considered (Chai et al., 2019; Yu et al., 2012), especially in the high-velocity impact cases. From the cavity expansion theory (Wang et al., 2007), the shank friction $f_{s}$ can be expressed as:

f_{s} = μ^{'} p l π d,

(23)

where

μ^{'}

is the friction coefficient between the projectile side wall and target, which is chosen to be 0.01. l represents the length of the projectile and d denotes the diameter of the projectile. p is the shank stress on the projectile, and it can be obtained by:

p = f_{c} a_{s},

(24)

where

f_{c}

is the compressive strength of the target, and

a_{s}

is a fitted cavity expansion coefficient.

In peridynamics, the contact force and shank friction force are applied as the force density and not imposed directly on the bonds. The peridynamic governing equation (1) can be reformulated as:

ρ \ddot{u} (x, t) = \int_{H_{x}} {\underline{T} [x, t] ⟨ x^{'} - x ⟩ - \underline{T} [x^{'}, t] ⟨ x - x^{'} ⟩} d V_{x^{'}} + b (x, t) + f_{c o n} (x, t) + f_{s} (x, t) .

(25)

Search strategy

For the search strategy in peridynamic simulations, there are mainly three typical choices in the numerical implementations brute-force search method, the grid search method, and the KD-tree search method. The brute-force search method solves problems by exhaustively enumerating all possible solutions. It is simple and easy to implement, but its efficiency is relatively low and causes significant computational consumption. The grid search method divides the search region into regular grids of points, and each point only searches other points in its ambient grid. It applies to both continuous and discrete search spaces. Nevertheless, its efficiency depends on the point distribution and even degrades to the brute-force search when the point distribution is considerably non-uniform. The core of the KD-tree search method is the binary tree structure. The points are partitioned into different regions, and the tree structure allows for the rapid exclusion of regions that cannot contain the nearest neighbor, thereby reducing the computational cost.

To have better calculation efficiency, numerical simulations are often implemented through parallel computing. It is evident that the grid search method and KD-tree search method are much more efficient than the brute-force search method. Whereas the grid search method and KD-tree search method still face restrictions when they are parallelized. On the one hand, if a shared-memory hash table is used in the grid search method, ensuring process safety incurs significant resource overhead. Additionally, hash collisions may occur when a large number of points are inserted simultaneously, thus leading to the parallel efficiency being lower than the single-process efficiency. On the other hand, there are two common approaches to parallelization for the KD-tree search method, single-process tree construction with multi-process search and multi-process tree construction with merging and multi-process search. The efficiency improvement of the former one is relatively limited due to the bottleneck in tree construction. The latter one introduces significant overhead to ensure the balance of the merged KD-tree, which may outweigh the benefits of parallelization.

In this study, an integration of the grid search method and the KD-tree search method is proposed to address these limitations and achieve satisfactory parallel performance. The process begins by dividing the search region into regular grids based on the number of parallel processes, following the grid search methodology. To account for edge effects, the boundary of each grid is extended by a certain distance, chosen as the horizon size in this study, to form an expansion space. Within this expansion space, a KD-tree is constructed to facilitate efficient query and search operations. This hybrid approach leverages the strengths of both methods, achieving a balance between computational efficiency and parallel scalability.

To evaluate the efficiency of the proposed strategy, a simple numerical experiment is conducted. A total of 30,000 points are randomly distributed within a three-dimensional space measuring 100 mm × 100 mm × 100 mm, with a search radius of 3 mm. Table 1 illustrates the search times for each method, comparing single-process execution and parallel execution using four processes. To further illustrate the performance of the proposed model and approach, benchmark tests at different problem scales (e.g., 30,000, 100,000, 300,000, and 900,000 material points) are conducted, as shown in Table 2. The proposed hybrid search strategy achieves a satisfactory acceleration on the search time using four parallel processes, especially when the number of material points is sufficient. Notably, the Grid search and KD-tree methods exhibit a great performance, while the proposed method further improves the search efficiency by combining these two methods. The greatest speedup ratio of the proposed method with four processes reaches approximately 100×, demonstrating good parallel scalability.

Table 1.

Comparison of search time across each method with randomly distributed material points.

	Brute-force search (s)	Grid search (s)	KD-tree search (s)	Integration of grid search and KD-tree search (s)
Single-process search	13.249	9.912	8.105
Parallel search	9.324	7.257	7.52/15.003	5.12

Table 2.

Comparison of search time across each method with uniformly distributed material points in different amounts.

Number of material points	Brute-force search (s)	Grid search (s)	KD-tree search (s)	Integration of grid search and KD-tree search (s)
30,000	3.768	2.707	2.496	2.032
100,000	16.605	4.177	3.689	3.175
300,000	218.681	8.716	8.125	4.94
900,000	1943.532	23.316	22.782	16.663

Numerical validation

The oblique impact tests on PMMA were performed in Dorogoy et al. (2010), and experimental data available from two cases are chosen here for numerical validation. The first case is conducted with a 0.3″ projectile impacting on a single PMMA plate, and the second test case involves a 0.5″ projectile impacting on the three-layer PMMA plates. The schematic diagram and geometric model for these two cases are shown in Figures 4 and 5. According to the experimental conditions, the incidence angle of the projectile in these two cases is 30°, and its initial velocity is 720 and 928 m/s, respectively.

In this work, the projectile is assumed to be a steel impactor. The material information in peridynamic simulations includes: for the projectile, we use Young's modulus $E = 210 GPa$ , Poisson's ratio $ν = 0.3$ , and density $ρ = 7800 kg / m^{3}$ . For the PMMA plate, Young's modulus $E = 5.76 GPa$ , Poisson's ratio $ν = 0.42$ , and density $ρ = 1190 kg / m^{3}$ . Values for $θ_{b 1}$ and $θ_{b 2}$ are chosen to be 1.33 and 4.05e⁻⁵ s. Besides, a horizon size of 0.6 mm is used for the peridynamic modeling here, and we use a uniform discretization with the horizon three times the grid spacing. For the current discretization, the resulting time step is approximately Δt = 2.0 × 10⁻⁸ s. The total simulation time is 500 μs for the single plate case and 600 μs for the three-layer plate case.

Figure 6 presents a comparison of the final failure patterns in two cases between experimental observations and peridynamic results. From Figure 6, we can see that the ricochet phenomenon is accurately predicted by the peridynamic simulations, both in terms of the direction and separated position of the projectile, whether in the single plate case or the three-layer plate configuration. In the first case, failure on the rear face of the plate can be evidently observed in both experimental and numerical results, despite the projectile not penetrating the plate. In the second case, the cracks in the bottom plate are effectively reproduced in the peridynamic simulations, although the crack propagation paths exhibit slight discrepancies compared to the experimental results. Notably, the experimental observations in Figure 6(b) reveal that the projectile exits the three-layer targets through the middle plate, resulting in severe failure of both the top and middle plates. The corresponding failure patterns for each of the three plates from peridynamic simulations are illustrated in Figure 7, providing a detailed depiction of the failure characteristics for each plate. Overall, one can conclude that the peridynamic results have a satisfactory agreement with the corresponding experimental findings.

To demonstrate the advantages of the proposed model, a comparative analysis is conducted using the single PMMA plate case between the proposed model and the basic short-range force contact model, as shown in Figure 8. The results indicate that the basic short-range force contact model with a constant friction coefficient fails to accurately capture the ricochet trajectory. The proposed model, which integrates both the dynamic friction coefficient and shank friction force, provides the best agreement with experimental observations in terms of failure patterns, projectile trajectory, and residual velocity. These findings confirm that both the dynamic friction coefficient and shank friction force are essential components for accurately simulating the oblique penetration process.

A convergence study is conducted to evaluate the robustness of the proposed model, where three different grid spacings are chosen here, Δx = 0.15, 0.2, and 0.25. Overly fine grid spacings may lead to prohibitively high computational cost, whereas overly coarse grid spacings may introduce discretization-induced uncertainty in local fracture details. The failure morphology of the PMMA plate and projectile trajectory are compared across the three grid spacings, as depicted in Figure 9. Numerical results with different grid spacings present a great convergence, which has satisfactory agreements compared with experimental measurements, further confirming the stability and reliability of the proposed model and approach.

To have a better understanding of projectile ricochet during oblique impacts, trajectories of the projectile over time for the two cases are displayed in Figures 10 and 11. Figure 10 provides a comprehensive view of the projectile's deflection and yaw throughout the process. In this scenario, the force acting on the projectile body is asymmetrical, thus leading to deflection during the penetration process. The position of the projectile is effectively predicted in the peridynamic simulations, as evidenced by the comparison between Figures 6(a) and 10(e). In Figure 6(a), the projectile has not completely exited the plate at 250 μs, and Figure 10(e) demonstrates a similar position of the projectile at 265 μs. Additionally, the projectile does not penetrate the plate, and its orientation shifts before reaching the bottom surface. Interestingly, the failure area on the rear face is nearly equivalent to that on the impact face.

For the three-layer PMMA plates, the projectile remains nearly horizontal upon initial contact with the middle plate. The projectile ultimately exits the middle plate, with its inclination angle evidently decreasing as it leaves the three-layer plates. In Figure 6(b), the times recorded for the three positions obtained from experimental tests are considerably greater than those in the peridynamic results, attributed to the initial distance between the project and target in the experimental setup. Whereas, the time interval between the projectile at position 1 and position 2 in Figure 6(b) is approximately 220 μs, closely matching the 230 μs interval obtained from the peridynamic results, as shown in Figures 11(a) and 11(d).

The curves depicting the projectile displacement, derived from experimental and numerical results, are plotted in Figure 12. Given the potential for slight deformation of the projectile, displacements at the centroid rather than the warhead are chosen here. Due to the restriction on the experimental conditions, directly recording the projectile displacement is exceedingly difficult. Consequently, the corresponding data presented in Figure 12 is gained from the flash X-ray pictures captured during the experimental tests. The numerical curves are in great coincidence with the experimental data for both cases, aligning closely in terms of both variation trends and actual values.

Based on the projectile displacement curves in Figure 12, we extract discrete data points and compute the relative errors and root-mean-square error. The relative errors for the single plate case at three main positions of the projectile are 2.53%, 3.25%, and 5.03%, respectively. The relative errors for the three-layer plates case at four main positions of the projectile are 1.49%, 1.77%, 2.71%, and 2.97%, respectively. The root-mean-square errors for these two cases are 3.691 and 3.683 mm, respectively. These quantitative indicators have been provided with more rigorous data support for the reliability of the proposed model.

Discussion

The influence of incidence angle on the penetration process

In the previous section, the ricochet phenomenon was observed for both the single plate and three-layer plates when the incidence angle was 30°. It is well established that the incidence angle is a critical factor in the penetration process, particularly in oblique impact scenarios, as it significantly influences the projectile's trajectory. To further investigate this influence, additional cases with incidence angles of 10°, 20°, 45°, and 60° are analyzed, while maintaining the initial velocity consistent with the previous section.

Trajectories of the projectile and failure patterns of the targets for these cases are displayed in Figures 13 and 14, representing the single plate and three-layer plates, respectively. For the single plate, the ricochet phenomenon is still evident at incidence angles of 10° and 20°, even if the projectile does not enter the plate at the incidence angle of 10°, so that the projectile just slides across the surface upon contact. Among the cases of 10°, 20°, and 30°, the failure area increases progressively, as does the time required for the projectile to ricochet off the plate. At higher incidence angles of 45° and 60°, the projectile directly penetrates the plate, with minimal deviation in its trajectory during the penetration process.

For the three-layer plates, as illustrated in Figure 14, the effect of incidence angle on the penetration process exhibits a similar trend to that observed for the single plate. The ricochet phenomenon occurs at incidence angles of 10° and 20° but is absent in the cases of 45° and 60°. Moreover, the failure area of the former two cases is significantly larger than that of the latter two cases. Notably, the projectile's trajectory undergoes a complete alteration after exiting the target in the case of 20°, as shown in Figure 14(b). In contrast, the projectile is even stuck in the three-layer plates in the case of 30°, as depicted in Figure 14(c). The reason for these differences can be further understood by analyzing the variation in the projectile's velocity, as discussed below.

The residual velocity of the projectile for these cases is recorded in Figure 15. As shown in Figure 15(a), we can see that the residual velocity of the projectile is the lowest among these cases when the incidence angle is 30°, corresponding to the longest duration of projectile interaction with the plate. Although the projectile does not penetrate the plate at incidence angles of 10°, 20°, and 30°, its residual velocity decreases dramatically as the incidence angle increases. This law is quite the opposite for the cases of 45° and 60°, where the residual velocity of the projectile gradually increases, reaching values much higher than that observed at 30°.

Additionally, for the three-layer plates, the influence of the incidence angle on the penetration process follows a similar trend to that of the single plate. It is worth noting that at the incidence angle of 45°, the residual velocity still decreases at the final stage. Based on the sectional view, it is inferred that the projectile can penetrate the three-layer plates and exit through the bottom plate with a thoroughly altered trajectory.

The influence of initial velocity on the penetration process

To investigate the effect of initial velocity on the penetration process, four representative velocities are chosen here, 500, 600, 800, and 900 m/s for the single PMMA plate, and 700, 800, 1000, and 1100 m/s for the three-layer PMMA plates. In all cases, the incidence angle is fixed to be 30°. Figures 16 and 17 show the trajectories of the projectile and corresponding failure patterns of the PMMA plates.

For the single plate, the projectile ricochet is observed at lower initial velocities of 500 and 600 m/s. The projectile does not totally enter the plate in the former case, and its trajectories and failure area closely resemble those observed when the initial velocity is 720 m/s with the incidence angle of 20° (as shown in Figure 13(b)). Besides, a comparison between Figures 6(a) and 16(b) reveals no significant differences in terms of the projectile trajectory and failure area, apart from slight variations in the crack propagation path on the top surface of the plate. At higher initial velocities of 800 and 900 m/s, the projectile directly penetrates the plate. Nevertheless, the projectile's direction is different between these two scenarios, as well as the failure area in the target.

For the three-layer plates, the ricochet phenomenon occurs across all analyzed initial velocities, as shown in Figure 17. However, the projectile trajectories exhibit significant variations with changes in the initial velocity. When the initial velocity is 700 and 800 m/s, the projectile moves in an opposite direction after exiting the target, moving back toward the incidence location. At the same time, the projectile ricochets from the three-layer plates at the boundary between the top and middle plates in Figure 17(a) and (b), while at the middle plate of 1000 and 1100 m/s in Figures 17(c) and 17(d). Additionally, the number of cracks in the three-layer plates increases dramatically with the rise of the initial velocity, indicating more severe damage at higher impact energies.

Figure 18 presents the residual velocity of the projectile for these cases, revealing distinct trends in Figure 18(a) and (b). In Figure 18(a), the residual velocity of the projectile is significantly higher when it directly penetrates the plate compared to when the ricochet phenomenon occurs. There is a substantial gap between the final residual velocities of the green and purple curves, despite the initial velocity differing by only 100 m/s. Furthermore, the residual velocity of the projectile when the initial velocity is 500 m/s is the highest among the three cases of 500, 600, and 720 m/s, where the projectile ricochets in both cases. In Figure 18(b), the ricochet phenomenon is observed in all these cases for the three-layer plates, but the variation in residual velocity across the curves is nearly consistent, with only slight differences in the final residual velocities for different initial velocities. Figure 18(b) illustrates that the projectile loses more kinetic energy at higher initial velocities during the impact process, which is a reasonable explanation for the substantial increase in the number of cracks in the three-layer plates as the initial velocity increases.

The influence of attack angle on the penetration process

In this section, the influence of attack angle, defined as the angle between the axis of the projectile and its velocity direction, on the penetration process is discussed through several numerical cases. A positive attack angle indicates that the projectile pitches upward, with its head slightly higher than the velocity direction, whereas a negative attack angle means the projectile pitches downward, with its head slightly lower than the velocity direction. The incidence angle is fixed at 30° for both the single PMMA plate and the three-layer PMMA plates, with initial velocities of 720 and 928 m/s, respectively. The attack angles are chosen to be +1°, +2°, +3°, +4°, +5°, and −1°, −2°, −3°, −4°, −5°.

Trajectories of the projectile and final failure patterns of the single plate with different attack angles are illustrated in Figures 19 and 20. For positive attack angles, we can see that the ricochet phenomenon exists in all these cases. Trajectories of the projectile and failure areas of the plate are generally consistent across these cases, except for the case of +5°, where an obvious deflection of the projectile direction is observed, and the failure area of the plate is smaller than in other cases as well. Compared to the results with no attack angle in Figure 6, the projectile ricochets much more easily here, even if the attack angle is +1°. In contrast, results are different when the attack angle is negative. As shown in Figure 20, the projectile can readily penetrate the plate when the attack angle is −2°, −3°, −4°, and −5°, just with variations in the deflection of the projectile after exiting the plate. Nevertheless, the result with the attack angle of −1° in Figure 20(a) is quite similar to that with the attack angle of 0° in Figure 6(b), in terms of both trajectories of the projectile and failure patterns of the plate.

Residual velocities of the projectile with different attack angles are shown in Figure 21, where the results for an attack angle of 0° are inserted to have a better comparison. The curves in Figure 21(a) are similar to those in Figure 21(b) (except for the case for attack angle of −1°) in terms of the variation trend and final residual velocity, despite differences in the projectile directions among these cases. Besides, the residual velocity of the projectile increases gradually with larger positive attack angles. Notably, the residual velocities of the projectile in cases of positive attack angles are consistently higher than those for the attack angle of 0°. Conversely, when the attack angles are negative, the residual velocities are much lower than those when the attack angle is 0°, particularly at the beginning of the impact process, where the projectile experiences a sharp velocity drop. Moreover, we can see that the black curve in Figure 21(b) is fairly close to the blue curve, which explains why the projectile direction in Figure 20(a) is contrary to other cases with negative attack angles.

For the three-layer plates, the final failure patterns and trajectories of the projectile with different attack angles are illustrated in Figures 22 and 23. Combining the results from these figures, it is evident that no projectile penetrates the plates in any of the cases, regardless of whether the attack angle is positive or negative. This contrasts with the results for the single plate, where penetration occurs for negative attack angles. The failure patterns of the plates remain nearly unchanged compared to the results with the attack angle of 0° in Figure 6(b). However, the deflection of the projectile is more pronounced for positive attack angles in Figure 22 than for negative attack angles in Figure 23, indicating that positive attack angles have a greater effect on the projectile's direction. In Figure 23(d), (c), and (f), the projectile is still stuck in the plates, particularly in the latter two cases, demonstrating that the deceleration effect of negative attack angles on the projectile is more evident.

Figure 24 shows the residual velocities of the projectile with different attack angles. There is no clear reduction in residual velocities of the projectile when the attack angle is 0°, compared with other curves in Figure 24, unlike the results in Figure 21, where the residual velocities of the projectile are smaller when the attack angle is 0° than others. Additionally, the residual velocity curves for positive attack angles are fairly close to each other, except for the special case of the attack angle of +4°, where the residual velocity is unexpectedly the lowest. The reason for this anomaly is not mentioned in this study, which will be thoroughly investigated in further work. For negative attack angles, the residual velocity curves for 0°, −1°, and −2° are almost identical, while a downward trend is observed for −3°, −4°, and −5°, which is why the projectile becomes lodged in the plates in Figure 23(d), (c), and (f).

Conclusion

In the present work, an improved peridynamic contact model has been proposed, incorporating the dynamic friction coefficient and shank friction force, along with an integration of the grid search method and KD-tree search method for the contact search. The proposed model and approach are validated through the oblique impact tests on PMMA, demonstrating satisfactory agreement between the peridynamic simulation results and corresponding experimental data. The ricochet phenomenon is accurately predicted in the peridynamic simulations, including the failure patterns of the plate and trajectories of the projectile, confirming the model's capability to effectively capture the failure characteristics of PMMA under oblique penetration.

The failure characteristics of PMMA under oblique penetration are comprehensively analyzed, with particular emphasis on the ricochet phenomenon. The effects of key parameters, including the incidence angle, initial velocity, and attack angle, on the oblique penetration process are systematically investigated. These effects are evaluated in terms of failure patterns, ballistic trajectories, and projectile deflection angles. Additionally, the residual velocity of the projectile is also analyzed, providing critical insights for the analysis of the oblique impact failure mechanism under various impact loading conditions.

Despite the progress achieved in this work, several issues remain unresolved, warranting further investigation. On the one hand, the deflection of the projectile is evident in oblique impacts, but no explicit patterns or laws governing this behavior under varying impact conditions have been established. On the other hand, the direction of the projectile is opposite to its initial direction after exiting the target in certain scenarios, potentially striking the target again. This introduces additional complexity to the impact failure mechanism, which requires further exploration to fully understand its implications and effects on the overall penetration process. Besides, the numerical consumption for such impact cases is still significant. To further improve the numerical efficiency, extending the proposed contact model into the dual-horizon PD-FEM platforms (Bie et al., 2020, 2025) will be considered in future work.

Footnotes

Acknowledgments

The authors would like to acknowledge the financial support from the National Natural Science Foundation of China (No. 12202136) and the Youth Funding of Rocket Force University of Engineering (No. 2023QN-S003).

ORCID iD

Liwei Wu

Credit statement

Junbin Guo: Methodology, Investigation, Software, Validation, Writing-original draft. Han Wang: Resources, Investigation, Funding acquisition, Writing-review & editing. Liwei Wu: Conceptualization, Methodology, Funding acquisition, Supervision, Writing-review & editing. Chuanqiang Yu: Project administration, Resources. Shuai Hou: Resources, Software.

Funding

The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work was supported by the Youth Funding of Rocket Force University of Engineering, National Natural Science Foundation of China, (grant number 2023QN-S003, 12202136).

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Data availability statement

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

References

Amani

Oterkus

Areias

, et al. (2016) A non-ordinary state-based peridynamics formulation for thermoplastic fracture. International Journal of Impact Engineering 87: 83–94.

Ben-Dor

Dubinsky

Elperin

(2007) Localized interaction models with non-constant friction for rigid penetrating impactors. International Journal of Solids and Structures 44(7-8): 2593–2607.

Bie

Liu

Yang

, et al. (2020) Abaqus implementation of dual peridynamics for brittle fracture. Computer Methods in Applied Mechanics and Engineering 372: 113398.

Bie

Ren

Bui

, et al. (2024a) Dual-horizon peridynamics modeling of coupled chemo-mechanical-damage for interface oxidation-induced cracking in thermal barrier coatings. Computer Methods in Applied Mechanics and Engineering 430: 117225.

Bie

Ren

Duan

, et al. (2025) Abaqus implementation of coupled thermo-mechanical dual-horizon bond-based peridynamics for thermally induced fractures in thermal barrier coatings. Theoretical and Applied Fracture Mechanics 141: 105298.

Bie

Ren

Rabczuk

, et al. (2024b) The fully coupled thermo-mechanical dual-horizon peridynamic correspondence damage model for homogeneous and heterogeneous materials. Computer Methods in Applied Mechanics and Engineering 420: 116730.

Bie

Wei

Rabczuk

, et al. (2024c) The implicit stabilized dual-horizon peridynamics-based strain gradient damage model. Applied Mathematical Modelling 128: 630–658.

Bobaru

(2012) Damage progression from impact in layered glass modeled with peridynamics. Central European Journal of Engineering 2: 551–561.

Børvik

Olovsson

Dey

, et al. (2011) Normal and oblique impact of small arms bullets on AA6082-T4 aluminium protective plates. International Journal of Impact Engineering 38(7): 577–589.

10.

Chai

, et al. (2019) On the friction effects on rigid-body penetration in concrete and aluminium-alloy targets. Defence Technology 15(4): 576–581.

11.

Chen

(1989) Penetration into dry porous rock: A numerical study on sliding friction simulation. Theoretical and Applied Fracture Mechanics 11(2): 135–141.

12.

Chen

(2022) A multiscale method coupling peridynamic and boundary element models for dynamic problems. Computer Methods in Applied Mechanics and Engineering 401: 115669.

13.

Dai

Xie

Jiang

(2023) A coupled peridynamics–smoothed particle hydrodynamics model for fracture analysis of fluid–structure interactions. Ocean Engineering 279: 114582.

14.

Dai

Xie

, et al. (2021) Numerical modeling on crack propagation based on a multi-grid bond-based dual-horizon peridynamics. Mathematics 9(22): 2848.

15.

Dorogoy

Rittel

Brill

(2010) A study of inclined impact in polymethylmethacrylate plates. International Journal of Impact Engineering 37(3): 285–294.

16.

Esfahlani

(2021) Ballistic performance of polycarbonate and polymethyl methacrylate under normal and inclined dynamic impacts. Heliyon 7: e06856.

17.

Faye

Parameswaran

Basu

(2016) Dynamic fracture initiation toughness of PMMA: A critical evaluation. Mechanics of Materials 94: 156–169.

18.

Forquin

Nasraoui

Rusinek

, et al. (2012) Experimental study of the confined behaviour of PMMA under quasi-static and dynamic loadings. International Journal of Impact Engineering 40: 46–57.

19.

Zhang

Xia

(2017) Voronoi-based peridynamics and cracking analysis with adaptive refinement. International Journal for Numerical Methods in Engineering 112(13): 2087–2109.

20.

Guan

Yan

Guo

(2023) An adaptive contact model involving friction based on peridynamics. European Journal of Mechanics - A/Solids 100: 104966.

21.

Guo

Zhang

, et al. (2020) Multi-scale peridynamic formulations for chloride diffusion in concrete. Engineering Analysis with Boundary Elements 120: 107–117.

22.

Hakan

Gunes

Apalak

, et al. (2023) Normal and oblique ballistic impact damage behaviour of functionally graded plates: Experimental and numerical. International Journal of Impact Engineering 181: 104756.

23.

Han

Lubineau

Azdoud

(2016) Adaptive coupling between damage mechanics and peridynamics: A route for objective simulation of material degradation up to complete failure. Journal of the Mechanics and Physics of Solids 94: 453–472.

24.

Huang

Jiang

(2021) Modeling and studies of fracture in functionally graded materials under thermal shock loading using peridynamics. Theoretical and Applied Fracture Mechanics 111: 102852.

25.

Chen

Spencer

, et al. (2018) Thermomechanical peridynamic analysis with irregular non-uniform domain discretization. Engineering Fracture Mechanics 197: 92–113.

26.

Huang

Zhang

L-W

Liew

(2023) A hybrid polymer–water peridynamics model for ballistic penetration damage of soft materials. Computer Methods in Applied Mechanics and Engineering 415: 116216.

27.

Kamensky

Behzadinasab

Foster

, et al. (2019) Peridynamic modeling of frictional contact. Journal of Peridynamics and Nonlocal Modeling 1(2): 107–121.

28.

Kazarinov

Bratov

Morozov

, et al. (2020) Experimental and numerical analysis of PMMA impact fracture. International Journal of Impact Engineering 143: 103597.

29.

Huang

, et al. (2021) Nonlocal steady-state thermoelastic analysis of functionally graded materials by using peridynamic differential operator. Applied Mathematical Modelling 93: 294–313.

30.

Liu

Kong

Fang

, et al. (2025) Peridynamics modelling of projectile penetration into concrete targets. International Journal of Impact Engineering 195: 105110–105110.

31.

Liu

Huang

(2009) Numerical simulations of oblique-angle penetration by deformable projectiles into concrete targets. International Journal of Impact Engineering 36(3): 438–446.

32.

Oterkus

, et al. (2021) Modelling of cracks with frictional contact based on peridynamics. Theoretical and Applied Fracture Mechanics 116: 103082.

33.

Huang

, et al. (2023) An extended peridynamic model for analyzing interfacial failure of composite materials with non-uniform discretization. Theoretical and Applied Fracture Mechanics 125: 103854.

34.

Macek

Silling

(2007) Peridynamics via finite element analysis. Finite Elements in Analysis and Design 43(15): 1169–1178.

35.

Montgomery

(1976) Friction and wear at high sliding speeds. Wear 36(3): 275–298.

36.

Pathrikar

Rahaman

Roy

(2019) A thermodynamically consistent peridynamics model for visco-plasticity and damage. Computer Methods in Applied Mechanics and Engineering 348: 29–63.

37.

Rosenberg

Dekel

(2012) Terminal ballistics. Springer, Cham: Switzerland.

38.

Silling

(2000) Reformulation of elasticity theory for discontinuities and long-range forces. Journal of the Mechanics and Physics of Solids 48(1): 175–209.

39.

Silling

Askari

(2005) A meshfree method based on the peridynamic model of solid mechanics. Computers & Structures 83(17-18): 1526–1535.

40.

Silling

Epton

Weckner

, et al. (2007) Peridynamic states and constitutive modeling. Journal of Elasticity 88(2): 151–184.

41.

Silling

Lehoucq

(2008) Convergence of peridynamics to classical elasticity theory. Journal of Elasticity 93(1): 13–37.

42.

Stronge

(2018) Impact mechanics. Cambridge: Cambridge University Press.

43.

Tian

Zhu

Sheng

, et al. (2024) A penetration correction contact model involving friction based on Peridynamics. Computers and Geotechnics 169: 106250.

44.

Wang

Huang

(2019) A non-ordinary state-based peridynamic formulation for thermo-visco-plastic deformation and impact fracture. International Journal of Mechanical Sciences 159: 336–344.

45.

Wang

Yin

Z-Y

Chen

(2024) Peridynamic modelling of dynamic damage and fragmentation of cracked solids during impact contact. Engineering Fracture Mechanics 295: 109784.

46.

Wang

Huang

, et al. (2007) Analysis of rigid motion for penetration of concrete. In: Proceedings of the 7th International Conference on Shock & Impact Loads on Structures, pp.627–638.

47.

Xia

(2004) Experimental study of tensile properties of PMMA at intermediate strain rate. Materials Letters 58(29): 3681–3685.

48.

Huang

(2022) Peridynamic modeling and simulations on concrete dynamic failure and penetration subjected to impact loadings. Engineering Fracture Mechanics 259: 108135.

49.

Huang

Bobaru

(2021) A reformulated rate-dependent visco-elastic model for dynamic deformation and fracture of PMMA with peridynamics. International Journal of Impact Engineering 149: 103791.

50.

Huang

, et al. (2022) Peridynamic modeling for impact failure of wet concrete considering the influence of saturation. International Journal of Damage Mechanics 31(9): 1448–1474.

51.

Huang

Wang

, et al. (2023) A comparison study on numerical analysis for concrete dynamic failure using intermediately homogenized peridynamic model and meso-scale peridynamic model. International Journal of Impact Engineering 179: 104657.

52.

Wang

Huang

, et al. (2020) An ordinary state-based peridynamic modeling for dynamic fracture of laminated glass under low-velocity impact. Composite Structures 234: 111722.

53.

Haijun

Fenglei

, et al. (2012) Study on the Behavior of Projectile Shank Friction Penetrating into Concrete. In: Transactions of Beijing institute of Technology, pp.12–17.

54.

Zhang

Dai

(2024) A fully coupled thermomechanical peridynamic model for rock fracturing under blast loading considering initial pore damage. Computers and Geotechnics 171: 106382.

55.

Zhang

Dai

(2025) A novel fracture model for composite laminates based on bond-based peridynamics. Engineering Analysis with Boundary Elements 176: 106229.

56.

Zhang

Liu

(2022) A peridynamic model for contact problems involving fracture. Engineering Fracture Mechanics 267: 108436.

57.

Zhong

Han

, et al. (2024) Multi-time-step coupling of peridynamics and classical continuum mechanics for dynamic brittle fracture. Engineering Fracture Mechanics 307: 110264.

58.

Zhu

Zhao

(2021) Peridynamic modelling of blasting induced rock fractures. Journal of the Mechanics and Physics of Solids 153: 104469.