An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading

Engineering Fracture Mechanics xxx (2017) xxx–xxx

Contents lists available at ScienceDirect

Engineering Fracture Mechanics journal homepage: www.elsevier.com/locate/engfracmech

An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading Heng Zhang a, Pizhong Qiao a,b,⇑ a State Key Laboratory of Ocean Engineering, Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration, School of Naval Architecture, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240, PR China b Department of Civil and Environmental Engineering, Washington State University, Sloan Hall 117, Pullman, WA 99164-2910, USA

a r t i c l e

i n f o

Article history: Received 18 May 2017 Received in revised form 16 September 2017 Accepted 16 September 2017 Available online xxxx Keywords: Peridynamics State-based model Thermal Interface crack propagation Damage growth Bimaterials

a b s t r a c t An extended ordinary state-based peridynamic model considering thermomechanical loading is presented to predict damage growth of bimaterial structures, such as cermet. In this new model, the three-dimensional (3D) and two-dimensional (2D) (both plane stress and strain) cases are all considered. As examples, 2D bimaterial beams and 3D thick plates are analyzed under thermal loading and three-point bending. m-convergence and dconvergence are discussed in the cases of 2D verification, and comparison of displacement with finite element model shows great accuracy of the extended model. Damage growth (in term of crack propagation) of bimaterial beams due to incremental thermal loading and three-point bending is investigated. The new model successfully captures interface crack propagation in bimaterial beams under thermal loading as well as crack growth within substrate material and at bimaterial interface under quasi-static and impact loading. Distribution of elastic strain energy density is analyzed during dynamic crack propagation under impact loading. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction Bimaterials, such as cermet, have been widely used because of their two different phase properties, coupling high hardness, and deforming capacity. Analysis, especially failure prediction, of cermet due to thermal and mechanical loading is necessary, and it has attracted much attention. Analytical solutions [1,2] for thermomechanical crack problems of bimaterials are available only for a few cases. The numerical methods, such as finite element method [3,4], extended finite element method [5], and boundary element method [6], were used to analyze behaviors of interfacial fracture in bimaterials under effect of thermomechanical loads. However, these available methods, which are based on classical local theory with assumption of displacement continuity, are naturally unsuitable to failure analysis of cermet. To ease inadequacies of classical local theory, theory of peridynamics [7] was formulated to handle problems involving discontinuities. Essentially, the theory of peridynamics is a reformulation of continuum mechanics, in which the integral-differential equations are established to replace the partial differential equations. This

⇑ Corresponding author at: Department of Engineering Mechanics, School of Naval Architecture, Ocean and Civil Engineering, Shanghai Jiao Tong University, Shanghai 200240, PR China. E-mail addresses: [email protected], [email protected] (P. Qiao). https://doi.org/10.1016/j.engfracmech.2017.09.023 0013-7944/Ó 2017 Elsevier Ltd. All rights reserved.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

2

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Nomenclature b d ed e½t 0 hni E G0 Hx KIC m s0 T t u v W

q l

X dV=V

edij b

H f ðHÞ h u

x d Dx 0 k ; a; b0

applied body force density damage history scalar state deviatoric extension state extension value of bond n at time t 0 Young’s modulus critical fracture energy release rate neighborhood of point x critical stress intensity factor ratio between the horizon size and the grid spacing critical stretch value force vector state scalar force state as the magnitude of T displacement of material point x Poisson’s ratio peridynamic strain energy density mass density shear moduli classical elastic strain energy density volume dilatation deviatoric strain tensor coefficient of thermal expansion temperature variation energy variation due to temperature change peridynamic volume dilatation volume weight of broken bonds influence function horizon value uniform grid spacing positive peridynamic constants

means that compared to displacement derivatives used in classical mechanics equations, which are not defined at discontinuities, the peridynamics-based formulation is applicable for fracture and discontinuity analysis. Thus, initiation and propagation of crack can be modeled by peridynamics without any special techniques. The original formulation of peridynamics, so called ‘‘bond-based”, assumes that points are connected with bonds through spring-like interactions and response in a bond is independent of other bonds. However, this formulation has a restriction on material properties, i.e., the Poisson’s ratio requires to be 1/3 for the 2D plane stress case and 1/4 for both cases of plane strain and 3D [7,8]. Even the proven model [9] and computational techniques [10] were proposed to overcome the constraint limitation of Poisson’s ratio, the bond-based model still cannot capture general characteristics of materials. To reduce the constraint, a more general framework, called ‘‘state-based” peridynamics [11], was proposed, and the bond force density between points depends on deformations of points of the whole family. The state-based peridynamic model can be classified into ‘‘ordinary state-based peridynamics” and ‘‘non-ordinary statebased peridynamics” [11]. The ordinary state-based model employs explicit dependence of volumetric and distortional deformations for material response, and it can reproduce material behaviors in conventional theory of solid mechanics, for not only linear elastic material state of 3D [12] and 2D cases [13,14], but also plastic [15,16], viscoelastic [17], and viscoplastic [18] cases. Unlike the ordinary one, the non-ordinary state-based model does not necessarily require force state parallel to deformed position of connected bond, and it was used to model peridynamic beam [19], plates and flat shells [20]. With the capability of handling discontinuities, the peridynamic theory has been successfully applied to damage growth prediction of various problems. Silling [21] presented numerical analysis of Kailthoff-Winkler experiments by peridynamics. Gerstle [9] analyzed plain and reinforced concrete structures through ‘‘micropolar peridynamics model”. Askari et al. [22] investigated failure modes of laminated composites with a large center notch under tension or shear-tension loads. Xu et al. [23] performed simulation of cruciform composite specimens under biaxial loads. Colavito et al. [24] studied nanoclay-epoxy nanocomposites subjected to lower-velocity impact. Ha et al. [25,26] evaluated dynamics crack propagation and crack branching with peridynamics. Oterkus et al. [27] analyzed fracture properties of stiffened composite curved panels with a pre-crack under combined axial tension and internal pressure. Agwai et al. [28] presented crack propagation in multilayer thin-film structures of electronic packages by peridynamics. Zhou et al. [29,30] investigated crack propagation in rock materials by extended non-ordinary state-based peridynamic model.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

3

