Bridging the gap between local and nonlocal numerical methods—A unified variational framework for non-ordinary state-based peridynamics

Highlights

• A unified variational framework is proposed to bridge the gap between PD and CCM.

• A new force state vector is introduced to ensure the consistency between PD and CCM.

• The proposed variational framework unifies boundary conditions in PD and CCM.

• A fully implicit algorithm is developed to simulate general nonlinear problems.

• A penalty method is employed to obtain a stable numerical solutions in NOSB-PD.

Abstract

The paper aims to develop a unified variational framework to bridge the gap between the non-ordinary state-based peridynamics (NOSB-PD) and the classical continuum mechanics (CCM). First, a new force state vector is proposed by introducing the first Piola–Kirchhoff stress. This new force state vector enables the stress divergence of each material point to be expressed by averaging all the force state vectors in its support domain. The new force state vector also ensures the mathematical consistency between the strong form of PD and CCM when the horizon of a material point approaches to zero. Second, the displacement and traction boundaries in CCM are transformed into the non-local fictious boundary layers in PD, and a non-local Gauss’s formulation is presented by transforming the displacement and traction boundaries in CCM into the non-local fictious boundary layers in PD, and this formulation unifies the variational framework and boundary conditions of PD and CCM. Third, a fully implicit algorithm is developed to obtain the general nonlinear problems such as fracture and large deformation of solid materials. Further, a penalty method is employed to eliminate the zero-energy mode oscillation inherently observed in NOSB-PD, and the penalty force and penalty stiffness matrix are derived for the proposed implicit algorithm and numerical implementation. Numerical results demonstrate that the proposed method is accurate and can well capture the fracture and large deformation of solid materials. Results also indicate that the method can effectively prevent the zero-mode oscillations inherently observed in the original NOSB-PD, and thus ensures the computational stability.

Introduction

Discontinuous problems, such as cracks and voids, are fundamental and important issues in the field of solid mechanics. Up to now, the crack simulation remains to be challenging in computational mechanics, and one of the critical issues is how to accurately capture the stress field in discontinuities. Recently, various numerical methods have been developed to simulate the fracture process of solid materials, such as Extended Finite Element Method (XFEM) [1], [2], Generalized Finite Element Method (GFEM) [3], Partition of Unity Method (PUM) [4], Diffuse Element Method (DEM) [5], Element-Free Galerkin (EFG) Method [6] and Material Point Method (MPM) [7], [8]. These methods either require additional variables and equations to capture cracks, or need complex shape functions to interpolate displacement and stress fields, which would be time-consuming for three-dimensional problems.

Peridynamics (PD) first proposed by Silling [9] represents a new theory for the simulation of fracture mechanics, where a continuum is modeled with material points connected by long range interactions. PD can be regarded as a nonlocal reformulation of the classical continuum mechanics (CCM). It formulates mechanical problems based on integral equations instead of partial differential equations, making it more suitable for discontinuous problems. Generally, two kinds of peridynamic models are developed, i.e. bond-based and state-based peridynamic models. The bond-based peridynamic (BB-PD) model was the first PD model proposed by Silling, in which two material points located in a same material horizon will interact with each other in the form of a pairwise force function. However, one of the defects of the bond-based PD is that the response of a bond is independent of other bonds, which makes the Poisson’s ratio restricted to $1∕4$ for the plane strain problem and $1∕3$ for the plane stress problem [9], [10]. Later, the ordinary state-based (OSB) PD and the non-ordinary state-based (NOSB) PD were introduced by Silling to overcome the limitations of the bond-based PD [11], [12]. The term of “ordinary” denotes that the force vector state is parallel to the deformation vector state. Compared to the BB-PD and the OSB-PD, one distinguished feature of NOSB-PD is that it incorporates constitutive model from the CCM such as plasticity [13], viscoelasticity [14] and visco-plasticity [15] into the force state. However, the zero-energy mode oscillation induced by the “PD correspondence material models”[16] in NOSB-PD is still a challenging issue. A detail review of NOSB-PD and OSB-PD models can be referred to [11]. Note, however, that, one defect of PD is that the long-range force between each two material points is non-physical at meso and macro scales [17]. The other one is that, the simulation of standard PD is based on the strong form of balance equation of linear momentum, which is different from the Bubnov–Galerkin method and weighted residual method in the CCM, and thus may result in difficulties in imposing zero traction conditions and truncated error in the interaction domain of a material point near the boundary surface [18]. An alternative way to address such issues is to unify or couple the PD and the CCM from the view of theory and numerical implementation.

