Home Contact Sitemap

DTD Handbook

Handbook for Damage Tolerant Design

  • DTDHandbook
    • About
    • Contact
    • Contributors
    • PDF Versions
    • Related Links
    • Sections
    • Examples
      • 1. FAC Problems
        • 0. Computing Stress Intensity Factor Histories
        • 1. Predicting Residual Strength - Fuselage Section
        • 2. Crack Interaction and Multi-Site Damage
        • 3. Predicting Ductile Tearing and Residual Strength - Flat Sheet
      • 2. Merc Problems
      • 3. NRC Problems
      • 4. SIE Problems
      • 5. UDRI Problems

Problem FAC-1

Title:  Computing Stress Intensity Factor Histories with the Finite Element Method and Interfacing with AFGROW and NASGRO

Objective:

To illustrate the process of using a finite element program, such as FRANC2D/L, to compute stress intensity factor histories accurately, so that they can be transmitted to AFGROW or NASGRO for fatigue crack growth rate (FCGR) and life predictions.

General Description:

This problem details the processes of using the finite element method to compute accurate relationships between stress intensity factors and crack length, and of communicating these relationships to AFGROW or NASGRO for further processing for FCGR and life predictions. Two problems are used for illustration: a typical detail in the standard library of AFGROW or NASGRO to confirm stress intensity factor accuracy, and a more complex problem for which no standard library solution is available.

 

Topics Covered:          Finite element analysis, stress intensity factor computation, non-standard geometry and/or boundary conditions, geometry factors for AFGROW or NASGRO

Type of Structure:       any single or multiple layer planar structure amenable to 2D or axisymmetric structural idealization

Relevant Sections of Handbook:  Sections 2, 5, 11

Author:                     Dr. A. R. Ingraffea

Company Name:   Fracture Analysis Consultants, Inc.

                                    121 Eastern Heights Drive
                                    Ithaca, NY  14850
                                    607-257-4970

Contact Point:       Dr. Paul Wawrzynek

Phone:       607-257-4970

e-Mail:      wash@fracanalysis.com

 

 


Overview of Problem Description

An essential task in predicting fatigue crack growth rate and remaining fatigue life is the accurate computation of a relationship between crack length and stress intensity factor.  In Example 1 of this problem a finite element code is used to compute such a relationship on a relatively simple structural detail, Figure FAC-1.1.  The pin load, P, is distributed as a normal pressure following a cosine function on the lower half-circumference of the hole. To assess the accuracy of the finite element calculations, they are compared to standard solutions stored in NASGRO for this detail.


Figure FAC-1.1.  Example 1 of this problem: a structural detail from standard library in NASGRO (Problem TC03) and AFGROW (Single Through Crack at Hole-Standard Solution).  Here W =10, D = 1, B = 5, t = 1, So = 10, and P = 10, and c/D is varied from 0.09 to 1.35.  All in kips and inches. Figure taken from NASGRO 3.

Frequently, the structural detail of interest is not contained in standard libraries.  In this case, either the actual detail must be idealized to resemble a standard solution, or a finite or boundary element model must be constructed. The finite element method can be used effectively in the role of general stress intensity factor calculator.  The method can accommodate virtually any actual situation involving non-standard geometry and boundary conditions, and include multiple cracks, multiple materials, material anisotropy, and initial stresses.

As an example of this situation, a variant of the detail shown in Figure FAC-1.1 is studied in Example 2 of this problem, Figure FAC-1.2. This figure shows a simple lug under non-symmetric loading. A contact-fit pin is inserted in the hole, and the pin load, P, is distributed to the hole by way of an elastic contact analysis.  A crack is initiated from a location of high tensile stress concentration along the lug bore, and then allowed to propagate in mixed mode.


Figure FAC-1.2.  Example 2 of this problem: a lug under unsymmetrical load. Lug and pin are both steel, E= 29,000ksi, n=0.30, and frictionless, contact-fit of the pin in the lug is assumed.  Thickness is 1.00 in. 

 

Computational Models: Example 1

Two finite element models shown in Figures FAC-1.3 and FAC-1.4 were used for Example 1 of this problem.  The purpose of using two meshes is to show the degree of mesh refinement that is needed to achieve accurate values of stress intensity factors using at least quadratic order triangular and quadrilateral elements. There are about 2300 degrees-of-freedom in the initial, crack-free mesh shown in Figure FAC-1.3, hereafter called the coarse mesh.

There are about 9800 degrees-of-freedom in the initial, crack-free mesh shown in Figure FAC-1.4, hereafter called the refined mesh.  The number of degrees-of-freedom in each model increases as the mesh is updated to accommodate crack growth. 