From the above, peridynamics was utilized to analyze and predict damage growth (e.g., crack initiation and propagation) in various materials. Additionally, peridynamics was considered to solve problems involving thermal effects. Kilic et al. [31] investigated crack paths in quenched glass plate with a bond-based thermomechanical response function. Bobaru et al. [32,33] proposed a peridynamic model for transient heat conduction when damage involves. Oterkus et al. [34] obtained fully coupled peridynamic thermomechanical equations using thermodynamics. However, the state-based peridynamic model considering thermomechanical loading and coupling effect has not been fully studied yet, and it is thus necessary to develop such a model for damage growth prediction of bimaterial structures (e.g., cermet) under thermomechanical loading. In this paper, an extended ordinary state-based peridynamic model for damage growth prediction of bimaterials under thermomechanical loading is presented, and damage growth and thermomechanical coupling effect of bimaterial cermet due to thermal loading or three-point bending are evaluated. First, the force state of the new model is obtained by taking Frechet derivative of peridynamic energy density. 3D and 2D (including both plane stress and strain) cases are all considered. Then, in thermal analysis, the new model is validated with a bimaterial beam subjected to uniform temperature loading in the 2D case as well as a thick plate subjected to gradient thermal loading in the 3D case. The damage growth (in term of crack propagation) prediction of bimaterial beams due to incremental thermal loading is also presented. Three-point bending analysis is conducted to validate the proposed peridynamic model and predict crack paths of bimaterial structures due to quasistatic or impact loading. Both m-convergence and d-convergence [35] are evaluated in 2D verifications, and finite element simulations by ABAQUS are performed to compare with predictions from the peridynamic model. Finally, adaptive dynamic relaxation (ADR) is introduced for the case of quasi-static analysis, while explicit integration method is utilized for the case analysis of bimaterial structures under impact loading. 2. State-based peridynamic model for thermomechanical problems A brief review of state-based peridynamic model is first introduced, followed by the proposed extended ordinary statebased peridynamic model under thermomechanical loading and damage growth prediction model considering thermal effect. Since the elastic strain energy density of point in peridynamics and classical mechanics is equal, the unknown peridynamic constants in peridynamic energy density expression are obtained by equalizing two energies. Then, the peridynamic force state is obtained by taking Frechet derivative of peridynamic energy density. The three-dimensional (3D) and twodimensional (2D) (both plane stress and strain) cases are all considered. 2.1. General state-based peridynamic theory The state-based peridynamic theory was first proposed by Silling [11], as a reformulation of classical continuum mechanics to deal with problems evolving discontinuities and to complement deficiencies of bond-based peridynamic model [7]. In the state-based peridynamic model, a physical system is discretized into finite material points with scalar volumes, and material points interact with all other points within a distance d, called ‘horizon’. The equation of motion for material point x in state-based peridynamics is given as:

Z

2

qðxÞ

d uðx; tÞ dt

2

¼ Hx

fT½x; thx0  xi  T½x0 ; thx  x0 igdV x0 þ bðx; tÞ

ð1Þ

where Hx is the neighborhood of point x, q is the density, and u is the displacement of material point x at time t. x0 is another material point in the neighborhood of x, and bðx; tÞ is the body force density of point x. As shown in Fig. 1, for the ordinary state-based peridynamic model [11], n ¼ x0  x is the vector state, as the bond between point x and x0 . T½x; t and T½x0 ; t are the force vector states given in the constitutive model of points x and x0 , respectively. The internal force density of point x acted by point x0 is comprised of T½x; thx0  xi  T½x0 ; thx  x0 i. Thus, the force density vector of point x is obtained through the integrand of T½x; thx0  xi  T½x0 ; thx  x0 i in Hx . Therefore, the form of T½x; t is the core of peridynamic model. As constitutive model, it establishes connection of kinematic field to mechanics field. In this study, the state-based peridynamic model is extended to consider thermomechanical loading, in form of T½x; t including temperature variation. 2.2. Peridynamic model 2.2.1. Classical elastic strain energy density The linear elastic strain energy density in classical continuum mechanics is [36]:

X¼

 2   X k dV dV þ f ðHÞ þl edij edij  3kbH 2 V V i;j¼1;2;3

ð2Þ

where k and l are the bulk and shear moduli, respectively. dV=V is the volume dilatation, edij is the deviatoric strain tensor, b is the coefficient of thermal expansion, and H is the temperature variation between instantaneous and reference Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

4

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Fig. 1. Ordinary state-based peridynamic model.

temperatures of point x. f ðHÞ is the energy variation only due to temperature change, which is not related to deformation. Therefore, the strain energy density is decomposed into four parts: the volumetric energy density, the distortional energy density, the thermal dilatational energy density, and the temperature-associated energy variation f ðHÞ. In 2D cases, the strain in z direction can be expressed by the strains in x-y coordinate. Thus, the volume dilatation and elastic energy density can be rewritten in the following. In the plane stress model, the stress and strain tensor components are:

2

3

2

r11 r12 0 e11 e12 7 6 r¼6 4 r12 r22 0 5 e ¼ 4 e21 e22 0

0

0

0

0

0

3

7 0 5

ð3Þ

e33

The relationship of stress vs. strain using Hooke’s law and considering thermal loading can be described as:

8 1 > < e11 ¼ E ðr11  v r22 Þ þ bH e22 ¼ 1E ðr22  v r11 Þ þ bH > : e33 ¼ Ev ðr11 þ r22 Þ þ bH

ð4Þ

where E and v are the Young’s modulus and Poisson’s ratio, respectively. The volume dilatation can be expressed as a function of H and the strain in x-y coordinate as:

dV 2v  1 v þ1 ¼ e11 þ e22 þ e33 ¼ ðe þ e Þ  bH V v  1 11 22 v  1 Using Eqs. (4) and (5), the deviatoric component

ed33

ð5Þ

ed33 is then given by

  1 dV v þ 1 dV v þ 1 ¼  ¼ e33  bH 3 V 3ð2v  1Þ V 2v  1

ð6Þ

Substituting Eq. (6) into Eq. (2), the strain energy density in the case of plane stress can be rewritten as:

"