Generally, there are three commonly approaches commonly used to unify or couple the CCM and the PD. The first is to divide the continuum into a PD subdomain and a CCM subdomain, and connect the two subdomains using a transition region or at an interface. The PD subdomain is employed to describe the meso-mechanical or micro-mechanical response of materials and local failure processes, whereas the CCM subdomain captures macro-scale response characteristics and further reduces the computational cost. Lubineau and Han et al. [19], [20] proposed a morphing strategy to couple non-local and local continuum mechanics, and introduced a unified model to encompass both local and non-local continuum representations. Galvanetto and Zaccariotto et al. [21], [22] presented a grid-mesh switching technique to couple PD material points and FEM meshes for static and dynamic problems. In this technique the FEM subdomain applies forces at the PD-FEM interface to achieve a direct coupling of FEM and PD. Bie et al. [23] developed a novel approach for the coupling of OSB-PD and node-based smoothed FEM, in which forces and displacements are directly transmitted from the CCM subdomain to the PD subdomain at the interface. Wang et al. [24] introduced a concurrent Arlequin method to couple the PD and the CCM for dynamic fracture problems. Note, however, that, in these coupling strategies, the transition zone or interface has to precedently defined in the model, and thus they are not applicable for scenarios where potential fracture or damage areas are unknown. Moreover, high-frequency waves would be spuriously reflected at locations where there is a transition zone or an interface [25]. The second approach is to develop new force state vectors of NOSB-PD to ensure the consistency of the force state vector in PD and the stress tensor in CCM when the horizon of a material point vanishes. Significant contributions have been made by a number of researchers to the development of this kind of approach. Gu et al. [26] proposed a Refined NOSB-PD by introducing an alternative form of the force density vector, which can be obtained from balance equations of linear momentum in CCM. Later, Gu et al. [27] introduced a bond-associated and higher-order stabilized force state vector to obtain a unified integral formulation for NOSB-PD, Smoothed Particle Hydrodynamics (SPH), Corrected Smooth Particle Hydrodynamics (CSPH), Reproducing kernel Particle Method (RKPM) and Gradient Reproducing kernel Particle Method (G-RKPM). However, a strong form of NOSB-PD is used in these methods, and thus, it is distinct from the Bubnov–Galerkin method and weighted residual methods in CCM. The last approach aims to interpret the long-range force in NOSB-PD from a mathematical view, and take the long-range force as a non-local approximation to describe the divergence of the stress in CCM. Bergel and Li [28] redefined the differential operator in CCM by introducing the integration of force state vectors, and proved the consistency of the new definition of differential operator and the classical differential operator in CCM. Madenci et al. [29], [30] also presented the differentiation of a scalar or a vector function in CCM with the integration of “orthogonal PD functions” and this new technique is referred as the Peridynamic Differential Operator (PDDO). Similarly, Rabczuk et al. [31] and Ren et al. [32], [33], [34] proposed the Nonlocal Operator Method (NOM) by introducing the idea of non-local approximation in dual-horizon Peridynamics [35]. Note that in the NOM, the differentiation of a scalar or a vector field function at a material point can be obtained by averaging the function values in its support domain. Similar to the dual-horizon Peridynamics [35], the NOM also allows for the h-adaptive refinement and provides a novel approach to solve partial differential equations (PDEs). Note that, the linear momentum balance equation in CCM is expressed in the form of a PDE related to the stress tensor. The critical issue is that, in essence, how to unify the connection between the force state vector in PD and the stress tensor in CCM, as well as the weak form and the numerical implementation between them.

The goal of this paper is to propose a new NOSB-PD model named by Unified Variational NOSB-PD (UVPD) that provides a unified variational framework and boundary conditions for PD and CCM at both continuous and discrete levels. First, a new force state vector is developed by introducing the first Piola–Kirchhoff stress, and the divergence of stress at each material point in CCM is obtained by averaging all the force state vectors in support domains. The new force state vector ensures that the strong form of PD is theoretically consistent with that of CCM when the horizon of each material point approaches to zero. Second, converting displacement and traction boundaries in CCM into “non-local fictious boundary layers” in PD, a non-local Gauss’s formulation is presented in the framework of PD to unify the variational framework and boundary conditions of NOSB-PD and CCM. Similar to the nonlocal numerical methods, such as PDDO [30], [18] and NOM [32], [33] in nonlocal math communities, the proposed unified variational formulation enables displacement and traction boundary conditions easily incorporated in numerical implementation, which is difficult or even impossible for the variable horizon Peridynamics. More importantly, there is no need of the surface correction technique [10], [36], [37] in the numerical implementation of the proposed method. Third, a proposed implicit PD formulation is derived from the unified variational framework developed in the paper, which is consistent with the Galerkin formulation widely used in the Finite element Method and the Meshfree Method. Nevertheless, the variable horizon Peridynamics is based on the strong form of the linear momentum equation, and thus it is difficult for both priori and posterior error analyses. A fully implicit algorithm is also developed to solve the general nonlinear problems such as large deformation and fracture of solid materials. Moreover, a penalty method is employed to eliminate zero-mode oscillation of NSOB-PD and ensure the computational stability. The penalty force and penalty stiffness matrix are also derived for the proposed implicit algorithm and numerical implementation. Finally, several numerical examples are presented to demonstrate the validation of the proposed method by comparing its predictions with those of analytical solutions and experimental results.

The remainder of this paper is organized as follows. In Section 2, a new force state vector is introduced and the theoretical consistency of the strong form between PD and CCM is proved. Then, a unified variational framework at continuous level is developed using the non-local Gauss’s theorem. In Section 3, two stress-based failure criteria, i.e. the maximum tensile stress failure criterion and the Mohr–Coulomb failure criterion, are presented and the penalty method for eliminating zero-energy mode oscillation is also introduced. In Section 4, a fully implicit algorithm is developed for the proposed UVPD model. Numerical examples are presented in Section 5 to validate the proposed UVPD model. Finally, we conclude in Section 6.