DOCTORA L T H E S I S
Department of Applied Physics and Mechanical Engineering Division of Polymer Engineering
Numerical Analysis of the Single Fiber Fragmentation Test Including
the Effect of Interfacial Friction
Enrique Graciani
ISSN: 1402-1544 ISBN 978-91-7439-165-7 Luleå University of Technology 2010
Enrique Graciani Numerical Analysis of the Single Fiber Fragmentation Test Including the Effect of Interfacial Friction
ISSN: 1402-1544 ISBN 978-91-7439-XXX-X Se i listan och fyll i siffror där kryssen är
Numerical Analysis of the Single Fiber Fragmentation Test Including
the Effect of Interfacial Friction
by
Enrique Graciani
Division of Polymer Engineering
Department of Applied Physics and Mechanical Engineering Luleå University of Technology
S-971 87 Luleå, Sweden December 2010
Printed by Universitetstryckeriet, Luleå 2010 ISSN: 1402-1544
ISBN 978-91-7439-165-7 Luleå 2010
www.ltu.se
To Carmen, my wife, and our son Enrique
Abstract
A numerical analysis of the single fiber fragmentation test is presented in this thesis. The aim of the numerical study is elucidating the influence of friction in the growth of the debond cracks along the fiber-matrix interface.
Two distinct approaches and numerical techniques have been employed in the analysis. Firstly, debond growth is studied using Linear Elastic Fracture Mechanics for interface cracks. Thus, an accurate solution of the near-tip stresses and displacements is needed, which is obtained using a self-developed numerical tool, based on the Boundary Element Method.
Energy Release Rate associated to debond growth is evaluated using a modified version of Irwin’s Virtual Crack Closure Technique, which takes into account the presence of frictional stresses along the crack faces.
Numerical results show that friction opposes to debond growth and, consequently, it has to be taken into account if quantitative values of the interfacial fracture toughness are going to be derived from test results.
In the abovementioned approach, numerical solutions are obtained neglecting the influence of load history on the frictional contact solution.
For that reason, a second approach has been employed in which debond onset and growth is simulated using cohesive elements, thus following the actual load history of the sample.
A commercial Finite Elements code has been employed for the numerical simulation in this approach, since there is no need for such a high refinement in the mesh as in the previous analyses.
Numerical results confirm that the influence of load history is negligible in the analysis, since analogous results are obtained with both approaches.
Preface
This thesis is a product of the research that I have carried out both in the Group of Elasticity and Strength of Materials of the University of Seville and the Division of Polymer Engineering of the Luleå University of Technology.
First of all, I wish to thank Professor Janis Varna for giving me the opportunity of carrying out these studies and for all the help and advice provided during this period. This gratitude is also extended to Professor Lars Berglund, Professor Federico París and Professor Vladislav Mantiþ who’s supervising have been highly instructive.
I whish to express my gratitude also to all the members of both research groups for their friendship and support, that has made the work much easy and enjoyable. My relationship with the Division of Polymer Engineering started with my stays in Luleå during the years 1996-1997, and I must make a special mention to Bobs, Kristofer, Anders, Johan, Fredrik, Mats, Lennart, Kristina, Christian, Conny, Marta-Lena, Leif, Erik… and all the people I met in those visits for being so kind and for making my stay so warm despite of the outdoor temperatures.
I’m glad to see that collaboration and contact have been maintained during these years and I hope it will continue growing in the future.
Seville, September 2010.
Enrique Graciani
Contents
1 Introduction ... 1
1.1 Micromechanical tests employed for interface characterization ... 2
1.2 Description of the Single Fiber Fragmentation Test ... 3
1.3 Theoretical aspects of the Single Fiber Fragmentation Test ... 4
1.4 Summary and Objectives ... 5
2 Description of the BEM formulation ... 6
2.1 Axisymmetric BEM formulation of the elastic problem ... 8
2.2 Weak boundary formulation of the coupling conditions ... 9
2.3 Solution of the final system of equations ... 10
3 Preliminary BEM analyses ... 11
3.1 Analysis of the distinct failure mechanisms ... 11
3.2 Analysis of the debond growth in absence of friction ... 16
3.3 Accuracy of the numerical results ... 20
4 BEM analysis of the effect of friction on debond growth ... 21
4.1 Near tip solution ... 22
4.2 Energy Release Rate and debond propagation ... 24
4.3 Axial stresses in the fiber ... 25
4.4 Concluding remarks ... 26
5 Cohesive zone analysis of the effect of friction on debond growth ... 27
5.1 Debond propagation in absence of friction ... 29
5.2 Debond propagation in presence of friction ... 31
6 Concluding remarks ... 33
7 Ongoing and future work ... 34
7.1 Evaluating ERR from near tip solution in presence of friction ... 34
7.2 Comparison of different approaches ... 35
7.3 Analysis of samples with other constituents ... 36
7.4 Comparison with experimental results ... 37
8 Acknowledgements ... 38
9 References ... 39
10 List of relevant own publications ... 42
10.1 Journal papers ... 42
10.2 Book chapters ... 43
10.3 Conference papers ... 43
Annex 1:
Graciani E, Mantiþ V, París F and Blázquez A. “Weak formulation of axi- symmetric frictionless contact problems with boundary elements.
Application to interface cracks”. Computers and Structures, 83:836-855, 2005.
Annex 2:
Graciani E, Mantiþ V, París F and Varna J. “Single fiber fragmentation test.
A BEM analysis”. Collection of Technical Papers -
AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics and Materials Conference, 2:988-997, Norfolk (Virginia), United States, 2003.
Annex 3:
Graciani E, Mantiþ V, París F and Varna J. “A BEM analysis of the single fibre fragmentation test. Influence of curing temperature in crack
propagation”. Proceedings of the 11th European Conference on Composite Materials (CD, File: A030.pdf), Rodas, Greece, 2004.
Annex 4:
Mantiþ V, Graciani E, París F and Varna J. “An Axisymmetric Boundary Element Analysis of Interface Cracks in Fiber Reinforced Composites”, in Advances in Boundary Element Techniques VI. Eds: AP Selvadurai, CL Tan ana MH Aliabadi, 21-26, Montreal, Canada, 2005.
Annex 5:
Graciani E, Mantiþ V, París F and Varna J. “Numerical analysis of debond propagation in the single fibre fragmentation test”. Composites Science and Technology, 69:2514-2520, 2009.
Annex 6:
Graciani E, Blázquez A, París F and Varna J. “Numerical analysis of the single fibre fragmentation test using cohesive elements”. Proceedings of the 14th European Conference on Composite Materials (Paper ID: 227- ECCM14), Budapest, Hungary, 2010.
Annex 7:
Erratum
1 Introduction
Among the wide range of fibrous composite materials that are used nowadays fiber reinforced composites are one of the most extended, especially for building primary structures subjected to high loads.
A typical composite of this kind is formed by a continuous matrix (usually constituted by a polymeric resin) and a fibrous reinforcement (mainly carbon or glass fibers) distributed within it. The resulting material combines a high stiffness and strength with a light weight, with fiber reinforcement being responsible for the high mechanical properties and the matrix providing consistency to the material.
Fiber reinforced composite materials are usually designed to work in the direction of the fibers, being their failure properties in this case mainly controlled by the strength properties of the fibers. However, the existence of differently oriented plies inside the laminates and/or the presence of impact loads, which produce stresses transferred along many directions within the laminate, are responsible for the appearance of failure mechanisms in which the failure properties of the matrix and the interface play a fundamental role.
Many experimental techniques have been developed to study the failure mechanisms of the interface. In most of them, the sample is constituted by a limited amount of fibers surrounded by polymeric matrix. These tests are generally known as micromechanical tests due to the small size of the samples.
1.1 Micromechanical tests employed for interface characterization
Among the micromechanical tests in which the sample contains a single fiber to which the load is directly applied, the following can be cited: pull- out test, push-out test, microbond test and microindentation test (see Figure 1).
Other widespread micromechanical tests, in which the load is applied to the sample, are the single fiber fragmentation test and the transverse single fiber test, depicted in Figure 2.
The main reason for developing such micromechanical tests is to obtain a clear failure of the interface within the sample to gain in understanding on interface failure mechanisms and properties.
All the above mentioned micromechanical tests present a common disadvantage associated to the small size of the samples, which makes it extremely difficult to obtain accurate measurements of the elastic state (displacements, deformations or stresses) within the samples during the
(a) (b) (c) (d) Figure 1: Micromechanical tests for studying interfacial failure properties:
(a) pull-out, (b) push-out, (c) microbond, (c) microindentation.
(a) (b) Figure 2: (a) single fiber fragmentation test, (b) transverse single fiber test.
tests. To overcome these difficulties, analytical or numerical solutions are employed during data reduction to obtain quantitative results of interfacial failure properties. Therefore, the validity of the results is clearly influenced by the assumptions made to obtain the elastic solution employed for the data reduction.
Single fiber fragmentation test has been chosen as the object of the analysis for being one of the most frequently used test for characterizing the interface of composite materials.
1.2 Description of the Single Fiber Fragmentation Test
Since it was introduced by Kelly and Tyson[1], single fiber fragmentation test has been extensively used for fiber-matrix interface characterization. Fragmentation samples consist of a sufficiently long fiber embedded in a resin matrix, subjected to tensile load. As maximum allowable strain is much lower in the fiber than in the matrix, and the fiber strength is characterized by a two-parameter Weibull distribution, an increase in the applied load, after reaching a certain value, results in successive breaking of the fiber. After fragmentation starts, the tensile load, applied by the grips at the ends of the sample, is transferred to the fiber fragments through shear stresses at the fiber-matrix interface.
If the failure properties of the interface and the matrix are high, fragmentation continues until a critical fragment length is reached where the fragments are not long enough to transfer through the interface an axial load sufficient to cause subsequent fiber breakages. This fact is usually referred to as saturation. In this case, in which the only appreciable failure in the sample takes place in the fiber, the final average length is mainly dependent on fiber failure properties.
Nevertheless, in most cases the above described saturation is not reached during the test because of the appearance of a series of cracks which arise from the ends of fiber fragments and grow through the matrix or the interface. There are three different crack paths observed in the experiments (see, for example, Ohsawa et al.[2] and Sjögren et al.[3]): a penny-shape crack and/or a bi-conical crack growing through the matrix, or a bi-cylindrical crack, growing along the interface, which is often referred to as debond crack. The appearance of these cracks implies a relaxation of the axial stress in the fiber, which prevents the appearance of new fiber fragments, and therefore the final average length of the fragments also depends on the interface or matrix failure properties.
1.3 Theoretical aspects of the Single Fiber Fragmentation Test
Traditionally, two distinct approaches (based on measuring interface strength or interface fracture toughness) have been employed for characterizing interface failure properties from the results of the single fiber fragmentation test.
The interface strength based approach relies on the initial shear lag stress distribution model proposed by Cox[4]. The correlation between the interface strength and the fiber strength was initially obtained by a simple equilibrium of the final fragment assuming constant interfacial stresses. It has been corrected afterwards to take into account several features like more realistic shear stress distributions on the interface, thermal residual stresses or statistical distribution of fragment lengths (see, for example, Ohsawa et al.[2], Whitney and Drzal[5], Netravali et al.[6], Henstenburg and Phoenix[7], Lacroix et al.[8-9], and Tripathi and Jones[10]). Interesting reviews in the early developments of this approach have been carried out, among
others, by Verpoest et al.[11], Herrera-Franco and Drzal[12], and Tripathi and Jones[13].
According to the shear lag model, the shear stress at the interface is finite and hence the interface failure takes place when its maximum value equals the interface strength. This criterion can be expressed as a certain relationship between the interface strength and the fiber strength and length at saturation.
The interface fracture toughness based approach, assumes that debond propagation would take place when the total energy released per unit length during the propagation (i.e., the energy release rate) equals the work needed to create the new unit debond crack surfaces (see, for example, the works of Nairn[14], Varna et al.[15], Nairn and Liu[16] and Wu et al.[17-20]). As the propagation takes place under pure mode II, failure of the interface is characterized by the mode II interface fracture toughness, which is the critical value of the mode II energy release rate (usually denoted as G ). IIc
Finally, the works carried out by Galiotis[21-22] and Huang and Young[23], measuring fiber strain by means of laser Raman spectroscopy, and the phase-stepping photoelasticity analysis made by Zhao et al.[24] can be cited as an example of the few existing experimental techniques which are able to provide direct stress or strain measurements during the execution of the test.
1.4 Summary and Objectives
The main objective of the present study is to carry out a numerical simulation of the single fiber fragmentation test, to show the influence of interfacial friction in the debond propagation. For that reason, a self- developed tool based on the Boundary Element Method has been
employed, which provides an accurate solution of the displacement and stresses within the sample.
The methods employed for developing the numerical tool and the details of the analyses are reported in Annexes 1-6. As some of the graphs included in the original publications contained misprints, corrected versions of them can be found in Annex 7.
In Annex 1, the axisymmetric BEM formulation and the weak boundary formulation of the coupling conditions along the interface developed for the numerical analyses are explained in detail. The results obtained in the numerical analyses of the single fiber fragmentation test are presented in Annexes 2-5.
Finally, in Annex 6, a numerical analysis of the problem using cohesive elements (carried out with a commercial Finite Elements code) is presented and compared with the previous BEM solution.
2 Description of the BEM formulation
In order to carry out the numerical analyses, collocation Boundary Element Method (BEM) for elastic problems with axial symmetry has been implemented in a self-developed numerical code. The fundamentals of BEM can be found, among others, in the books of Bakr[25], Baláš et al.[26]
and París and Cañas[27]. Coupling conditions along the bonded part and contact conditions along the debonded part of the fiber-matrix interface have been imposed using a weak boundary formulation, analogous to the one developed in Blázquez et al.[28] for plane elastic problems, which allows non-conforming discretizations to be used along the interface. This approach has been chosen for the numerical analysis of the single fiber fragmentation test since, in this particular case, it offers clear advantages if
compared with other alternative numerical methods (like, for example, the Finite Element Method).
Firstly, as the elastic solution of the single fiber fragmentation test presents axial symmetry with respect to the fiber axis, only the radial section of the model has to be analyzed. As BEM requires a mesh of the boundary of the domain under study, only radial section of the boundaries and interfaces of the sample need to be meshed. Therefore, a 1D BEM mesh will be sufficient for the analysis, being very easy to refine the mesh in the vicinity of the crack tips.
Secondly, as it is well known, crack propagation can be analyzed employing the near-tip elastic solution along crack faces and ahead of the crack tip. Since the debond crack in which the present analysis is focused grows along the fiber-matrix interface, the displacements and tractions needed for the study of crack propagation are the actual primary unknowns in the nodes of the BEM mesh.
Finally, the use of a BEM approach along with a weak boundary formulation of the coupling conditions permits the use of non conforming meshes in the interfaces. As the aspect ratio of the fiber is very high (since it is much longer than its radius) a fine mesh is required to obtain accurate results within the fiber. The use of a non conforming mesh along the fiber- matrix interface yields an important reduction on the number of elements needed for the boundary mesh in the matrix, where larger elements can be employed since sample radius is much larger than fiber radius.
The main and novel aspects of the numerical formulations presented in Annex 1 are briefly described in the following subsections.
2.1 Axisymmetric BEM formulation of the elastic problem
The axisymmetric collocation BEM formulation of the elastic problem[25-26] is based on the application of the axisymmetric Boundary Integral Equation (BIE) of the displacements for a point at the boundary (* ) of the domain (:) under study.
For a collocation point outside the symmetry axis this BIE can be written as:
³
* ED DJ J ED D
D
ED 'u U O 't T 'u d*
C* (x) (x) [ * (x,y) (y) (y) * (x,y) (y)] (1)
where x* is a fixed collocation point and y* is the integration point which moves along the boundary. 'uD( y), with D r,z, represent the global components of the displacements in cylindrical coordinates, while
)
J( y
't , with J n,s, are respectively the normal and tangential components of the tractions associated to a plane tangent to the boundary at
y . )ODJ( y is a rotation matrix relating global ( zr, ) and intrinsic ( sn, ) coordinates. )CED* (x are the components of the free term, which depend upon the local geometry of the boundary at x . Finally, U*ED(x,y) and
) ,
* (x y
TED are, respectively, the components of the axisymmetric fundamental solution and their associated tractions along the boundary. As
) ,
* (x y
UED and TED* (x,y) are singular functions when yo , the boundary x integral in (1) only exists in the sense of Cauchy Principal Value.
The main difficulty in the development of the axisymmetric BEM formulation is associated to the fact that U*ED(x,y) and TED* (x,y) are long non-explicit singular expressions written in terms of the complete elliptic integrals of the first and second kind. Series expansions of U*ED(x,y) and
) ,
* (x y
TED in the vicinity of the collocation point are presented in Annex 1, in which all non regular terms are identified. The novel analytical treatment of these non regular terms yields a high accuracy in the results.
For a collocation point xˆ at the symmetry axis it holds that 'ur(xˆ) 0, and the BIE of the axial component of the displacements can be written as:
³
* D DJ 'J D ' D *'u U O t T u d
Cˆzz*(xˆ) z(xˆ) [ˆz* (xˆ,y) (y) (y) ˆz* (xˆ,y) (y)] (2)
where )Uˆ*zD(xˆ,y and Tˆz*D(xˆ,y) can be obtained taking the limit of )
,
*D(x y
Uz and Tz*D(x,y) when the collocation point x tends to the symmetry axis. The free term Cˆ*zz(xˆ) is usually evaluated indirectly using numerical techniques. A novel analytical expression of Cˆzz*(xˆ) is presented in Annex 1, and its use has resulted in an increased accuracy of the results.
The collocation BEM employed for the numerical analysis is based in the application of equations (1) and (2) in a series of nodes in the boundary of the domain under study. Boundary integrals are approximated by the sum of the integrals along linear elements defined between the nodes.
Displacements and tractions within the elements are approximated by linear interpolation of their nodal values. Therefore, primary unknowns of the collocation equations are the components of the displacements and tractions at the nodes placed along the boundary of the solids under study.
2.2 Weak boundary formulation of the coupling conditions
Coupling of the collocation equations for different solids is carried out employing a novel weak boundary formulation of the equilibrium and compatibility equations at the interfaces and contact zones.
Let us consider domains :A and :B which are in contact along a part of their boundaries denoted by *CA and *CB. Equilibrium, that is
0 ) ( )
(y t y tB ' A
' , is guaranteed by the fulfillment of the following weak equation extended to *CB:
0 )
( )]
( ) (
[
³
*CB tB tA uB yrd CB D\D
D ' ' *
' y y y (3)
for any field of compatible displacements 'uB\( y).
Analogously, compatibility of the displacements along the contact region, that is 'įA(y) 'uA(y)'uB(y) with 'įA( y) being the relative displacements (opening and sliding) between opposite points of the contacting solids, is guaranteed by the fulfillment of the following weak equation extended to *CA:
0 )]
( )
( ) ( [ )
(
³
*CA tA\ uDA uDB DA yrd CAD ' ' 'G *
' y y y y (4)
for any tractions field in equilibrium 'tA\( y).
Employing the same boundary elements described in the previous section, with linear interpolations of the unknowns, to discretize equations (3) and (4) yields a system of equations which allow the coupling of the collocation equations corresponding to domains :A and :B to be made, since the primary unknowns of the equilibrium and contact equations (displacements and tractions along the contacting boundaries) are in fact a subset of the primary unknowns of the collocation equations.
The primary unknowns of the BEM formulation along the contact zone are the components of the displacements, the relative opening and sliding and the interfacial normal and tangential stresses. Therefore, both frictionless (see Annex 1) and frictional[29] contact conditions can be easily introduced in the formulation.
2.3 Solution of the final system of equations
It has to be emphasized that the final system of equations is constituted by the collocation equations described in section 2.1, the equilibrium and compatibility equations described in section 2.2 and the boundary and coupling conditions along the external boundaries, the interfaces and the contact zones. Accordingly, the final unknowns are the non prescribed
components of the displacements and tractions along the whole boundaries and the opening and sliding along the contact zones.
Since in many situations the sizes of the contact zones (and the adhesion or sliding regions inside them) are not known a priori, and usually depend upon the applied load, the corresponding contact boundary conditions are expressed by inequalities which make the problem to become nonlinear. The final solution of the system of equations is achieved following an incremental scheme which assures that all contact conditions are fulfilled during the whole load history. For this reason, in some cases, iterative stages are needed until a correct definition of the sizes of the contact zones is found which allows the solution procedure to be continued.
3 Preliminary BEM analyses
In this section, the results obtained in Annexes 2-4 are briefly described. In the majority of these preliminary results the effect of interfacial friction and/or the effect of thermal residual stresses resulting from the manufacturing of the sample are neglected.
3.1 Analysis of the distinct failure mechanisms
Comparison of the different kinds of failure of the sample in the single fiber fragmentation test is carried out in Annex 2, in terms of the energy release rate (ERR). Therefore, results of the computed ERR associated to the growth of a crack in the different orientations described (i.e., a penny- shape matrix crack, a bi-conical matrix crack or a bi-cylindrical interface crack) are presented for a range of crack lengths and orientations and fragment sizes. These crack configurations are depicted in Figure 3.
In this analysis, the effect of the thermal residual stresses and the friction between crack faces in closed cracks has been neglected. Analyses of the matrix cracks are carried out using a classical Fracture Mechanics approach, while the interface crack is studied using two distinct approaches: the open crack model introduced by Williams[30] (in which traction free crack faces are assumed) and the frictionless contact model developed by Comninou[31] (in which frictionless contact conditions are imposed along crack faces).
ERR is evaluated in all cases with an appropriate formulation of the virtual crack closure technique (VCCT). Material properties employed correspond to a typical E-glass fiber with an epoxy matrix.
Due to the large aspect ratio of the fiber no interaction between the cracks appearing in the sample is assumed. Therefore, in view of the existing symmetries, the model consists in the radial section of the portion of sample corresponding to one half of a typical fragment, with Lf being half of the fragment length and rf being the fiber radius.
a x y U T
a x y T U
x a U y T
(a) (b) (c) Figure 3: Configurations considered: (a) penny-shape matrix crack, (b) bi-
conical matrix crack and (c) bi-cylindrical interface crack (debond crack).
Results of the ERR, corresponding to the growth of a penny-shape matrix crack, are plotted versus crack length ( a ) in Figure 4 for selected fragment lengths.
For a conical crack, an initial straight growth oriented forming a certain angle D respect to the fiber axis is assumed. Results of the total ERR and its components are plotted versus crack orientation in Figure 5(a) for a crack with a and rf Lf 40rf.
Evolution of the ERR with crack length ( a ) is shown in Figure 5(b) for the particular case in which D 20º. Curves like the ones shown in Figure 5, along with the corresponding toughness of the matrix and the elastic solution in the vicinity of the crack tip can be used to predict the path followed by the conical crack.
As no quantitative information about the interface failure properties can be obtained when the previously described matrix failures occurs in the sample, the focus of this study has been set in the analysis of the case in which a bi-cylindrical debond crack grows along fiber-matrix interface.
0 10 20 30 40 50 60 70
0 2 4 6 8
Crack size (a/rf)
] 80
[ f f
I L r
G
] 40
[ f f
I L r
G
] 10
[ f f
I L r
G
] 20
[ f f
I L r
G
] 80
[ f f
I L r
G
] 40
[ f f
I L r
G
] 10
[ f f
I L r
G
] 20
[ f f
I L r
G
Energy release rate [GII,J/m2]
Figure 4: ERR in the penny-shape matrix crack (H 1%):
f
f a r
r 8
2 .
0 , 10rf Lf 80rf.
In this case, ERR has been evaluated using two different approaches. In Figure 6(a), results of the total ERR obtained with the open model is plotted as a function of the debond length ( a ), while Figure 6(b) shows the mode II ERR obtained using the frictionless contact model. It has to be emphasized that, although the asymptotic solution in both cases are clearly different, as described in Annex 2, results of the ERR are almost coincident in both approaches.
(a) 0 2 4 6 8 10 12 14
0 15 30 45 60 75 90
Crack orientation (D,degrees) GI
G GII
GI
G GII
Energy release rate [GII,J/m2]
(b) 0 5 10 15 20 25 30 35
0 1 2 3 4
Crack size (a/rf)
] 80 [Lf rf G
] 40 [Lf rf
G
] 10 [Lf rf
G
] 20 [Lf rf G
] 80 [Lf rf G
] 40 [Lf rf
G
] 10 [Lf rf
G
] 20 [Lf rf G
Energy release rate [GII,J/m2]
Figure 5: ERR in the bi-conical matrix crack (H 1%): (a) a , rf
f
f r
L 40 , º0º D 90 ; (b) 0.2rf a 4rf, 10rf Lf 80rf, ºD 20 .
As can be seen, for short debond lengths, ERR is decreasing steeply while, for larger debond lengths, it tends to a constant value, at least for longs fragments (having Lf t40rf). Therefore, debond growth will be stable in a first stage and will become nearly unstable when the plateau value of the ERR is being reached.
The plateau value can be estimated analytically studying the steady- state self-similar growth of a semi-infinite debond crack growing along the interface of a fiber fragment with infinite length. For the present case, the ERR evaluated numerically yields GII 8.85J/m2 for a 50Pm and
(a) 0 5 10 15 20 25 30
0 2 4 6 8
Crack size (a/rf)
] 80
[ f f
II L r
G
] 40
[ f f
II L r
G
] 10
[ f f
II L r
G
] 20
[ f f
II L r
G
] 80
[ f f
II L r
G
] 40
[ f f
II L r
G
] 10
[ f f
II L r
G
] 20
[ f f
II L r
G
Energy release rate [GII,J/m2]
(b) 0 5 10 15 20 25 30
0 2 4 6 8
Crack size (a/rf)
] 80
[ f f
II L r
G
] 40
[ f f
II L r
G
] 10
[ f f
II L r
G
] 20
[ f f
II L r
G
] 80
[ f f
II L r
G
] 40
[ f f
II L r
G
] 10
[ f f
II L r
G
] 20
[ f f
II L r
G
Energy release rate [GII,J/m2]
Figure 6: ERR in the debond crack (H 1%) with 0.2rf a 8rf ,
f f
f L r
r 80
10 : (a) open crack model; (b) contact crack model.
m 400P
Lf , which is in a very good agreement with the analytical plateau value GII 8.69J/m2 (for aof and Lf of)[32].
Although not included here for the sake of conciseness, the effect of debobd propagation on the axial stress along the fiber axis is also presented for all configurations, showing that, in all cases, the debobd propagation results in a relaxation of axial stresses along the fiber.
3.2 Analysis of the debond growth in absence of friction
The debond crack propagation is studied with greater depth in Annex 3 using the frictionless contact approach[31], taking also into consideration the effect of residual stresses due to the cooling stage after solidification of the sample.
As, both during the cooling stage and during the loading stage, the matrix radial shrinkage is larger than the fiber radial shrinkage, the debond crack is closed along its whole length. Moreover, if friction between crack faces is neglected, the solution can be obtained linearly combining two load cases: a first case in which only the decrease of temperature, 'T 80K, and its corresponding axial shrinkage, H'T |0.4%, are applied, and a second case in which a unit average axial strain, H 1%, is prescribed at constant temperature.
Therefore, the solution for a certain average mechanically applied axial strain, H, is given by
% 1 , 0 ,
80
% 1 , 0 ,
80
) , ( )
, ( ) , (
) , ( )
, ( ) , (
H ' H
H '
H ' H
H '
H H
' '
T K
T
T K
T
z r z
r z r
z r z
r z r
T T
ı ı
ı
u u
u (5)
Since notation is not uniform in the Annexes, a clear distinction is made in this document between the total axial strain of the sample, H0, which is equal to the average axial strain of the fiber, and the mechanically applied
axial strain, H. If thermal residual stresses are neglected both coincide, while if they are considered it holds that: H0 H'T H .
Firstly, the fragmentation of the sample is studied. In Figure 7(a) the axial stress along the fiber axis of a fragment with Lf 80rf , corresponding to a 1% average mechanically applied strain (H 1%), is shown in conjunction with the solution obtained after successively splitting of the original and the subsequent fragments by their middle points until
f
f r
L 10 is reached.
(a) 0
0.1 0.2 0.3 0.4 0.5
0 20
40 60
80
Axial stress [GPa]
Axial position [z/rf]
Szz(Lf=80rf)
Szz(Lf=40rf)
Szz(Lf=20rf)
Szz(Lf=10rf)
Lf= 10 rf Lf= 80 rf
H= 1%
Lf= 40 rf Lf= 20 rf
(b) 0
0.1 0.2 0.3 0.4 0.5
0 20 40 60 80 100 120
Max. axial stress [GPa]
Fragment size [Lf /rf]
H= 1%
Figure 7: (a) axial stress along the fiber during fragmentation; (b) maximum axial stress in the fragment.
As can be seen, when Lf t40rf the axial stress in the fragment reaches a plateau value in the middle of the fragment, in which the stress is equal to the nominal stress corresponding to a fragment having infinite length. On the contrary, shorter fragments do not reach this nominal stress, and therefore are subjected to lower stresses than the ones that caused the previous breakages.
In Figure 7(b) the maximum axial stress within the fragment is plotted against the fragment length. As can be seen, for Lf t40rf the maximum axial stress becomes constant. Therefore, it can be concluded that fragments having Lf t40rf have an equal transition zone in which the axial load in the fiber rises from cero to the nominal stress, and after this zone an uniform solution is obtained, which is almost equal to the one corresponding to the intact fiber.
To study the debond crack propagation, a set of configurations, in which the debond length is varied in the range 0.2rf a 8rf, has been analyzed (in fragments with Lf t40rf). The ERR obtained from the solution of these models, corresponding to a 1% average mechanically applied strain (H 1%), are shown in Figure 8(a).
As can be seen, comparing Figure 8(a) with Figure 6(b), the thermal residual stresses diminish the ERR, and therefore delay debond growth.
As the crack is closed and propagation occurs under pure mode II, a simple propagation criterion GII(H,a) GIIc has been employed for the analysis, where H is the average mechanically applied axial strain, a is the debond length and GIIc is the mode II interfacial fracture toughness. As shown in (5), the solution of displacement and stresses depend linearly upon H and, therefore, solution of GII( aH, ) is a second order polynomial in H whose coefficients can be obtained from the solution of the two load cases described at the beginning of this section.
Solving equation GII(H,a) GIIc, assuming a value of GIIc 50J/m2 yields the results shown in Figure 8(b), where the average strain needed to propagate a crack with a certain length is shown for two different cases:
considering and neglecting the thermal residual stresses.
As can be seen, crack growth is stable in an initial stage, followed by a nearly unstable (a vertical asymptote) crack propagation at higher strains.
Furthermore, these results clearly show that the presence of residual stresses delays crack propagation, as mentioned before.
(a) 0
2 4 6 8 10 12
0 1 2 3 4 5 6 7 8
] 80
)[
, (
f f a
II L r
GH
] 40 [
) , (
f f a
II L r
GH
] 80
)[
, (
f f a
II L r
GH
] 40 [
) , (
f f a
II L r
GH
% H 1
Energy release rate [J/m2]
Crack size, a/rf
(b) 0
1 2 3 4 5 6 7 8 9
0 0.5 1 1.5 2 2.5 3
] 0 , 80
) [
,
( G L r T K
G IIc f f
a
IIH '
] 0 , 40
) [
,
( G L r T K
G IIc f f
a
IIH '
] 80 , 80 [
) ,
( G L r T K
GIIHa IIc f f ' ] 80 , 40
) [
,
( G L r T K
GIIHa IIc f f ' ] 0 , 80
) [
,
( G L r T K
G IIc f f
a
IIH '
] 0 , 40
) [
,
( G L r T K
G IIc f f
a
IIH '
] 80 , 80 [
) ,
( G L r T K
GIIHa IIc f f ' ] 80 , 40
) [
,
( G L r T K
GIIHa IIc f f '
Crack size, a/rf
Applied strain, H [%]
J/m2 IIc 50 G
Figure 8: (a) ERR during debond crack propagation; (b) average strain needed to cause debond crack propagation.
Obtaining a series of curves like the one shown in Figure 8(b), with a parametric variation of GIIc, would permit the calculation of the actual interfacial fracture toughness by best fit of experimental debond propagation measurements. However, prior to the comparison with test results, the effect of friction between crack faces on debond propagation needs to be elucidated.
3.3 Accuracy of the numerical results
In order to show the accuracy of the results, a comparison of several aspects of the numerical solution of the BEM models with existing analytical or semianalytical results obtained from the literature is shown in Annex 4.
Firstly, the behavior of the BEM solution in the vicinity of the tip of the fiber crack terminating at the interface is compared with the behavior predicted by the semianalytical singularity analysis of multimaterial corners presented by Barroso et al.[33] According to this semianalytical analysis, when Uo0, with U being the distance to the crack tip, the following asymptotic behavior should be obtained: uD O(UO) and
) ( O1
DE U
V O , with D,E r,z and O 0.187. Numerical evaluation of the singularity exponent O carried out from the BEM near tip solution yields values in the range 0.190 O 0.198.
Secondly, the behavior of the BEM solution in the vicinity of the debond crack tip is compared with the asymptotical behavior corresponding to an interface crack. In particular, using the frictionless contact model, a singularity exponent in the range 0.501 O 0.511 is evaluated from the near tip stresses obtained with BEM, which is in very good agreement with the theoretical result O 0.5 predicted by Comninou[31].
Finally, according to the theoretical analysis presented by Leguillon and Sanchez-Palencia[34], limit behavior of the ERR when ao0, with a being the debond crack length, should be in accordance with the previously mentioned singularity exponent of the solution corresponding to the fiber crack terminating at the interface. Using the frictionless contact model, the BEM solution yields an asymptotic behavior G O(a0.609) which is again in a very good agreement with the expected theoretical behavior:
) ( )
(a2 1 O a0.626 O
G O .
4 BEM analysis of the effect of friction on debond growth
The propagation of the debond crack is studied in Annex 5, using an approach based on Linear Elastic Fracture Mechanics (LEFM) approach.
Numerical solution is obtained with the BEM formulation described in Section 2. The effects of interfacial friction and the thermal residual stresses resulting from the manufacturing of the sample are included in the analysis.
As in previous analyses, the material properties employed correspond to a typical E-glass fiber with an epoxy matrix. The main assumptions made in the numerical analysis are: linear elastic isotropic behavior of the constituents, axial symmetry (with respect to the fiber axis), repetitive solutions in the vicinity of all fragment ends and local symmetry with respect to the plane containing the fiber crack.
The model employed is depicted in Figure 9(a), where the main dimensions and the boundary conditions are indicated. A detail of the mesh employed is shown in Figure 9(b), where the advantages of using a boundary mesh can be clearly appreciated.
As described in section 3.2, the solution of the problem is obtained assuming a preexisting debond crack of length a , as shown in Figure 9(a), and linearly combining two load cases: a first case in which the decrease of temperature is applied (along with the corresponding average axial shrinkage) and a second case in which a unit average mechanically applied axial strain is applied at constant temperature. Therefore, the solution for a certain average mechanically applied axial strain, H, is given by (5).
4.1 Near tip solution
When interfacial friction between debond crack faces is considered, singular shear stresses are obtained both ahead of the crack tip and along crack faces, as can be appreciated in Figure 10(a). A friction coefficient
3 .
P 0 has been used for calculations.
In Figure 10(b), where a logarithmic scale is employed for clarity, the singular shear stresses in presence of friction (curve {3}) can be compared with the shear stresses in absence of friction (curve {2}).
It can be appreciated that friction has three distinct effects: the appearance of shear stresses along crack faces, a reduction of the shear
Lf
u0 H r
z
Lf
rf
rs
contact bonded
a matrix
fiber
0 5 10 15
-5 0 5 10 15
0 crack tip
z
r
fiber matrix
, Pm
, Pm
(a) (b) Figure 9: (a) Sketch of the BEM model; (b) detail of the BEM mesh.
stresses ahead of the crack tip and, what is more significant, a change on the singularity order of stresses.
As predicted by Comninou[35], when Uo0, with U being the distance to the crack tip, the following asymptotic behavior should be obtained:
) ( O
D O U
u and VDE O(UO1), where, for this particular problem, it holds that 5O!0. for any friction coefficient.
(a)
-150 -100 -50 0
0 10
20 30
Vrz(rf,z) [MPa]
Axial coordinate z [Pm]
debonded bonded
H=1%
'T = 80 K, P= 0.3, a = 10 Pm
(b)
1 E+1 1 E+2 1 E+3 1 E+4
1 E-4 1 E-3 1 E-2 1 E-1 1 E+0 1 E+1 Vrz(rf,arU) [MPa]
Distance to the crack tip U[Pm]
H=1%, a = 10 Pm {1}: 'T = 0 K, P = 0
{2}: 'T = 80 K, P = 0
^` 'T = 80 K, P = 0.3 Vrz(rf,a+U)
Vrz(rf,aU) 104
102 103
10
10-4 10-3 10-2 10-1 1 10
Figure 10: near-tip shear stresses (a) in presence of friction; (b) in presence/absence of friction and/or thermal residual stresses.
4.2 Energy Release Rate and debond propagation
Energy Release Rate has been numerically calculated in Annex 5 using a modification of Irwin’s Virtual Crack Closure Technique (VCCT) proposed by Sun and Qian[36]:
³
³
a
z rz
a
z rz
II
d a a r u a
a a r
d a a r u a
a r a
G
' '
U U ' '
U ' ' V
U U ' '
U ' V
'
0 0
) ,
( ) ,
2 ( 1
) ,
( ) , 2 (
) 1 ˆ (
(6)
where a' is the virtual crack increment considered, Vrz are the interfacial shear stresses in the vicinity of the debond crack tip (at both sides of the crack) and 'uz is the relative sliding between crack faces in the vicinity of the debond crack tip.
Since the crack is closed, a pure fracture mode II is obtained. In presence of friction, the energy released during crack propagation is different from that of the non-frictional case for two reasons. First, because both singular shear stresses ahead of the crack tip and relative displacement between crack faces change in the presence of friction, thus changing the value of the first integral in (6). Second, because after debond propagation the shear stress in the new created zone is non-zero due to friction. The second integral in (6) accounts for the frictional stresses along crack faces.
If fiber-matrix interfacial friction is considered, GˆII('a) vanishes when o0
'a , as a result of the change of the singularity order in the near-tip elastic solution[36]. A value 'a 0.012Pm has been considered in the calculations.
Once the ERR is evaluated, the following propagation criterion is employed to characterize the debond crack growth:
IIc
II a G
Gˆ (' ) (7)