The mesh is refined in the region of crack growth to the right of the hole in each case.  A rule-of-thumb in creating the initial mesh in the region of expected crack growth is to make the average element size there smaller than (say one-half to one-quarter) the specified increment in crack length.  Examples of this practice are shown in Figures FAC-1.3(b) and FAC-1.4(b) that show a detail of the mesh around a crack tip.  FRANC2D/L automatically remeshes after each increment of crack growth. Figure FAC-1.5 shows a detail of the refined mesh after a few steps of crack growth.

Figure FAC-1.3.  Example 1. (a) Coarse finite element model for structural detail shown in Figure FAC-1.1.  (b) Detail of meshing near initial crack tip, c/D=0.09.

 

Figures FAC-1.3(b) and FAC-1.4(b) show some other important details of meshing around a crack tip.  These include the relative size, type, and distribution of elements immediately surrounding the crack tip.  It is well known that mixed-mode stress intensity factors can be computed with excellent accuracy using the equivalent domain formulation of the elastic J-Integral (Nikishkov 1987). FRANC2D/L can automatically compute stress intensity factors with this technique.  Other techniques, such as displacement correlation (Chan 1970) and modified crack closure (Rybicki 1977), are also available, but these have proven to be generally less accurate for a given discretization. This formulation of the elastic J-Integral technique does not require extreme levels of mesh refinement in the vicinity of a crack tip. 

 

Figure FAC-1.4.  Example 1. (a) Refined finite element model for structural detail shown in Figure FAC-1.1. (b) Detail of meshing near initial crack tip, c/D=0.09.

 

Experience has shown that the treatment shown in this example, a rosette of eight six-noded triangular elements immediately surrounding the tip, quarter-point versions of these elements, and element size ranging from a few to as much as 25 percent of crack length, will produce very accurate values of stress intensity factors.  The use of quarter-point versions of these elements produces the singular stress field known to exist in linear elastic fracture mechanics, is automatically done in FRANC2D/L, and is well documented in the literature (Barsoum 1976). 

Computational models like these require little effort to generate.  An experienced finite element user can produce one in less than one-half hour of person-time using codes such as PATRANtm or ProEtm.

 

 

Figure FAC-1.5.  Example 1.  Detail of refined finite element model for structural detail shown in Figure FAC-1.1 after a few steps of crack propagation, showing automatic remeshing.

 

Computational Results: Example 1

Table FAC-1.1 shows a comparison between the stress intensity factor histories from NASGRO and finite element analysis for the coarse mesh.  Each step of analysis took less than one second of total computer time on a PC running Windows 2000 on a 300MHz Pentium II processor. Shown are stress intensity factors due to the pin load, S3, and the bypass load, S0, and the total of these. The average difference in total stress intensity factor across the range of crack lengths analyzed was less than 1 percent.

Table FAC-1.2 shows a comparison between the stress intensity factor histories from NASGRO and finite element analysis for the refined mesh.  Each step of analysis took about 10 seconds of total computer time on a PC running Windows 2000 on a 300MHz Pentium II processor.  The average difference in total stress intensity factor across the range of crack lengths analyzed remained less than 1 percent, but did not improve over the results from the coarse mesh. This implies that the results from finite element analysis have converged to values acceptably close to those produced from the boundary element method used to generate the values in NASGRO. It also implies that relatively coarse meshes, such as those shown in Figure FAC-1.3, can be used to compute stress intensity factors accurately when used with the equivalent domain formulation of the elastic J-Integral. Finally, the person-effort and computer time needed to produce an accurate stress intensity factor history for such 2D details are not overwhelming.

 

Table FAC-1.1.  Example 1. Comparison between NASGRO and finite element analysis stress intensity factor histories for Example 1, using the coarse mesh shown in Figure FAC-1.3.

Table FAC-1.2.  Example 1. Comparison between NASGRO and finite element analysis stress intensity factor histories for Example 1, using the refined mesh shown in Figure FAC-1.4.

Computational Models: Example 2

Figure FAC-1.6 shows the initial finite element model for the lug problem shown in Figure FAC-1.2. The lug and its frictionless, contact-fit pin are both steel. The elastic contact problem between the pin and the lug is solved in FRANC2D/L.  Six-noded, zero-thickness interface elements with non-linear constitutive capability are inserted around the contact surface.  This process is automated: the lug and pin are assigned different material property set numbers (but the same properties), and a feature in FRANC2D/L finds the closed curve defining the boundary between these two sets and automatically inserts interface elements along the entire curve.

The constitutive model used to represent the normal contact conditions is shown in Figure FAC-1.7.  No tension is allowed across the contact, and a high compressive normal stiffness is assigned to minimize intrusion of the pin into the lug.