X¼

  #  X k v þ 1 2 dV 2 þl þl edij edij 2 3ð2v  1Þ V i;j¼1;2 "   #    2 v þ1 2 dV v þ1 þl bH  3 k þ 2l bH þ f ðHÞ 3ð2v  1Þ V 2v  1

ð7Þ

While in the plane strain model, the stress and strain tensor components are:

2

r11 r12 6 r ¼ 4 r12 r22 0

0

3

2

3

e11 e12 0 6 7 7 0 5 e ¼ 4 e21 e22 0 5 0 0 0 r33 0

ð8Þ

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

5

The volume dilatation is:

dV ¼ e11 þ e22 þ e33 ¼ e11 þ e22 V

ed33 is given by

The deviatoric component



ed33 ¼ e33 

1 dV 3 V

ð9Þ



¼

  1 dV 3 V

ð10Þ

Substituting Eq. (10) into Eq. (2), the strain energy density in the case of plane stress is reformulated as:

   2 X k l dV dV þ þl edij edij  3kbH 2 9 V V i;j¼1;2



X¼

ð11Þ

ed33 and X are not only a function of the strain in x-y coordinate but also that of the temperature change H. While in the plane strain case, dV ; ed33 are not a function of temperature change H, but the strain energy denV In the plane stress case,

dV V

;

sity X is still a function of H.

2.2.2. Peridynamic strain energy density In the peridynamic model, h is defined as the nonlocal dilatation which is equal to the local volume dilatation [11]. Similar to the definition for 3D in [11] and according to expressions of volume dilatation in Eqs. (5) and (9) for 2D, the expressions of h for 3D, plane stress and plane strain cases can be, respectively, given as:

h¼

8 > > <

3

xxe

> :

2

3D

q

2ð2v 1Þ xxe > ðv 1Þ q

 vv þ1 bH 1

xxe

Plane Stress

ð12Þ

Plane Strain

q

where v is the Poisson’s ratio, x is the influence function, e is the extension scalar, and q is the weighted volume defined by q ¼ ðxxÞ  x. The dot product ðÞ is defined in [11]. Similar to the case in [11], considering thermal effect, the general peridynamic energy density is supposed to take form of: 0

Wðh; ed ; HÞ ¼

k h2 a 0 þ ðxed Þ  ed  3b0 Hh þ f ðHÞ 2 2

ð13Þ

0

where k , a and b0 are the positive constants, h is the volume dilatation defined in Eq. (12), ed is the deviatoric extension state 0 defined as ed ¼ e  hx=3, and f ðHÞ is the energy variation only due to temperature change in peridynamic model. 2.2.3. Peridynamic force state In ordinary model [11], the force vector state T ¼ tM, where t is called the ‘scalar force state’ as the magnitude of T, and M is the unit vector parallel to deformed bond direction. t can be obtained by taking the Frechet derivative of energy density function W, as t ¼ re W. Similar to [13], implementing Frechet derivative of function in Eq. (13) with respect to e, the form of scalar force state t is obtained as

  a 0 t ¼ k h  3b0 H  ðxed Þ  x re h þ axed 3

ð14Þ

where re h is the Frechet derivative of h with respect to e. According to the definition of Frechet derivative in [11] and the form of h in Eq. (12), re h can be written as:

re h ¼

8 xx 3 > > < q

3D

2ð2v 1Þ xx > ðv 1Þ q > : xx 2 q

Plane Stress

ð15Þ

Plane Strain

Now, the general form of force state of peridynamic model considering thermal loading is obtained in Eq. (14), where h 0 and re h are defined in Eqs. (12) and (15), respectively. Next, the peridynamic constant k , a and b0 in Eq. (14) can be obtained by equalizing the strain energy densities in classical continuum mechanics and peridynamics in the three respective cases as follows. In 3D case, according to [11]:

a 2

ðxed Þ  ed ¼

aq X d d ee 15 i;j¼1;2;3 ij ij

ð16Þ

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

6

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Equalizing the energy densities in classical mechanics of Eq. (2) and peridynamics of Eq. (13) and using Eq. (16), the relationships between peridynamic and classical mechanics constants are obtained as: 0

k ¼ k;

a¼

15l ; q

b0 ¼ kb;

0

f ðHÞ ¼ f ðHÞ

ð17Þ

In 2D cases, the peridynamic energy density of Eq. (13) needs to be reformatted in correspondence to classical strain energy density expressed in Eqs. (7) and (11). According to [13], in 2D cases:

a

ðxed Þ  ed ¼

2

aq d 2 aq X d d ee ðe Þ þ 16 33 8 i;j¼1;2 ij ij

ð18Þ

Substituting Eq. (18) into Eq. (13) and using Eq. (6), the peridynamic energy density in the plane stress case takes form of:

Wðh; ed ; HÞ ¼

"   #  0 k aq v þ 1 2 dV 2 aq X d d þ ee þ V 8 i;j¼1;2 ij ij 2 144 2v  1 !   !  2 aq v þ 1 2 dV aq v þ 1 0 b þ 3b H þ  b2 H2 þ f ðHÞ 24 2v  1 V 16 2v  1

ð19Þ

0

Comparing the energy density of plane stress model in Eqs. (7) and (19), the peridynamic constant parameters k ; a; b0 in the case of plane stress are expressed as





v þ1 2 8l k ¼kþl ; a¼ ; 3ð2v  1Þ q   l v þ1 2 2 2 0 f ðHÞ ¼ f ðHÞ þ b H 2 2v  1 0



v þ1 kþl 3ð2v  1Þ

0

b ¼

2 ! b; ð20Þ

Similarly, substituting Eq. (18) into Eq. (13) and using Eq. (10), the peridynamic energy density in the plane strain case is:

Wðh; ed ; HÞ ¼

 0   k aq dV 2 aq X d d dV þ e e  3b0 H þ f ðHÞ þ V 8 i;j¼1;2 ij ij V 2 144

ð21Þ 0

Then, comparing the energy densities in Eqs. (11) and (21), the constant parameters k ; a; b0 for the case of plane strain are given as: 0

k ¼kþ

l 9

;

a¼

8l ; q

b0 ¼ kb;

0

f ðHÞ ¼ f ðHÞ

ð22Þ

The scalar force state t of state-based peridynamic model considering thermomechanical loading is given as:

