- Email: [email protected]
Stress intensity factors for cracks in panels with uniformly spaced holes
0013-7944/x7 33.caf.00 Pergamon Journals Ltd.
Enp,emn~ Fracrure Mrchanw., Vol. 26, No. 2, pp. 261-271, 1987 Printed in Great Britain
STRESS INTENSITY FACTORS PANELS WITH UNIFORMLY Department
FOR CRACKS IN SPACED HOLES
ROBERT H. DODD& Jr? of Civil Engineering, University of Kansas, Lawrence,
KS 66045, U.S.A.
and
Fracture
and Deformation
Division,
DAVID T. READ National Bureau of Standards, CO 80302, U.S.A.
325 S. Broadway,
Boulder,
Abstract-Stress intensity factors are presented for cracks at the edges of the central holes in a panel containing a uniformly spaced array of circular holes. The finite element method is employed to obtain numerical values of the J-integral and the near-tip displacements from which the stress intensity factor, K,, is inferred. Six representative ratios of hole spacing to hole radius are studied: 1.2, 1.25, 1.4, 1.6, 2.0 and 4.0. For each hole spacing, the stress intensity factor is computed for a number of crack lengths for each of two crack configurations. Uniaxial loading in two orthogonal directions is considered. Computed stress intensity factors are presented in graphical and functional form for all cases. The functional form should prove convenient for fatigue crack growth computations.
INTRODUCTION FOR THE practical fracture and fatigue analysis of complex 3-D structures, simpler 2-D models are often constructed to determine approximate stress intensity factors. For example, the flat panel shown in Fig. 1 contains a uniformly spaced array of circular holes. The center holes have throughthickness edge cracks in two different configurations, denoted Case A and Case B. This 2-D model may be employed to obtain approximate stress intensity factors for cylindrical pressure vessels containing many closely spaced bore holes for the passage of water tubes, e.g. superheater headers, steam drums and water drums. As illustrated in Fig. 2, semi-elliptical surface cracks have been observed to form perpendicular to the bore holes due to the high stress concentrations, the cyclic thermal and pressure loading and the weldments between the tubes and the walls. The computation of exact stress intensity factors for this geometry is complicated by the wall curvature, adjacent holes and the crack shape. The 2-D model in Fig. 1 ignores the effects of wall curvature and treats the long but shallow surface crack as a through-thickness crack. This model permits the effects of hole spacing and adjacent cracks to be determined by numerical methods using reasonable amounts of computer resources[l]. Several solutions for cracks at the edges of holes are readily available. Tada et al.[2] present stress intensity factors for flat plates containing one or two holes with a through-thickness crack placed at a simplifying location. In addition to these 2-D solutions, Newman and Raju[3] reported the stress intensity factors for surface cracks in a bore hole of a flat plate containing a single penetration. Newman and Raju’s results are useful to verify the approximation of a long but shallow surface crack by a through-thickness crack. However, neither of these sets of solutions address the effects of hole spacing and adjacent cracks. This paper presents the results of a study to determine accurate stress intensity factors for the 2-D crack configurations (Cases A and B) shown in Fig. 1. For each case, six representative holespacing to hole-radius ratios (b/R) are studied: 1.2, 1.25, 1.4, 1.6, 2.0 and 4.0. Axial loading in the X and Y directions is considered separately to permit superposition of the stress intensity factors. Finite element models are used to obtain numerical values for the J-integral and the near-tip displacements from which the stress intensity factor is inferred. The computed stress intensity factors are presented in graphical and functional form to facilitate both fracture and fatigue analyses.
7 To whom all correspondence
should
be addressed.
267
R. H. DODDS
268
I
-
Jr and D. T. READ
-ly-T--r--i CASE A
CASE Fig. 1. 2-D models
B
with through-thickness
cracks
FINITE ELEMENT MODELS Finite element models were generated for the two cases, A and B, as illustrated in Fig. 3. The element grid shown in Fig. 3 modeled the hole-spacing to hole-radius ratio (b/R) of 1.4. The X = 0 and Y = 0 planes were taken as planes of symmetry in the analysis. Uniform tractions were applied over remaining ligaments on the right and top boundaries. The intensity of these tractions was specified to model an equivalent unit stress applied over the gross section. The ligaments on the top boundary (Y = 4b) were constrained to displace uniformly in the Y-direction. Similarly, ligaments on the right boundary (X = 4h) were constrained to displace uniformly in the X direction. These constraints provided a model with cyclic symmetry in both X and Y directions. Eight-node isoparametric elements were employed for the modeling. These elements provide a quadratic variation of displacement with a corresponding linear variation of strain and stress. Reduced (2 x 2) Gauss quadrature was used to generate the element stiffness matrices. The element grid was analyzed in plane stress with a Young’s modulus of 1000 and a Poisson’s ratio of 0.3.‘The adoption of specific material properties and the use of the plane stress or plane strain assumption are required for the finite element computations but do not affect the computed stress intensity factors. In the typical element grid, there were 1024 of the eight-node isoparametric elements. Another 20 or 22 elements were used to define the crack-tip vicinity. Substructuring and static condensation were applied at several levels to take advantage of the repeated geometry. The numbers shown on the element grid in Fig. 3 are the key nodes prior to the final condensation, which eliminated all nodes except 17 l- 187. Analyses were conducted for six hole-spacing to hole-radius ratios that bracket the range of ratios commonly found in steam generation equipment. For modeling, the hole radius was fixed at 2.54 cm (1.0in.) and the spacing between holes varied to model the following h/R ratios: 1.2,
in panels
Outlet
Cross-Section
Header
Tube
Fig. 2. Superheater
269
SIFs for cracks
header
Hole
with longitudinal
bore hole crack (l.Oin. = 2.54cm).
Finite
,f tr~r$l
-~Stress Area
on Gross t tl = 1.0
Elements
1,024, B-Node Isooarametrics. Plane Stress, Lfj8000*
“=“.3
i
168
Fig. 3. Typical
element grid pattern
used for all analyses.
This particular
grid is for h/R = 1.4.
270
R. H. DODDS
Jr and D. T. READ
1.25, 1.4, 1.6, 2.0 and 4.0. The same grid pattern shown in Fig. 3 was employed for each hole spacing. For each b/R ratio, solutions were generated for Case A and Case B crack configurations. Cracks were incorporated into the element grid across the ligament between nodes 171 and 187. To model each crack length, two of the eight elements across this ligament were replaced by the crack-tip substructure shown in Fig. 4. All nodes except 57-65 of the crack-tip substructure were statically condensed to produce an equivalent “element” for connection to the adjacent eight-node elements. Full displacement compatibility of the model was maintained with this procedure. At the crack tip in both Cases A and B, four eight-node isoparametric elements (denoted elements 1-4 in Fig. 4) were degenerated into triangles by coalescing nodes common to one edge. The adjacent mid-side nodes were displaced to the quarter point. This forced the strain field on all rays focusing at the crack tip in these elements to exhibit a l/J r singularity[4]. Crack-face nodes 65,. . ,37,. ,23 and 14 remained free to displace. Coincident crack-tip nodes I 9 and nodes 10, 15,. . ,57 on the crack plane were constrained to have no Y displacement. The coincident-tip nodes of the singularity elements were also constrained to displace equally in the X direction. With this modeling strategy, minimal computer time was required to generate the solution for each crack length. The stiffness matrix of the region shown in Fig. 3 and the crack-tip region shown in Fig. 4 were generated and condensed only once to provide equivalent models with fewer nodes. The grid in Fig. 3 containing 1024 right-node elements, was reduced to an equivalent model containing 17 nodes. For each crack length, this 17-node equivalent model was combined with one or two copies of the condensed crack tip model and five or six additional eight-node elements. The final model for each crack length had only 39 nodes for Case A and 35 nodes for Case B. The analyses were performed with the POLO-FINITEIS, 61 structural mechanics system. The substructuring, static condensation and multiple restarts needed for the analyses are standard features of this system. Stress intensity factors were computed with two independent methods. The first method used the J-integral; the second method used nodal displacements of the degenerated crack tip element that lies on the crack face. The J-integral was evaluated over paths passing through Gauss points of elements in the crack tip substructure. A pair of paths was defined to pass through each of three rings of elements: elements 5-8, 9 - 12 and I3316 (see Fig. 4). The variation between J-values computed over each pair of paths was found to be less than 3%. The average value for the six paths was doubled to obtain the full J-value. The stress intensity factor was obtained from the J-integral using the relation K, = \/EJ The stress intensity
factor was also obtained
Nodes
57-65
Retained
After
(plane stress). by substituting
(1) the computed
Condensation
Quarter
Fig. 4. Crack-tip
substructure
showing
nodal displacements
degenerated
isoparametric
triangles
Point
Elements
used at the crack
tip.
SIFs for cracks
in panels
271
on the crack face into the near-tip solution for displacement. The stress intensity factor becomes the only unknown and is easily computed. The procedure has been specialized for application to eight-node, quarter-point elements by Shih et al.[7]. The stress intensity factor is given by K, = {2GJ27~/(ti
+ 1)}(4VB - I/c)/&
(2)
where G is the shear modulus, K is the bulk modulus, L is the dimension of the crack-tip element and V,, Vc are the opening displacements of the crack-tip-element nodes. These quantities are illustrated in Fig. 5. The K, values calculated using the J-integral and nodal displacement methods agreed to within 1% for all hole spacings and crack lengths.
RESULTS AND DISCUSSION Stress concentration factors The adequacy of the modeling procedure was verified through an analysis with no cracks. An equivalent tensile stress of unit value on the gross section was imposed first over the top boundary and then over the right boundary. Stress concentration factors (SCFs) were determined at nodes 171, 187, 188, 204, 86, 102, etc. of the element grid shown in Fig. 3. The d,, stress at the edge of each hole was obtained by extrapolating stresses computed directly at the Gauss point locations in adjacent elements. Two stress concentration factors, denoted K,, and K,,, were computed for each location. Both factors refer to the resulting g‘ystress, but caused by unit tensile stresses applied in the X direction (KTX) and the Y direction (K,,). The computed stress concentration factors were identical at each of the nodes listed above. Numerical values are listed in Table 1. The K,, values agree closely with the solutions given in graphical form in [S]. The K,, values also agree well for b/R > 1.6. The solutions for closer hole-spacing in [S] do not exhibit as rapid a decrease toward zero as do the finite element results in Table 1.
Singularity
Fig. 5. Crack-tip
elements
Elements
showing
-
I
nodal displacements
Table
used to compute
1. Stress concentration factors
h,lR
KTX
K TY
I .20
-0.010 - 0.029 -0.170 -0.350 - 0.580 - 1.000
6.8 6.1 4.4 3.7 3.3 3.1
I .25 1.40 1.60 2.00 4.00
the stress intensity
factor
212
R. H. DODDS
Jr and D. T. READ
Stress intensity factors
The computed
stress intensity
factors
are presented
in the form
where e is the stress on the gross area, a is the crack length, c is the ligament and b/R is the hole-spacing ratio. These dimensions are detailed in Fig. 1. intensity coefficients, F, for Case A are shown in Figs 6 and 8 and for Case Nonlinear curve fitting procedures were employed to obtain functional forms corresponding to the curves shown in these figures. Stress intensity coeficients
forCase
length between holes Graphs of the stress B in Figs 9 and 10. for the F coefficients
A
The stress intensity coefficients for uniaxial Y loading loading in Fig. g. The form of the curves for both loading
are shown in Fig. 6 and for uniaxial states is given by
X
where a =
1 - a/c,
/I1 -+ /I6 are curve-fitting
parameters,
6 = b/R,
Numerical
values for each /I parameter
are given in Tables
\_
b/R
2 and 3. The F value for zero crack
= 1.2 =-.
\
= 1.25
0-1 0
0.2
0.4
0.6
’
’
’
0.8
--J
1.0
a/c Fig. 6. Stress intensity
factor
coefficients
for Case A, unit tensile stress Y on gross sectvm
SIFs for cracks
4
I
in panels
I
I
I
A
Ref. With
.
3
F > ," ; L
273
This b=4.
I
[21,Infinite
Plate
Single
Hole
Study R=l
Case
A,
2
1
0 0
I
I
I
I
I
1
2
3
4
5
Crack Fig. 7. Comparison
Length,
a
of Case A (b/R = 4) with solution
for a single hole in an infinite plate
. blR = . b/R = .
b/R
1
4.0 2.0
= 1.6
i 0.2
0
0.6
0.4
0.8
1.0
a/c Fig. 8. Stress intensity
length, a/c = 0, was obtained and the relation
factor coefficients
for Case A, unit tensile stress X on gross section.
from the computed
stress concentration
factor, K, = K,,
Table 2. Stress intensile factor b parameters for Case A, unit tensile stress Y on gross section h/R =6 1.20 1.25 1.40 J.60 2.00 4.00
/I1 8.835 8.476 8.029 7.273 4.079 2.683
B,
I%
84
Ps
- 6.694 - 6.960 - 4.677 -3.419 -2.142 -2.160
1.544 1.806 3.036 4.195 3.442 3.172
- 1.261 - 1.584 -2.371 - 3.265 -2.910 -3.016
- 1.424 - 1.265 -0.515 -0.323 -0.477 -0.507
86 2.704 2.113 2.782 3.48 1 3.081 0.732
or K,,,
214
R. H. DODDS Jr arId D. T. READ Table 3. Stress intensity factor /I parameters for Case A, unit tensile stress X on gross section hJR=6
1.60 2.00 4.00
[I,
B2
83
a4
8,
I&
8.458 6.284 3.950
-9.834 -8.152 - 5.001
11.706 18.034 18.807
-8.196 - 15.055 -17.869
-0.730 -0.347 -0.190
0.4610.487 0.051
K, = 1.12K,JG
(5)
F = 1.12K,
(6)
where 1.12 is the free edge correction factor. For uniaxial Y loading, all of the curves for Case A exhibit the same general shape. The stress intensity factor first decreases with an increasing crack length as the crack tip extends beyond the influence of the stress concentration at the edge of the hole. Once the crack, extends beyond midligament (a/c = 0.5) the increased stress due to the adjacent hole begins to raise the stress intensity factor. At a crack length of a/c = 0.9, the stress intensity factor rapidly increases and approaches infinity for a/c = 1.0. The trends of the F functions for a/c > 0.5 are masked by the 4’1 -- u/c scale factor, which was selected to obtain a clear separation of the curves for the larger a/c ratios. The influence of hole spacing is recognized by the much higher stress intensity coefficients for closely spaced holes at each crack length. The curves show that the coefficients decrease very rapidly between b/R of 1.2 and 1.6 and then more slowly for the larger hole spacings. The stress on the remaining ligament is considerably reduced for the larger hole-spacings. Figure 7 shows, in a non-parameterized form, the effects of an adjacent hole on the stress intensity factor for a single case. A hole radius of 1.0 unit and a hole spacing (center to center) of 8.0 units (b/R = 4.0) is considered. The computed results for Case A are compared with the solution for an infinite plate containing a single hole with cracks located as in Case A. For this larger holespacing, the single hole and Case A solutions match very well until the crack length exceeds 4.0 units (the SCF for Case A is very near 3.0). For longer crack lengths the stress concentration due to the adjacent hole rapidly raises the stress intensity factor which approaches DZat a crack length of 6.0 units. The single hole solution approaches the solution (F = 1.0) for an infinite plate with a center crack of length 2(R + a). Several checks can be made on the numerical results. For short cracks (a/R << l.O), the approximate stress intensity factor is given by K, = 1.12&J&
(7)
.For longer cracks (a/c = 0.4-0.5), an approximate upper bound for the stress intensity factor is given by the solution for an infinite panel with a periodic array of cracks along a line[2] with the crack length taken as (R + a). The addition of holes above and below the line of cracks reduces the actual stress intensity factor below this solution. An approximate lower bound for the stress intensity factor is given by the solution for a finite crack with a hole centered above and below the crack[2]. The addition of holes to the left and right of the crack increases the actual stress intensity above this solution. These checks verify that the finite element stress intensity factors are properly bounded. For uniaxial X loading, Fig. 8 shows that the stress intensity coefficients decrease rapidly to zero as a/c --f 0.5. The coefficients are negative due to the tensile sense of the applied X loading. The curves also show that the maximum stress intensity coefficient, which occurs for a/c = 0.0, decreases very rapidly with decreasing hole-spacing. These results reflect the very small stress concentrations, KTX, for the smaller hole spacings. The X loading has the most influence for the largest hole-spacing, b/R = 4.0. The maximum stress intensity coefficient is 33% of the magnitude obtained for the Y loading. Rigorous error estimates cannot be stated for the computed stress intensity factors. However, experience has shown that finite element models having much less grid refinement for test geometries (CCP, SEN and three-point bend) predict stress intensity factors within 2-3% of exact values. The high degree of grid refinement employed in this study should provide stress intensity factors having
SIFs for cracks
in panels
275
less than 3% error. Stress intensity coefficients for crack lengths beyond the last computed point (u/c = 0.875) are subject to larger curve-fitting errors. coeficients for Case B The stress intensity coefficients for uniaxial Y loading are shown in Fig. 9 and for uniaxial X loading in Fig, 10. The form of the curves for both loading states is given by
Stress intensity
eSlaa
F(a,‘c, b/R) = -‘_ Ja
84 w + 0.256 -i- a 46 + 4&E) + 131ac
(8)
where a = 1 - a/c, PI -b fi4 are curve-fitting parameters,
6 = b/R,and a/c < 0.375.
Values for the fi parameters are given in Tables 4 through 6 for the different crack tips and loadings. For uniaxial Y loading, Fig. 9 provides the stress intensity coefficients for crack-tip A located as shown in Fig. l(b). The stress intensity factors for tip B are slightly smaller than those for tip A due to the absence of a crack on the side of the hole opposite tip B. Although a graph is not provided for tip B, Table 5 lists the parameters for use in eq. (8) to compute the stress intensity factors. The curves reveal that a much larger degree of interaction between cracks growing from adjacent holes exists for b/R < 2.0. The stress intensity coefficient increases very rapidly once a/c exceeds 0.3. For the larger hole-spacings, the cracks grow into a rapidly diminishing stress field as the effects of the hole-induced stress concentrations decrease. The crack tips must grow very close for significant interaction to occur. Again, the scale factor applied to F in Fig. 9 somewhat masks this trend.
Fig. 9. Stress intensity
factor coefficients
for Case 8, crack tip A, unit tensile stress Y on gross section.
276
R. H. DODDS Jr and D. T. READ
. . l
b/R = 4.0 b/R = 2.0 b/R = 1.6
‘J
i -1.21
'I' 0
; 0.1
j
"1
c 0.2
I""", 0.3
" 0.4
0.5
a/c Fig. 10. Stress intensity factor coefficients for Case B, crack tips A and B, unit tensile stress X on gross section. Table 4. Stress intensity factor fl parameters for Case B, crack tip A, unit tensile stress Y on gross section bJR = 6 1.20 1.25 1.40 1.60 2.00 4.00
fi, 2.188 4.007 2.941 2.841 -0.695 0.389
82 0.201 6.065 2.293 3.955 2.124 -0.915
P3 - 3.966 -0.858 - 1.033 -0.628 3.614 1.027
P4 2.183 -5.160 - 1.890 - 4.468 2.257 I.536
Table 5. Stress intensity factor fi parameters for Case B, crack tip B, unit tensile stress Y on gross section b/R = 6 1.20 1.25 1.40 1.60 2.00 4.00
/I, 2.097
Pz 0.150
P, -4.365
3.737 2.718 2.521 - 0.689 0.616
3.925 1.455 1.884 2.254 2.382
- 1.023 - 1.297 -0,905 3.558 3.358
84 2.174 - 2.807 - 1.164 -- 2.244 2.065 1.522
Table 6. Stress intensity factor B parameters for case B, crack tip A, unit tensile stress X on gross section b/R =6 1.60 2.00 4.00
P* B, -26.926 -0.396 -24.344 -0.637 -61.262 - 1.062
8, 84 3.604 -0.00493 3.717 0.01489 3.196 -0.11420
Figure 10 shows the stress intensity coefficients at crack tip A for the uniaxial X loading. Coefficients for crack tip B are slightly smaller and are not shown. The functional form of the curves is the same as that for the Y loading (eq. 8), with values of the parameters given in Table 6. The variation of the coefficients with b/R and a/c is identical to that observed for Case A for the X loading. The stress intensity coefficients decrease rapidly to zero as a/c -+ 0.5 for each hole spacing. The values are negative due to the tensile sense of the applied X loading.
SIFs for cracks
in panels
211
Error estimates of the stress intensity factors for Case B are the same as the error estimates for Case A. The approximate upper and lower bound solutions for Case A also apply for Case B and demonstrate that the finite element solutions are properly bounded. Acknowledgements-Support for this study was provided by the Fracture Control Technology Program, Naval Sea Systems Command (SEA-05 R25). Mr John P. Gudas was program monitor. Numerical computations were performed on the Harris 500 and 800 computers in the Computer Aided Engineering Laboratory at the University of Kansas. Mr. Shahin Neghaban performed the nonlinear curve-fitting procedures.
REFERENCES engineering and 111 R. H. Dodds, Stress intensity factors for cracks in panels with uniformly spaced holes. Structural engineering materials, SL Report 85-1, University of Kansas Center for Research Inc., Lawrence, Kansas, pp.35 (April 1985). PI H. Tada, P. C. Paris and G. R. Irwin (Eds), Stress Analysis of Cracks Handbook. Del Research Corp., Hellertown, Pennsylvania (1973). factor equations for cracks in three-dimensional finite bodies. NASA c31 J. C. Newman and I. S. Raju, Stress-intensity Technical Memorandum 83200 (August 1981). quarter-point elements as elastic and perfectly-plastic crack tip elements. Inr. J. numer. Meth. c41 R. Barsoum, Triangular Engng 11, 85-98 (1977). software for nonlinear analysis. Znt. J. Adtr. Engng Software 2, 161-168 I51 R. H. Dodds and L. A. Lopez, Generalized (1981). in linear and nonlinear analysis. Int. J. numer. Meth. Engng 15, 583-597 WI R. H. Dodds and L. A. Lopez, Substructuring (1980). isoparametric c71 C. F. Shih, H. G. deLorenzi and M. D. German, Crack extension modeling with singular quadratic elements. Inc. J. Fracture 12, 647-651 (1976). PI W. Griffel, Handbook of Formulas for Stress and Strain. (Received 6 June 1986)
Enp,emn~ Fracrure Mrchanw., Vol. 26, No. 2, pp. 261-271, 1987 Printed in Great Britain
STRESS INTENSITY FACTORS PANELS WITH UNIFORMLY Department
FOR CRACKS IN SPACED HOLES
ROBERT H. DODD& Jr? of Civil Engineering, University of Kansas, Lawrence,
KS 66045, U.S.A.
and
Fracture
and Deformation
Division,
DAVID T. READ National Bureau of Standards, CO 80302, U.S.A.
325 S. Broadway,
Boulder,
Abstract-Stress intensity factors are presented for cracks at the edges of the central holes in a panel containing a uniformly spaced array of circular holes. The finite element method is employed to obtain numerical values of the J-integral and the near-tip displacements from which the stress intensity factor, K,, is inferred. Six representative ratios of hole spacing to hole radius are studied: 1.2, 1.25, 1.4, 1.6, 2.0 and 4.0. For each hole spacing, the stress intensity factor is computed for a number of crack lengths for each of two crack configurations. Uniaxial loading in two orthogonal directions is considered. Computed stress intensity factors are presented in graphical and functional form for all cases. The functional form should prove convenient for fatigue crack growth computations.
INTRODUCTION FOR THE practical fracture and fatigue analysis of complex 3-D structures, simpler 2-D models are often constructed to determine approximate stress intensity factors. For example, the flat panel shown in Fig. 1 contains a uniformly spaced array of circular holes. The center holes have throughthickness edge cracks in two different configurations, denoted Case A and Case B. This 2-D model may be employed to obtain approximate stress intensity factors for cylindrical pressure vessels containing many closely spaced bore holes for the passage of water tubes, e.g. superheater headers, steam drums and water drums. As illustrated in Fig. 2, semi-elliptical surface cracks have been observed to form perpendicular to the bore holes due to the high stress concentrations, the cyclic thermal and pressure loading and the weldments between the tubes and the walls. The computation of exact stress intensity factors for this geometry is complicated by the wall curvature, adjacent holes and the crack shape. The 2-D model in Fig. 1 ignores the effects of wall curvature and treats the long but shallow surface crack as a through-thickness crack. This model permits the effects of hole spacing and adjacent cracks to be determined by numerical methods using reasonable amounts of computer resources[l]. Several solutions for cracks at the edges of holes are readily available. Tada et al.[2] present stress intensity factors for flat plates containing one or two holes with a through-thickness crack placed at a simplifying location. In addition to these 2-D solutions, Newman and Raju[3] reported the stress intensity factors for surface cracks in a bore hole of a flat plate containing a single penetration. Newman and Raju’s results are useful to verify the approximation of a long but shallow surface crack by a through-thickness crack. However, neither of these sets of solutions address the effects of hole spacing and adjacent cracks. This paper presents the results of a study to determine accurate stress intensity factors for the 2-D crack configurations (Cases A and B) shown in Fig. 1. For each case, six representative holespacing to hole-radius ratios (b/R) are studied: 1.2, 1.25, 1.4, 1.6, 2.0 and 4.0. Axial loading in the X and Y directions is considered separately to permit superposition of the stress intensity factors. Finite element models are used to obtain numerical values for the J-integral and the near-tip displacements from which the stress intensity factor is inferred. The computed stress intensity factors are presented in graphical and functional form to facilitate both fracture and fatigue analyses.
7 To whom all correspondence
should
be addressed.
267
R. H. DODDS
268
I
-
Jr and D. T. READ
-ly-T--r--i CASE A
CASE Fig. 1. 2-D models
B
with through-thickness
cracks
FINITE ELEMENT MODELS Finite element models were generated for the two cases, A and B, as illustrated in Fig. 3. The element grid shown in Fig. 3 modeled the hole-spacing to hole-radius ratio (b/R) of 1.4. The X = 0 and Y = 0 planes were taken as planes of symmetry in the analysis. Uniform tractions were applied over remaining ligaments on the right and top boundaries. The intensity of these tractions was specified to model an equivalent unit stress applied over the gross section. The ligaments on the top boundary (Y = 4b) were constrained to displace uniformly in the Y-direction. Similarly, ligaments on the right boundary (X = 4h) were constrained to displace uniformly in the X direction. These constraints provided a model with cyclic symmetry in both X and Y directions. Eight-node isoparametric elements were employed for the modeling. These elements provide a quadratic variation of displacement with a corresponding linear variation of strain and stress. Reduced (2 x 2) Gauss quadrature was used to generate the element stiffness matrices. The element grid was analyzed in plane stress with a Young’s modulus of 1000 and a Poisson’s ratio of 0.3.‘The adoption of specific material properties and the use of the plane stress or plane strain assumption are required for the finite element computations but do not affect the computed stress intensity factors. In the typical element grid, there were 1024 of the eight-node isoparametric elements. Another 20 or 22 elements were used to define the crack-tip vicinity. Substructuring and static condensation were applied at several levels to take advantage of the repeated geometry. The numbers shown on the element grid in Fig. 3 are the key nodes prior to the final condensation, which eliminated all nodes except 17 l- 187. Analyses were conducted for six hole-spacing to hole-radius ratios that bracket the range of ratios commonly found in steam generation equipment. For modeling, the hole radius was fixed at 2.54 cm (1.0in.) and the spacing between holes varied to model the following h/R ratios: 1.2,
in panels
Outlet
Cross-Section
Header
Tube
Fig. 2. Superheater
269
SIFs for cracks
header
Hole
with longitudinal
bore hole crack (l.Oin. = 2.54cm).
Finite
,f tr~r$l
-~Stress Area
on Gross t tl = 1.0
Elements
1,024, B-Node Isooarametrics. Plane Stress, Lfj8000*
“=“.3
i
168
Fig. 3. Typical
element grid pattern
used for all analyses.
This particular
grid is for h/R = 1.4.
270
R. H. DODDS
Jr and D. T. READ
1.25, 1.4, 1.6, 2.0 and 4.0. The same grid pattern shown in Fig. 3 was employed for each hole spacing. For each b/R ratio, solutions were generated for Case A and Case B crack configurations. Cracks were incorporated into the element grid across the ligament between nodes 171 and 187. To model each crack length, two of the eight elements across this ligament were replaced by the crack-tip substructure shown in Fig. 4. All nodes except 57-65 of the crack-tip substructure were statically condensed to produce an equivalent “element” for connection to the adjacent eight-node elements. Full displacement compatibility of the model was maintained with this procedure. At the crack tip in both Cases A and B, four eight-node isoparametric elements (denoted elements 1-4 in Fig. 4) were degenerated into triangles by coalescing nodes common to one edge. The adjacent mid-side nodes were displaced to the quarter point. This forced the strain field on all rays focusing at the crack tip in these elements to exhibit a l/J r singularity[4]. Crack-face nodes 65,. . ,37,. ,23 and 14 remained free to displace. Coincident crack-tip nodes I 9 and nodes 10, 15,. . ,57 on the crack plane were constrained to have no Y displacement. The coincident-tip nodes of the singularity elements were also constrained to displace equally in the X direction. With this modeling strategy, minimal computer time was required to generate the solution for each crack length. The stiffness matrix of the region shown in Fig. 3 and the crack-tip region shown in Fig. 4 were generated and condensed only once to provide equivalent models with fewer nodes. The grid in Fig. 3 containing 1024 right-node elements, was reduced to an equivalent model containing 17 nodes. For each crack length, this 17-node equivalent model was combined with one or two copies of the condensed crack tip model and five or six additional eight-node elements. The final model for each crack length had only 39 nodes for Case A and 35 nodes for Case B. The analyses were performed with the POLO-FINITEIS, 61 structural mechanics system. The substructuring, static condensation and multiple restarts needed for the analyses are standard features of this system. Stress intensity factors were computed with two independent methods. The first method used the J-integral; the second method used nodal displacements of the degenerated crack tip element that lies on the crack face. The J-integral was evaluated over paths passing through Gauss points of elements in the crack tip substructure. A pair of paths was defined to pass through each of three rings of elements: elements 5-8, 9 - 12 and I3316 (see Fig. 4). The variation between J-values computed over each pair of paths was found to be less than 3%. The average value for the six paths was doubled to obtain the full J-value. The stress intensity factor was obtained from the J-integral using the relation K, = \/EJ The stress intensity
factor was also obtained
Nodes
57-65
Retained
After
(plane stress). by substituting
(1) the computed
Condensation
Quarter
Fig. 4. Crack-tip
substructure
showing
nodal displacements
degenerated
isoparametric
triangles
Point
Elements
used at the crack
tip.
SIFs for cracks
in panels
271
on the crack face into the near-tip solution for displacement. The stress intensity factor becomes the only unknown and is easily computed. The procedure has been specialized for application to eight-node, quarter-point elements by Shih et al.[7]. The stress intensity factor is given by K, = {2GJ27~/(ti
+ 1)}(4VB - I/c)/&
(2)
where G is the shear modulus, K is the bulk modulus, L is the dimension of the crack-tip element and V,, Vc are the opening displacements of the crack-tip-element nodes. These quantities are illustrated in Fig. 5. The K, values calculated using the J-integral and nodal displacement methods agreed to within 1% for all hole spacings and crack lengths.
RESULTS AND DISCUSSION Stress concentration factors The adequacy of the modeling procedure was verified through an analysis with no cracks. An equivalent tensile stress of unit value on the gross section was imposed first over the top boundary and then over the right boundary. Stress concentration factors (SCFs) were determined at nodes 171, 187, 188, 204, 86, 102, etc. of the element grid shown in Fig. 3. The d,, stress at the edge of each hole was obtained by extrapolating stresses computed directly at the Gauss point locations in adjacent elements. Two stress concentration factors, denoted K,, and K,,, were computed for each location. Both factors refer to the resulting g‘ystress, but caused by unit tensile stresses applied in the X direction (KTX) and the Y direction (K,,). The computed stress concentration factors were identical at each of the nodes listed above. Numerical values are listed in Table 1. The K,, values agree closely with the solutions given in graphical form in [S]. The K,, values also agree well for b/R > 1.6. The solutions for closer hole-spacing in [S] do not exhibit as rapid a decrease toward zero as do the finite element results in Table 1.
Singularity
Fig. 5. Crack-tip
elements
Elements
showing
-
I
nodal displacements
Table
used to compute
1. Stress concentration factors
h,lR
KTX
K TY
I .20
-0.010 - 0.029 -0.170 -0.350 - 0.580 - 1.000
6.8 6.1 4.4 3.7 3.3 3.1
I .25 1.40 1.60 2.00 4.00
the stress intensity
factor
212
R. H. DODDS
Jr and D. T. READ
Stress intensity factors
The computed
stress intensity
factors
are presented
in the form
where e is the stress on the gross area, a is the crack length, c is the ligament and b/R is the hole-spacing ratio. These dimensions are detailed in Fig. 1. intensity coefficients, F, for Case A are shown in Figs 6 and 8 and for Case Nonlinear curve fitting procedures were employed to obtain functional forms corresponding to the curves shown in these figures. Stress intensity coeficients
forCase
length between holes Graphs of the stress B in Figs 9 and 10. for the F coefficients
A
The stress intensity coefficients for uniaxial Y loading loading in Fig. g. The form of the curves for both loading
are shown in Fig. 6 and for uniaxial states is given by
X
where a =
1 - a/c,
/I1 -+ /I6 are curve-fitting
parameters,
6 = b/R,
Numerical
values for each /I parameter
are given in Tables
\_
b/R
2 and 3. The F value for zero crack
= 1.2 =-.
\
= 1.25
0-1 0
0.2
0.4
0.6
’
’
’
0.8
--J
1.0
a/c Fig. 6. Stress intensity
factor
coefficients
for Case A, unit tensile stress Y on gross sectvm
SIFs for cracks
4
I
in panels
I
I
I
A
Ref. With
.
3
F > ," ; L
273
This b=4.
I
[21,Infinite
Plate
Single
Hole
Study R=l
Case
A,
2
1
0 0
I
I
I
I
I
1
2
3
4
5
Crack Fig. 7. Comparison
Length,
a
of Case A (b/R = 4) with solution
for a single hole in an infinite plate
. blR = . b/R = .
b/R
1
4.0 2.0
= 1.6
i 0.2
0
0.6
0.4
0.8
1.0
a/c Fig. 8. Stress intensity
length, a/c = 0, was obtained and the relation
factor coefficients
for Case A, unit tensile stress X on gross section.
from the computed
stress concentration
factor, K, = K,,
Table 2. Stress intensile factor b parameters for Case A, unit tensile stress Y on gross section h/R =6 1.20 1.25 1.40 J.60 2.00 4.00
/I1 8.835 8.476 8.029 7.273 4.079 2.683
B,
I%
84
Ps
- 6.694 - 6.960 - 4.677 -3.419 -2.142 -2.160
1.544 1.806 3.036 4.195 3.442 3.172
- 1.261 - 1.584 -2.371 - 3.265 -2.910 -3.016
- 1.424 - 1.265 -0.515 -0.323 -0.477 -0.507
86 2.704 2.113 2.782 3.48 1 3.081 0.732
or K,,,
214
R. H. DODDS Jr arId D. T. READ Table 3. Stress intensity factor /I parameters for Case A, unit tensile stress X on gross section hJR=6
1.60 2.00 4.00
[I,
B2
83
a4
8,
I&
8.458 6.284 3.950
-9.834 -8.152 - 5.001
11.706 18.034 18.807
-8.196 - 15.055 -17.869
-0.730 -0.347 -0.190
0.4610.487 0.051
K, = 1.12K,JG
(5)
F = 1.12K,
(6)
where 1.12 is the free edge correction factor. For uniaxial Y loading, all of the curves for Case A exhibit the same general shape. The stress intensity factor first decreases with an increasing crack length as the crack tip extends beyond the influence of the stress concentration at the edge of the hole. Once the crack, extends beyond midligament (a/c = 0.5) the increased stress due to the adjacent hole begins to raise the stress intensity factor. At a crack length of a/c = 0.9, the stress intensity factor rapidly increases and approaches infinity for a/c = 1.0. The trends of the F functions for a/c > 0.5 are masked by the 4’1 -- u/c scale factor, which was selected to obtain a clear separation of the curves for the larger a/c ratios. The influence of hole spacing is recognized by the much higher stress intensity coefficients for closely spaced holes at each crack length. The curves show that the coefficients decrease very rapidly between b/R of 1.2 and 1.6 and then more slowly for the larger hole spacings. The stress on the remaining ligament is considerably reduced for the larger hole-spacings. Figure 7 shows, in a non-parameterized form, the effects of an adjacent hole on the stress intensity factor for a single case. A hole radius of 1.0 unit and a hole spacing (center to center) of 8.0 units (b/R = 4.0) is considered. The computed results for Case A are compared with the solution for an infinite plate containing a single hole with cracks located as in Case A. For this larger holespacing, the single hole and Case A solutions match very well until the crack length exceeds 4.0 units (the SCF for Case A is very near 3.0). For longer crack lengths the stress concentration due to the adjacent hole rapidly raises the stress intensity factor which approaches DZat a crack length of 6.0 units. The single hole solution approaches the solution (F = 1.0) for an infinite plate with a center crack of length 2(R + a). Several checks can be made on the numerical results. For short cracks (a/R << l.O), the approximate stress intensity factor is given by K, = 1.12&J&
(7)
.For longer cracks (a/c = 0.4-0.5), an approximate upper bound for the stress intensity factor is given by the solution for an infinite panel with a periodic array of cracks along a line[2] with the crack length taken as (R + a). The addition of holes above and below the line of cracks reduces the actual stress intensity factor below this solution. An approximate lower bound for the stress intensity factor is given by the solution for a finite crack with a hole centered above and below the crack[2]. The addition of holes to the left and right of the crack increases the actual stress intensity above this solution. These checks verify that the finite element stress intensity factors are properly bounded. For uniaxial X loading, Fig. 8 shows that the stress intensity coefficients decrease rapidly to zero as a/c --f 0.5. The coefficients are negative due to the tensile sense of the applied X loading. The curves also show that the maximum stress intensity coefficient, which occurs for a/c = 0.0, decreases very rapidly with decreasing hole-spacing. These results reflect the very small stress concentrations, KTX, for the smaller hole spacings. The X loading has the most influence for the largest hole-spacing, b/R = 4.0. The maximum stress intensity coefficient is 33% of the magnitude obtained for the Y loading. Rigorous error estimates cannot be stated for the computed stress intensity factors. However, experience has shown that finite element models having much less grid refinement for test geometries (CCP, SEN and three-point bend) predict stress intensity factors within 2-3% of exact values. The high degree of grid refinement employed in this study should provide stress intensity factors having
SIFs for cracks
in panels
275
less than 3% error. Stress intensity coefficients for crack lengths beyond the last computed point (u/c = 0.875) are subject to larger curve-fitting errors. coeficients for Case B The stress intensity coefficients for uniaxial Y loading are shown in Fig. 9 and for uniaxial X loading in Fig, 10. The form of the curves for both loading states is given by
Stress intensity
eSlaa
F(a,‘c, b/R) = -‘_ Ja
84 w + 0.256 -i- a 46 + 4&E) + 131ac
(8)
where a = 1 - a/c, PI -b fi4 are curve-fitting parameters,
6 = b/R,and a/c < 0.375.
Values for the fi parameters are given in Tables 4 through 6 for the different crack tips and loadings. For uniaxial Y loading, Fig. 9 provides the stress intensity coefficients for crack-tip A located as shown in Fig. l(b). The stress intensity factors for tip B are slightly smaller than those for tip A due to the absence of a crack on the side of the hole opposite tip B. Although a graph is not provided for tip B, Table 5 lists the parameters for use in eq. (8) to compute the stress intensity factors. The curves reveal that a much larger degree of interaction between cracks growing from adjacent holes exists for b/R < 2.0. The stress intensity coefficient increases very rapidly once a/c exceeds 0.3. For the larger hole-spacings, the cracks grow into a rapidly diminishing stress field as the effects of the hole-induced stress concentrations decrease. The crack tips must grow very close for significant interaction to occur. Again, the scale factor applied to F in Fig. 9 somewhat masks this trend.
Fig. 9. Stress intensity
factor coefficients
for Case 8, crack tip A, unit tensile stress Y on gross section.
276
R. H. DODDS Jr and D. T. READ
. . l
b/R = 4.0 b/R = 2.0 b/R = 1.6
‘J
i -1.21
'I' 0
; 0.1
j
"1
c 0.2
I""", 0.3
" 0.4
0.5
a/c Fig. 10. Stress intensity factor coefficients for Case B, crack tips A and B, unit tensile stress X on gross section. Table 4. Stress intensity factor fl parameters for Case B, crack tip A, unit tensile stress Y on gross section bJR = 6 1.20 1.25 1.40 1.60 2.00 4.00
fi, 2.188 4.007 2.941 2.841 -0.695 0.389
82 0.201 6.065 2.293 3.955 2.124 -0.915
P3 - 3.966 -0.858 - 1.033 -0.628 3.614 1.027
P4 2.183 -5.160 - 1.890 - 4.468 2.257 I.536
Table 5. Stress intensity factor fi parameters for Case B, crack tip B, unit tensile stress Y on gross section b/R = 6 1.20 1.25 1.40 1.60 2.00 4.00
/I, 2.097
Pz 0.150
P, -4.365
3.737 2.718 2.521 - 0.689 0.616
3.925 1.455 1.884 2.254 2.382
- 1.023 - 1.297 -0,905 3.558 3.358
84 2.174 - 2.807 - 1.164 -- 2.244 2.065 1.522
Table 6. Stress intensity factor B parameters for case B, crack tip A, unit tensile stress X on gross section b/R =6 1.60 2.00 4.00
P* B, -26.926 -0.396 -24.344 -0.637 -61.262 - 1.062
8, 84 3.604 -0.00493 3.717 0.01489 3.196 -0.11420
Figure 10 shows the stress intensity coefficients at crack tip A for the uniaxial X loading. Coefficients for crack tip B are slightly smaller and are not shown. The functional form of the curves is the same as that for the Y loading (eq. 8), with values of the parameters given in Table 6. The variation of the coefficients with b/R and a/c is identical to that observed for Case A for the X loading. The stress intensity coefficients decrease rapidly to zero as a/c -+ 0.5 for each hole spacing. The values are negative due to the tensile sense of the applied X loading.
SIFs for cracks
in panels
211
Error estimates of the stress intensity factors for Case B are the same as the error estimates for Case A. The approximate upper and lower bound solutions for Case A also apply for Case B and demonstrate that the finite element solutions are properly bounded. Acknowledgements-Support for this study was provided by the Fracture Control Technology Program, Naval Sea Systems Command (SEA-05 R25). Mr John P. Gudas was program monitor. Numerical computations were performed on the Harris 500 and 800 computers in the Computer Aided Engineering Laboratory at the University of Kansas. Mr. Shahin Neghaban performed the nonlinear curve-fitting procedures.
REFERENCES engineering and 111 R. H. Dodds, Stress intensity factors for cracks in panels with uniformly spaced holes. Structural engineering materials, SL Report 85-1, University of Kansas Center for Research Inc., Lawrence, Kansas, pp.35 (April 1985). PI H. Tada, P. C. Paris and G. R. Irwin (Eds), Stress Analysis of Cracks Handbook. Del Research Corp., Hellertown, Pennsylvania (1973). factor equations for cracks in three-dimensional finite bodies. NASA c31 J. C. Newman and I. S. Raju, Stress-intensity Technical Memorandum 83200 (August 1981). quarter-point elements as elastic and perfectly-plastic crack tip elements. Inr. J. numer. Meth. c41 R. Barsoum, Triangular Engng 11, 85-98 (1977). software for nonlinear analysis. Znt. J. Adtr. Engng Software 2, 161-168 I51 R. H. Dodds and L. A. Lopez, Generalized (1981). in linear and nonlinear analysis. Int. J. numer. Meth. Engng 15, 583-597 WI R. H. Dodds and L. A. Lopez, Substructuring (1980). isoparametric c71 C. F. Shih, H. G. deLorenzi and M. D. German, Crack extension modeling with singular quadratic elements. Inc. J. Fracture 12, 647-651 (1976). PI W. Griffel, Handbook of Formulas for Stress and Strain. (Received 6 June 1986)