Figure FAC-1.7. Example 2.  Constitutive model for normal stress-displacement on the pin/lug contact surface.

 

A non-linear solver is now needed, and FRANC2D/L uses a dynamic relaxation procedure (Underwood 1983).  This is a relatively slow but extremely robust algorithm for solving nonlinear equilibrium equations. Figure FAC-1.8 presents the uncracked, deformed shape and shows that separation of the pin from the lug has occurred from about the 8 o’clock to the 2 o’clock positions, while non-overlapping contact occurs along the rest of the contact.  Figure FAC-1.9 shows contours of major principal stress for the initial, uncracked configuration.  This solution involved about 9200 DOF, and, with an error tolerance of 0.0005 on both equilibrium and displacement change between time steps, required about 5200 time steps and 2.5 minutes on a PC running Windows 2000 on a 1GHz Pentium III processor.

 

Figure FAC-1.6. Example 2.  Initial finite element mesh.

Figure FAC-1.8. Example 2.  Deformed shape of uncracked configuration.  Amplification factor is 225.


 


Figure FAC-1.9. Example 2.  Contours of major principal stress in uncracked condition.

Computational Results: Example 2

Figure FAC-1.8 indicates that there are two locations of high stress concentration around the lug bore, as expected.  A slightly higher concentration occurs in the lower left quadrant, at about the 8 o’clock position, and a short crack (0.067 in.) is initiated at the location of highest circumferential tensile stress near this point. Growth of this crack is then simulated by:

1.      Computing KI and KII using the equivalent domain formulation of the elastic J-Integral;

2.      Calculating the direction of crack growth using the maximum circumferential tensile stress theory (Erdogan 1963);

3.      Extending the crack by 0.075 inch in this direction;

4.      Resolving the finite element problem, including the non-linear contact between lug and pin;

5.      Repeating steps 1-4 until the crack has extended about 2.5 inch at which point fatigue life or residual strength limits are likely to have been reached.

Figure FAC-1.10 shows the corresponding amplified displaced shape, while Figure FAC-1.11 shows the resulting stress intensity factor histories. Figure FAC-1.10 shows that the crack trajectory is not quite radial, and is responding to the asymmetrical loading and geometry. The trajectory is dictated here by the maximum circumferential tensile stress theory that requires that KII remain zero along the crack path.  With a finite element model that discretizes the trajectory into finite, straight segments, there will always be residual, non-zero values of KII computed at each crack tip location. If the segments are short enough, these residuals should be small compared to the KI values.  Figure FAC-1.12 shows that the values of KII are oscillating around zero, and are indeed small, in this case, never reaching more than 3.5% of KI.  The highest values usually occur early in the trajectory while the finite element model is adjusting to the stress field that is evolving as a result of crack growth, as shown in Figure FAC-1.12.  The key practical issue suggested here is the length of crack growth increment.  This length should be sufficiently short to accurately discretize a curvilinear trajectory and provide enough data points for the accurate integration required for FCGR calculations, while not being so short that excessive computation times accrue.  Here 32 increments were used.  The number of DOF’s grew to nearly 14,000 at the last increment, and a total of about 2 hours of computing time was required.

Figure FAC-1.10. Example 2.  Final deformed shape of single cracked configuration.  Amplification factor is 150.

 

As previously mentioned, a finite element code with fracture mechanics features can be thought of as a general stress intensity factor calculator.  As such it can be used to attack practically interesting variants of problems.  For example, it is possible that two fatigue cracks might initiate in this lug problem, one from each of the locations of initially high stress concentration. This possibility is also simulated here, under the assumptions that

 

Figure FAC-1.11. Example 2.  Computed stress intensity factor histories. Single crack case.

Figure FAC-1.12. Example 2.  History of computed ratio of stress intensity factors.  Single crack case.

 

the cracks initiate simultaneously, and that they have equal rates of growth.  The resulting trajectories under these simplifying assumptions are shown in Figure FAC-1.13. The corresponding mode I stress intensity factor histories are given in Figure FAC-1.14.  This figure shows that, even if initiation were simultaneous, rate of growth would not be equal; the left crack would have higher growth rate.  However, even under these simplifying assumptions, Figure FAC-1.14 also shows that the growth rate of the left crack would be higher than it would be if it were the only crack to occur.

 

 

Figure FAC-1.13. Example 2.  Final deformed shape of multiply cracked configuration.  Amplification factor is 100.

 

Figure FAC-1.14. Example 2.  Comparison of stress intensity factor histories for single and multiple crack cases.
It is currently possible to couple the results of a stress intensity factor history analysis from finite element analyses for a single crack to FCGR simulators like NASGRO and AFGROW. This process is described next.

 

