Microsoft Word - numero_46_art_16 V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 158 Developments in the fracture and fatigue assessment of materials and structures Analysis of cylindrical delamination cracks in multilayered functionally graded non-linear elastic circular shafts under combined loads Victor Rizov Department of Technical Mechanics, University of Architecture, Civil Engineering and Geodesy, 1 Chr. Smirnensky blvd., 1046 – Sofia, Bulgaria v_rizov_fhe@uacg.bg ABSTRACT. This paper is focused on delamination fracture analyses of a multilayered functionally graded circular shaft under two loading combinations (centric tension and torsion, and bending and torsion) assuming non-linear elastic mechanical behavior of the material. The loading combinations under consideration generate mixed-mode II/III delamination crack loading conditions (the centric tension and bending generate mode II crack loading, while the torsion is responsible for mode III crack loading). The shaft is made by concentric longitudinal layers. The number of layers is arbitrary. Besides, each layer has individual thickness and material properties. The material in each layer is functionally graded in radial direction. Hyperbolic laws are used to describe the continuous variation of material properties in radial direction. A cylindrical delamination crack (the crack front is a circle) is located arbitrary between layers. The delamination fracture is studied in terms of the strain energy release rate by analyzing the energy balance. In order to verify the solution obtained, the strain energy release rate is derived also by differentiating the complementary strain energy with respect to the delamination crack area. Parametric investigations of the behavior of the cylindrical delamination crack are carried-out. The present paper is a contribution in the fracture mechanics of multilayered functionally graded non-linear elastic circular shafts under combined loads. KEYWORDS. Multilayered circular shaft; Functionally graded material; Cylindrical delamination crack; Material non-linearity; Combined loads. Citation: Rizov, V., Analysis of cylindrical delamination cracks in multilayered functionally graded non-linear elastic circular shafts under combined loads, Frattura ed Integrità Strutturale, 46 (2018) 158-177. Received: 04.05.2018 Accepted: 03.06.2018 Published: 01.10.2018 Copyright: © 2018 This is an open access article under the terms of the CC-BY 4.0, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. INTRODUCTION ne of the main advantages of the functionally graded materials over the traditional structural materials is the material property graded distribution [1 - 6]. Very often, fracture is the critical failure mode for the functionally graded structures to lose their structural capacity [7 - 10]. Therefore, the study of fracture mechanics of O http://www.gruppofrattura.it/VA/46/16.mp4 V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 159 functionally graded materials plays an important role in the design of various structural members and devices made by these novel un-homogeneous materials. Understanding the fracture behavior of functionally graded structures under various loading conditions is vital for the further development of the methods for safety design. The present paper deals with analyses of a cylindrical delamination crack in a multilayered functionally graded non-linear elastic circular shaft under combined loads. It should be mentioned that in one of the previous works of the author, non- linear analyses of cylindrical delamination cracks in circular shafts have been developed assuming that the shafts are loaded in pure torsion only [10]. However, in reality, the circular shafts usually are under various load combinations which include torsion (this fact is the basic motive for writing the present paper). FRACTURE ANALYSIS IN TERMS OF THE STRAIN ENERGY RELEASE RATE Shaft under centric tension and torsion he multilayered functionally graded circular shaft, shown schematically in Fig. 1, is under consideration. The shaft is made of adhesively bonded concentric longitudinal layers. In each layer, the material is functionally graded in radial direction. Besides, the functionally graded material exhibits non-linear mechanical behavior. Figure 1: Multilayered functionally graded circular shaft loaded in centric tension and torsion. The number of layers is arbitrary. Also, each layer has individual thickness and material properties. The shaft cross-section is a circle of radius, R . The length of the shaft is 2l . The shaft is loaded in centric tension and torsion, respectively, by longitudinal forces, F , and torsion moments, T , applied at the end sections of the shaft as shown in Fig. 1. A circular notch is cut-out in the middle of the shaft in order to generate conditions for delamination fracture. It is assumed that a cylindrical delamination crack of length, 2a , is located symmetrically with respect to the middle of the shaft. The delamination crack represents a cylindrical surface (the crack front is a circle of radius, br ). Thus, the internal crack arm is a shaft of length, 2a , and circular cross-section of radius, br . The external crack arm is a shaft of length, 2a , and ring- shaped cross-section of internal radius, br , and external radius, R . The delamination crack is located arbitrary between layers. The circular notch divides the external crack arm in two symmetric segments of length, a , each. Apparently, the two segments of the external crack arm are free of stresses (Fig. 1). Due to the symmetry, only half of the shaft, 2l x l  , is analyzed. In the present paper, the delamination fracture is studied in terms of the strain energy release rate. It is obvious that the longitudinal force, F , induces mode II crack loading conditions. The mode II component of the strain energy release rate, IIG , is determined by analyzing the energy balance. By assuming an increase of the crack length, a , the energy balance is written as F II C U F u a G l a a        (1) T V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 160 where u is the increase of the longitudinal displacement of the end section of the shaft, FU is the strain energy cumulated in half of the shaft as a result of the centric tension by F , Cl is the length of the crack front. By substituting of 2C bl r in (1), IIG is obtained as 1 2 2 2 F II b b UF u G r a r a         (2) The expression in brackets in (2) is doubled due to the symmetry (Fig. 1). It should be specified that the present delamination fracture analysis is valid for non-linear elastic behavior of the material. The analysis can also be applied for elastic-plastic behavior if the shaft undergoes active deformation, i.e. if the external loading increases only [11, 12]. It should also be mentioned that the present analysis is carried-out assuming validity of the small strains assumption. By using methods of Mechanics of materials, one obtains  L Hu a l a    (3) where L and H are, respectively, the longitudinal strains in the internal crack arm and the un-cracked shaft portion, 2l a x l   , induced by the longitudinal force, F . Figure 2: Cross-section of the internal crack arm loaded in centric tension and torsion. The longitudinal strain in the internal crack arm is determined from the following equation for equilibrium of the cross- section of the internal crack arm: 1 1 i i n i i A F dA     (4) where 1n is the number of layers in the internal crack arm, i is the distribution of the longitudinal normal stresses in the i-th layer, iA is the area of the cross-section of the same layer. In the present paper, the mechanical behavior of the functionally graded material in the i-th layer is described by the following non-linear stress-strain relation [13]: i i is p      (5) where is and ip are the distributions of the material properties in the same layer,  is the longitudinal strain. V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 161 The properties, is and ip , vary continuously in the radial direction of the i-th layer according to the following hyperbolic laws: 1 1 i i B i i D i i s s r r s r r     (6) 1 1 i i B i i D i i p p r r p r r     (7) where iB s , iD s , iB p and iD p are material properties ( iD s and iD p control the material gradient of is and ip , respectively), ir and 1ir  are shown in Fig. 2. In (6) and (7), the radius, r , varies in the interval  1;i ir r  . It should be mentioned that the distribution of the longitudinal strains is analyzed assuming validity of the hypothesis for plane sections, since the length to diameter of the cross-section ratio of the shaft under consideration is large. Thus, L is distributed uniformly in the cross-section of the internal crack arm. Hence, by substituting of (5), (6) and (7) in (4), one derives     1 2 2 3 3 1 1 1 2 2 3 i n i i L i i i i i F r r r r                (8) where 1 i i i B L B i i s p       (9) 2 2 2 i i i i i i B D L B D i i i i i B L B i i s s p p s p                  (10) 1i i ir r   (11) 1 i D i i i s r     (12) 1 i D i i i p r     (13) It should be noted that by substituting of 0 iD s  and 0 iB p  in (8), one obtains   1 2 2 1 1 1 i i n L i i i B F r r s       (14) V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 162 which is exact match of the equation for equilibrium of multilayered circular shaft made by homogeneous linear-elastic layers loaded in centric tension [14]. This fact is an indication for consistency of Eqn. (8) since at 0 iD s  and 0 iB p  the non-linear stress-strain relation (5) transforms in the Hooke’s law assuming that 1/ iB s is the modulus of elasticity in the i-th layer. Eqn. (8) should be solved with respect to L by using the MatLab computer program. Eqn. (8) is applied also to determine H . For this purpose, 1n and L are replaced, respectively, with n and H in (8), (9) and (10). Here, n is the number of layers in the un-cracked shaft portion. Figure 3: Non-linear   diagram. Since the external crack arm is free of stresses (Fig. 1), the strain energy cumulated in half of the shaft as a result of the centric tension is written as F FL FHU U U  (15) where FLU and FHU are the strain energies in the internal crack arm and the un-cracked shaft portion, respectively. The strain energy in the internal crack arm is obtained by addition of strain energies cumulated in the layers 1 0 1 i i i n FL FL i A U a u dA     (16) where 0 iFLu is the strain energy density in the i-th layer. The strain energy density is equal to the area, OPQ, enclosed by the stress-strain curve (Fig. 3). Thus, 0 iFLu is written as 0 0 L iFL i u d     (17) By substituting of (5) in (17), one derives 0 1 ln ln i i i i i FL L L i i i i i s s s s u p p p p p                (18) The strain energy cumulated in the un-cracked shaft portion is expressed as V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 163   0 1 i i i n FH FH i A U l a u dA      (19) where the strain energy density in the i-th layer, 0 iFHu , is obtained by formula (18). For this purpose, L is replaced with H . By substituting of (3), (15), (16) and (19) in (2), one obtains   1 0 0 1 1 1 i i i i i n i n II L H FL FH i ib b A A F G u dA u dA r r                       (20) Apparently, the torsion moment, T , induces mode III crack loading conditions (Fig. 1). By analyzing the balance of the energy, the mode III component of the strain energy release rate, IIIG , is written as 1 2 2 2 T III b b UT G r a r a           (21) where  is the angle of twist of the end section of the shaft, TU is the strain energy cumulated in half of the shaft as a result of the torsion. In (21), the expression in the brackets is doubled in view of the symmetry (Fig. 1). By applying methods of Mechanics of materials, one obtains  qm b a l a r R      (22) where m and q are the shear strains at the periphery of the cross-sections of the internal crack arm and the un-cracked shaft portion, respectively. The shear strain at the periphery of the cross-section of the internal crack arm is determined by using the following equation for equilibrium of the cross-section of the internal crack arm: 1 1 i i n i i A T rdA     (23) where i is the distribution of the shear stresses in the i-th layer induced by the torsion. In the present paper, the mechanical behavior of the functionally graded material in torsion is described by the following non-linear stress-strain relation [13]: i i if g      (24) where  is the shear strain, if and ig are the distributions of the material properties in the i-th layer. The continuous variation of if and ig in the radial direction of the i-th layer is described by the following hyperbolic laws: 1 1 i i B i i D i i f f r r f r r     (25) V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 164 1 1 i i B i i D i i g g r r g r r     (26) where 1i ir r r   (27) The radiuses, ir and 1ir  , are shown in Fig. 2. In (25) and (26), iBf , iDf , iBg and iDg are material properties ( iDf and iD f govern the gradient of if and ig , respectively). The distribution of shear strains in radial direction is written as m b r r   (28) By substituting of (24), (25), (26) and (28) in (23), one derives     1 4 4 5 1 1 1 2 1 1 5 4 5 i n m i i i i i i i b T r r r r r                 (29) where 1 i i    (30)  22 21 1 i i i B i B m i D ii i i i b i f g g rr r                (31) iD i i f    (32) 1 iB i i i f r     (33) It should be mentioned that at 0 iB g  and 0 iD f  Eqn. (29) transforms in   1 4 4 1 1 1 2 i i n m i i i b B T r r r f        (34) The fact that (34) is exact match of the equation for equilibrium of a multilayered circular shaft made by linear-elastic homogeneous layers loaded in torsion [14] is an indication for consistency of (34) since at 0 iB g  and 0 iD f  the non- linear stress-strain relation (24) transforms in the Hooke’s law assuming that 1/ iB f is the shear modulus in the i-th layer. Eqn. (29) should be solved with respect to m by using the MatLab computer program. Eqn. (29) is used also to determine the shear strain at the periphery of the cross-section of the un-cracked shaft portion. For this purpose, 1n , br and m are replaced, respectively, with n , R and q in (29) and (31). V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 165 Figure 4: Non-linear   diagram. The strain energy cumulated in half of the shaft as a result of the torsion is obtained as T TL THU U U  (35) where TLU and THU are the strain energies in the internal crack arm and the un-cracked shaft portion, respectively. The strain energy in the internal crack arm is written as 1 0 1 i i i n TL TL i A U a u dA     (36) where 0 iTLu is the strain energy density in the i-th layer as a result of the torsion. In principle, the strain energy density is equal to the area, OPQ, enclosed by stress-strain curve in Fig. 4. Thus, formula (18) can be used to obtain 0 iTLu . For this purpose, 0 iFLu , L , is and ip are replaced, respectively, with 0 iTLu ,  , if and ig , where  is expressed by (28). The strain energy cumulated in the un-cracked shaft portion as a result of the torsion is expressed as   0 1 i i i n HL TH i A U l a u dA      (37) where the strain energy density in the i-th layer, 0 iTHu , is obtained by (18). For this purpose, 0 iFLu , L , is and ip are replaced, respectively, with 0 iTHu , H , if and ig . Here, the distribution of the shear strains is written as H q r R   (38) By substituting of (22), (35), (36) and (37) in (21), one obtains 1 0 0 1 1 1 i i i i i n i n qm III TL TH i ib b b A A T G u dA u dA r r R r                         (39) The total strain energy release rate, G , is written as II IIIG G G  (40) V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 166 By substituting of (20) and (39) in (40), one arrives at   1 0 0 1 1 1 i i i i i n i n L H FL FH i ib b A A F G u dA u dA r r                        1 0 0 1 1 1 i i i i i n i n qm TL TH i ib b b A A T u dA u dA r r R r                         (41) The integration in (41) should be performed by the MatLab computer program. In order to verify (41), the strain energy release rate is derived also by differentiating the complementary strain energy with respect to the crack area. The total strain energy release rate is written as [15] *dU G dA  (42) where *dU is the change of the complementary strain energy, dA is an elementary increase of the crack area. For the cylindrical delamination crack (Fig. 1), dA is expressed as 2 bdA r da (43) where da is an elementary increase of the crack length. By substituting of (43) in (42), one arrives at * 2 b dU G r da  (44) The complementary strain energy cumulated in half of the shaft as a result of the centric tension and torsion is obtained as * * *L HU U U  (45) where *LU and * HU are the complementary strain energies in the internal crack arm and the un-cracked shaft portion, respectively. The complementary strain energy in the internal crack arm is expressed as 1 * * 0 1 i i i n L L i A U a u dA     (46) where *0 iLu is the complementary strain energy density in the i-th layer. The complementary strain energy density is equal to the area, OQR, that supplements the area enclosed by the stress-strain curve to a rectangle (Fig. 3 and Fig. 4). Thus, * 0 iL u is written as *0 0 0iL i L FL i TLu u u       (47) By substituting of (5), (18), (24) and (28) in (47), one obtains V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 167 2 * 0 1 ln ln i i i i iL L L L i i L i i i i i s s s s u s p p p p p p                     2 1 ln lni i i i i i i i i i i f f f f f g g g g g g                    (48) where  is determined by (28). The complementary strain energy in the un-cracked shaft portion as a result of centric tension and torsion is written as  * *0 1 i i i n H H i A U l a u dA      (49) where the complementary strain energy density, *0 iHu , is obtained by (48). For this purpose, * 0 iL u , L and  are replaced, respectively, with *0 iHu , H and H , where H is determined by (38). The expression, obtained by substituting of (45), (46) and (49) in (44), is doubled in view of the symmetry (Fig. 1). The result is 1 * * 0 0 1 1 1 i i i i n i n L H i ib A A G u dA u dA r                (50) Integration in (50) should be carried-out by the MatLab computer program. It should be noted that the strain energy release rate calculated by (50) is exact match of the strain energy release rate determined by (41). This fact verifies the analysis of the cylindrical delamination crack in the multilayered functionally graded circular shaft loaded in centric tension and torsion (Fig. 1). Shaft under bending and torsion The cylindrical delamination crack is analyzed also when the external loading consists of bending moments, M , and torsion moments, T , applied at the two ends of the multilayered functionally graded circular shaft (Fig. 5). Figure 5: Multilayered functionally graded circular shaft loaded in bending and torsion. Obviously, the bending moments induce mode II crack loading. By considering the balance of the energy, the mode II component of the strain energy release rate is derived as V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 168 1 2 2 2 M MII b b UM G r a r a           (51) where  is the angle of rotation of the end section of the shaft due to the bending, MU is the strain energy cumulated in half of the shaft as a result of the bending. It should be noted that the bending induces stresses not only in the un-cracked shaft portion and the internal crack arm, but also in the external crack arm. By using methods of Mechanics of materials,  is obtained as  L Ha l a     (52) where L and H are the curvatures of the crack arms and the un-cracked shaft portion, respectively. Since the bending generates mode II crack loading conditions, the two crack arms deform with the same curvature. Therefore, L is determined in the following way. First, the equation for equilibrium of the cross-section of the internal crack arm is used 1 1 1 i i n d i i A M z dA     (53) where dM is the bending moment in the internal crack arm. The distribution of the longitudinal normal stress, i , in the i-th layer, induced by the bending of the shaft, are expressed by (5). The distribution of the longitudinal strains,  , is written as 1L z  (54) By substituting of (5), (6) and (7) in (23), one derives     1 4 4 5 5 1 1 1 1 1 4 5 i n d L i i i L i i i i M r r r r                 (55) where i i i Bs    (56) i i D i B i s s    (57) The radiuses, ir and 1ir  , in (55) are shown in Fig. 6. It should be noted that at 0 iB p  and 0 iD s  Eqn. (55) transforms in   1 4 4 1 1 1 4 i i n L d i i i B M r r s      (58) which is exact match of the equation for equilibrium of multilayered circular shaft made of homogeneous linear-elastic layers loaded in bending [14] assuming that 1/ iB s is the modulus of elasticity in the i-th layer. V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 169 Figure 6: Cross-section of the internal crack arm loaded in bending and torsion. Figure 7: Two three-layered functionally graded circular shafts loaded in centric tension and torsion with cylindrical delamination crack located between (a) layers 2 and 3 and (b) layers 1 and 2. In (55), there are two unknowns, dM and L . One more equation with unknowns dM and L is derived by considering the equilibrium of the cross-section of the external crack arm. Obviously, (55) can be used as equation for equilibrium of the cross-section of the external crack arm. For this purpose, 1n has to be replaced with 2n ( 2n is the number of layers in the external crack arm). Besides, dM has to be replaced with dM M (this follows from the fact that the sum of the bending moments in the two crack arms is equal to M ). Thus, the equation for equilibrium of the cross-section of the external crack arm is written as     2 4 4 5 5 1 1 1 1 1 4 5 i n d L i i i L i i i i M M r r r r                  (59) Eqns. (55) and (59) should be solved with respect to L and dM by using the MatLab computer program. V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 170 Eqn. (55) is used also to determine H . For this purpose, dM , 1n and L are replaced, respectively, with M , n and H . The strain energy cumulated in half of the shaft as a result of the bending is obtained as M ML MQ MHU U U U   (60) where MLU , MQU and MHU are the strain energies in the internal crack arm, the external crack arm and the un- cracked shaft portion, respectively. Formula (16) is applied to determine MLU . For this purpose, FLU and 0 iFLu are replaced with MLU and 0 iMLu , respectively. The strain energy density, 0 iMLu , in the i-th layer of the internal crack arm as a result of the bending is obtained by formula (18). For this purpose, 0 iFLu and L are replaced with 0 iMLu and  , respectively (  is expressed by formula (54)). Figure 8: Two three-layered functionally graded circular shafts loaded in bending and torsion with cylindrical delamination crack located between (a) layers 2 and 3 and (b) layers 1 and 2. Formula (16) is used also to determine MQU by replacing of FLU , 1n and 0 iFLu with MQU , 2n and 0 iMQu , respectively. 0 iMQu is determined by replacing of L with  . MHU is found by formula (19). For this purpose, FHU and 0 iFHu are replaced with MHU and 0 iMHu , respectively. The strain energy density, 0 iMHu , in the i-th layer of the un-cracked beam portion as a result of bending is obtained by formula (18). For this purpose, 0 iFLu and L are replaced with 0 iMHu and b , respectively. The distribution of the longitudinal strains, b , in the cross-section of the un-cracked shaft portion is found by (54). For this purpose,  , L V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 171 and 1z are replaced with b , H and 2z , respectively ( 2z is the vertical centroidal axis of the cross-section of the un- cracked portion of the shaft). By substituting of (52), MLU , MQU , MHU and (60) in (51), one derives the following expression for the mode II component of the strain energy release rate as a result of the shaft bending:   1 2 0 0 0 1 1 1 1 i i i i i i i n i n i n MII L H ML MQ MH i i ib b A A A M G u dA u dA u dA r r                            (61) The mode III component of the strain energy release rate induced by the shaft torsion is obtained by formula (39). The total strain energy release rate is found by addition of (39) and (61). The result is   1 2 0 0 0 1 1 1 1 i i i i i i i n i n i n L H ML MQ MH i i ib b A A A M G u dA u dA u dA r r                             1 0 0 1 1 1 i i i i i n i n qm TL TH i ib b b A A T u dA u dA r r R r                         (62) Integration in (62) should be carried-out by the MatLab computer program. Figure 9: The strain energy release rate in non-dimensional form plotted against 1D s property (curve 1 – for the shaft configuration shown in Fig. 7a, curve 2 – for the shaft configuration shown in Fig. 8a, curve 3 – for the shaft configuration shown in Fig. 7b, curve 4 – for the shaft configuration shown in Fig. 8b). Formula (62) is verified by obtaining of G with the help of (44). The complementary strain energy cumulated in half of the shat as a result of bending and torsion is found as * * * *MTL MQ MTHU U U U   (63) where *MTLU , * MQU and * MTHU are the complementary strain energies in the internal crack arm, the external crack arm and the un-cracked shaft portion, respectively. It should be specified that *MQU is due to the bending only, since the external crack arm is not loaded in torsion. V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 172 Figure 10: The strain energy release rate in non-dimensional form plotted against 1 3 /B Bs s ratio: curve 1 – for the shaft configuration shown in Fig. 7a (linear-elastic solution), curve 2 – for the shaft configuration shown in Fig. 7a (non-linear solution), curve 3 – for the shaft configuration shown in Fig. 8a (linear-elastic solution), curve 4 – for the shaft configuration shown in Fig. 8a (non-linear solution). Formula (46) is applied to determine *MTLU . For this purpose, * LU and * 0 i Lu are replaced with * MTLU and * 0 i MTLu , respectively. Formula (48) is used to obtain the complementary strain energy density, *0 i MTLu , in the i-th layer of the internal crack arm. For this purpose, *0 i Lu and L are replaced, respectively, with * 0 i MTLu and  , where  is expressed by (54). * MQU is obtained by replacing of * LU , 1n and * 0 i Lu , respectively, with * MQU , 2n and * 0 i MQu in (46). * 0 i MQu is determined by replacing of *0 i Lu and L , respectively, with * 0 i MQu and  in (48) and taking into account the fact that the external crack arm is loaded in bending only. Thus, *0 i MQu is written as 2 * 0 1 ln ln i i i i i MQ i i i i i i i s s s s u s p p p p p p                    (64) The complementary strain energy cumulated in the un-cracked shaft portion as a result of bending and torsion is calculated by (49). For this purpose, *HU and * 0 i Hu are replaced with * MTHU and * 0 i MTHu , respectively. Formula (48) is used to obtain the complementary strain energy density, *0 i MTHu . For this purpose, * 0 i Lu , L and  are replaced, respectively, with *0 i MTHu , b and H , where H is found by (38). Finally, *MTLU , * MQU and * MTHU are added-up and substituted in (44). The result is 1 2 * * * 0 0 0 1 1 1 1 i i i i i i n i n i n MTL MQ MTH i i ib A A A G u dA u dA u dA r                     (65) The integration in (65) should be performed by the MatLab computer program. It should be noted that the strain energy release rate obtained by (65) is exact match of the strain energy release rate calculated by (62), which is a verification of the V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 173 delamination fracture analysis of the multilayered functionally graded non-linear elastic circular shaft loaded in bending and torsion. Figure 11: The strain energy release rate in non-dimensional form plotted against 1D p property at three 1 3 /B Bp p ratios for the shaft configuration shown in Fig. 7a (curve 1 – at 1 3 / 0.5B Bp p  , curve 2 – at 1 3/ 1B Bp p  and curve 3 – at 1 3/ 1.5B Bp p  ). PARAMETRIC INVESTIGATIONS arametric investigations of delamination fracture in the multilayered functionally graded non-linear elastic circular shaft are performed in order to elucidate the effects of material gradients, cylindrical delamination crack location, non-linear mechanical behavior of the material and load combinations. For this purpose, calculations of the strain energy release rate are carried-out by formulae (41) and (62). The results obtained are presented in non-dimension form by using the formula 3 /N BG G s R . Two three-layered functionally graded circular shafts loaded in centric tension and torsion are analyzed in order to elucidate the influence of the cylindrical delamination crack location on the fracture behavior (Fig. 7). A cylindrical delamination crack is located between layers 2 and 3 in the shaft configuration shown in Fig. 7a. A shaft with cylindrical delamination crack located between layers 1 and 2 is also under consideration (Fig. 7b). In both shaft configurations, the thickness of the layers is t (Fig. 7). It is assumed that 50T  Nm, 300F  N and 0.01t  m. In order to elucidate the influence of the load combination on the fracture behavior, two three-layered functionally graded circular shafts loaded in bending and torsion are also analyzed (Fig. 8). In the shaft shown in Fig. 8a, a cylindrical delamination crack is located between layers 2 and 3. A cylindrical delamination crack is located between layers 1 and 2 in the shaft configuration in Fig. 8b. The thickness of each layer is 0.01t  m in both shafts (Fig. 8). The loading is 50T  Nm and 40M  Nm. The strain energy release rate in non-dimensional form is presented as a function of 1D s material property ( 1D s controls the gradient of 1s property in layer 1) in Fig. 9 for the four shaft configurations shown in Fig. 7 and Fig. 8. It is assumed that 2 3 / 0.3D Ds s  , 1 3/ 2B Bs s  , 2 3/ 1.8B Bs s  , 1 3/ 0.3D Dp p  , 2 3/ 0.2D Dp p  , 3 3/ 0.4D Dp s  , 1 3/ 1.7B Bp p  , 2 3 / 1.9B Bp p  , 3 3/ 0.9B Dp p  , 1 3/ 0.4D Df f  , 2 3/ 0.3D Df f  , 3 3/ 0.5D Df s  , 1 3/ 1.7B Bf f  , 2 3/ 1.6B Bf f  , 3 3 / 1.1B Bf s  , 1 3/ 0.2D Dg g  , 2 3/ 0.5D Dg g  , 3 3/ 0.4D Dg s  , 1 3/ 1.5B Bg g  , 2 3/ 1.9B Bg g  and 3 3/ 0.8B Bg s  . Curves in Fig. 9 indicate that the strain energy release rate decreases with increase of 1D s . This behavior is due to the decrease of the shaft stiffness. It can also be observed in Fig. 9 that the strain energy release rate increases when the cylindrical delamination crack position is changed from this shown in Fig. 7a and Fig. 8a to that shown in Fig. 7b and Fig. P V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 174 8b. This finding is attributed to the increase of the stiffness of the internal crack arm. Fig. 9 shows also that the loading combination “bending and torsionˮ generates higher strain energy release rate in comparison with the loading combination “centric tension and torsionˮ for the considered values of F , T and M . Figure 12: The strain energy release rate in non-dimensional form plotted against 1D f property at three 1 3 /B Bf f ratios for the shaft configuration shown in Fig. 7a (curve 1 – at 1 3 / 0.7B Bf f  , curve 2 – at 1 3/ 1.2B Bf f  and curve 3 – at 1 3/ 1.7B Bf f  ). Figure 13: The strain energy release rate in non-dimensional form plotted against 1D g property at three 1 3 /B Bg g ratios for the shaft configuration shown in Fig. 7a (curve 1 – at 1 3 / 0.8B Bg g  , curve 2 – at 1 3/ 1.3B Bg g  and curve 3 – at 1 3/ 1.8B Bg g  ). The strain energy release rate in non-dimensional form is presented as a function of 1 3 /B Bs s ratio in Fig. 10 for the shaft configurations shown in Fig. 7a and Fig. 8a. One can observe in Fig. 10 that the strain energy release rate increases with increasing of 1 3 /B Bs s ratio (this due to the decrease of the shaft stiffness). Curves in Fig. 10 confirm the finding that when the shat is loaded in bending and torsion the strain energy release is higher in comparison with the case when the V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 175 shat is loaded in centric tension and torsion for the considered values of F , T and M . In order to evaluate the effect of the material non-linearity on the delamination fracture, the strain energy release rate obtained assuming linear-elastic behavior of the three-layered functionally graded shafts is plotted in non-dimensional form against 1 3 /B Bs s ratio in Fig. 10 for comparison with the strain energy release rate generated by the non-linear solution. It should be mentioned that that the linear-elastic solution to the strain energy release rate is derived by substituting of 0 iB p  and 0 iB g  , where 1, 2, 3i  , in formulae (41) and (62). It can be observed in Fig. 10 that the material non-linearity leads to increase of the strain energy release rate. The strain energy release rate in non-dimensional form is presented as a function of 1D p material property in Fig. 11 at three 1 3 /B Bp p ratios for the shaft configuration shown in Fig. 7a. The curves in Fig. 11 indicate that the strain energy release rate decreases with increase of 1D p . It can be observed also that the strain energy release rate increases with increasing of 1 3 /B Bp p ratio (Fig. 11). The influence of 1D f material property on the delamination fracture behavior is shown in Fig. 12. The shaft configuration in Fig. 7a is considered. One can observe that the strain energy release rate decreases with increasing of 1D f (Fig. 12). The curves in Fig. 12 show that the increase of 1 3 /B Bf f ratio leads to increase of the strain energy release rate. Figure 14: The /III IIG G ratio plotted against /T F ratio for the shaft configuration shown in Fig. 7a. The effect of 1D g material property and 1 3 /B Bg g ratio on the delamination fracture is illustrated in Fig. 13. The shaft configuration shown in Fig. 7a is analyzed. Fig. 13 shows that the strain energy release rate decreases when 1D g increases. It can be observed in Fig. 13 that the strain energy release rate increases with increasing of 1 3 /B Bg g ratio. The influence of the torsion moment - to - longitudinal force, /T F , ratio on the mode III component of the strain energy release rate - to - mode II component of the strain energy release rate, /III IIG G , ratio is shown in Fig. 14. The shaft configuration in Fig. 7a is considered. One can observe in Fig. 14 that /III IIG G ratio increases with increasing of /T F ratio. V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 176 Finally, the effect of the bending moment - to - torsion moment, /M T , ratio on the strain energy release rate is illustrated in Fig. 15. The shaft configuration shown in Fig. 8a is analyzed. The curve in Fig. 15 indicates that the total strain energy release rate increases with increasing of /M T ratio. Figure 15: The strain energy release rate in non-dimensional form plotted against /M T ratio for the shaft configuration shown in Fig. 8a. CONCLUSIONS delamination fracture analysis of multilayered functionally graded circular shaft is developed in terms of the strain energy release rate. The shaft is made of an arbitrary number of adhesively bonded concentric longitudinal layers which have different thicknesses and material properties. The material in each layer is functionally graded in radial direction. Besides, the material exhibits non-linear mechanical behavior in each layer. The continuous variation of material properties in radial direction is described by hyperbolic laws. A cylindrical delamination crack is located arbitrary between layers (the internal crack arm is a shaft of circular cross-section; the external crack arm is a shaft of ring-shaped cross-section). Two load combinations (centric tension and torsion, and bending and torsion) are investigated. These load combinations generate mixed mode II/III delamination fracture (the centric tension and bending generate mode II crack loading conditions, the torsion generates mode III crack loading conditions). The strain energy release rate is derived by analyzing the balance of the energy. In order to verify the solution obtained, the strain energy release rate is determined also by differentiating the complementary strain energy with respect to the crack area. Parametric investigations of the delamination fracture are carried-out in order to evaluate the effects of material gradients, the crack location, the material non-linearity and the load combinations. It is found that the strain energy release rate increases with increasing of 1 3 /B Bs s , 1 3/B Bp p , 1 3/B Bf f and 1 3/B Bg g ratios. The increase of 1Ds , 1Dp , 1Df and 1Dg leads to decrease of the strain energy release rate. Concerning the influence of the delamination crack location on the fracture behavior, it is found that the strain energy release rate decreases when the diameter of the cross-section of the internal crack arm increases. The analysis reveals that the non-linear mechanical behavior of the material leads to increase of the strain energy release rate. The comparison between the strain energy release rates generated by the two loading combinations shows that the strain energy release rate is higher when the shaft is loaded in bending and torsion for the considered values of the longitudinal force, bending and torsion moments. The results obtained in the present paper show that the strain energy release rate in A V. Rizov, Frattura ed Integrità Strutturale, 46 (2018) 158-177; DOI: 10.3221/IGF-ESIS.46.16 177 multilayered functionally graded circular shafts can be controlled by using appropriate material gradients in radial direction. The present paper is a contribution in delamination fracture mechanics of multilayered functionally graded circular shafts exhibiting non-linear mechanical behavior of the material under combined loads. REFERENCES [1] Hirai, T. and Chen, L. (1999). Recent and prospective development of functionally graded materials in Japan, Material Science Forum, 308-311, pp. 509-514. [2] Mortensen, A. and Suresh, S. (1995). Functionally graded metals and metal-ceramic composites: Part 1 Processing, International Materials Review, 40, pp. 239-265. [3] Nemat-Allal, M.M., Ata, M.H., Bayoumi, M.R. and Khair-Eldeen, W. (2011). Powder metallurgical fabrication and microstructural investigations of Aluminum/Steel functionally graded material, Materials Sciences and Applications, 2, pp. 1708-1718. [4] Neubrand, A. and Rödel, J. (1997). Gradient materials: An overview of a novel concept, Zeit f Met, 88, pp. 358-371. [5] Suresh, S. and Mortensen, A. (1998). Fundamentals of functionally graded materials, IOM Communications Ltd, London. [6] Bohidar, S.K., Sharma, R. and Mishra, P.R. (2014). Functionally graded materials: A critical review, International Journal of Research, 1, pp. 289-301. [7] Paulino, G.C. (2002). Fracture in functionally graded materials, Engineering Fracture Mechanics, 69, pp. 1519-1530. [8] Pan, S.-D., Feng, J.-C., Zhou, Z.-G. and Zhi, W.-L. (2009). Four parallel non-symmetric Mode –III cracks with different lengths in a functionally graded material plane, Strength, Fracture and Complexity: an International Journal, 5, pp. 143-166. [9] Tilbrook, M.T., Moon, R.J. and Hoffman, M. (2005). Crack propagation in graded composites, Composite Science and Technology, 65, pp. 201-220. [10] Rizov, V.I. (2016). Elastic-plastic fracture of functionally graded circular shafts in torsion, Advances in Materials Research, 5, pp. 299-318. [11] Lubliner, J. (2006). Plasticity theory (Revised edition), University of California, Berkeley, CA. [12] Chakrabarty, J. (2006). Theory of plasticity, Elsevier Butterworth-Heinemann, Oxford. [13] Lukash, P.A. (1998). Fundamentals of non-linear structural mechanics, Stroizdat. [14] Varvak, P.M. (1997). New methods in strength of materials, Vishta shkola. [15] Rizov, V.I. (2017). Analysis of longitudinal cracked two-dimensional functionally graded beams exhibiting material non-linearity, Frattura ed Integrità Strutturale, 41, pp. 498-510. << /ASCII85EncodePages false /AllowTransparency false /AutoPositionEPSFiles true /AutoRotatePages /None /Binding /Left /CalGrayProfile (Dot Gain 20%) /CalRGBProfile (sRGB IEC61966-2.1) /CalCMYKProfile (U.S. Web Coated \050SWOP\051 v2) /sRGBProfile (sRGB IEC61966-2.1) /CannotEmbedFontPolicy /Error /CompatibilityLevel 1.4 /CompressObjects /Tags /CompressPages true /ConvertImagesToIndexed true /PassThroughJPEGImages true /CreateJobTicket false /DefaultRenderingIntent /Default /DetectBlends true /DetectCurves 0.0000 /ColorConversionStrategy /CMYK /DoThumbnails false /EmbedAllFonts true /EmbedOpenType false /ParseICCProfilesInComments true /EmbedJobOptions true /DSCReportingLevel 0 /EmitDSCWarnings false /EndPage -1 /ImageMemory 1048576 /LockDistillerParams false /MaxSubsetPct 100 /Optimize true /OPM 1 /ParseDSCComments true /ParseDSCCommentsForDocInfo true /PreserveCopyPage true /PreserveDICMYKValues true /PreserveEPSInfo true /PreserveFlatness true /PreserveHalftoneInfo false /PreserveOPIComments true /PreserveOverprintSettings true /StartPage 1 /SubsetFonts true /TransferFunctionInfo /Apply /UCRandBGInfo /Preserve /UsePrologue false /ColorSettingsFile () /AlwaysEmbed [ true ] /NeverEmbed [ true ] /AntiAliasColorImages false /CropColorImages true /ColorImageMinResolution 300 /ColorImageMinResolutionPolicy /OK /DownsampleColorImages true /ColorImageDownsampleType /Bicubic /ColorImageResolution 300 /ColorImageDepth -1 /ColorImageMinDownsampleDepth 1 /ColorImageDownsampleThreshold 1.50000 /EncodeColorImages true /ColorImageFilter /DCTEncode /AutoFilterColorImages true /ColorImageAutoFilterStrategy /JPEG /ColorACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /ColorImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000ColorACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000ColorImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasGrayImages false /CropGrayImages true /GrayImageMinResolution 300 /GrayImageMinResolutionPolicy /OK /DownsampleGrayImages true /GrayImageDownsampleType /Bicubic /GrayImageResolution 300 /GrayImageDepth -1 /GrayImageMinDownsampleDepth 2 /GrayImageDownsampleThreshold 1.50000 /EncodeGrayImages true /GrayImageFilter /DCTEncode /AutoFilterGrayImages true /GrayImageAutoFilterStrategy /JPEG /GrayACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /GrayImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000GrayACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000GrayImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasMonoImages false /CropMonoImages true /MonoImageMinResolution 1200 /MonoImageMinResolutionPolicy /OK /DownsampleMonoImages true /MonoImageDownsampleType /Bicubic /MonoImageResolution 1200 /MonoImageDepth -1 /MonoImageDownsampleThreshold 1.50000 /EncodeMonoImages true /MonoImageFilter /CCITTFaxEncode /MonoImageDict << /K -1 >> /AllowPSXObjects false /CheckCompliance [ /None ] /PDFX1aCheck false /PDFX3Check false /PDFXCompliantPDFOnly false /PDFXNoTrimBoxError true /PDFXTrimBoxToMediaBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXSetBleedBoxToMediaBox true /PDFXBleedBoxToTrimBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXOutputIntentProfile () /PDFXOutputConditionIdentifier () /PDFXOutputCondition () /PDFXRegistryName () /PDFXTrapped /False /CreateJDFFile false /Description << /ARA /BGR /CHS /CHT /CZE /DAN /DEU /ESP /ETI /FRA /GRE /HEB /HRV (Za stvaranje Adobe PDF dokumenata najpogodnijih za visokokvalitetni ispis prije tiskanja koristite ove postavke. Stvoreni PDF dokumenti mogu se otvoriti Acrobat i Adobe Reader 5.0 i kasnijim verzijama.) /HUN /ITA /JPN /KOR /LTH /LVI /NLD (Gebruik deze instellingen om Adobe PDF-documenten te maken die zijn geoptimaliseerd voor prepress-afdrukken van hoge kwaliteit. De gemaakte PDF-documenten kunnen worden geopend met Acrobat en Adobe Reader 5.0 en hoger.) /NOR /POL /PTB /RUM /RUS /SKY /SLV /SUO /SVE /TUR /UKR /ENU (Use these settings to create Adobe PDF documents best suited for high-quality prepress printing. Created PDF documents can be opened with Acrobat and Adobe Reader 5.0 and later.) >> /Namespace [ (Adobe) (Common) (1.0) ] /OtherNamespaces [ << /AsReaderSpreads false /CropImagesToFrames true /ErrorControl /WarnAndContinue /FlattenerIgnoreSpreadOverrides false /IncludeGuidesGrids false /IncludeNonPrinting false /IncludeSlug false /Namespace [ (Adobe) (InDesign) (4.0) ] /OmitPlacedBitmaps false /OmitPlacedEPS false /OmitPlacedPDF false /SimulateOverprint /Legacy >> << /AddBleedMarks false /AddColorBars false /AddCropMarks false /AddPageInfo false /AddRegMarks false /ConvertColors /ConvertToCMYK /DestinationProfileName () /DestinationProfileSelector /DocumentCMYK /Downsample16BitImages true /FlattenerPreset << /PresetSelector /MediumResolution >> /FormElements false /GenerateStructure false /IncludeBookmarks false /IncludeHyperlinks false /IncludeInteractive false /IncludeLayers false /IncludeProfiles false /MultimediaHandling /UseObjectSettings /Namespace [ (Adobe) (CreativeSuite) (2.0) ] /PDFXOutputIntentProfileSelector /DocumentCMYK /PreserveEditing true /UntaggedCMYKHandling /LeaveUntagged /UntaggedRGBHandling /UseDocumentProfile /UseDocumentBleed false >> ] >> setdistillerparams << /HWResolution [2400 2400] /PageSize [612.000 792.000] >> setpagedevice