CC BY
62
DOI: 10.17277/amt.2017.01.pp.020-035
Innovative Modeling of the Crack Path and Stress Intensity Factor for Arbitrary Shaft Configurations
R. Dimitri1, Y. Li2, N. Fantuzzi3, F. Tornabene3*
1 Department of Innovation Engineering, University of Salento, Lecce, Italy
2 State Key Laboratory of Coal Mine Disaster Dynamics and Control, Chongqing University,
Chongqing, China 3 DICAM-Department University of Bologna, 2, Via del Risorgimento, Bologna, 40136, Italy
* Corresponding author: Tel.: +39 051 2090312; E-mail: [email protected]
Abstract
The lifetime of most engineering structures and components is known to depend on the presence of defects, such as holes, cracks or voids usually introduced during a manufacturing process. In many cases, the crack growth, extension and propagation within a body, still remains a challenging problem in fracture mechanics. The present paper proposes an extended analytical model based on a section method to predict the fracture direction and compute the stress intensity factors for a cracked shaft under mixed-mode loading conditions. The advantage of the present formulation is mainly related to its capability of predicting the direction of crack propagation within a shaft under coupled longitudinal, flexural and torsional loading conditions. The analytical results are straightforwardly compared with the theoretical expressions available from the handbooks and the numerical solutions found with the extended finite element method. The present approach agrees quite well with the theoretical and numerical results already proposed in the literature, thus confirming its potentials for accurate computations of the crack propagation and stress intensity factors for arbitrary configurations.
Keywords
Crack propagation; cracked shafts; extended finite element method; stress intensity factor.
© R. Dimitri, Y. Li, Ni. Fantuzzi, F. Tornabene, 2017
Introduction
The accurate modeling of 3D cracks in continuum bodies remains a challenging problem in computational mechanics. The importance of computing the fracture parameters and simulating the crack growth in engineering components, stems from the widespread use of numerical fracture mechanics for static and fatigue problems. Fatigue or static failure usually occurs due to the initiation and propagation of single or multiple cracks with different shapes of multiphysical nature along surfaces or nearby them [1, 2]. Despite the large effort by the scientific community for the study of the main factors controlling the crack paths within materials, the results are still not satisfactory. The works by Pook [3, 4] provide an interesting overview
of the active research on several crack path topics of current interest, e.g. fatigue cracking in aircraft structures, fracture toughness testing, fractography, mixed mode fatigue thresholds, crack path stability, failure of laminated and generally jointed structures. A linear elastic study of cracked bodies, crack paths and a collection of some realistic cases can be also found in Pook [5] to study the fatigue crack paths within metallic materials. A large variety of crack path topics, however, makes a useful generalization of the problem within standards and codes quite difficult.
Within a linear elastic fracture mechanics (LEFM) framework, the propagation of a crack is mainly governed by the singular behavior of the stress field in the vicinity of its crack tip. For a 2D quasi-static crack, this behavior is generally governed by
at] (r, 9) = -^ Tj (9 )+K= j)^0), (1) V2nr V2nr
where £ (9) and £^(9) are universal functions
describing the angular variation of the stress field, and KI and Kn are the stress intensity factors (SIFs). Thus, except for the cases where higher order terms are required to characterize the crack paths, SIFs usually play an important role for a correct characterization of the mechanical properties of cracked materials. The state of art for the stress analysis of cracks, based on the fracture mechanics concept of the SIF, was already reviewed in 1965 by Paris and Sih [6]. Many 2D theoretical solutions of the SIF available from the literature still represent basic formulations for many engineering applications, design procedures, standards and failure assessment codes [7-10].
However, the simplified analyses based on 2D theories, are largely recognized to ignore the 3D effects, including the effects of the plate thickness [11], vertex singularities [12], and coupling of fracture modes I, II and III on the stress and strain fields near the crack front [13-15]. Thus, an accurate computation of the SIFs for 3D surface cracks would be required for practical applications, which is based on the implicit assumption of continuity along a crack front. The understanding of 3D effects has long been recognized as the most important and challenging problem in fracture mechanics, and the increasingly powerful computers have made the numerical investigation possible, especially in the vicinity of corner points. In the last decades, a large amount of theoretical and experimental information about the 3D effects has been discussed by several researchers, as it can be seen in the overview provided by Pook [16], or other recent works from the literature [17-21].
There is an important class of crack problems, which involve the surface damage in the form of part-through cracks in round bars or beams, where a 3D analysis is involved. Miao et al. [22] analyzed numerically 3D mixed mode central cracked plates under biaxial loading conditions, by including the elastic and elastic-plastic ./-integrals. Valiente [23] obtained numerically the solutions for the sample compliance and the energy release rate based on the finite element approach. The stiffness derivative technique was applied by the same author to determine the energy release rate, based on the virtual crack extension. Levan and Royer [24] and Couroneau and Royer [25] applied the boundary integral equation method to determine the SIFs for round bars with transverse circular cracks. The finite element approach was then applied by Daoud and Cartwright [26, 27] to compute the SIFs for a circular bar with a straight crack under tension and bending. Ng and Fenner [28]
built a 3D finite element model to treat fracture mechanics and to determine the SIF for a cracked circular bar with a straight and circular crack fronts on one side, when subjected to an axial and flexural loading.
Shin and Cai [29] compared some experimental and FEM-based results about the SIFs for an elliptical crack within a circular shaft under tension and bending. Carpinteri [30] computed numerically the dimensionless SIFs for a straight-front and semi-elliptical crack by using a 3D finite element modeling. Erjian [31] described the normalized SIFs as a function of the crack depth for round bars under axial forces and bending moments. Naik and Maiti [32] presented the formulation for analyzing the triply coupled free vibration of beams based on a compliance approach in the open edge crack in an arbitrary angular orientation. Some reviews of the SIF solutions for surface cracks in round bars subjected to tension loading can be found in the literature [33, 34]. Due to the complexity of cracked problems, a large variety of numerical techniques have been employed in the literature to determine 3D SIFs nearby a crack. In this context, the finite element method (FEM) is the most commonly used procedure for computing the SIFs due to the progressive increase in the computational power and FEM-based codes available for computational studies. Since 1969, Dixon and Pook [35] announced the adoption of the FEM as a useful tool to predict the SIFs in the study of opening mode crack growth problems where the nominal stress level is known to be lower than the yield stress of a material. Chen and Wang [36] determined the mixed mode geometry factors for a finite plate with an inclined edge crack by means of a weight function approach. This was also adopted by the same authors to calculate the mixed-mode SIFs and T-stress for an offset double edge-cracked plate. Jogdand and Murthy [37] proposed an innovative method for estimating the mixed-mode SIFs based on the formation of over-determined system of equations, which was then implemented in existing finite element codes. A FEM-based dislocation model was also proposed by Feng et al.
[38] in order to determine the SIFs for an arbitrary number of through-thickness cracks. Whitcomb et al.
[39] and Chaudhuri [40] commented on the accuracy of FEM results, based on a displacement potential energy approach, in the neighborhood of stress discontinuities and stress singularities. As mentioned in [39, 40], the numerical distribution of the stress and/or displacement as obtained with FEM for cracked specimens usually loses the accuracy in the closest region around the discontinuity or singularity.
The SIFs determined from local displacements or stresses, require some unconventional crack tip elements, also known as quarter point elements, which
account for the i/Vr singular stress behavior of the crack tip. A quarter-point element is obtained from the standard quadratic element by simply moving the mid-node coordinates three-fourths of the side towards the crack tip. It is well known from the literature that the solution can converge quite slowly for cracked problems, when adopting the only conventional elements, which do not include stress singularities. It seems that the most common approach for modeling 3D cracks using the FEM is to introduce the singularity stress behavior in the crack tip elements by relocating the element mid-side nodes to new locations. This introduces a singularity in the Jacobian inverse of the geometric transformation. Closed-form solutions for the SIFs are available for simple crack geometries in 3D. However, if the crack surface is not aligned with the element boundaries, the displacement discontinuity and the traction conditions on the crack surface cannot be treated by standard FEMs. In the actual cracked structures, the numerical prediction of the fracturing process is very hard to perform, especially for arbitrarily-shaped cracks in finite specimens. Based on the above mentioned numerical limitations of FEM when treating physical models in presence of cracks, some alternative numerical strategies have been recently developed by the scientific community to analyze the fracture process of elastic structures. To this end, the Cell Method [41, 42] the Embedded Mesh formulation [43, 44], the Generalized Differential Quadrature (GDQ) Method [45-48], the Isogeometric Method [49-52], or the extended FEM (XFEM) [53-55] have been proposed as computationally efficient techniques to treat the problem. In this study, we have employed the XFEM to treat stress singularities in fracture mechanics for comparative purposes, as introduced by Belytschko and Black [56]. This approach enables a local enrichment of the discrete space, based on the concept of partition of unity. This results to be very useful for the approximation of solutions with a pronounced non-smooth behavior in some parts of the computational domain (i.e. near discontinuities or singularities), where a standard FEM could be insufficient to offer good quality results. The XFEM offers significant advantages by enabling optimal convergence for these applications, since it alleviates the main shortcomings of a classical FEM approach associated with meshing of the cracked surface. The FEM would represent the essence of the XFEM, and hence a lot of the theoretical and numerical developments in finite elements can be readily extended and applied. The XFEM-based approach splits a modeling problem in two distinct parts: the mesh generation for the geometric domain (without considering any possible crack), and the enrichment step of the discrete
approximation by means of additional functions that model a geometric discontinuity. Anyway, the fracture criterion is the primary basis for the safe evaluation in cracked materials. A theoretically integrated fracture criterion should be capable to predict the fracture direction and loading of cracked structures under varying loading conditions.
In this sense, different fracture criteria were developed in the past decades to describe a mixed-mode fracturing. Erdogan and Sih [57] proposed the maximum circumferential stress criterion, which assumed that fracture occurs in the direction where the circumferential stress surrounding the crack tip reaches the maximum value. Smith et al. [58] developed a generalized maximum tensile stress criterion based on the T-stress. A generalized maximum tangential stress criterion was employed in Ayatollahi et al. [59-61] to predict the fracture initiation angles, which is strongly dependent on the elastic T-stress, under mixed mode loading conditions in different brittle or composite materials. Sih [62, 63] introduced the minimum strain energy density factor criterion which assumes that fracture occurs in the direction where the strain energy density reaches the minimum value. A similar fracture criterion was also proposed by Richard et al. [64] to treat 3D-mixed mode problems, which revealed to be suitable for the description of the crack growth, in agreement with [62].
Additional bases and applications of the averaged strain energy density criterion for the computation of notch SIFs for coarse meshes, are detailed in Lazzarin et al. [65, 66]. The T-stress effects in the application of a conventional minimum strain energy density criterion has been also recently discussed by Ayatollahi et al. [67], by considering higher order terms in the William's series expansion. A clear influence of the magnitude and sign of the T-stress on the crack trajectory and the fracture strength of different specimens is highlighted by the authors in their work. Hussain et al. [68] proposed the maximum energy release rate criterion based on Griffith's theory [69]. Cornetti et al. [70] introduced a fracture criterion based on a double enforcement of the stress and energy balance, named as finite fracture mechanics, to describe the finite crack propagation and advancement, recently applied for mode-I loaded specimens in Dimitri et al. [71].
A criterion for 3D cracking under mixed-mode loading conditions, was additionally described by Schollmann et al. [72]. The criterion allows for the prediction of the 3D crack advancement with the help of two different deflection angles. Similarly, Pook [73] proposed a 3D criterion based on the definition of two crack deflection angles. Chang et al. [74] presented
a general mixed-mode brittle fracture criterion based on the concept of the maximum potential energy release rate. Additional criteria of fracturing are reviewed and discussed in [75, 76] together with their advantages and drawbacks for different applications. For an engineering estimation of the maximum strength for cracked bodies, the use of efficient procedures involving a limited computational effort can be suitable, also for lower levels of accuracy. The present paper introduces a simple engineering method to calculate the SIF, based on the section method, as introduced in [77] and modified by Nobile [78, 79]. The section method considers the stress singularity at the tip of an elastic crack, and is based on an elementary beam theory equilibrium condition for internal forces evaluated in the cross-section passing through the crack tip. The derived simple formulation for SIFs of transversely cracked circular section beams subjected to bending and torsional moment, as well as to the shear force, are in a reasonable agreement with the results from the available literature. Moreover, three different criteria as proposed by Richard [64], by Sih [62, 63], and by Schollmann et al. [72] are used here to determine the crack direction of propagation. The results, both in terms of SIFs and crack angles of propagation, are compared with the available results from the literature as well as with results obtained with the XFEM as implemented in the ABAQUS code.
Analytical Modeling
Problem Statement. Let us consider a 3D structure under mixed-mode loading condition, which is made of two perpendicular circular beams subjected to a vertical concentration load F at the free ending. The lengths of the two beams are L and L/2, respectively, while the crack is located in the middle of the first beam of length L. The structure has a circular section with radius R and a crack length a (see Fig. 1). The cracked section is subjected to a shear force V = F, a bending moment M = Fl/2, and a torsional moment T = Fl/2. For this structure, we have implemented the section method as proposed
Fig. 2. Rectangular cross section and stress distribution
by Parton and Morozov [77] in order to compute the SIF. Based on this method, the crack tip is subjected to an additional force related to an increase in the stress expressed as
rt
Act = f Jo
ct
dr,
(2)
Fig. 1. The structural problem
where b is the length, cs is the nominal stress, and Ac is an additional force per unit length, yielding to the equivalence between the areas Aabdc and Acfik/ (see Fig. 2, a).
Nobile [78, 79] introduced a slight variation in the section method, i.e. the additional stress Ac due to the cracked section was computed through a balance between the equivalent force of the areas ACFIG and ACEIG (see Fig. 2, a). In other words, b is the distance where the singular stress is equal to the normal stress when the neutral axis moves along CG (see Fig. 2, a). As clearly seen from Fig. 2, a, however, the force associated with the area ACEIG is lower than the force associated with the area ACFIG. This would make the assumption by [78, 79] almost inaccurate. Thus, in this work we propose an improved
version of the section method, starting from the original section method by Parton and Morozov [77] and modified by Nobile [78, 79]. Based on our improved version, the additional force due to the stress concentration at the crack tip is computed by comparing the equivalent forces associated to the
areas
Ac
D E 'F 'I'H'G '
and A
A 'BD HG C
in
Fig. 2, b.
The distance b can be determined from the following balance condition (see Fig. 2, b)
(3)
G
s I r=b
-h = G,
y '=y'
where y is the distance of the reduced cross section at the crack tip from the neutral axis.
Rectangular Beam under a Bending Moment. As reference example, we start the work considering a beam with a rectangular cross section of thickness t = 1, height h, and a straight through transverse crack of length a. The cross section is subjected to a bending moment M, as seen in Fig. 2, b, with a linear distribution of the normal stress expressed as
M
G =-y,
I
(4)
where I = th3/12 is the second moment of inertia for the uncracked cross section.
Similarly, the reduced portion of the cross section due to the presence of the crack, features the following linear distribution
M ,
Gc = — y , Ic
(5)
where Ic = t(h - a')3/12 is the corresponding second moment of inertia, y' = y + e is the general distance of a point in the system of coordinates O' x' y', e = OO' = a '/2 is the distance of the neutral axis from the x-axis. It is well known from the literature that the stress component as in mode-I condition near a singular corner point or discontinuity, is related to the SIF, defined as
KT
s V2nr
(6)
where r is the general distance from the crack tip.
The equilibrium of the forces in the horizontal direction reads
h/2 b
JGdy =|g sdr,
(7)
h ' b —a -b 2
where the exact value of the distance b can be found by comparing Eqs. (5) and (6), i.e.
Kr
M-
•J2%b Icy'
(8)
After some algebraic manipulations, the mode-I SIF can be expressed as follows
6M
Kl = th32
F (n),
(9)
where n = a'/h is the dimensionless crack depth, whereas the geometry function is defined as
V2np
F (n) =
(1 -n)2
(10)
with
ß=
1
1
—n— i
2 (i-n)2 H(i-n)4 (i-n)2 4
1 ■ 2n 1 ■1 for 0.5-n-ß> 0
1 - __
2 n+(1-n)2 V(1-n)4 (1-n)2 4
2n-1 ■ I elsewhere.
(11)
The geometry function F (n) defined by the Eqs. (10) and (11) can be compared with the corresponding expressions given by Brown [80], Tada [8] and Nobile [78]:
- least square fitting method, precision of 0.2 % for n< 0.6 [80]
F(n) = -\/nn(1.122- 1.40n + 7.33n2 - 13.08n3 + 14.0n4);
(12)
- precision better than 0.5 % for any n [8]
F (n)^ 2 tan ^
0.923 + 0.1991 1 - sin^
cos
nn
2
-; (13)
- modified section method [78] 0.482
F (n) =
VÓ-ñ3'
(14)
The main analytical results in terms of geometry functions are plotted comparatively in Fig. 3, where it is worth noticing the good agreement between our
0.2 0.3 0.4 0.5
C-a: k dep:h ratio T]
Fig. 3. Geometric function for rectangular beam
4
2
proposed formulation and the results by Brown [80], Tada [8] and Nobile [78], especially for n- 0.6. A relative difference of few units percent can be generally noticed with respect to the predictions by Brown [80] or Tada [8], with larger differences with respect to the predictions by Nobile [78].
Circular Beam. The previous analysis is now repeated for a circular beam with radius R and precrack length a (see Fig. 4, a) subjected to a bending moment, shear force, and torsional moment, alternatively. Due to the symmetry of the problem, we analyze only half of the geometry in what follows.
Bending Moment. The cracked circular section is here assumed to be divided into independent rectangular strips of width t, which do not exchange any internal stress along their interface. This approximation would be less appropriate near the ends of the crack tip. Each strip has a general height h and crack depth a' which depend on the general distance s of the strip from the y-axis, as follows:
h
= 2>/R 2 -,!
' = a = — (R - a).
(15)
2
The distribution of the normal stresses for the uncracked cross section is defined by Eq. (4), where I is the second moment of inertia for the uncracked circular section. Similarly, the distribution of the normal stresses for the reduced cross section passing through the crack tip is defined by Eq. (5), where Ic is the second moment of inertia for the reduced portion of the circular cross section.
The singular stress component is related to the mode-I SIF as in the Eq. (6), whereas the equilibrium of the forces along the longitudinal direction of the beam (z-axis) is expressed as
h/2 b Jcdy =jcsdr. 0
(16)
h u —a—b 2
The value of b in the Eq. (16) is found from the equilibrium
Kj = m .
V2nb Ic
y'=y=h/ 2—a+e'
(17)
while the combination of Eqs. (16) and (17) yields to the following expression for the SIF
M
Kj = R 2.5
F (n),
(18)
where n = a/h is the crack depth ratio, and F(n) is defined as
f (n)=^^P^er R f (1 -n),
2Ic
with
n = 1 + V 2 D — 1
1 —
i \2 ( s
V h1 J
1—11—2 D
2
(19)
(20)
and
R
h
2S- -1
h
21 2 D -1
(21)
The geometry function F(n)is plotted in Fig. 4, b vs. the dimensionless ratio s / h1 for different lengths
of the precrack a/D. As predicted in the literature [8], the value of F (n) decreases for increasing distances from the y-axis and decreasing values of the dimensionless ratio a/D. The analytical results based on our proposed formulation agree quite well with
a/D 0.2-Tada
afll=0.4-Tada
ajl) ().6-'l'ada
a/Tl=0.2-Prcscrit
aJD=0.4-Preseiit
a/r>=fl.6-Prcscnt
a)
Fig. 4. Circular cross section (a) and geometry function vs. the dimensionless ratio s/h1 for K (b)
1
2
2
a
the predictions by Tada [8] adapted for a circular section, given by
^Tada=
4 h
R
1-
R
2Л
2tan
-0.923+ 0.199 1-sin^
nn
cos— 2
(22)
For KI , the relative difference between results seems to increase for higher values of the ratios a/D and s / h1.
Shear Force. Consider now the circular cross section in Fig. 5, a subjected to a shear force V . In this case, the actual direction of the transverse shear stress would depend on the general distance (R - a), such that we could no longer assume that the shear stress acts parallel to the y-axis. Considering, first, the only contribution of the shear stress along the y - direction, we have
4V cos S
4V (r2
3nR
2
3nR
4
(23)
where cos$ = V2Ra - a2 / R and h = VR 2 - y2.
For y = 0, the vertical shear stress tzy reaches the maximum value equal to
4V
zy _max 2'
(24)
3nR
The vertical shear stress contributes to the mode-II fracturing. As described in [79], we approximate the average shear stress acting on the reduced cross-section throughout the crack tip with the average shear stress acting on the un-cracked cross section. Thus, we apply the same method by Parton and Morozov [77] to compute the SIF.
The singular stress due to crack is related to the mode-II SIF as follows
KjI V (25)
T
zy _ s
V2nr
The equilibrium of the forces along the longitudinal direction of the beam can be written as
h/2
i Tzydy = j
= I Tzy _ sdr>
(26)
h —a 2
where b can be determined by enforcing the following condition
Tzy s \r=b Tzy \ h
y=—a
(27)
Upon substitution of Eqs. (23) and (25) into Eqs. (26) and (27), it is
Tzy_max
0.4 0.6
s/h,
C)
Fig. 5. Approximated shear stress distribution within the circular cross section (a), and geometry functions for Kn V (b) Кш V (c) due to the shear stress vs.
the dimensionless ratio s/h1
K
II V
= V2np
4V h
3nR1 V R
(
2
Л
1 -IR J(i-n)2
(28)
with
{ 13 12 fR? 1 ^ -n +-n + I — I n — 3l 2 ' I h J 4
2
^ R
1
—--
n
2
(29)
h J \ 2
By substituting Eq. (29) into Eq. (28), the mode-II SIF reads
4
2
s
2
V
KII _V =1T5 F(n), R
where the geometry function is
F (n) = ^ (l-( R T(1 -n)2 A
(30)
(31)
Figure 5, b plots the geometry function F(n) vs. the dimensionless ratio s /h1 for different lengths of the precrack a/D. A comparative assessment of the proposed formulation is performed with respect to previsions by Tada [8], adapted for a circular section, given by
v( ) 1 /h 1.3-0.65n + 0.37n2 + 0.28n3
F(n)Tada =-Vonn- l -L (32)
rcV R V1 -n
where a good matching between the two approaches can be noticed especially for lower values of a/D. An increasing difference between the proposed approach and predictions by Tada [8], instead, can be observed for higher values of a/D, whereas F(n) never gets a null value in the present formulation even when s / h1 = 1.
Considering now, the additional shear stress component Tzx along the x -direction, expressed as
5 „ 4sV sin e 4Vsy
Tzx =~Tzv tane=-— =--
hi zy 3nR 3nR4
(33)
whose singular value is related to the mode-III SIF as
KIII V
42%r
(34)
the equilibrium of the forces along the x -direction of the beam is expressed as
h/2
f Tzxdy = fTzx _ A, 0
h —a 2
(35)
where b can be determined by enforcing the following condition
zx s lr=b
lr=b Tzx I
(36)
y=—a 2
The substitution of Eqs. (33) and (34) into Eqs. (35) and (36) gives
Kr
4Vs ( h Y.5 ( 1
and
3nR ( R ) ( 2
ß= (1 -n)n
-n]V2nß
(37)
2(1 - 2n)
The substitution of Eq. (29) into Eq. (28), yields
(38)
with
kiii _ v =~t7 F (n),
R
f (n)=3 R ( R f ( 2-n^.
(39)
(40)
The geometry function F (n) related to the shear stress Tzx is plotted in Fig. 5, c vs. the dimensionless ratio s /h1 for varying lengths of the precrack a/D. A comparative assessment of the proposed formulation is performed with respect to previsions by Tada [8], adapted for a circular section, given by
F(n)Tada = - SJR -iПnl-, n Ri R y sin nn
(41)
where a quite good agreement is reached between the two different formulations, especially for lower values of s /h1 and a/D.
From Fig. 5, c, it can be noticed, once again, that the increasing values of a/D always yield to an increasing geometry function F(n) whose value
increases from zero up to a maximum value, for s / h1 ranging between zero and 0.7, thus decreasing for higher values of s/h1 up to the unit value.
Torsional Moment. Consider now the circular cross section in Fig. 6, a subjected to a torsional moment T . In this case, the torsional shear stress t'zx parallel to the crack tip is related to the mode-III fracturing as follows
T T
T'zx = — P cos 9 = — У, (42)
Jp J p
where JP is the polar moment of inertia of the circular section.
If we consider that the torsional stress distribution in the reduced cross section is equal to the corresponding distribution in the whole circular section, except for the fact that the reference axis is shifted by a distance e from the original coordinate, the stress distribution on the reduced cross section passing through the crack tip can be expressed as
■' = L '
' zx _ c = ~j~ y , Jc
(43)
where Jc is the polar moment of inertia for the reduced portion of the circular cross section defined in dimensionless form as
Jc n 1 f 2a \ —r =---arccosl 1--1 +
R4 2 2
D )
+
(1 - 2a / D)(3 - 8a / D + 8(a / D)2 )a / D-(a / D)2
(44)
T
zx s
h
3
0.4 0.6
s/h,
c)
Fig. 6. Geometry function for Km T (a) and Kjj t (b) due to a torsional moment vs.
the dimensionless ratio s/h1
where b can determined by enforcing the following condition
T zx s \r=b = T zx c I ' - h • (47)
- y = y=—a+e 2
This yields
K
T
111 ^Rï
F(n), (48)
where
F (n)=
f R4 ^
V y
1.5
and
ß 1 n
ß = — n — 2 2
R
fR4 \
(1 -n) (49)
V Jc y
(1 -n)
(
f R4 \
2n-1 + n — (1 -n) +n(l-n).
V Jc y
(50)
A comparative evaluation of the results in terms of F(n) as given by Eq. (49) and from the literature by Tada [8] is given in Fig. 6, b. Figure 6, b shows a comparative evaluation of the proposed formulation with respect to solution proposed by Tada [8], and adapted for the a circular section as follows
F(n)Tada =h{R ^H^ n)(51)
where a quite good agreement is reached between the two different formulations.
The vertical shear stress induced by the torsional moment is related to mode II, and it is expressed as follows
T
T
Tzy = — pSin^^ — X
J
P
J
(52)
P
The singular stress component related to the mode-III SIF is
kiii _ T
V5nr
(45)
whose value is independent of the crack depth a, whereas it reaches a minimum value x'zy min = 0 in the middle of the
crack tip x = 0 and a maximum value =—vr2 - y2 at x = -Rcos9 = -vR2 - y2 .
zy _max
The equilibrium of the forces along the x - direction of the beam is expressed as
h/2 b
J
P
i Tzxdy = |
= I T zx s dr'
(46)
—a-b
Thus, the tensional stress of the beam under the torsional moment is based on the double definition as in the following
T'y c = Ts (53)
' zy _ c
Jc
h
2
t
zx. s
h
2
and its singular form
zy _ s
KII _ T л/ 2nr
(54)
By equating the expressions (53) and (54) for ' - h
y = y = ^ - a + e we obtain the following expression for K
II T
T
KII t =—W 2nb,
Jc
(55)
where b can be determined with the section method by considering the following relation
h/2 b
J T'zydy = JT zy _ sdr-0
h b —a—b 2
(56)
Thus, based on the Eqs. (55) and (56), b is computed in dimensionless form as
ß = n Jc
2 Jp — Jc
(57)
while the substitution of Eq. (57) in the Eq. (55) yields to the following expression for the mode-II SIF
-T^F(n), (58)
K
II T
with
f (n)=^ №.
V V Jc R\ R
(59)
Figure 6, c plots the main results in terms of F(n) obtained with our proposed formulation,
compared to the corresponding predictions from the literature [8], adapted for the a circular section as in the following
2 s h 1.122-0.561n + 0.085^ + 0.18n
F (n)lada =- R<JR -r==-.
n RVR -n
(60)
A good agreement is reached between the two different formulations only for lower values of a/D, with an increasing discrepancy for higher ratios of a/D and s/h1. This is, mainly, related to the 2D approximated formulation provided by Tada [8] instead of the more appropriate 3D theory herein proposed for the purpose. For a/D ratios lower than 0.5, the geometry function F(n) seems to increase from zero up to a maximum value for s/h1 ranging between zero and 0.7, thus decreasing up to the null value for s/h1 approaching the unit value. Differently, F(n) seems to increase monotonically for a/D ratios higher than 0.5, in accordance to the predictions by Tada [8].
Comparison with Numerical Results Based on XFEM
Basis Notions of XFEM and Computation of the Mixed-Mode SIF. A short description of the SIF in mixed-mode conditions is provided in this section, remembering that the complete SIF for a 3D structure can be easily obtained as an algebraic summation of the components KI, Ku and Km defined as:
Ki =
M5 FI (n);
R
kii = kii к + KII T
R
Ts fii (n);
KIII = KIII V + KIII T
T
= t27 fiii (n). (61)
R
All the geometric functions determined analytically are now compared with those ones proposed in the literature [8] or obtained numerically with the commercial software ABAQUS by applying the XFEM approach. For the sake of brevity, the detailed implementation of the XFEM is not reported since it can be found in some specific references within the ABAQUS code.
The ABAQUS software provides the following methods to study fracture mechanics problems. In other words, quasi-statics can be applied to study the crack onset, growth and propagation along predefined paths, including the low-cycle fatigue, through the adoption of contour integrals. Part-through cracks in shells can be modeled inexpensively by using some line spring elements within a static procedure. An XFEM-based procedure models a crack as an enriched feature by adding degrees of freedom in elements with special displacement functions. The XFEM can be applied to simulate initiation and propagation of a discrete crack along an arbitrary, solution-dependent path without any requirement of remeshing. In addition, XFEM can be adopted to evaluate the contour integral without the need of refining the mesh around the crack tip (see Fig. 7). This approach also removes the requirement
Fig. 7. Coordinate configuration for the crack front enrichment functions
of defining explicitly the crack front or specifying the virtual crack extension direction when evaluating the contour integral. The data required for the contour integral are determined automatically based on the level set signed distance functions at nodes within an element.
As mentioned before, SIFs KI, KII and KIII play a key role in LEFM, since they characterize the influence of load or deformation on the magnitude of the crack-tip stress and strain fields and measure the entity of crack propagation or crack driving forces. Furthermore, the stress intensity can be related to the energy release rate, i.e. the ./-integral, computed directly in ABAQUS for a linear elastic material, as follows
/ =1 (K + Kl)) Kjn, (62)
E 2|
with E = E for plane stress and E = E/(l -v2) for plane strain conditions. However, for a crack under a mixed-mode loading, the integral method [82] is implemented in the ABAQUS code to compute all the SIFs. This method is well known to be applicable to
cracks in isotropic and anisotropic linear materials. Figures 8, a - i show the XFEM-based results in terms of geometry functions FI(n), FII(n) and Fm(r\), for a single beam and the whole structure, compared to the theoretical predictions provided by our formulation and proposed in the literature by Tada [8]. All these plots clearly show the quite good agreement between results, especially for small values of s/h1. Increasing deviations between theoretical and numerical results can be observed for larger values of s/h1 together with some visible discontinuities in the X-FEM trend for s/h1 ratios approaching the unit value, where the solution could be affected by some boundary effects.
Computation of the Crack Propagation. For completeness of the analysis, we finally apply the main values of the SIFs as given by our proposed model to determine the crack propagation angle 0O (see Fig. 9) both numerically and analytically. In the last case, we apply three different fracture criteria, i.e. the criterion by Richard [64] (here labeled as CR1), the minimum strain energy density factor criterion by Sih [62, 63] (here labeled as CR2), and the criterion by Schollmann et al. [72] (here labeled as CR3). The analytical results
-XFtM pStruclLrei 0 XFEM (Single beam) a Hr«ent Tilda
n R o A □ V « . „ A- 1
>00000000$
02 01 06 OB
slh,
g)
Fig. 8. Geometry functions for K,: a/D = 0.2 (a), a/D = 0.4 (b), a/D = 0.6 (c); Kn: a/D = 0.2 (d), a/D = 0.4 (e), a/D = 0.6 (f); Kln: a/D = 0.2 (g); a/D = 0.4 (h); a/D = 0.6 (i), vs. the dimensionless ratio s/Aj
a) b)
Fig. 9. Schematic representation of the negative tilting angle 0 (a) and twisting angle y (b)
are compared with the corresponding predictions provided by XFEM.
Additional details about the three selected criteria are given in Appendix A, B and C, respectively. Figures 10, a - c show the tilt angle 90 as given by the
Fig. 10. Propagation angle 00 vs. the dimensionless ratio s/h1 for a/D = 0.2 (a), a/D = 0.4 (b), a/D = 0.6 (c)
three cracking criteria, which clearly show a very good agreement of the results with those ones obtained numerically with the ABAQUS code. This confirms once more the ability of our proposed theoretical formulation combined with the selected criteria of fracturing to capture well the cracking development under mixed mode loading conditions. As also seen from Figs. 10, a - c, the propagation angle increases for increasing values of the dimensionless ratio s/h1 and it seems to be quite unaffected by the dimensionless crack depth a/D, especially for low values of the ratio s/h1.
Conclusion
An extended analytical model based on a section method has been presented in this work for a correct estimation of the SIFs and fracture direction of propagation within a cracked shaft under a coupled longitudinal, flexural and torsional loading conditions. The extended version of the section method here presented, considers the stress singularity at the tip of an elastic crack, and is based on an elementary beam theory equilibrium condition for internal forces evaluated in the cross-section throughout the crack tip. The quality of the analytical geometry functions in modes I, II, and III is evaluated by comparison with the theoretical predictions available from the handbooks. Based on a comparative assessment, and an error estimation of the geometry functions for Kj, Kjj, and Kjjj, a reasonable quality of the results is demonstrated in most cases yielding to a percentage error within 20 %. The only exception is represented by the torsional case and the associated computation of Kjj T, whose quality is reached only for low values
of the geometry cracking ratios a/D and s/h1. This is mainly related to the 2D approximated formulation provided in the handbooks [8] instead of the 3D analytical theory herein proposed for the purpose. A similar parametric investigation is also performed numerically via the XFEM-based approach as implemented in the ABAQUS code, for the computation of the geometry functions and propagation angle of cracking for 3D circular beams and structures. In the last case, the numerical results are compared to the analytical ones as given by three different fracture criteria selected from the literature, and here combined with the SIFs provided by our formulation. The cracking propagation angle under mixed mode loading condition is demonstrated to be accurately captured due to the perfect matching of the numerical and analytical curves. The mixed-mode geometry functions Fj(n), Fjj(n), and Fjjj(n), instead, are of acceptable accuracy especially for small values
0
0
of s/h1 when the error maintains lower than 20 %. The present study is not applicable only to rectangular and circular shafts because it has been already shown in the literature that the applied methodologies can be employed for the evaluation of the SIFs in general shape beams such as T and I beam in [83]. It will be the scope of a future contribution to study general shape beams within the present theoretical framework. Such analyses can lead to a numerical modeling where the crack behavior is modeled through line-springs for approximating the structural behavior [84, 85].
Acknowledgments
The research topic is one of the subjects of the Center of Study and Research for the Identification of Materials and Structures (CIMEST) "M. Capurso" of the University of Bologna (Italy), supported by the National Natural Science Foundation Project of China (51474039, 51404046, U1361205), by the Scientific Research Foundation of State Key Laboratory of Coal Mine Disaster Dynamics and Control (2011DA105287-ZD201302, 2011DA105287-MS201403), by the Fundamental Research Funds for the Central Universities (106112015CDJXY240003) and by the Basic Research of Frontier and Application of Chongqing (cstc2015jcyjA90019).
Appendix A. The Criterion by Richard et al. [64]
According to the generalized fracture criterion by Richard et al. [64], an equivalent stress intensity factor Keq is defined as follows
=Kl ~eq 2 2
As a function of the SIFs KI, KII, and KIII, with ai = KIc/KIIc = 1.155 and a2 = KIc/KIIIc = 1.0. Consequently, unstable crack growth occurs when Keq exceeds the fracture toughness KIc for mode I. For the determination of the crack deflection angle 90 there does exist a simple relation, which has been verified by a large number of experiments
( ik1 "
140--- - 70 -u
where y0 < 0° for KIir > 0 and y0 > 0° for KIir < 0
Keq = ^ + -JK2 + 4(a]K// )2 + 4(a2Km)2 < KIC.(A1)
00 =+
Ki +|kii| +|kiii| [Ki +|kii| + {Km
Ill
,(A2)
where 00 < 0 ° for Kn > 0 and 00 > 0 ° for Ku < 0 and Ki > 0.
The same approach could be also applied to compute the twisting angle y0 as follows
y0 = 3 78-N , , -3si KJ
Ki +|Kii|+\Krn\ [KI +\Ku\ +K
ii
VIII y
(A3)
and KI > 0.
Appendix B. The Strain Energy Density Factor Criterion by Sih [62, 63]
As far as brittle fracture is concerned, the strain energy density (S) theory, as formulated by Sih [62, 63] provides a more general treatment of fracture mechanics problems due to its ability to describe the multi-scale features of the material damage or the crack growth and propagation under mixed mode conditions. From the linear theory of isotropic and homogeneous elasticity, the factor S is found to feature a singularity of the order 1/r near the crack front as in the following
dW = S dVr
where W is the strain energy, and V is the volume. The core of the S - criterion is the parameter of strain energy density factor S which is a function of the SIFs for linear elasticity, defined as (see Figs. 9, a, b)
S (0, V) = T7Z-L— (a
16n| cos y
\anKj + 2a12 KjKn +
+ Ü22 Kir + Ü33KItt ).
^22 kii^U33kiii). (B2)
In the previous equation, aij (with i, j = 1, 2, 3) are some auxiliary functions which depend on the elastic properties (i.e. the Young modulus E, Poisson ratio v, and shear modulus of the material, as well as on the polar angle 9 of the crack tip, as follows:
a11 = (3 - 4v - cos 9)(1 + cos 9);
a12 = 2 sin 9(cos 9 - (1 - 2v));
a22 = 4(1 - v)(1 - cos 9) + (1 + cos 9)(3 cos 9 -1);
a33 = 4. (B3)
Based on the assumption of the S - criterion, the prediction of the crack initiation and direction is based on the following relations:
da]
11
1
50 16n|cosy
dan 1
50 16n|cosy
daii 1
2sin 0(cos 0 + 2v-1);
1-(cos 0(cos 0 + 2v -1) - 2sin2 0);
50 16n|cosy and
3S_ 50
2sin 0(1 - 2v - 3cos 0)
= sin 0(cos 0- 0.4) KI +
(B4)
0=60
+
(cos 0(cos0- 0.4) - sin2 0)lKIKII
+ sin 0(0.4 - 3cos 0) KI = 0,
(B5)
where v is herein set to 0.3. Similarly, the crack angle
y = y0 is derived by minimizing S with respect to y,
i.e. dS / dyl = 0. T|y=yo
Appendix C. The Maximum Principle Stress Criterion by Schollmann et al. [72]
The 3D criterion as proposed by Schollmann et al. [72], considers all the fracture modes I, II and III, and is otherwise called as maximum principle stress Gi'-criterion, where al is the maximum principle stress on a virtual cylindrical surface around the crack front. The present criterion reduces to the well-known maximum tangential stress-criterion in case of in-plane mixed-mode loading conditions (i.e. modes I and II) along the crack front, and thus can be considered as its generalization in 3D. Based on this criterion, the crack kinks in mode II with an angle 90 (see Fig. 9, a) which is determined through the following equation (for more details, see [72])
- 6Kj r - Kjj (6 -12r2)+
+ (4Kj - 12Kjjr)(- 6Kjr - Kn (6 - 12r2)) - 32Kl//r(l + r yl(4Kj -12KUr)2 + 64K]jj (1 + r2 )2
= 0
(C1)
where r = tan(90/2). At the same time, the crack front in mode III loading conditions grows into a new plane, which twists with an angle y0 (see Fig. 9, b) determined as
f TF F A
1 -% = 2tan
_KjL_
VKI (rcos!(e0 /2) - 2ur)-Kn (3+3sin2(90 /2) - 2u),
(C2)
An unstable crack growth occurs if the comparative SIF Kv reaches the fracture toughness KIC, where Kv is expressed as
Kv = lcosfe0
v 2 I 2
J Kj
cos2 (e° 1- 3Kjj sin(e0 )-
+ .
(C3)
Kjcos2I- 3Kjj sin(e0)
+ 4 K
jjj
References
1. Sih, G.C. Short and long crack data for fatigue of 2024-T3 Al sheets: binariness of scale segmentation in space and time. Fatigue Fract. Engng. Mater. Struct., 2014, vol. 37, no. 5, pp. 484-493.
2. Sih, G.C. Redemption of the formalism of
segmented linearity: multiscaling of non-equilibrium and non-homogeneity applied to fatigue crack growth. Fatigue Fract. Engng. Mater. Struct., 2014, vol. 38, no. 5, pp. 621-628.
3. Pook, L.P. An aircraft company, a government laboratory, and a university. In: Rossmanith HP editor. Fracture research in retrospect. Rotterdam: A.A. Balkema, 1997, pp. 493-505.
4. Pook, L.P. Five decades of crack path research. Eng. Fract. Mech, 2010, vol. 77, pp. 1619-1630.
5. Pook, L.P. The linear elastic analysis of cracked bodies, crack paths and some practical crack path examples. Eng. Fract. Mech., 2016, vol. 167, pp. 2-19.
6. Paris, P.C. and Sih, G.C. Stress analysis of cracks. In: Fracture toughness testing and its applications. ASTM STP 381. American Society for Testing and Materials, Philadelphia, USA, 1965, pp. 30-81.
7. Sih, G.C. Handbook of stress intensity factors for researchers and engineers. Lehigh University, Bethlehem: Institute for Fracture and Solid Mechanics, 1973.
8. Tada, H., Paris, P.C. and Irwin, G.R. The stress analysis of crack handbook. Paris Production Incorporated and Del Research Corporation, St. Louis, 1985.
9. Murakami, Y. Handbook of stress intensity factors.Oxford: Pergamon Press, 1987.
10. Pook, L.P. Approximate stress intensity factors obtained from simple plate bending theory. Eng. Fract. Mech., 1979, vol. 12, no. 4, pp. 505-522.
11. Hartranft, R. J. and Sih, G.C. An approximate three-dimensional theory of plates with application to crack problems. Int. /. Eng. Sci., 1970, vol. 8, pp. 711-729.
12. Benthem, J. P. State of stress at the vertex of a quarter-infinite crack in a half space. Int. /. Solids Struct., 1977, vol. 13, pp. 479-492.
13. Nakamura, T. and Parks, D.M. Antisymmetric 3D field near the crack front of a thin plate. Int. /. Solids Struct., 1989, vol. 25, pp. 1411-1426.
14. Berto, F., Lazzarin, P., Kotousov, A. and Pook, L.P. Induced out-of-plane mode at the tip of blunt lateral notches and holes under in-plane shear loading. Fatigue Fract. Engng. Mater. Struct., 2012, vol. 35, pp. 538-555.
15. Kotousov, A., Lazzarin, P., Berto, F. and Pook, L.P. Three-dimensional stress states at crack tip induced by shear and anti-plane loading. Eng. Fract. Mech., 2013, vol. 108, pp. 65-74.
16. Pook, L.P. A 50-year retrospective review of three-dimensional effects at cracks and sharp notches. Fatigue Fract. Engng. Mater. Struct., 2013, vol. 36, pp. 699-723.
17. Zhang, Z.N., Wang, D.Y., Zheng, H. and Ge, X.R. Interactions of 3D embedded parallel vertically inclined cracks subjected to uniaxial compression. Theor. Appl. Fract. Mech, 2012, vol. 61, pp. 1-11.
18. Garcia-Manrique, J., Camas, D. and Lopez-Crespo, P. and Gonzalez-Herrera, A. Stress intensity factor analysis of throug thickness effects. Int. /. Fatigue, 2013, vol. 46, pp. 58-66.
19. Aygul, M., Al-Emrani, M., Barsoum, Z. and Lenader, J. An investigation of distortion-induced fatigue cracking under variable amplitude loading using 3D crack propagation analysis. Eng. Fail. Anal., 2014, vol. 45, pp. 151-163.
20. Berto, F., Campagnolo, A. and Lazzarin, P. Fatigue strength of severely nothced specimens made of Ti-6Al-4V
2
under multiaxial loading. Fatigue Fract. Engng. Mater. Struct., 2015, vol. 38, pp. 503-517.
21. Storgards, E., Simonsson,, K. and Sjostrom S. Three-dimensional crack growth modeling of a Ni-based superalloy at elevated temperature and sustained loading. Theor. Appl. Fract. Mech, 2016, vol. 81, pp. 2-10.
22. Miao, X.T., Zhou, C.Y., Li, J. and He X.H. Studies of elastic and elastic-plastic J integral for mixed mode cracked plate under biaxial loading. Fatigue Fract. Engng. Mater. Struct., 2016, vol. 39, pp. 536-550.
23. Valiente, A. Criterios de fractura para alambres. PhD thesis, Polytechnic University of Madrid, Spain. <http://oa.upm.es/745/01/04198001.pdf>, 1980.
24. Levan, A. and Royer, J. Part-circular surface cracks in round bars under tension, bending and twisting. Int. J. Fract, 1993, vol. 61, pp. 71-99.
25. Couroneau, N. and Royer, J. Simplified model for the fatigue growth analysis of surface cracks in round bars under mode I. Int. J. Fract, 1998, vol. 20, pp. 711-718.
26. Daoud, O.E.K., Cartwright, D.J. and Carney, M. Strain-energy release rate for a single-edge-cracked circular bar in tension. J. Strain Anal. Eng. Des., 1978, vol. 13, no. 2, pp. 83-89.
27. Daoud, O.E.K. and Cartwright, D.J. Strain energy release rate for a straight-fronted edge crack in a circular bar subject to bending. Eng. Fract. Mech., 1984, vol. 19, pp. 701-707.
28. Ng, C.K. and Fenner, D.N. Stress intensity factors for an edge cracked bar in tension and bending. Int. J. Fract., 1988, vol. 36, pp. 291-303.
29. Shin, C.S. and Cai, C.Q. Experimental and finite element analyses on stress intensity factors of elliptical surface crack in a circular shaft under tension and bending. Int. J. Fract, 2004, vol. 129, pp. 239-264.
30. Carpinteri, A. Stress intensity factors for straight-fronted edge cracks in round bars. Eng. Fract. Mech., 1992, vol. 42, no. 6, pp. 1035-1040.
31. Eijian,S. Stress intensity factors for edge cracks in round bars. Eng. Fract. Mech, 1990, vol. 37, pp. 805-812.
32. Naik, S.S. and Maiti, S.K. Triply coupled bendingtorsion vibration of Timoshenko end Euler-Bernoulli shaft beams with arbitrarily oriented open crack. J. Sound Vib., 2009, vol. 324, pp. 1067-1085.
33. James, L.A. and Mills, W.J. Review and synthesis of stress intensity factor solutions applicable to cracks in bolts. Eng. Fract. Mech., 1988, vol. 30, pp. 641-654.
34. Toribio, J., Alvarez, N., González, B., Matos, J.C. A critical review of stress intensity factor solution for surface cracks in round bars subjected to tension bending. Eng. Fail. Anal., 2009, vol. 16, pp. 794-809.
35. Dixon, J.R. and Pook, L.P. Stress intensity factors calculated generally by the finite element technique. Nature, 1969, vol. 224, no. 5215, pp. 166-167.
36. Chen, C.H. and Wang. C.L. Stress intensity factors and T-stresses for offset double edge-cracked plates under mixed-mode loadings. Int. J. Fract., 2008, vol. 152, pp. 149-162.
37. Jogdand, P.V. and Murthy, K.S.R.K. A finite element based interior collocation method for the computation of stress intensity factors and T-stresses. Eng. Fract. Mech., 2010, vol. 77, pp. 1116-1127.
38. Feng, X.Q., Shi, Y.F., Li, B., Yu, S.W. and Yang, Q. Dislocation-based semi-analytical method for calculating stress intensity factors of cracks: two-dimensional cases. Eng. Fract. Mech., 2010, vol. 77, pp. 3521-3531.
39. Whitcomb, J.D., Raju, I.S. and Goree, J.G. Reliability of the finite element method for calculating free edge stresses in composite laminates. Comput. Struct., 1982, vol. 15, pp. 23-27.
40. Chauduri, R.A. A combined FEM and three-dimensional asymptotic solution based analysis on stress concentration/intensity around elliptical/circular cylinder shaped surface flaws in stretched plates. Appl. Math. Model., 2015, vol. 39, pp. 5341-5353.
41. Ferretti, E. Crack propagation modeling by remeshing using the Cell Method (CM). CMES-Comput. Model Eng., 2003, vol. 4, pp. 51-72.
42. Ferretti, E. A discrete nonlocal formulation using local constitutive laws. Int. J. Fract., 2004, vol. 130, no. 3, pp. 175-182.
43. Dolbow, J. and Harari, I. An efficient finite element method for embedded interface problems. Int. J. Numer. Methods Eng., 2009, vol. 78, no. 2, pp. 229-252.
44. Fang, X. and Charalambides, P.G. The fracture mechanics of cantilever beams with an embedded sharp crack under end force loading. Eng. Fract. Mech, 2015, vol. 149, pp. 1-17.
45. Viola, E., Tornabene, F., Ferretti, E. and Fantuzzi, N. GDQFEM numerical simulations of continuous media with cracks and discontinuities. CMES-Comput. Model Eng., 2013, vol. 94, no. 4, pp. 331-369.
46. Tornabene, F., Fantuzzi, N., Ubertini, F., Viola, E. Strong formulation finite element method based on differential quadrature: a survey. Appl. Mech. Rev., 2015, vol. 67, pp. 1-55.
47. Fantuzzi, N., Dimitri, R. and Tornabene, F. A SFEM-based evaluation of mode-I Stress Intensity Factor in composite structures. Comput. Struct., 2016, vol. 145, pp. 162-185.
48. Dimitri, R., Fantuzzi, N., Tornabene, F., Zavarise, G. Innovative numerical methods based on SFEM and IGA for computing stress concentrations in plates with discontinuities. Int. J. Mech. Sci., 2016, vol. 118, pp. 166-187.
49. Verhoseel, C.V., Scott, M.A., de Borst, R. and Hughes, T.J.R. An Isogeometric approach to cohesive zone modeling. Int. J. Numer. Methods Eng, 2011, vol. 87, no. 1-5, pp. 336-360.
50. Nguyen, V.P. and Nguyen-Xuan, H. High-order B-splines based finite elements for delamination analysis of laminated composites. Compos. Struct., 2013, vol. 102, pp. 261-275.
51. Dimitri, R., De Lorenzis, L., Wriggers, P. and Zavarise, G. NURBS- and T-spline-based isogeometric cohesive zone modeling of interface debonding. Comput. Mech., 2014, vol. 54, pp. 369-388.
52. Dimitri, R., Isogeometric treatment of large deformation contact and debonding problems with T-splines: a review. Curved Layer. Struct., 2015, vol. 2, pp. 59-90.
53. Campilho, R.D.S.G., Banea, M.D., Pinto, A.M.G., da Silva, L.F.M. and de Jesus, A.M.P. Strength prediction of single- and double-lap joints by standard and extended finite element modeling. Int. /. Adhes., 2011, vol. 31, pp. 363-372.
54. Benvenuti, E. and Tralli, A. Simulation of finite-width process zone in concrete-like materials by means of a regularized extended finite element model. Comput. Mech., 2012, vol. 50, no. 4, pp. 479-497.
55. Sadeghirad, A., Chopp, D.L., Ren, X., Fang, E. and Lua, J. A novel hybrid approach for level set characterization and tracking of non-planar 3D cracks in the extended finite element method. Eng. Fract. Mech., 2016, vol. 160, pp. 1-14.
56. Belytschko, T. and Black, T. Elastic crack growth in finite elements with minimal remeshing. Int. /. Num. Methods Eng., 1999, vol. 45, pp. 601-620.
57. Erdogan, F. and Sih, G.C. On the crack extension in plates under plane loading and transverse shear. /. Basic Eng., 1963, vol. 85, pp. 519-527.
58. Smith, D.J., Ayatollahi, M.R. and Pavier, M.J. The role of T-stress in brittle fracture for linear elastic materials under mixed-mode loading. Fatigue Fract. Engng. Mater. Struct., 2001, vol. 24, no. 5, pp. 137-150.
59. Ayatollahi, M.R., Aliha, M.R.M., Hassani, M.M. Mixed mode brittle fracture in PMMA—An experimental study using SCB specimens. Mater Sci Eng., 2006, vol. 417, no. 1-2, pp. 348-356.
60. Ayatollahi, M.R. and Aliha, M.R.M. On the use of Brazilian disc specimen for calculating mixed mode I-II fracture toughness of rock materials. Eng. Fract. Mech., 2008, vol. 75, no. 16, pp. 4631-4641.
61. Aliha, M.R.M. and Ayatollahi, M.R. Analysis of fracture initiation angle in some cracked ceramics using the generalized maximum tangential stress criterion. Int. /. Solids Struct., 2012, vol. 49, no. 13, pp. 1877-1883.
62. Sih, G.C. Mechanics of fracture 1: a special theory of crack propagation. Leyden: Noordhoff International Publishing, 1973.
63. Sih, G.C. Strain energy-density factor applied to mixed mode crack problem. Int. /. Fract., 1974, vol. 10, no. 3, pp. 305-321.
64. Richard, H.A., Fulland, M. and Sander, M. Theoretical crack path prediction. Fatigue Fract. Engng. Mater. Struct., 2005, vol. 28, no. 1, pp. 3-12.
65. Lazzarin, P., Berto, F. and Zappalorto, M. Rapid calculations of notch stress intensity factors based on averaged strain energy density from coarse meshes: Theoretical bases and applications. Int. /. Fatigue, 2010, vol. 32, pp. 1559-1567.
66. Berto, F. and Lazzarin, P. A review of the volume-based strain energy density approach applied to V-notches and welded structures. Theor. Appl. Fract. Mech., 2009, vol. 52, no. 3, pp. 183-194.
67. Ayatollahi, M.R., Rashidi Moghaddam, M., Razavi, S.M.J. and Berto, F. Geometry effects on fracture trajectory of PMMA samples under pure mode-I loading. Eng. Fract. Mech., 2016, vol. 163, pp. 449-461.
68. Hussain, M.A. and Pu, S.L. Underwood J. Strain energy release rate for a crack under combined mode I and mode II. In: Paris PC, Irwin GR, editors. Fracture analysis, ASTM STP 560, Philadelphia, 1974, pp. 2-28.
69. Griffith, A.A. The phenomena of rupture and flow in solids. Philos. Trans. R. Soc. London, 1920, vol. 221, pp. 163-198.
70. Cornetti, P., Pugno, N., Carpinteri, A., Taylor, D. Finite fracture mechanics: A coupled stress and Energy failure criterion. Eng. Fract. Mech., 2006, vol. 73, pp. 2021-2033.
71. Dimitri, R., Cornetti, P., Mantic, V., Trullo M. and De Lorenzis, L. Mode-I debonding of a double cantilever beam: a comparison between cohesive crack modeling and finite fracture mechanics. Submitted.
72. Schollmann, M., Richard, H.A., Kullmer, G. and Fulland, M. A new criterion for the prediction of crack development in multiaxially loaded structures. Int. J. Fract., 2002, vol. 117, no. 2, pp. 129-141.
73. Pook, L.P. Linear Elastic Fracture Mechanics for Engineers: Theory and Application. WIT Press, Southampton, 2000.
74. Chang, J., Xu, J.Q. and Mutoh, Y. A general mixed-mode brittle fracture criterion for cracked materials. Eng. Fract. Mech, 2006, vol. 73, pp. 1249-1263.
75. Sih, G.C. and Cha, C.K. A fracture criterion for three-dimensional crack problems. Eng. Fract. Mech, 1974, vol. 6, pp. 699-723.
76. Kong, X.M., Schluter, N. and Dahl, W. Effect of triaxial stress on mixed-mode fracture. Eng. Fract. Mech., 1995, vol. 52, pp. 379-388.
77. Parton, V.Z. and Morozov, E.M. Elastic-plastic fracture mechanics. Moscow: Mir Publishers, 1978.
78. Nobile, L. Mixed mode crack initiation and direction in beams with edge crack. Theor. Appl. Fract. Mech., 2000, vol. 33, pp. 107-116.
79. Nobile, L. Mixed mode crack growth in curved beams with radial edge crack. Theor. Appl. Fract. Mech., 2001, vol. 36, pp. 61-72.
80. Brown, W.F. and Srawley, J.E. Plane strain crack toughness testing of high strength metallic materials. ASTM STP, 1966, vol. 410, pp. 63-65.
81. Shin, C.S. and Cai, C.Q. Experimental and finite element analyses on stress intensity factors of an elliptical surface crack in a circular shaft under tension and bending. Int. J. Fract, 2004, vol. 129, pp. 239-264.
82. Shih, C.F. and Asaro, R.J. Elastic-Plastic Analysis of Cracks on Bimaterial Interfaces: Part I, Small Scale Yielding. J. Appl. Mech, 1988, vol. 8, pp. 537-545.
83. Ricci, P., Viola, E. Stress intensity factors for cracked T-sections and dynamic behaviour of T-beams. Eng. Fract. Mech, 2006, vol. 73, pp. 91-111.
84. Viola, E., Fantuzzi, N. and Marzani, A. Cracks interaction effect on the dynamic stability of beams under conservative and nonconservative forces. Key Eng. Mat., 2012, vol. 488-489, pp. 383-386.
85. Viola, E., Marzani, A. and Fantuzzi, N. Interaction effect of cracks on flutter and divergence instabilities of cracked beams under subtangential forces. Eng. Fract. Mech., 2016, vol. 151, pp. 109-129.