Interfacing with NASGRO and AFGROW

This section describes the process of interfacing the stress intensity factor history from finite element analyses with NASGRO (NASGRO 2000) or AFGROW (AFGROW 2000) to do FCGR and remaining life assessments.

 

Interfacing with NASGRO

In addition to its built in library of stress intensity factor solutions for many standard structural geometries, NASGRO offers the capability of performing crack growth analysis for non-standard geometries through its user defined data table (DT) approach.  In this approach NASGRO requests stress intensity factor input in the form of a one-dimensional table (NASGRO:Reference Manual, Appendix C):

 

a/D

F3

0.067

0.52

0.142

0.96

0.217

1.34

0.292

1.66

0.367

1.97

0.442

2.23

 

Table FAC-1.3. Portion of NASGRO stress intensity factor data table for a non-NASGRO-standard geometry. D= 4.0 in. and S3 = 7.5 ksi.

 

NASGRO’s general expression for stress intensity factors is,

NASGRO Eqtn. (2.2)

where the stress quantities S0, S1 and S2, S3, and S4 are for applied tension/compression, bending in the thickness and width directions, pin bearing pressures, and special cases, respectively. For Example 2 only S3 is non-zero, so NASGRO expects a relationship of the form

           

where F3 is the geometry correction factor for the present lug problem.  The values of F3 needed for the NASGRO data table are contained within the stress intensity factor values computed with FRANC2D/L, and can be explicitly obtained through

 

           

 

Table FAC-1.3 contains a portion of the geometry correction factor history obtained from FRANC2D/L for Example 2 with a single crack.  Herein D is the pin diameter, 4.0 in., and S3 is the pin bearing stress

 

            S3 = P/Dt = 30 kips/((4 in.)(1.0 in)) = 7.5 ksi

 

Once stress intensity factor history data in this form has been supplied, all standard NASGRO fatigue crack growth capabilities become available for any cracked structural model for which FRANC2D/L is an appropriate modeling system.

 

Interfacing with AFGROW

In addition to its built in library of stress intensity factor solutions for many standard structural geometries, AFGROW offers the capability of performing crack growth analysis for non-standard geometries through its user-defined Beta Table approach (AFGROW: Users Guide).  In this approach for a single through-crack, AFGROW requests stress intensity factor input in the form of a one-dimensional table of crack length versus stress intensity factor geometry correction, or beta, values, where beta is defined by

           

and where s is the relevant applied stress and c is crack length.  The beta values are obtained from finite element results in the same way as the equivalent Fi values are obtained for NASGRO. Figure FAC-1.15 shows the User-Input Beta Table dialogue box in AFGROW with the first few beta values obtained for Example 2 with a single crack, and with s set to unity.

Once stress intensity factor history data in this form has been supplied, all standard AFGROW fatigue crack growth capabilities become available for any cracked structural model for which FRANC2D/L is an appropriate modeling system.


 

Figure FAC-1.15. AFGROW dialogue box for creating user-input beta table, showing portion of table for Example 2 with a single crack.  Here s has been set to unity.

 

References

AFGROW: Users Guide and Technical Manual, AFGROW for Win95/98/NT4/2K, Version 4.0001011.8, AFRL-VA-WP-TR-2000-XXXX, May 2000.

Barsoum, R.S. (1976), "On the use of Isoparametric Finite Elements in Linear Fracture Mechanics," International Journal for Numerical Methods in Engineering, Volume 10 (1976) pp 25-38.

Chan, S.K, Tuba, I.S., and Wilson, W.K., "On the Finite Element Method in Linear Fracture Mechanics," Engineering Fracture Mechanics, Volume 2 (1970) pp 1-17.

Erdogan, F., and Sih, G.C., "On the Crack Extension in Plates Under Plane Loading and Transverse Shear," Journal of Basic Engineering, Volume 85 (1963) pp 519-525.

NASGRO: Reference Manual, Fatigue Crack Growth Computer Program “NASGRO” Version 3.0, JSC-22267B, March 2000.

Nikishkov, G.P. and Atluri, S.N. (1987), "Calculation of Fracture Mechanics Parameters from an Arbitrary Three-Dimensional Crack by the Equivalent Domain Integration Method," International Journal for Numerical Methods in Engineering, Volume 24 (1987) pp 1801-1821.

Rybicki, E.R., and Kanninen, M. (1977), "A Finite Element Calculation of Stress Intensity Factors by a Modified Crack Closure Integral," AIAA Journal Volume 2, 1977.

Underwood, P., (1983)"Dynamic Relaxation", in Computational Methods for Transient Analysis, Volume 1. Belytchko T., and Hughes, J.R. Eds. North-Holland, Amsterdam, 1983.