8 0 xx 3ðk h  3b0 H  a3 ðxed Þ  xÞ q þ axed > > < xx v 1Þ 0 ðk h  3b0 H  a3 ðxed Þ  xÞ q þ axed t ¼ 2ð2 ðv 1Þ > > : 0 xx 2ðk h  3b0 H  a3 ðxed Þ  xÞ q þ axed

3D Plane Stress

ð23Þ

Plane Strain

0

where h is defined in Eqs. (12), and k ; a; b0 are given in Eqs. (17), (20) and (22) for the cases of 3D, plane stress, and plane strain, respectively. In summary, the strain energy density and scaler force state of state-based peridynamics considering thermomechanical loading are given in Eqs. (13) and (23), respectively. Similar to classical thermoelasticity [36], the uniform temperature variation H directly causes its volume dilatation variation, and the volume variation in turn leads to the changes of strain energy density and scaler force state. Specifically, in 3D model, the third part of t is equal to zero according to [11], and substituting Eqs. (15) and (17) into Eq. (14), the scalar force state t in Eq. (23) can be rewritten as:

t ¼ 3kðh  3bHÞ

xx q

þ

15l d xe q

ð24Þ

Similar to classical mechanics model, the scalar force state of 3D peridynamic model is decomposed into co-isotropic, codeviatoric, and thermal-extensional parts. Unlike that given in [11], the thermal dilatation needs to be deducted from the volume dilatation due to thermal loading. Obviously, the scalar force state of present new model can be reduced to that in [11] for 3D and [13] for 2D when the temperature change H is equal to 0. Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

7

2.3. Peridynamic damage growth model with thermomechanical effect Peridynamics is primarily proposed to deal with discontinuous problems, and it shows great potential in failure pattern and damage growth prediction. However, quantitative analysis of fracture problems is still difficult, and suitable damage models are needed to build the relationship between classical mechanical parameters and peridynamic damage constants. In [21], the bond stretch criterion of damage model for prototype microelastic brittle materials was proposed by building relationship of critical stretch value s0 with classical fracture parameter energy release rate G0. In the present study, the scalar state d is defined to track damage history of pair point vector n ¼ x0  x, as the ‘bond’ in the bond-based model. Similar to [21], the damage model can take the form:

