- Email: [email protected]
A note on numerically computed eigenfunctions and generalized stress intensity factors associated with singular points
Pergamon
Engineering Fracture Mechanics Vol. 54, No. 4, pp. 593-595, 1996 Copyright © 1996 ElsevierScienceLtd. Printed in Great Britain. All rights reserved 0013-7944(95)00090-9 0013-7944/96 $15.00+ 0.00
TECHNICAL NOTE
A NOTE ON NUMERICALLY COMPUTED EIGENFUNCTIONS AND GENERALIZED STRESS INTENSITY FACTORS ASSOCIATED WITH SINGULAR POINTS ZOHAR YOSIBASHI" and BARNA A. SZABO Center for Computational Mechanics, Washington University, Campus Box 1129, St. Louis, MO 63130, U.S.A. Abstract--The solution of linear elastostatic problems in the neighborhood of singular points is characterizedi by a sequence of eigenpairs and their coefficients, called generalized stress intensity factors (GSIFs). For general singular points, as cracks in anisotropic multimaterial interfaces, and V-notches in composite materials, only numerical approximations of the eigenfunctions are available. This note addresses the: question: how should the GSIFs computed using the numerically determined eigenfunctions be interpreted? Copyright © 1996 Elsevier Science Ltd.
DISCUSSION THE EXACTsolutions of linear elliptic partial differential equations which represent heat conduction or elastic deformation typically have a number of singular points. Each singular point is characterized by a sequence of eigenpairs (see refs [I-3]). For example, the solution of linear elasticity problems in two dimensions in the vicinity of a singular point, can be expanded in the form: M
u = ~.
~ C~,.r~J.lnmr.fjm(O) + u*,
(1)
jffil m ~ 0
where u* is smoother than the terms in the sum, Cjmare the coefficients of the asymptotic expansion (called the generalized stress intensity factors - - GSIFs), ~ and f~(0) are eigenpairs which depend on the differential operator and boundary conditions. The detelTnination of these eigenpairs is an essential step towards reliable computation of the GSIFs, which are of great practical importance because failure theories directly or indirectly involve these coefficients. Most of the research performed in the past concentrated on solutions of problems corresponding to isotropic materials. In this case, the eigenpairs can be computed analytically [4]. For example, in ref. [5] the eigenpairs for cracks along the interface of two dissimilar isotropic materials are explicitly given and in ref. [6] explicit eigenfunctions for a crack along the interface of anisotropic materials are provided as well. It is easier to obtain explicit eigenvalues than eigenfunctions (see, for example, refs [7, 8] for cases of up to three sub-domains). Complete exact solutions are restricted to rather simple geometries. Recently, the behavior of the solution for singular problems in general domains containing anisotropic materials has become of great interest because this is related to composite materials and electronic devices. For these general singular points, i.e. singularit:~es associated with corners, non-isotropic multi-material interfaces and abrupt changes in boundary conditions, the eige~£unctions, and very often the eigenvaiues as well, cannot he computed explicitly and have to be computed by numerical methods. Several numerical methods have been presented for the computation of the eigenvaiues during the last eight yrs, but only scant attention has been given to the computation of the eigenfunctions and their connection to the GSIFs. Leguillon and Sanchez-Palenci,'L [9] proposed the matrix and determinant methods based on the finite element approach to compute the eigenvalues and corresponding eigenfunctions. In reL [10] the eigenvalues are determined using the matrix method followed by an iterative approach and the "shooting method", while the eigenfunctions are evaluated using numerical integration based on the Runge-Kutta-Fehlberg method. This enables the calculation of the eigenpairs for anisotropic non-homogeneous materials with several types of boundary conditions. A detailed discussion on the evaluation of the eigenvalues, which involves the solution of a nonlinear eigenvaiue problem, is presented. Barsoum [11, 12] proposed the iterative finite eleme:at method for the computation of the first smallest eigenpair. This method, however, is much less efficient, robust and general than the matrix and the determinant methods. Gu and Belytschko [13] proposed a method based on stress function and the interpolation of displacements (based on the finite element method) for the computation of the eigenvalues, without discussing the computation of the eigenfunctions. This method, however, is not general and efficient. A new numerical method, called the modified Steklov method, for the reliable computation of the eigenpairs resulting from singularities due to corners, abrupt changes in material properties and boundary conditions was presented in refs [14] and [15]. The modified Steklov method, used in conjunction with the p-version of the finite element method, is robust, efficient and the computed eigenpairs are shown to converge strongly and accurately. The connection between ?Author to whom all correspondence should be addressed. Present address: Pearlstone Center for Aeronautical Engineering Studies, Dept of Meehanicai Engineering, Ben-Gurion University of the Negev, P.O. Box 653, Beer-Sheva 84105, Israel. 593
594
Technical Note
the numerical eigenpairs and the coefficients of the asymptotic expansion is presented in refs [16, 17], where high accuracy and efficiency are demonstrated for general singular points. This note addresses the question: how should the GSIFs computed using the numerically determined eigenfunctions be interpreted? Consider eq. (1) and assume for simplicity that M = 0 in eq. (1), i.e. the eigenfunctions do not contain log r terms. The eigenfunctions fj(0) are expressed as combinations of trigonometric functions such that max0 [fj(0) [ 4= 1. All numerical methods involve the computation of eigenvectors, obtained by solving an eigenproblem. These eigenvectors, defined up to a multiplying constant, are used for the computation of the eigenfunctions. The numerically computed eigenfunctions, denoted by gi(0), are related to [(0) by: fj(0) = ~jgj(0),
(2)
~,j being scalars. If one uses these gj(0) for the computation of the GSIFs, the obtained values depend, of course, on yj. Using g/instead of $ one obtains c~ = Q/?j, where we denote by cj the GSIFs corresponding to the numerical eigenfunctions. The physical quantity u is substantially unchanged whether one uses the numerically computed eigenfunctions or analytic ones. Correlation between the c/s and experimental data is achievable once u is measured and correlated with the computed data. The situation is more complicated when the j-th and j + 1-th eigenvalues are complex conjugates. In this case eq. (2) becomes: ~"'(0) _+ i~')(0) = yje'a,[g}R'(0)_+ i~"(0)],
(3)
where i = ~ - 1 and the superscript ~R)(respectively ")) represents the real (respectively imaginary) part of the function. The j-th and j + 1-th GSIFs are not independent, and are given as the real and complex parts associated with g~R>(0)_+ ig~J~(0).Therefore, the numerical eigenfunction involves a phase shift in addition to scaling when compared with the analytic one. Nevertheless, the physical quantity u is unchanged, while the GSIFs are related to the analytic ones by a complex number. If the analytic eigenfunctions are known, the scaling factors and phase shift between the "classical" and numerically computed GSIFs can be easily obtained. One has to apply on a model problem, called calibration model, tractions corresponding to the analytic eigenfunctions together with unit GSIFs, and compute the GSIFs using the numerical eigenfunctions. If the eigenfunctions are real, then ~j = l/c~. If the j-th eigenfunction is complex, both cj and cj + 1 have to be considered (see, for example, ref. [17]). Correlation with experimental data, for obtaining a critical value for cj at which failure initiates, does not require the analytic eigenfunctions. It suffices to correlate cj and gj(0) with experimental observations, provided that these numerical computed values are free from numerical errors. This can be assured by an adaptive numerical method in which the discretization errors can be controlled and estimated with high accuracy. These numerical methods are presented in ref. [I 7].
SUMMARY Classical fracture mechanics, applied to crack tip singularities, addresses mainly the computation of the stress intensity factors, the associated eigenfunctions usually are not addressed. This is because they are known explicitly and given as analytic functions. For general singular points, however, only numerical approximations of the eigenfunctions are available. An intrinsic property of all numerical approximations is that the eigenfunctions can be computed up to a scaling factor (sometimes this is a complex scaling factor). As a consequence, the generalized stress intensity factors are not independent of the numerical eigenfunctions and any reported value should include both. Furthermore, failure theories for general singular points have to be formulated such that explicit dependency on the eigenpairs and GSIFs is reflected.
REFERENCES [1] v. A. Kondratiev, Boundary value problems for elliptic equations in domains with conical or angular points. Trans. Moscow Math. Soc. 16, 227-313 (1967). [2] V. A. Kondratiev, The smoothness of a solution of Dirichlet's problem for second order elliptic equations in a region with picewise smooth boundary. Differ. Equat. 6, 1392-1401 (1970). [3] P. Grisvard, Elliptic Problems in Nonsmooth Domains. Pitman Publishing, England, U.K. (1985). [4] M. L. Williams, Stress singularities resulting from various boundary conditions in angular corners of plates in extension. Trans. ASME, J. appl. Mech. 19, 526-528 (1952). [5] J. R. Rice and G. C. Sih, Plane problems of cracks in dissimilar media. Trans. ASME, J. appl. Mech. 32, 418--423 (1965). [6] Z. Suo, Mechanics of interface fracture. Ph.D. Thesis, Harvard University, Cambridge, Massachusetts, U.S.A. (May 1989). [7] X. Ying, A reliable root solver for automatic computation with application to stress analysis of a composite plane wedge. Ph.D. Thesis, Washington University, St. Louis, Missouri, U.S.A. (December 1986). [8] J. P. Dempsey and G. B. Sinclair, On the singular behavior at the vertex ofa bi-material wedge. J. Elasticity 11, 317-327 (1981). [9] D. Leguillon and E. Sanchez-Palencia, Computation of Singular Solutions in Elliptic Problems and Elasticity. John Wiley, New York, NY, U.S.A. (1987). [10] P. Papadakis, Computational aspects oftbe determination oftbe stress intensity factors for two-dimensional elasticity. Ph.D. Thesis, University of Maryland at College Park, College Park, Maryland, U.S.A. (December 1988). [11] R. Barsoum, Theoretical basis of the finite element iterative method for the eigenvalue problem in stationary cracks. Int. J. numer. Meth. Engng 26, 531-539 (1988). [12] R. Barsoum, Application of the finite element iterative method to the eigenvalue problem of a crack between dissimilar media. Int. J. numer. Meth. Engng 26, 541-554 (1988). [13] L. Gu and T. Belytschko, A numerical study of stress singularities in a two-material wedge. Int. J. Solids Structures 31, 865-889 (1994).
Technical Note
595
[14] Z. Yosibash, Nmnerical analysis of singularities and first derivatives for elliptic boundary value problems in two-dimensionas. D.Sc. Thesis, Sever Institute of Technology, Washington University, St. Louis, Missouri, U.S.A. (August 1994). [15] Z. Yosibash and B. A. Szab6, Numerical analysis of singularities in two-dimensions. Part 1: computation of eigenpairs. Int. J. nurner. Meth. Engng 38, 2055-2082 (1995). [16] B. A. Szab6 and Z. Yosibash, Numerical analysis of singularities in two-dimensions. Part 2: computation of the generalized flux/stress intensity factors. Int. J. numer. Meth. Engng 39, 409-434 (1996). [17] Z. Yosibash and B. A. Szab6, Generalized stress intensity factors in linear elastostatics. Int. J. Fracture 72, 223-240 (1995). (Received 27 October 1994)
Engineering Fracture Mechanics Vol. 54, No. 4, pp. 593-595, 1996 Copyright © 1996 ElsevierScienceLtd. Printed in Great Britain. All rights reserved 0013-7944(95)00090-9 0013-7944/96 $15.00+ 0.00
TECHNICAL NOTE
A NOTE ON NUMERICALLY COMPUTED EIGENFUNCTIONS AND GENERALIZED STRESS INTENSITY FACTORS ASSOCIATED WITH SINGULAR POINTS ZOHAR YOSIBASHI" and BARNA A. SZABO Center for Computational Mechanics, Washington University, Campus Box 1129, St. Louis, MO 63130, U.S.A. Abstract--The solution of linear elastostatic problems in the neighborhood of singular points is characterizedi by a sequence of eigenpairs and their coefficients, called generalized stress intensity factors (GSIFs). For general singular points, as cracks in anisotropic multimaterial interfaces, and V-notches in composite materials, only numerical approximations of the eigenfunctions are available. This note addresses the: question: how should the GSIFs computed using the numerically determined eigenfunctions be interpreted? Copyright © 1996 Elsevier Science Ltd.
DISCUSSION THE EXACTsolutions of linear elliptic partial differential equations which represent heat conduction or elastic deformation typically have a number of singular points. Each singular point is characterized by a sequence of eigenpairs (see refs [I-3]). For example, the solution of linear elasticity problems in two dimensions in the vicinity of a singular point, can be expanded in the form: M
u = ~.
~ C~,.r~J.lnmr.fjm(O) + u*,
(1)
jffil m ~ 0
where u* is smoother than the terms in the sum, Cjmare the coefficients of the asymptotic expansion (called the generalized stress intensity factors - - GSIFs), ~ and f~(0) are eigenpairs which depend on the differential operator and boundary conditions. The detelTnination of these eigenpairs is an essential step towards reliable computation of the GSIFs, which are of great practical importance because failure theories directly or indirectly involve these coefficients. Most of the research performed in the past concentrated on solutions of problems corresponding to isotropic materials. In this case, the eigenpairs can be computed analytically [4]. For example, in ref. [5] the eigenpairs for cracks along the interface of two dissimilar isotropic materials are explicitly given and in ref. [6] explicit eigenfunctions for a crack along the interface of anisotropic materials are provided as well. It is easier to obtain explicit eigenvalues than eigenfunctions (see, for example, refs [7, 8] for cases of up to three sub-domains). Complete exact solutions are restricted to rather simple geometries. Recently, the behavior of the solution for singular problems in general domains containing anisotropic materials has become of great interest because this is related to composite materials and electronic devices. For these general singular points, i.e. singularit:~es associated with corners, non-isotropic multi-material interfaces and abrupt changes in boundary conditions, the eige~£unctions, and very often the eigenvaiues as well, cannot he computed explicitly and have to be computed by numerical methods. Several numerical methods have been presented for the computation of the eigenvaiues during the last eight yrs, but only scant attention has been given to the computation of the eigenfunctions and their connection to the GSIFs. Leguillon and Sanchez-Palenci,'L [9] proposed the matrix and determinant methods based on the finite element approach to compute the eigenvalues and corresponding eigenfunctions. In reL [10] the eigenvalues are determined using the matrix method followed by an iterative approach and the "shooting method", while the eigenfunctions are evaluated using numerical integration based on the Runge-Kutta-Fehlberg method. This enables the calculation of the eigenpairs for anisotropic non-homogeneous materials with several types of boundary conditions. A detailed discussion on the evaluation of the eigenvalues, which involves the solution of a nonlinear eigenvaiue problem, is presented. Barsoum [11, 12] proposed the iterative finite eleme:at method for the computation of the first smallest eigenpair. This method, however, is much less efficient, robust and general than the matrix and the determinant methods. Gu and Belytschko [13] proposed a method based on stress function and the interpolation of displacements (based on the finite element method) for the computation of the eigenvalues, without discussing the computation of the eigenfunctions. This method, however, is not general and efficient. A new numerical method, called the modified Steklov method, for the reliable computation of the eigenpairs resulting from singularities due to corners, abrupt changes in material properties and boundary conditions was presented in refs [14] and [15]. The modified Steklov method, used in conjunction with the p-version of the finite element method, is robust, efficient and the computed eigenpairs are shown to converge strongly and accurately. The connection between ?Author to whom all correspondence should be addressed. Present address: Pearlstone Center for Aeronautical Engineering Studies, Dept of Meehanicai Engineering, Ben-Gurion University of the Negev, P.O. Box 653, Beer-Sheva 84105, Israel. 593
594
Technical Note
the numerical eigenpairs and the coefficients of the asymptotic expansion is presented in refs [16, 17], where high accuracy and efficiency are demonstrated for general singular points. This note addresses the question: how should the GSIFs computed using the numerically determined eigenfunctions be interpreted? Consider eq. (1) and assume for simplicity that M = 0 in eq. (1), i.e. the eigenfunctions do not contain log r terms. The eigenfunctions fj(0) are expressed as combinations of trigonometric functions such that max0 [fj(0) [ 4= 1. All numerical methods involve the computation of eigenvectors, obtained by solving an eigenproblem. These eigenvectors, defined up to a multiplying constant, are used for the computation of the eigenfunctions. The numerically computed eigenfunctions, denoted by gi(0), are related to [(0) by: fj(0) = ~jgj(0),
(2)
~,j being scalars. If one uses these gj(0) for the computation of the GSIFs, the obtained values depend, of course, on yj. Using g/instead of $ one obtains c~ = Q/?j, where we denote by cj the GSIFs corresponding to the numerical eigenfunctions. The physical quantity u is substantially unchanged whether one uses the numerically computed eigenfunctions or analytic ones. Correlation between the c/s and experimental data is achievable once u is measured and correlated with the computed data. The situation is more complicated when the j-th and j + 1-th eigenvalues are complex conjugates. In this case eq. (2) becomes: ~"'(0) _+ i~')(0) = yje'a,[g}R'(0)_+ i~"(0)],
(3)
where i = ~ - 1 and the superscript ~R)(respectively ")) represents the real (respectively imaginary) part of the function. The j-th and j + 1-th GSIFs are not independent, and are given as the real and complex parts associated with g~R>(0)_+ ig~J~(0).Therefore, the numerical eigenfunction involves a phase shift in addition to scaling when compared with the analytic one. Nevertheless, the physical quantity u is unchanged, while the GSIFs are related to the analytic ones by a complex number. If the analytic eigenfunctions are known, the scaling factors and phase shift between the "classical" and numerically computed GSIFs can be easily obtained. One has to apply on a model problem, called calibration model, tractions corresponding to the analytic eigenfunctions together with unit GSIFs, and compute the GSIFs using the numerical eigenfunctions. If the eigenfunctions are real, then ~j = l/c~. If the j-th eigenfunction is complex, both cj and cj + 1 have to be considered (see, for example, ref. [17]). Correlation with experimental data, for obtaining a critical value for cj at which failure initiates, does not require the analytic eigenfunctions. It suffices to correlate cj and gj(0) with experimental observations, provided that these numerical computed values are free from numerical errors. This can be assured by an adaptive numerical method in which the discretization errors can be controlled and estimated with high accuracy. These numerical methods are presented in ref. [I 7].
SUMMARY Classical fracture mechanics, applied to crack tip singularities, addresses mainly the computation of the stress intensity factors, the associated eigenfunctions usually are not addressed. This is because they are known explicitly and given as analytic functions. For general singular points, however, only numerical approximations of the eigenfunctions are available. An intrinsic property of all numerical approximations is that the eigenfunctions can be computed up to a scaling factor (sometimes this is a complex scaling factor). As a consequence, the generalized stress intensity factors are not independent of the numerical eigenfunctions and any reported value should include both. Furthermore, failure theories for general singular points have to be formulated such that explicit dependency on the eigenpairs and GSIFs is reflected.
REFERENCES [1] v. A. Kondratiev, Boundary value problems for elliptic equations in domains with conical or angular points. Trans. Moscow Math. Soc. 16, 227-313 (1967). [2] V. A. Kondratiev, The smoothness of a solution of Dirichlet's problem for second order elliptic equations in a region with picewise smooth boundary. Differ. Equat. 6, 1392-1401 (1970). [3] P. Grisvard, Elliptic Problems in Nonsmooth Domains. Pitman Publishing, England, U.K. (1985). [4] M. L. Williams, Stress singularities resulting from various boundary conditions in angular corners of plates in extension. Trans. ASME, J. appl. Mech. 19, 526-528 (1952). [5] J. R. Rice and G. C. Sih, Plane problems of cracks in dissimilar media. Trans. ASME, J. appl. Mech. 32, 418--423 (1965). [6] Z. Suo, Mechanics of interface fracture. Ph.D. Thesis, Harvard University, Cambridge, Massachusetts, U.S.A. (May 1989). [7] X. Ying, A reliable root solver for automatic computation with application to stress analysis of a composite plane wedge. Ph.D. Thesis, Washington University, St. Louis, Missouri, U.S.A. (December 1986). [8] J. P. Dempsey and G. B. Sinclair, On the singular behavior at the vertex ofa bi-material wedge. J. Elasticity 11, 317-327 (1981). [9] D. Leguillon and E. Sanchez-Palencia, Computation of Singular Solutions in Elliptic Problems and Elasticity. John Wiley, New York, NY, U.S.A. (1987). [10] P. Papadakis, Computational aspects oftbe determination oftbe stress intensity factors for two-dimensional elasticity. Ph.D. Thesis, University of Maryland at College Park, College Park, Maryland, U.S.A. (December 1988). [11] R. Barsoum, Theoretical basis of the finite element iterative method for the eigenvalue problem in stationary cracks. Int. J. numer. Meth. Engng 26, 531-539 (1988). [12] R. Barsoum, Application of the finite element iterative method to the eigenvalue problem of a crack between dissimilar media. Int. J. numer. Meth. Engng 26, 541-554 (1988). [13] L. Gu and T. Belytschko, A numerical study of stress singularities in a two-material wedge. Int. J. Solids Structures 31, 865-889 (1994).
Technical Note
595
[14] Z. Yosibash, Nmnerical analysis of singularities and first derivatives for elliptic boundary value problems in two-dimensionas. D.Sc. Thesis, Sever Institute of Technology, Washington University, St. Louis, Missouri, U.S.A. (August 1994). [15] Z. Yosibash and B. A. Szab6, Numerical analysis of singularities in two-dimensions. Part 1: computation of eigenpairs. Int. J. nurner. Meth. Engng 38, 2055-2082 (1995). [16] B. A. Szab6 and Z. Yosibash, Numerical analysis of singularities in two-dimensions. Part 2: computation of the generalized flux/stress intensity factors. Int. J. numer. Meth. Engng 39, 409-434 (1996). [17] Z. Yosibash and B. A. Szab6, Generalized stress intensity factors in linear elastostatics. Int. J. Fracture 72, 223-240 (1995). (Received 27 October 1994)