( dhni ¼

1 if 0

e½t 0 hni jnj

 bH < s0 for all t0 6 t;

otherwise

ð25Þ

where jnj is the undeformed length of bond n, and e½t 0 hni is the extension value of bond n at time t0 . Unlike the one given in [21], the elastic stretch e½t 0 hni=jnj  bH rather than bond stretch e½t 0 hni=jnj is used to check if the bond is deformed large enough to break because of thermal expansion. At every time step, the elastic stretch of each pair of bonded material points is computed, and bond breaks when the elastic stretch exceeds the critical stretch value s0 . If bond x0  x is broken, the scalar force thx0  xi of point x is set to zero. The internal force density of point x acted by point x0 is reduced to T½x0 ; thx  x0 i, while the internal force density of point x0 acted by point x is reduced to T½x0 ; thx  x0 i. While the critical stretch value s0 is related to the fracture mechanics constants [37] as:

8 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3D > 0 =9kd < p5G ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s0 ¼ 2ð1  v ÞpG0 =3Ed Plane Stress > ffi : pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2ð1  v  v ÞpG0 =3Ed Plane Strain

ð26Þ



where G0 is the energy release rate of material, which can be calculated from G0 ¼ K 2IC =E and G0 ¼ K 2IC 1  v 2 =E, for the cases of plane stress and plane strain, respectively, and KIC is the critical stress intensity factor. As shown in Eq. (26), for 3D, plane stress and plane strain cases, the forms of critical stretch value s0 are different. Moreover, as for multi-material systems, there is no special interface layer points in discretized peridynamic model, and every point belongs to unique material and has unique material parameters. For pair points across interface, the smaller value of two materials is chosen to calculate the critical stretch value s0 of bond. For example, considering two material points x and x0 connected by a peridynamic bond with different fracture parameters G0 and G00 , leading to different respective critical stretch values s0 and s00 , the smaller value of s0 and s00 is chosen as the critical stretch value of this bond. In addition, the point variable parameter damage u as the volume weight of broken bonds is defined as:

Z

uðx; tÞ ¼ 1  ZHx Hx

xddV n ð27Þ

xdV n

where x is the influence function, meaning that bonds connected to the point have different weight for its damage. Damage parameter u is equal to 0 when bond connected to the point is not broken, while u is equal to 1 when all bonds to the point are broken. For crack propagation problems, crack grows with broken bonds, and the damage u of points along crack is kept changing. Thus, the variable parameter damage u can be used to track crack path. Unlike the bond-based peridynamic model, in the state-based peridynamic model, there is no real bond between material points. However, the relationship between pair points can be loosely regarded as the ‘bond’, and the bond stretch criterion still works. Furthermore, the energy-based damage model was introduced by Foster et al. [38]. The study to develop damage model with the peridynamic strain energy is needed in further research.

3. Results and discussions Validity of the proposed thermomechanical coupled peridynamic model is first demonstrated by considering a 2D bimaterial beam subjected to uniform temperature loading and a 3D plate subjected to a gradient temperature through thickness direction, both without any displacement constraints. The numerical results obtained by the finite element method (FEM) (which is only applied for the cases without discontinuity) via the commercial software ABAQUS are used to compare with and validate the peridynamic model prediction. Then, the damage growth of a 2D bimaterial beam due to incremental thermal loading is investigated. For quasi-static loading analysis, the adaptive dynamic relaxation (ADR) method [39] is utilized in this section to obtain the valid results. Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

8

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

3.1. 2D bimaterial beam under uniform temperature loading A 2D bimaterial beam under uniform temperature loading is considered to validate the proposed model for both the plane stress and plane strain problems. Similar to the model in [39], the geometrical sizes of the bimaterial beam are shown in Fig. 2, and the thickness of beam has the constant value of t = 0.1 mm. The beam without any displacement constraints and under the uniform temperature load of H ¼ 50  C is considered. Both the substrate materials are elastic and isotropic, and their material properties are given in Table 1. For numerical peridynamic simulation, the system is discretized into a finite number of nodes with uniform grid spacing Dx. m-convergence is considered to converge to the exact continuum peridynamic solution, where the horizon value d is fixed and the value m ¼ ðd=DxÞ ? 0. Meanwhile, d-convergence is considered to converge to the classical continuum mechanics solution, where m is fixed and the horizon d ? 0. The FEM results are obtained with the grid size value of 0.1 mm and same thermal loading condition. For the numerical FEM analysis, the symmetrical boundary condition is applied to avoid the rigid body deformation. The results captured from the peridynamic model and numerical finite element method are shown in Figs. 3–6. In the plane stress case, the profiles of displacement along interface of bimaterial with a fixed value of m = 4 and different values of horizon d decreasing from 0.8 mm to 0.2 mm are shown in Fig. 3; whereas the profiles of displacement along interface with a fixed value of d = 0.4 mm and different values of m ranging from 2 to 8 are illustrated in Fig. 4. Similarly, the profiles of displacement with the two convergences in the plane strain case are shown in Figs. 5 and 6, respectively. Since the coefficient of thermal expansion of material 1 is larger than that of material 2, the bimaterial bends as revealed by the displacement in y direction and the displacement in x direction increases linearly along interface as shown in Figs. 3–6. The maximum and minimum displacements in the plane strain case are larger than those in the plane stress case. For a fixed value of m = 4 and varying values of horizon d, the peridynamic results are very close to the FEM ones (with a maximum different of 2.0%), so does the m-convergence. The profiles almost coincide for different values of d and m, indicating that when only under thermal loading the deformations have little relations with respect to the grid size. Thus, the reasonable results can be obtained without very dense grid. As the thermomechanical coupled peridynamic model given above, the thermal load is applied directly as the variations of force state on material points. The effects of thermal loading can be obtained without gathering deformations in horizon. As a consequence, the horizon and grid sizes have minor effect on the peridynamic simulation results. The displacement contours of peridynamic model in the plane stress case are shown in Fig. 7. The values of m = 4 and d = 0.05 mm are used, and as shown in Fig. 7, the distribution trend of displacement is similar to that of [39], as well as the ones by the finite element analysis. As expected, the bimaterial bends symmetrically as revealed by the displacements contours in y direction as shown in Fig. 7(b).

3.2. 3D plate subjected to temperature gradient along thickness direction A 3D model of thick plate under temperature gradient is considered to verify the proposed peridynamic model under nonuniform thermal loading in 3D case. The 3D model in Fig. 8 is free of any displacement constraints, and the temperature increases linearly from H ¼ 0  C to

H ¼ 20  C along the z (thickness) direction, which can be expressed as: Hðx; y; zÞ ¼ 104 z. The material constants for ‘Material 1’ in Table 1 are considered. For 3D peridynamic simulation, the value of horizon d is set as 0.4 mm and m is 4. The same grid

Fig. 2. Bimaterial beam model.

Table 1 Material properties for cermet (metal-ceramic). Material number and type 1 (Al) 2 (Al2O3)

Elastic modulus (GPa) E

Poisson’s ratio v

Density (g/cm3)

71 340

0.22 0.31

q

Critical stress intensity factor (MPa) KIC

Thermal expansion coefficient (ppm/°C) b

2.70 3.92

30.0 4.2

8.5  106 23  106

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

9

10

5

Displacement, y component(× 10-3 mm)

Displacement, x component(× 10-3 mm)

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

δ = 0.2 mm δ = 0.4 mm δ = 0.05 mm FEM

0

-5

-10 -10

-5

0

5

10

8 4 0 -4

δ = 0.2 mm δ = 0.4 mm δ = 0.05 mm

-8

FEM

-12 -10

-5

0

x (mm)

x (mm)

(a)

(b)

5

10

10

5

Displacement, y component (×10-3 mm)

Displacement, x component (× 10-3 mm)

Fig. 3. Displacements along interface with different values of d and a fixed value of m = 4 in the plane stress case: (a) x component, and (b) y component.

m=2 m=4 m=8 FEM

0

-5

-10 -10

-5

0

5

10

8 4 0 -4 m=2 m=4 m=8 FEM

-8 -12 -10

-5

0

x (mm)

x (mm)

(a)

(b)

5

10

10

5

δ = 0.2 mm δ = 0.4 mm δ = 0.05 mm FEM

0

-5

-10 -10

-5

0

5

10

Displacement, y component (×10-3 mm)

Displacement, x component (×10-3 mm)

Fig. 4. Displacements along interface with different values of m and a fixed value of d = 0.4 mm in the plane stress case: (a) x component, and (b) y component.

10 5 0 -5 -10

δ = 0.2 mm δ = 0.4 mm δ = 0.05 mm

-15 -20 -10

FEM

-5

0

x (mm)

x (mm)

(a)

(b)

5

10

Fig. 5. Displacements along interface with different values of d and a fixed value of m = 4 in the plane strain case: (a) x component, and (b) y component.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

10

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Fig. 6. Displacements along interface with different values of m and a fixed value of d = 0.4 mm in the plane strain case: (a) x component, and (b) y component.

Fig. 7. Displacement distribution contours with m = 4 and d = 0.05 mm in the plane stress case: (a) x component, and (b) y component.

Fig. 8. 3D plate model.

size and thermal conditions are also used in the FEM model, and the symmetrical boundary conditions are applied to avoid rigid body deformation. The comparisons of displacement by the peridynamic and FEM models are shown in Fig. 9. As the coefficient of thermal expansion linearly increases along the z (thickness) direction, the plate bends as demonstrated by the z direction displace-

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

11

Fig. 9. Comparisons of displacement due to temperature gradient between the peridynamic model in (a) x, (c) y, (e) z directions, and FEM model in (b) x, (d) y, (f) z directions.

ment shown in Fig. 9(e). Also, the displacement by the peridynamic model shows consistent and close correlations with the finite element results. 3.3. Interface crack propagation of bimaterial beams under incremental thermal loading After validating accuracy of the proposed peridynamic model under thermal loading, the damage prediction of bimaterial beams subjected to incremental thermal loading is investigated. Based on the verification analysis, the value of m = 4 and d = 0.4 mm is chosen to obtain reliable results. The geometrical sizes and material parameters are the same as the ones shown in Fig. 2. But the temperature load increases 5.0  103 °C at each iteration rather than a constant temperature. The plane stress model is used. The initiation and growth of damage contours at several incremental iterations representing interface crack propagation are shown in Fig. 10. As shown in Fig. 10, symmetrical cracks initiate at interface at both the end edges, and they then grow along the bimaterial interface toward the middle of the beam abruptly. The failure of interface begins with small damage u which is defined in Eq. (27) and grows to full interface cracking (a complete separation) with a final damage value of 0.42. Unlike the model in [39], the material points interact with each other through ‘state’ rather than ‘bond’ in the state-based peridynamic model. Thus, in the given interface debonding analysis, the damage u of interface points increases progressively with respect to iterations rather than accomplished via the crack propagation along the bimaterial interface. 4. Three-point bending problem First, the convergence study of bimaterial beams under three-point bending is evaluate to validate effectiveness of the proposed peridynamic model for simulating bending problems. Then, the damage growth (interface crack propagation) prediction of bimaterial due to quasi-static and impact loading is presented. The plane stress model is considered in the analysis. Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

12

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Fig. 10. Interface crack propagation path of bimaterial, illustrating the damage growth, due to incremental thermal loading at iteration steps of (a) 25,000, (b) 27,000, (c) 28,000, (d) 30,000, and (e) 36,000.

4.1. Bimaterial beam under three-point bending The geometry of bimaterial beam is shown in Fig. 11, and the material constants are the same as listed in Table 1. Under three-point bending, the components of displacement in x and y directions of point (10 mm, 1 mm) in Fig. 11 are set to 0. Similarly, the displacement in y direction of point (10 mm, 1 mm) is set to 0, where those points are material points in peridynamics and the grid nodes in FEM as shown in Fig. 11. The constant normal stress rs = 10 MPa is applied over a length of 1 mm (see Fig. 11). In the peridynamic model, the line force is loaded as the body force equivalently on material points along the loading line. The boundary and loading areas are set as no-damage zones to avoid the undesired damage because of local effect due to boundary and loading conditions (i.e., the displacement constraints at the boundary and location of concentrated applied load). The ADR method is used for quasi-static analysis. The y-direction displacement profiles along the bimaterial interface with respect to m-convergence and d-convergence are shown in Figs. 12 and 13, respectively. As shown in Fig. 12, the displacement profile along bimaterial interface with a fixed value of m = 4 approaches that of FEM when d becomes smaller. While for a fixed value of d = 0.4 mm, the displacement profile along bimaterial interface converges to that of FEM as the m value increases (Fig. 13). The downward y-direction displacement profiles by the peridynamic model match relatively closely to the finite element results with the decreased d and increased m values with a maximum difference of 2.2%, thus indicating the convergence of peridynamic model to classical continuum mechanics. However, unlike the thermal loading case, the dense grid size is needed to obtain reasonably accurate displacement prediction under three-point bending problems. Since the peridynamics is nonlocal theory and external force is loaded as the body force on partial material points, the effect of localized force transmits among the whole system through interactions between material points as the force state evolves. However, the cost of computation increases rapidly with the decreasing d and increasing m values. Thus, the suitable grid sizes need to be chosen to obtain acceptable predictions. The displacement contours in x- and y-directions predicted by the peridynamic model with d = 0.05 mm and m = 4 are shown in Fig. 14. 4.2. Damage growth of bimaterial beam under incremental quasi-static loading Damage growth, in term of crack propagation, of bimaterial beams under quasi-static incremental loading in the plane stress case is presented. Based on the convergence study above, the values of d = 0.2 mm and m = 4 are chosen.

Fig. 11. Bimaterial beam model under three-point bending.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

Displacement, y component (×10-3 mm)

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

13

0 δ = 0.2 mm δ = 0.4 mm δ = 0.05 mm

-2

FEM

-4

-6 -10

-5

0

5

10

x (mm)

Displacement, y component (×10-3 mm)

Fig. 12. The y-direction displacement along the bimaterial interface under three-point bending with different d values and a fixed value of m = 4.

0 m=2 m=4 m=8 FEM

-2

-4

-6 -10

-5

0

x (mm)

5

10

Fig. 13. The y-direction displacement along the bimaterial interface under three-point bending with different m values and a fixed value of d = 0.4 mm.

Fig. 14. Displacement contours under three-point bending with m = 4 and d = 0.05 mm: (a) x component, and (b) y component.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

14

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

Fig. 15. Bimaterial beam model with a pre-crack under incremental quasi-static loading.

The same bimaterial beam model but with a transverse pre-crack in the substrate material (Material 2) at the central span is considered (see Fig. 15), and the normal stress load of 1.0  102 MPa is applied at each iteration rather than directly subjected to a constant value of stress load. For quasi-static analysis, the ADR method is used to induce loading to the beam gradually. The patterns of crack propagation at several typical iterations representing different extents of crack propagation within the substrate material and along the bimaterial interface are shown in Fig. 16. As expected, the crack initiates from the tip of pre-notch at the iteration step of 7000 when the normal stress load reaches 70 MPa, and it propagates upward and rapidly along the y-direction within the substrate material (Material 2). Then, crack propagation ceases at interface at the iteration step of 8000. When the applied load increases to 120 MPa, the bifurcating crack initiates to extend again along the bimaterial interface toward the supporting boundary symmetrically. Since the critical stretch of Material 1 is much larger than that of Material 2 due to material properties of cermet as defined in Table 1, the transverse crack ceases at the interface and then propagates along the bimaterial interface rather than keeping propagating into Material 1. Moreover, since the bending moment decreases from the center of beam toward two boundary ends, it is easier for crack to grow along vertical direction at the central span where there is a maximum bending moment at the early time, and it is much harder to propagate along the bimaterial interface as crack approaches two boundary ends. 4.3. Damage growth of bimaterial beams under impact loading The damage growth and crack propagation of bimaterial beams under impact loading is also analyzed by the proposed state-based peridynamic model. The values of peridynamic numerical parameters of m = 4 and d = 0.2 mm are again used. The sample model under impact loading is presented in Fig. 17, and the constant normal stress rs = 100 MPa is applied

Fig. 16. Crack propagation path of bimaterial beam due to incremental quasi-static loading at iterations of (a) 7000, (b) 8000, (c) 12,000, and (d) 26,000.

Fig. 17. Bimaterial model with a transverse pre-crack under impact loading.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

15

Fig. 18. Crack propagation path of bimaterial due to impact loading at time of (a) 9.5 ls, (b) 10 ls, (c) 11.25 ls, (d) 12.5 ls, and (e) 14 ls.

but at different location compared to the last example of quasi-static case in Fig. 15. Instead of the ADR method, the explicit time integration Velocity-Verlet algorithm [21] is used for dynamic analysis, and the stable time step of 2.5 ns (nanosecond) is considered. A transverse pre-crack in Material 2 is introduced at the central span as shown in Fig. 17. The boundary and loading areas are set to no-damage zones to avoid undesired damage due to local concentration. Crack propagations with progressive time steps are shown in Fig. 18. The crack initiates at the pre-crack tip at 9.5 ls (microsecond) and propagates within Material 2 toward the bimaterial interface, and it then grows with a curvilinear pattern as shown in Fig. 18(c) because of lack of symmetry of loading (see Fig. 17). Then, the crack ceases to extend upward at the bimaterial interface at 11.25 ls, and it propagates along the interface toward the loading zone. The elastic strain energy density distributions with the progressive time steps are shown in Fig. 19. The actual deformed shape and the crack pattern at 50 ls are illustrated in Fig. 20. The actual deformation of the upper layer material (Material 1) is quite large. However, the crack still propagates along the bimaterial interface because of large deformation capacity (ductility) of metal (Material 1).

Fig. 19. Strain energy density maps of bimaterial beam under impact loading at time of (a) 9.5 ls, (b) 10 ls, (c) 11.25 ls, (d) 12.5 ls, and (e) 14 ls (unit: J/m3).

Fig. 20. Crack propagation path and deformed shape of bimaterial beam under impact loading at time of 50 ls.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

16

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

As shown in Fig. 19(a), the strain energy density is first concentrated at the locations of impact loading zone, bimaterial interface (due to material mismatch), free boundary blow loading zone, and the crack tip. As the crack propagates, the strain energy density of crack tip becomes larger and the concentrated point moves with the crack tip (see Fig. 19(b) and (c)). Moreover, when the crack arrives at the bimaterial interface and then propagates along the interface, the concentrated energy at the free boundary blow impact zone is reduced (see Fig. 19(d)) and eventually fully released (see Fig. 19(e)); while the strain energy density of free boundary near the crack propagation path remains concentrated during crack propagation. At the bimaterial interface, though the displacement across the interface is continuous, the strain energy density is concentrated because of the mismatch of material properties. Thus, the damage can easily happen at the bimaterial interface. 5. Conclusions In this study, the ordinary state-based peridynamic model is extended to predict damage growth and crack propagation in bimaterial (e.g., cermet) structures under the thermomechanical loading, and the three-dimensional (3D) and twodimensional (2D) (plane stress and strain) cases are considered. A bimaterial beam in 2D and a thick plate in 3D under thermal loading and three-point bending are simulated by the proposed state-based peridynamic model, and the accuracy is demonstrated by the close correlations of displacement obtained by the present peridynamic model and the numerical finite element model. The damage growth of bimaterial beams under thermal loading and three-point bending is also investigated, and the interface crack propagation at different loading steps (iterations) is observed. The dense grid size is not needed in problems when only considering uniform thermal loading, but it is required when under three-point bending. In the case of interface cracking of bimaterial due to thermal loading, symmetrical cracks initiate at the interface near the end edges and propagate along the interface toward the central span. The damage variable u of interface points grows progressively along with iterations. For bimaterial interface crack under three point bending, in case of incremental quasi-static loading, the crack initiates from the tip of pre-notch and propagates upward within the substrate material along the y-direction abruptly, and it ceases at the bimaterial interface and then propagates along the bimaterial interface toward the end edges symmetrically. In case of impact loading, the crack initiates at the tip of the pre-crack and propagates with a curvilinear pattern because of lack of symmetry, and it ceases to extend upward at the interface and then propagates along the interface toward the impact loading area. Furthermore, the concentrated zones of strain energy density moves with the crack tip, and the strain energy density is concentrated at the bimaterial interface due to mismatch of material properties. In summary, the new ordinary state-based peridynamic model proposed is capable of successfully capturing the damage growth at bimaterial interface under thermal loading and predicting the crack propagation within substrate material and then along bimaterial interface under three-point bending. The crack grows naturally with respect to time steps without any special fracture criteria, demonstrating capability of the proposed model in predicting damage growth and crack propagation path in bimaterial structures under thermomechanical loading. Acknowledgements The authors would like to thank for the partial financial support from the National Natural Science Foundation of China (NSFC Grant Nos. 51478265 and 51679136) to this study. Appendix A. Supplementary material Supplementary data associated with this article can be found, in the online version, at https://doi.org/10.1016/j.engfracmech.2017.09.023. References [1] Ikeda T, Sun CT. Stress intensity factor analysis for an interface crack between dissimilar isotropic materials under thermal stress. Int J Fract 2001;111:229–49. [2] Khandelwal R, Chandra Kishen JM. The use of conservative integral in bi-material interface crack problems subjected to thermal loads. Int J Solids Struct 2008;45:2976–92. [3] Boutabout B, Chama M, Bouiadjra BAB, Serier B, Lousdad A. Effect of thermomechanical loads on the propagation of crack near the interface brittle/ductile. Comput Mater Sci 2009;46:906–11. [4] Drai A, Bouiadjra BB, Meddah M, Benguediab M. Analysis of interfacial fracture in ceramic-metal assemblies under effect of thermal residual stresses. Comput Mater Sci 2009;46:1119–23. [5] Pathak H, Singh A, Singh IV. Fatigue crack growth simulations of bi-material interfacial cracks under thermo-elastic loading by extended finite element method. Eur J Comput Mech / Rev Eur Mécanique Numérique 2013;22:79–104. [6] Keppas LK, Anifantis NK. Boundary element analysis of cracked homogeneous or bi-material structures under thermo-mechanical cycling. Comput Methods Appl Mech Eng 2010;199:3345–55. [7] Silling SA. Reformulation of elasticity theory for discontinuities and long-range forces. J Mech Phys Solids 2000;48:175–209. [8] Gerstle W, Sau N, Silling S. Peridynamic modeling of plain and reinforced concrete structures. In: SMiRT18. Int conf structural mech reactor technol; 2005. [9] Gerstle W, Sau N, Silling S. Peridynamic modeling of concrete structures. Nucl Eng Des 2007;237:1250–8. [10] Liu W, Hong JW. Discretized peridynamics for linear elastic solids. Comput Mech 2012;50:579–90.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023

H. Zhang, P. Qiao / Engineering Fracture Mechanics xxx (2017) xxx–xxx

17

[11] Silling SA, Epton M, Weckner O, Xu J, Askari E. Peridynamic states and constitutive modeling. J Elast 2007;88:151–84. [12] Silling SA. Linearized theory of peridynamic states. J Elast 2010;99:85–111. [13] Le QV, Chan WK, Schwartz J. A two-dimensional ordinary, state-based peridynamic model for linearly elastic solids. Int J Numer Methods Eng 2014;98:547–61. [14] Sarego G, Le QV, Bobaru F, Zaccariotto M, Galvanetto U. Linearized state-based peridynamics for 2-D problems. Int J Numer Methods Eng 2016;108:1174–97. [15] Madenci E, Oterkus S. Ordinary state-based peridynamics for plastic deformation according to von Mises yield criteria with isotropic hardening. J Mech Phys Solids 2016;86:192–219. [16] Mitchell JA. A nonlocal, ordinary, state-based plasticity model for peridynamics. Sandia National Laboratories, Albuquerque SAND2011-3166; 2011. [17] Mitchell JA. A non-local, ordinary, state-based viscoelasticity model for peridynamics. Sandia National Laboratories, Albuquerque SAND2011-8064; 2011. [18] Foster JT, Silling SA, Chen WW. Viscoplasticity using peridynamics. Int J Numer Methods Eng 2009;81:1242–58. [19] O’Grady J, Foster J. Peridynamic beams: A non-ordinary, state-based model. Int J Solids Struct 2014;51:3177–83. [20] O’Grady J, Foster J. Peridynamic plates and flat shells: A non-ordinary, state-based model. Int J Solids Struct 2014;51:4572–9. [21] Silling SA, Askari E. A meshfree method based on the peridynamic model of solid mechanics. Comput Struct 2005;83:1526–35. [22] Askari E, Xu J, Silling S. Peridynamic analysis of damage and failure in composites. In: 44th AIAA/ASME/ASCE/AHS/ASC aerospace sciences meeting and exhibit, 9–12 January 2006, Reno, NV. AIAA 2006–88. [23] Xu J, Askari A, Weckner O, Razi H, Silling S. Damage and failure analysis of composite laminates under biaxial loads. In: 48th AIAA/ASME/ASCE/AHS/ ASC structures, structural dynamics, and materials conference, 23–26 April 2007, Honolulu, HI. AIAA 2007–2315. [24] Colavito KW, Kilic B, Celik E, Madenci E, Askari E, Silling S. Effects of nanoparticles on stiffness and impact strength of composites. In: 48th AIAA/ ASME/ ASCE/AHS/ASC structures, structural dynamics, and materials conference, 23–26 April 2007, Honolulu, HI. AIAA 2007–2021. [25] Ha YD, Bobaru F. Characteristics of dynamic brittle fracture captured with peridynamics. Eng Fract Mech 2011;78:1156–68. [26] Ha YD, Bobaru F. Studies of dynamic crack propagation and crack branching with peridynamics. Int J Fract 2010;162:229–44. [27] Oterkus E, Madenci E, Weckner O, Silling S, Bogert P, Tessler A. Combined finite element and peridynamic analyses for predicting failure in a stiffened composite curved panel with a central slot. Compos Struct 2012;94:839–50. [28] Agwai A, Guven I, Madenci E. Crack propagation in multilayer thin-film structures of electronic packages using the peridynamic theory. Microelectron Reliab 2011;51:2298–305. [29] Zhou XP, Wang YT. Numerical simulation of crack propagation and coalescence in pre-cracked rock-like Brazilian disks using the non-ordinary statebased peridynamics. Int J Rock Mech Min Sci 2016;89:235–49. [30] Wang Y, Zhou X, Xu X. Numerical simulation of propagation and coalescence of flaws in rock materials under compressive loads using the extended non-ordinary state-based peridynamics. Eng Fract Mech 2016;163:248–73. [31] Kilic B, Madenci E. Prediction of crack paths in a quenched glass plate by using peridynamic theory. Int J Fract 2009;156:165–77. [32] Bobaru F, Duangpanya M. A peridynamic formulation for transient heat conduction in bodies with evolving discontinuities. J Comput Phys 2012;231:2764–85. [33] Bobaru F, Duangpanya M. The peridynamic formulation for transient heat conduction. Int J Heat Mass Transf 2010;53:4047–59. [34] Oterkus S, Madenci E, Agwai A. Fully coupled peridynamic thermomechanics. J Mech Phys Solids 2014;64:1–23. [35] Bobaru F, Yang M, Alves LF, Silling SA, Askari E, Xu J. Convergence, adaptive refinement, and scaling in 1D peridynamics. Int J Numer Methods Eng 2009;77:852–77. [36] Parkus H. Thermoelasticity. Wien New York: Springer-Verlag; 1976. [37] Bobaru F, Foster JT, Geubelle PH, Silling S. Handbook peridynamic modeling. CRC; 2017. [38] Foster JT, Silling SA, Chen W. An energy based failure criterion for use with peridynamic states. Int J Multiscale Comput Eng 2011;9:675–87. [39] Kilic B, Madenci E. An adaptive dynamic relaxation method for quasi-static simulations using the peridynamic theory. Theor Appl Fract Mech 2010;53:194–204.

Please cite this article in press as: Zhang H, Qiao P. An extended state-based peridynamic model for damage growth prediction of bimaterial structures under thermomechanical loading. Engng Fract Mech (2017), https://doi.org/10.1016/j.engfracmech.2017.09.023