ELASTOSTATIC STRESS ANALYSKS 057 ENE ANISOTROPIC PLATES WITH CENT RALLY LOCATED TRACTIONeFREE CRACKS Thesfs for the Degree of Pb D‘ MICHIGAN STATE umvzasm GEORGE STEPHEN cvm’bvts: 1972 unfit-"9 This is to certify that the thesis entitled ELASTOSTATIC STRESS ANALYSIS OF FINITE ANISOTROPIC PLATES WITH CENTRALLY LOCATED TRACTION-FREE CRACKS presented by GEORGE S. GYEKENYESI has been accepted towards fulfillment of the requirements for Ph . D o dggree in MEChaniQS or professor Date O I '72 0-7639 LIBRARY Michigan State Uni v entity ABSTRACT ELASTOSTATIC STRESS ANALYSIS OF FINITE ANISOTROPIC PLATES WITH CENTRALLY LOCATED TRACTION-FREE CRACKS By George Stephen Gyékényesi A mapping-collocation method is developed for the elaStostatic stress analysis of finite anisotropic plates with centrally located traction-free cracks. The essence of the method is as follows: 1. The crack is mapped into the unit circle. 2. The boundary conditions on the crack are satisfied exactly by expressing one of the stress potentials in the form of the other. 3. A form of representation is assumed for the remaining unknown stress potential in the parametric plane. 4. The boundary conditions on the outer boundary of the region are satisfied by the method of least squares boundary collocation. In order to demonstrate the feasibility of the method, rectan- gular orthotropic plates with centrally located traction-free cracks under constant tensile and shear loads are analyzed. A parametric study of the finite plate stress intensity factors is presented show— ing the effects of varying material properties, orientation angle, crack length to plate width and plate height to plate width ratios. In general, some of the more significant results can be summarized as follows: George Stephen Gyékényesi l. Rectangular orthotropic Plate with a centrally located traction- free crack under constant tensile load acting in the direction perpendicular to the crack. (a) The opening mode stress intensity factor increases with decreasing EZZ/Ell ratios, where E22 and E11 are Young's moduli of elasticity for an orthotropic mate— rial and 0 < E < l. 22/E11 (b) The opening mode stress intensity factor increases as the crack length to plate width ratio increases. (c) The opening mode stress intensity factor increases as the plate height to plate width ratio decreases. (d) Considering a square plate, the maximum opening mode stress intensity factor occurs at the value of the orientation angle of 900 while the minimum opening mode stress intensity factor is obtained at O0 for any constant crack length to plate width ratio. (e) The presence of the sliding mode stress intensity fac- tor is due to values of the orientation angle other than 00 or 900 and to the finite size of the plate. 2. Orthotropic square plate with a centrally located traction-free crack under constant shear loading. (a) The sliding mode stress intensity factor increases with decreasing EZZ/Ell ratios. (b) The sliding mode stress intensity factor increases as the crack length to plate width ratio increases. George Stephen GyékényeSi (c) The minimum sliding mode stress intensity factor always occurs at the zero value of the orientation angle while the maximum sliding mode Stress intensity factor can occur at various values of the orientation angle for conStant crack length to plate Width ratios. (d) The presence of the opening mode stress intensity fa;— tor is due to values of the orientation angle differ- ent from O0 or 900 and to the finite size of the plate. In addition to the development of the mapping-collocation method and to the parametric study of stress intensity factors, the infinite anisotrOpic plate solutions of Savin and Lekhnitskii are derived by considering a finite rectangular anisotropic plate with a centrally located traction-free crack and extending its dimensions to infinity. ELASTOSTATIC STRESS ANALYSIS OF FINITE ANISOTRDPIC PLATES WITH CENTRALLY LOCATED TRACTION-FREE CRACKS BY George Stephen.Gyékényesi A THESIS Submitted to Michigan State university in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Matallurgy, Mechanics and Materials Science 1972 To my parents Gvékénvesi Gyorgy Lészlé and Gy. Korcsmér Katalin ii ACKNOWLEDGMENTS The author wishes to express his Sincere gratitude to his major thesis advisor, Dr. Alexander Mendelson, for his generous assiStan3e and guidance which ultimately led to the completion of the work re- ported in this dissertation. The author is also grateful to the chairman of his doctoral com- mittee, Dr. W. N. Sharpe and the members thereof, Drs. R. W. Little, W. A. Bradley and N. L. Hills who with their helpful and understand— ing attitude greatly encouraged him in the course of his research. Many thanks are due to Mr. P. F. Michaelis who gave his most needed help in the programing of the IBM TSS/36O computer. The author wishes also to express his deep appreCiation to the National Aeronautics and Space Administration in general, under whose auspices this research was completed, and in particular, to Miss Gertrude R. Collins of the Training Office for her dedication in man- aging the author's program of study at NASA's Lewis Research Center. Last but not least, recognition is in due order to the author's wife, Sarolta and the children, Gyorgy, lstvén and Enikg who with their considerate patience most selflessly contributed to the author's efforts. iii TABLE OF CONTENTS Page LIST OF TABLES. . . . . . . . . . . . . . . . . . . . . . . . . . . vii LIST OF FIGURES . . . . . . . . . . . . . . . . . . . . . . . . . . viii LIST OF SYMBOLS . . . . . . . . . . . . . . . . . . . . . . . . . . x 1. INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . . . l 2. PLANE PROBLEM FORMULATIONS FOR ANISOTROPIC MATERIALS . . . . . . 5 2.1 INTRODUCTION. . . . . . . . . . . . . . . . . . . . . . 5 2.2 STRESS FUNCTION FORMULATION OF PLANE STRESS PROBLEMS. . . . . . . . . . . . . . . . . . . . . . . . 6 2.3 DISPLACEMENT FUNCTION FORMULATION OF PLANE STRESS PROBLEMS . . . . . . . . . . . . . . . . . . . . 9 3. GENERAL CONSIDERATIONS FOR PROBLEM SOLUTIONS CONCERNING MULTIPLY-CONNECTED REGIONS . . . . . . . . . . . . . . . . . . . 13 3.1 INTRODUCTION. . . . . . . . . . . . . . . . . . . . . . 13 3.2 SOLUTION OF THE DIFFERENTIAL EQUATION . . . . . . . . . 13 3.3 STRESS COMPONENTS . . . . . . . . . . . . . . . . . . . 15 3.4 DISPLACEMENTS . . . . . . . . . . . . . . . . . . . . . 17 3.5 BOUNDARY CONDITIONS . . . . . . . . . . . . . . . . . . 18 3.51 First Fundamental Problem . . . . . . . . . . 18 3.52 Second Fundamental Problem . . . . . . . . . 19 3.6 ON THE COMPLEX PARAMETERS OF LEKHNITSKII. . . . . . . . 20 3.7 RESULTANT FORCE AND RESULTANT MOMENT. . . . . . . . . . 23 3.8 FORMS OF THE STRESS POTENTIALS FOR FINITE MULTIPLY- CONNECTED REGIONS . . . . . . . . . . . . . . . . . . . 24 iv 4. AN APPROACH TO THE SOLUTION OF FINITE GEOMETRY TRACTION- FREE CRACK PROBLEMS. 4.1 INTRODUCTION. . . . . 4.2 THE MAPPING OF A CRACK AND ITS EXTERIOR INTO THE UNIT CIRCLE AND ITS EXTERIOR. 4.3 ON THE CONTINUATION AND FORM OF THE STRESS POTENTIALS. , . . . . . 4.4 BOUNDARY CONDITIONS FOR FINITE REGIONS CONTAINING A TRACTION-FREE CRACK o - o o o r: a c 4.5 METHOD OF LEAST SQUARES BOUNDARY COLLOCATION. 4.6 STRESS COMPONENTS 4.7 CRACK TIP STRESS COMPONENTS AND STRESS INTENSITY FACTORS O O O O O O .J- C O O O O O O O C 4.8 DISPLACEMENTS . . . . . . , . . . . . . . . . . . 5. THE ANISOTROPIC INFINITE PLATE CONTAINING A CENTRALLY LOCATED TRACTION-FREE CRACK UNDER CONSTANT LOAD. 5.1 INTRODUCTION. 5.2 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT LOAD IN THE y-DIRECTION . . . 5.3 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT SHEAR LOADING O O D O t O O O 5.4 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT LOAD IN THE x-DIRECTION 6. PROBLEM FORMULATION FOR A FINITE ORTHOTROPIC RECTANGULAR PLATE CONTAINING A CENTRALLY LOCATED TRACTION-FREE CRACK UNDER CONSTANT LOAD. . . . . . . . . . . . 6.1 INTRODUCTION. . . . . . . . . . . . . . 6.2 STRESS SYMMETRY . . . . . . . . . . . . . . . . 6.3 DETERMINATION OF INTEGRATION CONSTANTS. . . 6.4 BOUNDARY CONDITIONS . . . . . . . . . . . 6.5 SOME REMARKS ABOUT THE COMPUTER PROGRAM . 3O 30 31 38 43 46 49 51 59 59 59 64 67 7O 7O 7O 71 73 75 7. RESULTS AND DISCUSSION OF THE SOLUTIONS OF THE TENSION AND SHEAR PROBLEMS. . . . . . . . . . . . . . . . . . . . . . . 77 7.1 INTRODUCTION . . . . . . . . . . . . . . . . . . . . . 77 7.2 SPECIFICATIONS OF PARAMETERS TREATED . . . . . . . . . 77 7.3 EXAMINATION OF THE MAPPING-COLLOCATION METHOD. . . . . 79 7.4 SOLUTION OF THE TENSION PROBLEMS OF AN ORTHOTROPIC RECTANGULAR PLATE WITH A CENTRAL CRACK . . . . . . . . 84 7.41 Effects of Material Properties on the Stress Intensity Factors . . . . . . . . . . 84 7.42 Effects of Plate Width and Plate Size on the Stress Intensity Factors. . . . . . . 87 7.43 Effects of Orientation Angle and Plate Width on the Stress Intensity Factors . . . 89 7.44 Stress Distribution on the y = 0, x > a Line I O O O O O O O O O O O O O O O O 9 l 7.5 SOLUTION OF THE SHEAR PROBLEM OF AN ORTHOTROPIC SQUARE PLATE WITH A CENTRAL CRACK. . . . . . . . . . . 92 7.51 Effects of Material Properties on the Stress Intensity Factors . . . . . . . . . . 92 7.52 Effects of Orientation.Ang1e and Plate Width on the Stress Intensity Factors . . . 95 7.53 Stress Distribution on the y = O, x > a Line. . . . . . . . . . . . . . . . . 9b 8. SUMMARY AND RECOMMENDATIONS . . . . . . . . . . . . . . . . . . 122 9. LIST OF REFERENCES. . . . . . . . . . . . . . . . . . . . . . . 126 10. APPENDICES A. DISCUSSION OF ELASTIC CONSTANTS . . . . . . . . . . . . 128 B. COMPUTER PROGRAM. . . . . . . . . . . . . . . . . . . . 150 vi Table LIST OF TABLES Comparison of opening made stress intensity factors for rectangular isotropic plates with central cracks. Typical pattern of convergence of the stress intensity factorsofor tensile loading (h/c = l, a/c = 1/2, 6 = 45 , Tt = 1). Effects of material properties on the stress intensity factorsofor tensile loading (h/c = l, a/c = 2/3, 5 = 45 , Tt = 1). Effects of Eflate width and plate size on the stress intensity factors for tensile loading (6 = 45 , Tt = l). - . . . . . . . Effects of :plate Width and orientation angle on the stress intensity factors for tensile loading (h/c = l, Tt = 1). Stress distribution on the y = 0, x > a line 0 for tensile loading (h/c = l, a/c = 2/3, 5 = 45 , T = l). 0 C 3 C O x t O Q I' C O O O I Effects of material properties on the stress intensity factor gor shear loading (h/c = 1, a/c = 2/3, 6 = 45 , TS = 1). Effects of tflefie width and orientation angle on the stress intensity factors for shear loading (h/c = 1, TS = 1). Stress distribution on the y = O, > for shear loading (h/c = l, a/c = TS = I). . . . . . . . . - . . . x a ling 2/3, 5 = 45 , vii 98 101 104 107 110 113 116 119 Figure 10. 11. 12. LIST OF FIGURES Curve segment loaded by.the complex stress, X.n + iYn. . . The multiply connected region, L, in the z—, 21- and 22-p1anes. . . . . . . . The mapping of the unit circle into a crack . The mapping of the doubly connected regions L, L1 and L into their corresponding parametric planes 0 O O O O 0 I O O O O I O I O O O O The region LC as continued across the unit circle . The crack tip coordinate system . The tension problems of a rectangular plate with a central crack . . . . . . . . . . . . . . . The shear problem of a square plate with a central craCR O O C O O C C O O O O O O O O O I O I 0 Typical pattern of convergence of the opening mode .stress intensity factor for tensile loading (h/c = l, a/c = 1/2, 6 = 45°, T = 1) . . . . . . . . . t Typical pattern of convergence of the sliding mode stress intensity factor for tensile loading (h/c = l, a/c = 1/2, 6 = 45°, Tt = 1) . . . . . . . . . . . . . Effects of material properties on the opening mode stress intensity factor for tensile loading (h/c = l, a/c = 2/3, 6 - 45°, T: = 1) . . . . . . . . . . . . Effects of material.properties on the sliding mode stress intensity.gactor for tensile loading (h/c = l, a/C32/3,6=45 ’Tt=l) o o o o o o o o o 0 viii 29 55 56 57 57 58 58 99 100 102 103 Effects of plate width. and plate size on the opening mode stress intensity factor for tensile loading (5 = 45‘, T = 1). . . . . . . . . . . . . 105 Effects of plate width an sliding mode stress inter c tensile loading (5 = 45d, '). . . - . a . . . . . . . 106 Effects of orientation angle and plate width on the opening m:de stress untetsity rector ' ceding (h;c = I, T, 4 I). . . . 4 . . . . . . 108 k Effects of orientation angle and place width on the sliding mode stress intensity factor for tensile loading (h/c = 1, Tt = 1). , . . . . . . . . . 109 Normal Stress distribution on the y a O, x > a line for tensile loading (h/c = 1, a/C = 2/3, 6 a 4:)» , TE = l) v ' o o o a o o o o o o 11.1 Shear stress distribution on the y = 0, x > a line for tensrle loading (h/c = 1, 8/": = 2/3, 6 = ASU, II‘t 2' l) '. t -. a o o o o o o o o a o o 112 Effects of material properties on the opening mode Stress intensity factog for shear loading (hl'C - l, a/L 1 23i3, 5 = 45 , Ts = l) o o o a e o o o o o o 114 Effects of material properties on.the sliding mode stress intensity factog for shear loading (h/c = l, a/: = 2/3, 5 s 45 , T8 = 1). . . . . . . . . . . 115 Effects of orientation angle and plate width on the opening mode stress intensrty factor for shear loading (h/c = 1, TG = l). . . t . . . . . . . - . 117 Effects of orientation angle and plate width on the sliding mode stress intensity factor for shear loading (h/c = 1, T5 = l). . . . . . . . . . . . . . 118 Normal stress discribution on the y = 0, x > a line for shear loading (h/c = 1, a/c = 2/3, 5 = 45°, T5 = i) . . , . . . . . 120 Shear stress distribution on the y = 0, x > a line fgr shear loading (h/c = l, a/ 5 = 45 , T8 = 1) . . . . . . . . . . . . , . . . . . . . . 121 ix LIST OF SYMBOLS The following list contains the more commonly used symbols and their representative meaning unless otherwise defined in the text. In general, all symbols are defined when first introduced. a A1n 1n ‘4 %"41 E11’ E22’ E33 F F G G 1n’ 2n’ 1n’ 2n 613’ G12’ 23 half crack length th . . n complex constant in a Laurent series assocr- ated with positive exponent terms th . . . n complex constant in a Laurent series assoc1- ated with negative exponent terms half-width of a rectangular plate elastic compliance matrix ijth element of elastic compliance matrix real arbitrary integration constants nth modified complex constant associated with the positive exponent terms in a Laurent series th . . n modified complex constant associated Wlth the negative exponent terms in a Laurent series complex constant depending on material properties modified complex constants depending on the sum and difference of A and B , respectively 11 11 Young's.moduli of elasticity with respect to coor- dinate directions complex quantities utilized in the construction of the coefficient matrix in the least squares bound— ary collocation method shear moduli of an anisotropic material half-height of a rectangular plate Opening and sliding mode stress intensity factors, respectively total number of collocation points on the outer boundary of a region with a central crack truncation numbers on infinite sums associated with the negative and posrcive exponent terms in a Laurent series, respectively complex constants depending on material properties crack tip polar coordinates elastic stiffness matrix .,th . . . 13 element of elastic stiffness matrix constant shear stress applied to the boundary of a square plate constant tensile stress applied to the boundary of a rectangular plate displacements in the x- and y-direction, respectively rectangular coordinates applied boundary stresses acting in the x— and y—direction, respectively conventional complex variable complex variables modified by material properties reduced compliance matrix ,.th . . 1] element of reduced compliance matrix real parts of Lekhnitskii's material parameters reduced stiffness matrix ._th . 1] element of reduced stiffness matrix imaginary parts of Lekhnitskii's material parameters xi 8 8 xx YY XY C e ”12' ”21’ v23, V32’ ”13’ V31 “1, ”2 O ¢1(zl>. ¢2<22) a» a» material orientation angle components of the strain tensor, rectangular coordinates complex variable in parametric plane polar angle in parametric plane Poisson's ratios in an anisotropic material Lekhnitskii's material parameters complex variable on the unit circle components of the stress tensor, rectangular coordinates complex stress.potential functions Airy's stress function displacement function xii 1. INTRODUCTION The recently extended application of anisztrcpi: materials in vari- ous fields of structural design makes one aware :f a whole class of new problems in anisotropic elasticity theory. in order to study the frac- ture phenomena of these materials, the problem of a centrally cracked finite plate seems to be of immediate importanze. Methods of solution for centrally cracked finite isotropic plates were developed by A. S. Kobayashi (ref. 1) and O. L. Bowie and D. M. Neal (ref. 2) while the problem of the finite anisotropic plate remained fairly neglected. The plane problem of a centrally cracked finite anisotropic plate could be considered as a special case of the more general problem of the anisotr0pic.plate containing an elliptic hole. The limiting case of this problem, i. e., the problem of an elliptic hole in an infinite ani- sotrOpic plate, has extensively been investigated by various authors over the past thirty years. The complex variable approach of N. I. Muskhelishvili (ref. 3) as adopted by S. G. Lekhnitskii (ref. 4) to plane anisotropic theory was used by both Lekhnitskii and G. N. Savin to solve the problem of the infinite anisotropic plate bounded by an ellipse. Lekhnitskii (refs. 4, 5) employed the method of series ex- pansion while Savin (refs. 6, 7) applied.Schwarz‘s formula in order to obtain the complex Stress potentials of Lekhnitskii. Simultaneously with these achievements, A. E. Green, in a series of articles published.in the early 1940’s, developed his own method, 2 also solving the problem of an elliptic hole in an infinite anisotropic region. His results are condensed in his Theoretical Elasticity (ref. 8). Another approach to plane elastic problems is advanced in L. M. Milne-Thomson's Plane Elastic Systems (ref. 9) in handling both iso- trOpic and anisotropic plane problems by a "semi-unified" method. He also discusses the problem of the elliptic hole in an infinite aniso- tropic region. An entirely different method is due to D. D. Aug and M. L. Williams (ref. 10). They used a formulation in integral equations for a cen- trally cracked orthotropic plate. This method, however, applies only for infinite orthotropic plates with zero orientation angle. The semi-infinite specially orthotropic plate problem with a cen- trally located crack was first discussed and solved by A. Mendelson and S. W. Spero (ref. 11). The solution was obtained by the application of finite Fourier transforms resulting in an integral equation for the crack opening. The concept of the extension of Fracture Mechanics from isotropic to anisotropic media was proposed by P. C. Paris and G. C. Sih in 1961 (ref. 12). If their approach is adopted, the stress intensity factors can be determined directly from the anisotropic stress potentials in a manner similar to that introduced by Sih for isotropic materials. As a further development, E. M. Wu examined the conditions neces- sary for the application of Fracture Mechanics to anisotrOpic materials and also verified them experimentally for orthotropic plates (ref. 13). He concluded that the crack tip stress singularity is of the same order 3 as that of isotropic materials and the stress intenSity factors are sim- ilar to those of the isotropic case. The subject of this dissertation is a logical continuation of the above listed work as it is applied and extended to the problem of a fin— ite anisotropic region With a centrally io-a:ed traction—free crack. The proposed method of solution parallels Bow1e's modified mapping- collocation technique (ref. 2) which was successfully used in the case of isotropic materials. Since the presence of flaws in a material is of critical nature, the main objective of the work reported herein was the construction of a practical method toward defining and obtaining the stress intensity factors of finite anisotropic plates with central traction—free cracks. The method of solution can briefly be described as follows: 1. The crack 18 mapped into the unit Circle. 2. The complex stress potentials of Lekhnitskii are found by exactly satisfying the zero traction conditions on the crack with the application of Muskhelishvili's function extension concept a-ross the unit circle. 3. Upon expanding the remaining undefined stress potential in a Laurent series, the conditions on the outer boundary or the region are approximated by the method of least squares boundary collocation. The stress intensity factors for both the opening and sliding modes are then computed Since they are shown to be functions of the coefficients of the Laurent series already defined. In addition, expressrons for the 4 stress components and displacements are given in general for the case of the traction-free crack problem. As an illustration of the practicality of the method, the problem of an orthotrOpic rectangular plate with a central traction—free crack is solved under constant tension and shear loading. The effects of varying material properties, orientation angles, crack length to plate width and plate height to plate width ratios are investigated and dis— cussed at some length. It is shown for plates of various crack length to plate width ratios that in the tension problem, the maXimum opening mode stress intensity factor is obtained in each case when the "strong" axis of the material is perpendicular to the crack. In the shear prob- lem, the maximum sliding mode stress intensity factor can occur at var- ious values of the orientation angle. In addition to the finite plate results, the problem of infinite anisotropic plates with central cracks is also discussed for the cases of constant tension and shear loading. It is shown that the infinite plate solutions of Lekhnitskii and Savin (refs. 4, 6) can also be ob- tained by extending the dimensions of a rectangular plate to infinity. 2. PLANE PROBLEM FORMULATIONS FOR ANISOTROPIC MATERlALS 2.1 INTRODUCTION In linear anisotropic elasti;ity there are two plane problems which may be discussed from the classrcal point of view: 1. The state of plane stress, corresponding to a plane plate of constant thickness loaded by forces in the plane of the plate; 2. The state of plane strain, corresponding to a long (theo- retically infinite) body acted upon by loads which are uni- formly distributed in the infinite direction and have no components normal to the finite planes. These cases are similar to those in plane linear isotropic elasti- city (ref. 15)— It can easily be shown that the difference between the plane stress and the plane strain problems is inherent in the use of elastic constants within the stress-strain relations as given in Appendix A. Therefore, any method of problem formulation for the plane stress problem is also valid for the state of plane strain (or vice versa) With the proper interchange of the elastic constants. In the following, two methods of problem formulation will be demonstrated: 1. The conventional stress function formulation (ref. 4); 6 2. The displacement function formulation corresponding to Marguerre's displacement function approach (refs. 15, 16) in linear isotropic elasticity. The formulations developed on the following pages are given for plane stress problems with no body forces. 2.2 STRESS FUNCTION FORMULATION 0F PLANE STRESS PROBLEMS Assuming that the stress components are derivable from a function such that 32¢ 32¢ 32¢ Ox=—-2—;O =——2_;Ox=_3x8 2' X By YY 8x Y y one finds that the equilibrium equations 30 80 £+JX=O 8x 8y 2. 30 80 __§Z.+.__ZZ.= 0 3x 3y are identically satisfied. The displacements, u and v, satisfy the strain compatibility equation 825x azsy 328 + =2—-§Z 2. ay2 ax2 3x3y also identically. Therefore, the function, ¢, can be related to the displacements, u and v, through the use of the strain-displacement relations 1 2 3 Bu E2 = ”T" XX OX 6 = 3X. 2. yy 8y 26 = du +‘31 xy 8y 8x and Hooke's Law for the state of anisotropic plane scress o. = B..e\ (1,j = l, 2, 6) 2. in the following manner: 2 Li 3 _ .3. i 3y2 B11 dx + 816 a u I 816 3x + 812 ay V 322 a a a a 3x2 = (812 3§'+ 826 E?) u + (826 §§'+ B22 5;) V 2' __3_<.P_= _3 .9. .2. 3 8x3y 816 8x + 866 By u + 866 8x + 826 68y V Upon the elimination of u and v froanqs.2.6, the following fourth order partial differential equation is obtained: 4 _ 2 3 ¢ 3 2 (811866 816) 3X 4 I 2(8118 26 812816) a 3. X 6y 4 2 3 ¢ + (28168 26 2812866 + 811822 ‘ 812) 2 2 Bx 8y 349 2 342 + 2(8168 22 812826) 3 + (8228 66 826) 4 = O axay 8y 8 Using the relations 42 from Appendix A, this differential equation can be written as, 4 4 4 4 4 M _ ii. _3_1_ _ .12.. Li = C22 4 2C26 3 + (C66 + 2C12) 2 2 2C16 3 + C11 4 0 3x 3x 8y 8x 3y Bxay By 2.8 The displacements are found from the stress-strain relations 6 = C ,o (i,j = 1, 2, 6) 2.9 by substituting the strain-displacement relations and the stress func— tion, ¢. into the stress-strain relations, and integrating the result- ing expressions, I: I 2 - 33—2 ii _ .39. _ C11 / y2 dx + C12 8x Ci6 8y + “o “’63' X 2.10 < ll 2 Li 2i- 2i C22 0/“ ax2 dy + C12 8y C26 8x + vo + wox Y IhiEqs.2.10, the terms, uo, woy, v0, wox, designate rigid body motions. This formulation was first obtained by S. G. Lekhnitskii in the 1930's (ref. 4). The general solution of the differential equation (Eq. 2.8) sub- ject to displacement and/or stress boundary conditions would then be the solution of the linear anisotropic plane stress problem in the case of zero body forces. 9 2.3 DISPLACEMENT FUNCTION FORMULATION OF PLANE STRESS PROBLEMS The substitution of the strain-displacement relations xx 8x =22 yy 3y Bu 3v 28xy 3y 3x into Hooke's Law for the state of anisotropic plane stress 1 ' B 5 (131 = 19 2: 6) 12 By au 8v Bu 3v °xy 616 8x + 826 By + 866 l + > Using the equations of equilibrium 80 80 3x By 30 30 __£2.+.__22.. 0 8x 8y 2.11 2.12 2.13 2.14 in conjunction with the stresses given in 2.13, results in Navier's equ- ations for the state of anisotropic plane stress, 10 32 32 32 32 32 32 811 2'*2816 axay I 866 2 “‘I 816 2'*(812 B66) axay'*826 2 V=‘) 3x By L 3x 3y 32 32 32— 32 32 32 B16 2'I(812'I866) axay'I826 2‘ “'I 866 2'I2826 axay'I822 2 V“) 8x 3y 3x 8y u 2.15 If the displacements, u and v, are defined in terms of an unknown func- tion, w, such that 32 32 32 u I 816 2 I (812 I 866) axay I 826 2 w 3x 3y 2,16 < I 32 32 a2 ' ‘ B11 2 I 2816 axay I 866 2 w 8x 3y where w = W(x,y), then the first equation of 2.15 is identically satis- fied and the second equation becomes 4 4 2 a p _ 3 (866811 ' B16) 4 I 2(811866 816812) 3 3x 3x 3y + (B 8 + 28 B - 28 B - 82 ) ‘ififlL— ll 22 16 26 12 66 12 2 2 3x 3y 4 4 a 2 a w _ I 2(822816 ‘ 826812) 3 I (866822 I 826) . 4 ’ 0 3x3y 8y 2.17 This equation can be written with the elastic compliances as, 4 4 4 4 4 ill- _3_L __3_i’__ 3‘4 ii): C22 4 ZC26 3 I (C66 I 2C12) 2 2 2C16 3 I C11 4 0 8x 8x 8y 3x 3y axay 8y 11 It may be noted that Eqs. 2.18 and 2.8 are of the same form. The stress components can be obtained from 2.13 upon substituting for the displace- mnt'o a -(a a 43133 HS 6 B)—-L83 +(816 -e 31—9133 xx 11 66 16 2 11 26 8816 12 26 12 66 3 3x 3y Bxayz 3y 3 2 a3 H<812 16 '8115 26) a 3 I (812666 I512 ' 822811826816) 2 x 8x 8y + 2(8 -B e > 33 '+(82 -e B ) 313 12 B26 16 22 a 3 26 22 66 3 x3y 8y 6 -(62-ee)—‘P-3+(ee- e)—§-—-“’+(88 -66)33 xy 16 11 66 a 3 12 16 811826 2 12 866 16 26 2 x 8x 3y 8x3y 2.19 In order to compare the above formulation with Marguerre's displacement function approach, consider the case of isotropy 1 _ c11 ' C22 I E’ C16 C26 ' 0 2.20 v 2 l + v C12 I ' E C66 I E Then the displacements are written as, u . E 326 2(1 - v) axay 2.21 v _ _ E 32w __ E 82w 1 _ v2 ax2 2(1 + v) 8y2 With the definition . - ___EML___ ,4 2(1 + v) 2'22 12 the displacements, = _ l + v Bag? u l - v 3x3y 2.23 2 32% 32.11 V=l~v 2I 2 3x 3y are obtained, where .a1' is the displacement function which was proposed by Marguerre (ref. 16). Of course, in the case of isotropy, the fourth order partial differential equation with constant coefficients (Eq= 2.18) will simplify to the biharmonic form. In conclusion it can be stated that the general solution of the dif- ferential equation 2.18 subject to the apprOpriate boundary conditions would constitute the solution of the anisotropic plane stress problem with zero body forces. The displacement function method was presented only as an interest— ing point to note in plane anisotropic problem formulation paralleling that of Marguerre's proposed method in isotropic elasticity. 3. GENERAL CONSIDERATIONS FOR PROBLEM SOLUTIONS CONCERNING MULTIPLY-CONNECTED REGIONS 3.1 INTRODUCTION The problem of multiply connected regions in plane linear aniso- tropic elasticity is discussed in both Lekhnitskii's and Savin's works (refs. 4, 6). In the following, a brief treatment of the subject matter is given for reasons of completeness. First, the governing differential equation is solved in terms of two arbitrary complex functions (stress potentials) and expressions for the stress components and displacements are derived in terms of these two functions. The boundary conditions are established for both the First and the Second Fundamental Problems. Lekhnitskii's complex material parameters are examined and the resultant force and moment for a segment of a curve are also found. Finally, the forms of the stress potentials for a finite multiply connected region are defined. For a more comprehensive discussion, the reader is referred to the references given above. 3-2 SOLUTION OF THE DIFFERENTIAL EQUATION The fourth order partial differential equation with constant coefficients obtained from the stress function formulation .4 4 4 4 4 ii. .242. _§_L_ LL Li: C22 4 2C26 3 I (C66 I 2012) 2 2 2C16 3 I C11 4 0 8x 3x 3y 8x 3y Bxay 3y 3. l3 1 14 has been extensively studied by Lekhnitskii, Savin, Green and Milne- Thomson (refs. 4, 6, 8, 9). Following Lekhnitskii's method of solu— tion, it is observed that the above differential equation can be written in symbolic form as, D1D2D3D4¢ = 0 3.2 where 8 3 Dk — 3§-- “k 3;- (k - 1. 2. 3. 4) and uk are the roots of the characteristic equation given as, C - 2C u + (C + 2C )u2 - 2C u3 + C u4 = O 3 3 22 26 66 12 16 ll ' Integrating Eq. 3.2, the stress function, o, can be written as the sum of four arbitrary functions, 4 ¢ = E :Fk(x + uky) 3.4 k=l Considering the positive definiteness of strain energy, Lekhnitskii proved that all the roots, uk, must necessarily be complex. Since the roots, uk, are all complex and the coefficients of the characterisric equation defining the roots, uk, are all real; the roots, “k’ must be complex conjugates of each other. So, the roots can be designated as, pl, “2’ hi and B}. These roots, called complex parameters, depend entirely on the material constants. For real stresses, the stress function, ¢. must necessarily be real. Therefore, the general solution of Eq. 3.2 must be of the form 6 = 2Re[F1(zl) + F2(22)] 3.5 15 where 21 = x + uly z = x + and “I # “2 2 “2y If the roots are not distinct, i.e., “1 I L‘2’ “1 I u2 then the stress function becomes ¢ = 2Re[Fl(zl) + le2(zl)] 3.6 In particular, consider an isotropic material. The characteristic equa- tion becomes .1 2(1 + v) _ 2' 2 l_ 4 _ E + [}——7;——- 2 g] u + E u — O 3.7 Solving for the roots, one obtains U1 = U2 = 1 “1 = “2 = -i Hence 21 = 22 = x + 1y; 21 = 22 = x - iy and the well known expression 6 = 2Re[EF2(Z) + Fl(z)] 3-8 is obtained. 3.3 STRESS COMPONENTS It was shown in the preceding section that the stress function, ¢, can be expressed in terms of two arbitrary functions of the complex var- iables 21 and z2 when ul # “2’ in terms of 21 when ul - uz and in terms of 2 when isotropy is considered. Knowing the stress 16 function, one can readily find the stress components in terms of the above arbitrary functions by simple differentiation. Defining dFl dz1 I ¢1(21) sz 3.9 622 I ¢2(22) the stress components will be given by the following equations: 1. The complex parameters are distinct (ul # uz): 2 ' 2 ' xx I 2Re[%1¢1(z1) I “2¢2(22E] ' ' yy - 2R2 E1(Zl) + (112(22fl 3.10 °xy I ZRE[%1¢1(21) I “2¢2(22{} 2. The complex parameters are equal (u1 = uz): O 0' _ r'2 ' _. ._ 2 ' °xx I 233Lf1¢1(z1) I 2“‘1“1"’21¢1(zl)+Plzl¢2(zl)+ [33110111111 -ul) +012 IUICIJFZQI} + u - w o ~2Re (>‘(1E-23 V" qi¢1 z1 'Iqizi¢2 21. I pl + v + m 0 CV fi' . 1 _ - 171-) 4-C121J1 -C26] F2(zl} X 0 3.15 18 3. Isotr0pic case (111 = u2 = 1, ul = “2 = — i): 2Re{-1;v [411(2) + 2152(2)] + g—E-l F202} + no - woy 2Re{- i L%—Y- [¢1(z) + E¢2(z)] - i _3E_V F2(z)} + vO + (sex The fundamental displacement combination is written as, u + iv = 2 {- l—E-y- [Tb-1(a + 232(2):] + 11%;) F2(z)} + u0 - woy + i(vO + wox) 3.17 C II <2 II 3.5 BOUNDARY CONDITIONS 3.51 First Fundamental Problem Suppose that the boundary stresses, Xu and Yn, are given (Fig- ure 1). Then on the boundary, F, of the area, L, the following equa— tions must be satisfied: Xn = oxxcos(n,x) + Oxycos(n,y) 3.18 Yn = 0x cos(n,x) + oyycos(n,y) Taking into consideration that Q: dv cos(n,x) = 5:1 cos(n,y) = - ds and substituting for the stresses, the boundary conditions on the stress function can be written as, x—As Y--sa 3w 19 Upon integration, these two expressions become .. .33: S - ‘//) Ynds + Cl 8x 8 where C1 and C2 are arbitrary real constants of integration and s 3.20 L2: is measured from an arbitrarily chosen point on F. Since the stress function, o, can be expressed in terms of two arbitrary complex functions, Fl(z1) and F2(22), and consequently using the definitions, 3.9; the derivatives, 3%. and 333 can be written in terms of the complex functions, 61(zl) and ¢2(22). Hence, the equa— tions which must be satisfied by the functions, ¢l(zl) and ¢2(z2), on the boundary, P, of the area, L, are as follows: 41621) + 42oz) + Elczl) + $26,) = - / Yndx + c1 5 ul¢>l(zl) + u2¢2(22) + 316161) + "5262022) = f Xnds + 02 z in I‘ s 3.21 3.52 Second Fundamental Problem Suppose that the displacement components u(x,y) and' v(x,y) are given in the boundary, P, of the area, L. Then the boundary conditions will be given by Pl¢l(zl) + 92¢2(221)4+ 313151) + 32%(32) - woy + no = 110(8) qul(zl) + q2¢2(22) +‘Ei$i(§i) + EEEéCEé) + wox + v0 = vo(s) z in P 3.22 20 where u°(s) and v°(s) are specified values of the components of dis- placement on the boundary, P, which are functions of the arcs of the contour from an arbitrarily chosen reference point. Considering both the First and the Second Fundamental Problems, it is observed that both problems can be solved completely in terms of the two arbitrary complex functions, ¢1(zl) and ¢2(22). 3.6 ON THE COMPLEX PARAMETERS OF LEKHNITSKII Lekhnitskii has shown (ref. 4) that the complex parameters, pk, in a rotated coordinate system are related to the complex parameters, pg, in the original system by the equations u: cos 6 - sin 6 u = (k = 1, 2) 3.23 cos 6 + “a sin 6 where 6 is the angle between the x-axes of the original and the rota— ted systems. These transformation equations are of particular impor- tance when problems of orthotropic media are considered. For the case when the material and the reference axes coincide (special orthotropy) Eq. 3.3 reduces to 4 2 C22 + (C66 + 2C12)u + C11“ - 0 3.24 For the case when the material and reference axes are not aligned along the same directions (general orthotropy), the are found from o “k 3.24 with respect to the material axes first and then in accordance with Eq. 3.23, the values of “k can be determined with respect to the reference system. 21 Upon solving Eq. 3.24, one obtains the roots 1/2 2 1/2 u6 = 1 C66 I 2(312 + (C66 I 2C12> _ C22 1 2C11 2C11 Cll 1/2 2 1/2 76 = _ 1 C66 I 2C12 + (C66 I 2C12) _ C22 1 2011 2C11 C11 .2 C + 2C C + 2C 2 C 1/2 1/2 3 5 u6=1 66 12_<66 12)_22 2 2011 2Cll Cll 1/2 _ 2 1/2 T16 = _ 1 C66 I 2C12 _ (C66 I 2C312) _ C22 2 2Cll 2Cll Cll Since the complex parameters, pi, are of the form 0 o 0 uk — ok + 18k (k — l, 2) 3.26 one can immediately conclude that for the condition C + 2C 2 C 66 12 22 20 _. ——C _>_ O 3.27 11 11 the complex parameters, pi, become purely imaginary, i.e., o = o = ' uk in (k l, 2) 3.28 where o o 81 > O, 82 > 0 It also follows from Eq. 3.25 that o o 81 2.82 3.29 Note that condition 3.27 is equivalent to 22 E /E 22 < 1 3.30 2— (E11 _v) 2612 12 in terms of engineering constants. The rest of the range of material orthotropy is covered by the condition C + 2C 2 C 66 12 22 2C - Er-'< 0 3.31 11 11 resulting in ai = -dg and Bi = a; E-O. In terms of engineering con— stants, this condition is equivalent to E /E 11 22 > 1 3,32 2 (_IL1__,) 2G12 12 These dimensionless ratios are convenient to use and cover the whole range of material properties. The positive definiteness of B: guarantees that B in Eq. 3.23 k are positive definite. As a matter of fact, one can always construct such that B > 0, even in the case of complete anisotropy (ref. 6). “k k This fact is of importance if the following affine transformation is considered: xk = x + aky (k = l, 2) 3.33 yk = Sky where 2k is now defined as, zk = xk + iyk (k = 1, 2) 3.34 23 The Jacobians of the transformation become k Jk = = Bk > O (k = l, 2) 3.35 O Bk Since 81 > O and 82 > 0, the description of any two transformed con— tours, F1 and F2, in L1 and L2, respecrively, correspond to T in L such that when 2 describes I in L, zl describes II in L1 and 22 describes T2 in L2, all in the same sense (ref. 9, Figure 2). 3.7 RESULTANT FORCE AND RESULTANT MOMENT Consider the curve AB (Figure I) loaded by.the.complex stress, X.n + iYn. Suppose that one wants to find the resultant force and the resultant moment of the forces acting on the curve AB. The force vector acting on an element ds of the curve AB can be written as, (X +1Y)ds=d—a$_iai=_id_a_?i+i§$ 3.36 n n 3y 3x 3X 8y Hence, the total force acting on the curve AB becomes B B x+1r= (X+1Y)ds=-‘131+139i 3.37 n n 3x By 'A A where [...]i denotes the increase in the value of the bracketed ex- pression when moving along the curve AB from A to B. Upon substitution for 211,21 3x By the final form of the force acting on the curve AB is obtained 24 X + iY = - 1&1 + iu1)¢l(zl) + (l + iu2)¢2(zz) B + (l + iul)¢l(zl) + (1 + iu2)¢2(§2)}A 3.38 The moment, Mo’ with respect to the origin acting on the curve AB is given by M = (xYn - an)ds 3.39 A After substituting for Xnds and Ynds the moment, MO, becomes B M = - [xd(%§) + yd(%$—) 3.40 A Upon integration by parts, one obtains Mo=..[:x_gix+ygy]:+[6133'“ Then the moment of the forces acting on the curve AB can be written in terms of the complex stress potentials as, B Mo = 2Re[F1(zl) + F2(22flA - Reg [(1 - iu1)¢l(zl) + (1 - iu2)¢2(zz) B + (1 - iTIlWlGl) + (1 - 1F2)$2(?2fl}A 3.42 3.8 FORMS OF THE STRESS POTENTIALS FOR FINITE MULTIPLY CONNECTED REGIONS Savin has shown (ref. 6) that the functions, ¢l(zl) and ¢2(22) for a multiply connected region are expressed in general as, 25 L ()— Ar( +* I1 21 I 2 n 21 I 21,2) ¢1(zl) =1 L 3.43 - B 2 + * ¢2(22) - 2 u(z2 - 22,,) ¢2(22) i=1 where .21 i and 22 2 are arbitrarily chosen points inside the bound— ! 9 * * aries, Pl,£ and F2’£ (Figure 2); ¢l(zl) and ¢2(zz) are holomorphic functions in L1 and L2, respectively and A2 and Bi are complex constants which have to be determined from some known condition. In general, the ith boundary, F2’ is loaded by the force vector, Xg + iY£, as given in the preceding section. Due to force equilibrium, any closed curve C in L reconciliable with I will be subject to 2 2 a load given by - (X2, + 1Y1) = - i{(l + iul)¢1(zl) + (l + iu2)¢2(22) B=A + (1 + filwlfil) + (1 + 1U2)'6'2(72)} 3.44 A Considering the forms of ¢1(zl) and ¢2(22)(Eqs.3.43), it is noted * * that the holomorphic functions ¢l(zl) and ¢2(zz) do not increase in value over the complete circuit, C Therefore, only the logarithmic 2,. terms have to be investigated in accordance with Eq. 3.44. If one carries out this investigation, two equations will result which are complex conjugates of each other. These equations are given as, X2 + IYQ (1 + 1111M, - (1 + 1111M, + (1 + 1112)}31 — (1 + iu2)B£ a .. 2“ X9’ - 1Y2 2n (ref. 6) - (1 .. 161ml.» (1 - film, — (1 — iu2)B9v+ (1 - 1112», = _ 3.45 26 Obviously, Eqs. 3-45 are not sufficient in themselves to determine the values of the complex constants, Ag and 32' The other two equa- tions, necessary for the determination of Al and B2’ are based on the consideration of single valued displacements. Upon writing the fundamental displacement combination u + iv = Pl¢l(zl) + pl¢l(zl) + p2¢2(22) + p2¢2(22) + 1[ql¢l + qlrlol) + q2¢2<22> + q2$2(22)] 346 one requires that the condition B=A [p + iv] A II C be satisfied. The satisfaction of this condition will also result in two complex conjugate equations (Pl + iql)A£ - (31 + iql)A£ + (p2+ iq2)B£ - (p2 + 1q2)B£ = 0 - (pl -1q1)A,+ ('51 431%, - (p2 -iq,)B,+ (13’, 43213,: 0 (ref. 6) 3.47 For a solution of the system of equations consisting of Eqs. 3.45 and 3.47, the determinant of the coefficient matrix cannot vanish. In order to prove that this condition is satisfied, the determinant, D, of the coefficient matrix is obtained I16B182°11022 2 2 2 D = 2(a - a ) (B + B ) 2 2 2 2 {: 2 l 2 1 (ol + Bl)(a2 + 82) 2 2 4 + (82 - 81) (82 + 81) + (0.2 - 61)} 3.48 27 This determinant is the corrected form of that given by Savin which is apparently in error (ref. 6). Obviously, for B > 0, the determinant, k D, does not vanish. Hence, the system of equations constructed from A B Eqs. 3.45 and 3.47 can be solved for the complex constants, A 2’ 2 £9 and B2. There are two cases possible: 1. The resultant force vanishes on each boundary, i.e., X£+ 1Y2 = X2 - iY2 = 0 2. The resultant force does not vanish on each individual boundary. In the first case, the complex constants, A A B and Bk, must be 2’ 2’ 2 taken as zero since D # 0. Hence, the stress potentials, ¢l(zl) and ¢2(22), are holomorphic functions in their respective regions, L1 and £, BR and Bi, L . In the second case, the complex constants, A A 2 Q. will have values other than zero. .Hence, the complex stress potentials, ¢l(zl) and ¢2(z2), will become multiple-valued due to the presence of the logarithmic terms. 28 Figure 1. - Curve segment loaded by the complex stress, Xn + W”. 29 z-plane 21-plane Figure 2. - The multiply connected region, L, in the z-, 21-, and 22-planes. 4. AN APPROACH TO THE SOLUTION OF FINITE-GEOMETRY TRACTION-FREE CRACK PROBLEMS 4.1 INTRODUCTION It was established in the preceding section that two complex holo- morphic functions define the solution for problems for which the result- ant force vanishes on each individual boundary. Since the solution is defined in terms of two arbitrary complex holomorphic functions, it still remains to find the forms of these functions such that the appro- priate boundary conditions are satisfied on the boundaries of the region. In this section, the problem considered is that of a central, traction-free crack in a finite anisotropic region. The method of solu— tion proposed herein was called by Bowie "A.Modified Mapping-Collocation Technique" (ref. 2) and was successfully employed by him for isotropic problems. As a first step, the crack is mapped into the unit circle. Then the zero traction conditions on the crack are exactly satisfied by using Muskhelishvili's continuation concept (ref. 3) across the unit circle such that the form of either of the two functions is defined in the form of the other. These expressions, will for the first time, be derived in this section. Upon satisfying the zero traction conditions on the crack, a Laurent series form of representation is assumed for the remaining 30 31 arbitrary function. The boundary conditions for the outer boundary are then derived in terms of this Laurent series. In order to satisfy the conditions on the outer boundary, the least squares boundary collocation method is proposed. Bowie (ref. 2) also suggests the possible use of two stress and the moment conditions on the boundary in addition to the force boundary conditions for problems with local stress irregularities. In this dissertation, the force boundary conditions were found to be sufficient for the solution of the problems considered. In addition to the method of solution described above, the expres- sions for the stress components are derived and specialized to the neighborhood of the crack tip. It is then shown that the stress inten- sity factors for a finite plate for both the opening and the shear modes can be defined analogously to that of the infinite plate. Finally, expressions for the displacements are derived and in par- ticular, the displacements of the crack boundary are also given. 4.2 THE MAPPING OF A CRACK AND ITS EXTERIOR INTO THE UNIT CIRCLE AND ITS EXTERIOR Let the unit circle, y, in a c—plane be mapped into the crack boundary, E:,in the z-plane by w(o) 4.1 N II where i0 0 = e The unit circle, y, in the C-plane will then correspond to the crack, PC, in the z-plane with the same sense of description (Figure 3). 32 Since the transformation 2k = x + uky = xk + iyk (k = l, 2) 4.2 can be expressed in terms of the conventional complex variables as, l - iuk l + ink zk=--2——z+——2—-z (k=l,2) 4.3 it follows that z gflwfiy) +3.:im(l) : w (0) (k = l 2) 4 4 k 2 2 o I k ’ ‘ maps the unit circle, y, in the t-plane into the crack boundary, Fck’ in the zk-planes. Consequently 2k = wk(c) 4.5 will map the region, LC, determined by the unit circle into the regions, Lk’ in the 2 —planes. k The mapping of the unit circle, y, into the crack boundary, PC, is accomplished by the function z =% (o +—) (ref. 9) 4.6 where a is the half—length of the crack (Figure 3). Hence, the unit circle, y, is mapped into the crack boundary, Ick’ in the zk-planes by _a .1.) _ zk-2(o+0 (k-1,2) 4.7 Upon explicitly expressing the variable on the unit circle, one obtains o = (k = l, 2) 4.8 33 Therefore, the exterior region, LC, to the unit circle, y, is given by the mapping C = ' (k = l, 2) 4.9 Note that the positive sign preceding the radical guarantees the "exterior to exterior” transformation. 4.3 ON THE CONTINUATION AND FORM OF THE STRESS POTENTIALS In the following Muskhelishvili's extension concept (ref. 3) will be applied to the complex stress potentials, ¢l(zl) and ¢2(22), in the C—plane. Consider that a traction-free crack, PC, in the z-plane is mapped into the unit circle, y, in the c-plane by the functions C = , z in F 4.10 and C = , , zk 1n Fck (k = l, 2) 4.11 (Figure 4). Then on the unit circle, the following boundary conditions will prevail: 41(4) + ¢2(c) + ¢l(c) + ¢2(4) = o, r in y 4.12 ul¢l(c) + u2¢2(l;) + Fl$1® + U2$2® = 0. c in v where ¢k1 s ¢k(c) (k = 1. 2) 34 The functions, ¢1(C) and ¢2(C), are holomorphic in LC and are defined in LC. In order to define ¢1(C) and ¢2(t) in RC (inside the unit circle), substitute I/C for 2' into Eq. 4.12 and consider t in R . C Then the following expressions will result: 9 (C) + 9 (C) = - 3' I‘ - $l l‘ C in y 1. 2 1 C 2 C ’ 4.13 u 4 (C) + u ¢ (C) = - I'I' 1' - E'E' I) C in Y 1 l 2 2 l l t 2 2 t ’ Hence U ‘11 IT.“ - 1 2 1 2 l - 1 ¢ (‘) = ‘—‘:—=F’¢ (C)'+'—‘f:‘=’ ¢ (—), C in R 1 C “1 “1 2 “1 “1 2 C C 4.14 11 '11 L1 '11 -— l l 2 l 2 - 1 ¢(-)' _— ¢(C)+ _— ¢(—)9 C in R 2 C ”2 “2 l “2 “2 1 C C The substitution C - l/E’ will change these definitions to u - u 3' - u - 2 1 l 2 l-— $(C)--—-—_e ¢(=)+__— ¢(C). C in L l “1 "1 2 C “1 "1 2 C __ 4.15 n - u u - u -' - l 2 1 l 2-— -— .6. __ (z)._.__.....), ... L 2 ”2 ”2 1 g “2 “2 1 C Upon taking the complex conjugate expressions, one obtains fi' - i' 3' - u l 2 1 l 2 ¢(C)'——_'=—$(")+——:———¢(C). C in L 1 “1 "1 2 C ul “1 2 C 4.16 H' - 3' 3' - u 2 l —- l 2 1 9 (C) ' “‘1f1=' ¢ (') i"“1:1=‘ ¢ (C). C in L 35 Equations 4.16 show that either ¢1(t) or ¢2(C) can be expressed in terms of the other. The zero traction condition on the crack provided the extended definitions of either one of the two functions in terms of the other such that each is continuous across the unit circle in its C-plane. The inverse transformation, back to the z —p1anes, results in k _ _ 2 2 u - z - z '+ z -a 1 2— 1 1 1 . ¢l(zl) __.1¢2 a , 21 in L1 II - z - 22-a II - z -+ 22 -a2 ¢(z)_“2 “1,. 2 2 ,“2 “1¢ 2 2 z in L 2 2 u E 1 a In -H 1 a 2 2 2 2 2 4.17 Now, 61(21) and ¢2(zz) are defined in L1 and L2 such that using either one or the other in solving a finite geometry problem, will satisfy the zero traction condition on the crack. Furthermore, if ¢1(zl) or 62(22) is transformed to the t-plane, each can be consid- ered holomorphic within its doubly connected region enclosed by r and i' where r' is obtained by inversion of r with respect to the unit circle (Figure 5). The determination of either ¢l(§) or ¢2(C) depends on a form of representation and the satisfaction of boundary conditions on r cor- responding to I1 or P2, respectively. It will be assumed that ol(c) or ¢2(C) can be represented in the form of a Laurent series. This assumption appears to be reasonable for a certain class of problems although the boundaries T' and r are not circular (ref. 2). As 36 Bowie stated in reference 2: "There is no a priori reason to suspect (for a certain class of problems) that the region of convergence of such a series could not extend over the desired parameter range." In general, assume the form for ¢1(C) as n -n ¢1(C) - A10 + > :anc + BlnC 3 4.18 n=l considering that ¢2(22) will be found from the second equation of 4.17. Hence 2 2 n 2 2 n 21 + 21 - a z1 - zl — a I191) I A10 I A1n a I Bln a n: 4.19 z + Q/zz - a2 and ¢ 2 2 become 1 a n __ 22- 22-a _ K 22-3/221-8 a 10 In a n z + z - a +2 2 I 2 4.20 In a n= 37 n z + z2 - a2 z + 22 - a2 ¢ 2 2 = A + A 2 v 2 l a 10 In a n=l n z - 22 - a2 + B 2 I 2 4.21 In a Therefore, 62(22) is obtained as, 3' - 3' 3' - u 2 l —- 2 l ¢(z)=—-_—:-A +—T-:-A 2 2 “2 “2 10 “2 ”2 10 I _. _. _. + 2 2 “ U'IJ “'11 Z 2‘3 + 2 _ _1 A1 + 2 - _1 I} 2 “2 “2 n “2 “2 n 3 n=1 TJI-ILJI 17-11 2- + 2 _ _} -1 + 2 _ _} Bl 2 4.22 u2 2 n 112 “2 “ It is obvious that this expression is also a Laurent series and it is of the form so n — ¢2(C) - A20 + E :[AZnt + BZHC I] 4.23 n=l At this point, it is noted that upon assuming both ¢1(C) and ¢2(C) in the form of Laurent series' and satisfying the zero traction conditions on the crack, one obtains the very same form of ¢2(z2) as given by 4.22. This result is obtained by the elimination of the con- stants A20, A2n and an which is p0881ble upon the satisfaction of the zero traction conditions on the crack. 38 4.4 BOUNDARY CONDITIONS FOR FINITE REGIONS CONTAINING A TRACTION-FREE CRACK Considering the First Fundamental Problem, the boundary conditions 3.21 must be satisfied by the functions ¢1(zl) and 62(22) on the outer boundary of the region. Upon making the proper substitutions (using the results of the preceding section), the boundary conditions 3.21 are written as, 39 and 2 2 a a _ _ 2 2 z z a. ,a J a + _ 2 2 C G ) )q‘l l U 2 _U 2 U U _ . 2 - 2 _ U 2 U 2 /\ _U I\ _U 2 2 _U _U + + H Al/ .25 4 40 where Ci and c; are real arbitrary constants of integration cor- responding to the outer boundary of the doubly connected region. It may be noted, here, that the complex constants, A10 and 210’ remain undetermined. This fact, however, is of no importance since the stresses are not influenced by these constants and the rigid body dis— placements have already been expressed by other constants (see Eq. 3.14). Upon examination of the boundary conditions, 4.24 and 4.25, one finds that the n-l terms can be contracted due to the fact that and ‘ 4.26 111012 -F2)21+u2(fiz -u1)z2+'172(u1 -u2)?2 =- v(uz - 111M112 - u2)(u2 - “1) Hence, the boundary conditions are written as A " B 2 11 11 ., _. 2 2 2 2 _2 2 232 3012-172) [(112—112) zl—a +(1Iz—u1) zz-a +(u2-p1) V22_a] \ 4‘1 ' 2 2 I 2 2 I I n — I - - - + 2Re/‘ .- . 2322(1) -_ ) [(112 F2) (21 + 21 a ) +52 111) (22+ 22 a -- -4 2 H2 n-2 I +(u1 I“2)'22I (222 42)? z in P 4.27 41 and + B N 11 11 - '- 2Re{:(u2 _ F2) “#2 " H1)(l12 ‘ H2)(l12 ‘ #1)} — 3 A11 11 _ _ f 2 _ 2 _ _ 2 _ 2 + F2012 - "1) 22 - 82]} I \ n +2Re) A1“ u(u F)z+ zZ-a2 atl? ‘x. \ B I n + 2Re '. 1n ( _ —-) z _ z2 _ a2\ / n( __) “1 “2 “2 1 1 / / , a 112 112 n-2 1: H v G N l a N N N l m N v :3 I "2672 I + .. /_ + 2 n 112(111 - 112) \22 22 - a = fxnds+c‘2’, z in r 4.28 8 At this point, it is convenient to introduce the following definitions: ll 11 In Zn In 2n 42 A11 I B11 gs - 3(112 - F2) g = A11 I B d a(u2 - “2) g, . A1. 1n n( _-> = B1n (an n( _—> a ”2 ”2 <4 = (112 - 111)(u2 - 32 )(F, - ul) (uz-TIZM’Zi-azi-(IM-u )h/zg: -2ula+(u2- _ 2 2 _ ul(u2 - “2)4/21 - a + u2(u2- 2 n (112 - 5,)(21 + 21 - 82) + (3, MA, - 82 + u2(u2- 4.29 F :MFI—az + (u - 2)GC 2% n n (“2 I “2) (21 I 22 I 32) I (“2 “1)<22 I 22 I 82) + (u u )(2' + 32 1 2 2 2 n n u1(u2 u2)(l 21 - 82) + u2(u2 “1)<22 + z: - 82) 43 Now, observing that F e? + F G’ = 2[ReF 11 d 11 d 1R8%3 I ImF 11m%%] etc., 1 l the boundary conditions become ReFllReQE — ImFllImEE + E [ReFlnRegan - ImF n=2 lnlmffln] + E :[ReFZHRefian - ImFZnImegn] n=2 _ l [E [/3 Y ds + C5], 2 in F 4.30 — - n l 2 s yReQEReQ; - yImfiZImeg + ReGllRe%Z - ImGlllmQE + E : [ReGlnRe%in — ImGlnImein] +- E :[RernRegan - ImGZHImegn] =2 n=2 1 o _ E-[’L/1 Xnds + C%], 2 in P 4.31 s 4.5 METHOD OF LEAST SQUARES BOUNDARY COLLOCATION' Up to this point, it was shown that the solution of the governing differential equation, 3.1, can be constructed in terms of the stress potentials, ¢l(zl) and ¢2(zé);and in turn, ¢2(zz) can be defined in the form of ¢1(zl) such that the zero traction conditions on the crack boundary are exactly satisfied. It still remains to satisfy the 44 boundary conditions, 4.30 and 4.31, on the outer boundary of the region in order to complete the solution of the problem of a finite anisotropic region with a central crack. Perhaps the best method applicable to this end is the least squares boundary collocation method as given in references 17 and 18. The method is described briefly in the following: 1. 4. The boundary conditions are satisfied at discrete points on the boundary in the least squares sense. One takes an overdetermined system of equations, i.e., the number of equations must exceed the number of unknowns of the system. Then in indicial notation, one has the following system of equations: égqu = Ki 4.32 where K1 depend on the loads acting on the boundary, 431 are dependent on material properties and the boundary geometry and 43 are the unknowns of the system. Upon assuming an approximation to the unknowns of the system, one has, say, 435?; - Ki = R1 4.33 where R1 are the error terms due to the approximation, ‘6; to Then the square of the error, R1R1 I (61173 I K1)(駣?; I K1) 4'34 45 is minimized by taking d(R R ) _l_1_ = o 4.35 as" 2 ~which in turn results in * . éfléalkgk I 5.121(1 “'36 Or in matrix notation, one has 635%" = §TK 4.37 The applicability of the least squares boundary collocation method to the problem at hand is readily obvious. Considering the truncation of the infinite series' in the boundary conditions at NP correspond- ing to the series' with the unknown constants, flan, and at NN for the series' with g2n’ one has at a point on the outer boundary two equa- tions in 2(NN + NP + 2) unknowns on the left hand side and two arbi- trary integration constants on the right hand side. In general, these integration constants can be handled as unknowns in the collocation method; however, for certain orthotropic.cases they can readily be de- termined. Then at a point on the boundary for the anisotropic case, one has the folloWing two equations, r- 1 '- — -1- " o o ReF -ImF c — i Y dJ 2 ll. . . ZNN 2 n A S ‘ 1 ? 2 yReCc -yImCC ReGll...-ImG2NN C O HO O NO Ree: l—\/p X ds 2 n ImQ; L s - Reg’ IngNN L 4.38 46 Considering M points on the boundary, the coefficient matrix 6’ will be ofw 2M rows and 2(NN +INP + 3) columns subject to the condition that 2M > 2(NN +tNP + 3). 4.6 STRESS COMPONENTS Using expressions 3.10 in conjunction with Eqs. 4.19 and 4.22, the stress components are obtained. These expressions are somewhat modified by noting the contraction of the n=l terms and using some of the defi- nitions given in 4.29: oxx = Heiress] 2 _. z1 2 22 + 2Re % 111012 - 112) + 11252 - “1) z2 _ a2 22 _ 2 1 2 a —2 E2 + u2(u2 - “1) .2 _ 2 22 a NP 112012 - 712) 2 2 n + 2Re nfiin 1 (21 + 21 - a ) 22 - 82 n=2 1 +u2(u:- +142012 - ul) __ . _2 2 n 32 2 22 - 22 ‘ a IVE:- a NN A u2(u - 3') n _ 2Re ‘ nVZn 1' "2 2 (21 - 2i - a2) 22 - a2 - n=2 1 n n + u26:2 “1) 2 2 “2(“2 “1) 2 Z - z - a + + z - a z - a -— - 2 22 a +2Re - 2Re NP n=2 n=2 2 2 21 a (“2 - Hz) 2 Z -a (n2 “2) 2 Z ’8 47 2 - “1) 2 Z - a 2' 2 -ul) _2 2 Z ‘ a n z2_2 l a (112 ul) ( 22+ 2 2_2 z2 a (U2 - U1) (E _2 22 ‘ Z 2_2 22 a n 22_2) l a (uz-u1)(z- z 22 - 32 2 2 (U2 - “1) (_ _2 22+ 2 Z - 82 2 48 21 22 o = - 2R2 g? 11(u — E') + u (3' - u ) xy (1 1 2 2 2- 2 2 2 l 22- 2 ’21 a 2 a _2 _ z2 + u2(u2 - “1) 2 22 - a NP u (u - H") n - 2Re ng’ l 2 2 z + 22 - a2) 1n 1 1 z2 2 n=2 1 ”2(32 “1) 2 2)n + z + z - a 2 2 (2 2 z2 - a + 6%2 - “’22 - a 22-112 2 NN . - n u1(u2 “2) 2 2 + 2Re 2n 21 - zl — a 22 _ 2 n=2 1 4.41 49 4.7 CRACK TIP STRESS COMPONENTS AND STRESS INTENSITY FACTORS The stress components were expressed, in general, for any traction— free central crack problem in series form as given above. One can Specialize these expressions for the neighborhood of the crack tip by considering the transformation, 2 = a + r(cos 8 + i sin 8) (Fig. 6) 4.42 The substitution of 4.42 into 4.3 will result in z = a + r(cos 6 + u sin 6) (k = l, 2) 4.43 k k Subjecting the expressions for the stress components to these transform— ations and observing the condition for the neighborhood of the crack tip, (i.e., r/a << 1), the following expressions are obtained: U U U U OX = 4 Re [(1.12 - Ill-)2] Re{ 1. 2 [ 2 _ l 1 x V2ra , U1 U2 /cos 0-+u sin 6 /cos 6-+ulsin 6 2 2 2 4 - 1 [' u2 L11 .1 {a Re[:uz(u2 -ul)E] Re 11 - -uL q.' 2 l 2 /cos 6-+uzsin 6 Vhos 6-+t131n 9] 4.44 U 11 o = 4 Re[(u2 -ul)§:/] Re{ i l _ 2 yy /2ra H1 H2 /cos 8-+uzsin 6 /cos 9'+p151n 6 4 Rafi-2012 -ul)Z:] Re u Eu 1 __ 1 V2ra 1 2 Vhos 6-+u sin 8 /cos 9 +1 sin 6 2 l 4.45 50 [1.249314% 1 - 1 xy “3:2“ l '-“2 Vhos 6-+ulsin 6 /cos 6-+uzsin 6 _4__ “1 “2 - 32(112-111)2Re -u - V2raRe 2 /COs 6-+ulsin 6 Vhos 6-+uzsin 6 4.46 where K \ NP L]: a? + Enangln -Enan‘€ =2 n=2 :3 The forms of the crack tip stress components are the same as that of given by Paris and Sih (ref. 12) for an infinite anisotrOpic region with a central crack. Hence, the stress intensity factors for a finite anisotropic region with a central crack are automatically defined as 4- K1 “IE Re[(u2 - 1192] 4.47 and K2 = - énefi'zwz - ”92] 4.48 where K1 and K2 are the opening mode and sliding mode stress inten- sity factors, respectively. A further examination of the stress inten- sity factors reveals that 2 K+uK=--—(u -u) m¢(c) 4-49 2 2 1 J; l 2 CI’l l which form was first published in reference 12. 51 vIn addition, one may note that the stress intensity factors depend on the material properties and the geometry of the region, while the stress components at the crack tip have a singularity of the order of J?. This singularity is similar to that obtained for isotrOpic materials. 4.8 DISPLACEMENTS The substitution of Eqs. 4.19 and 4.22 into 3.14 results in the ex- pressions for the displacements. These expressions are somewhat modi- fied by using the definitions 4.29: u=2(xcll - yclémefgscuz - u2)(u2 - 111M112 — 111)} + 2 c Re <6 3(u -E)+u3 y 11 3 H1 2 2 2 2 1 “2 1 “2 _ [—2 "‘2 f——2 2 _ [22—7 + ZREEIEl-(uZ-uz) zl-a +P2m2-Ul) 22 - a +p2(u2 -111) 22—8} NP ' n < \ n + 2R3 i {g’lnEIWZ '32) ("‘1 + 41% '32) +p2(1-2 ”“1)(22 + V23 '82) n=2 p ' __ 2 n +132(“1 ”2) éz ' ‘22 "a )D NN 2 2 n 2 n , _ .. 2 + 2R3 {’EnEfluz "“2) (7‘1 " 21 “a > +p2(“2 "“1)(22 " 22 'a ) n=2 2 n +32(“1 “2) Z2 + 22 ”a 52 and v=2yclzne{gas(u2 - U2)(u2 - ul)(U2 - 111)) 2xC22 “1“1”:“zne n .. + 2C11Re [E g1n(p2 - u2)(u2 - “1)C—2 l)(cos n O + i sin n 2% n=2 NN n -— -— v o + ZCllRe E [a g3n(“2 - u2)(u2 - Ll)(U2 - ul)(cos n G - 1 Sln n 6% n=2 4.53 2C =—-—-2—2=-Re{a‘€sul (11 -?)(u - 11 )(17 - 11 )cos 0 u u u u 2 2 2 l 2 l l l 2 2 + iagu 111-012 WZHUZ - u1)(u2 - 111)sin 9) 2C2 2Re NP 3 p u > , a “(Z n“1(‘“2 112M211 -u l)(?2-ul)(cos n 9 + 1 Sin n 0:] “1 1 2 2 _ n—2 2C2 2R2 . “fin—17i- 2 :1: [an WUIF‘1(112- u2)(u2- 111) (Hz -Ul)(cos n O - 1 8111 n O]; 4.54 With the expressions of the crack displacements, the formulation of the method of solution of problems of central cracks in finite aniso- trOpic regions can be considered complete.. As a matter of fact, the 54 method described above could also be used for the Second Fundamental Problem by employing Eqs. 4.50 and 4.51 for the prescribed displace- ments on the outer boundary. }< \s V// 7 z-plane N \\ g-‘w \ \ &\ >3 W “\\\ .\\\‘ g-plane Figure 3. - The mapping of the unit circle into a crack. iii 56 22-plane $€ g-plane Figure 4. - The mapping of the doubly-connected regions, L, L1, and L2 into their corresponding parametric planes. 57 +3 $15 g-plane Figure 5. - The region, Lg, as continued across the unit circle. 9 (29% ix _> a l.— z-plane Figure 6. - The crack tip coordinate system. 58 T I I n 1 l i if ’1 Ti T h <- a» h ‘=.r/§9 J \px : I' BE X -- a+~ a .. al- .. 6 *C—. Ell-axis : c—-_,,, Ell-axis Ll l i i HT T T (a) t (b) Figure 7. - The tension problems of a rectangular plate with a central crack. Figure 8. - The shear problem of a square plate with a central crack. 5, THE ANISOTROPIC INFINITE PLATE CONTAINING A CENTRALLY LOCATED TRACTION-FREE CRACK UNDER CONSTANT LOAD 5.1 INTRODUCTION The conventional method of solving the problem of a centrally cracked infinite anisotropic plate under constant loading is due to Lekhnitskii (ref. 4) and Savin (ref. 6). They consider the superposi— tion of a constant state of stress in an infinite plate and that of the stress disturbance caused by the presence of the crack. In this section it will be demonstrated that the solutions of Lekhnitskii and Savin can also be obtained by a direct limiting process applied to the dimensions of a finite rectangular plate under constant loading. In the following, three problems will be considered: 1. Central crack in an infinite plate under constant load in the y—direction. 2. Central crack in an infinite plate under constant shear loading. 3. Central crack in an infinite plate under constant load in the.x-direction. 5.2 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT LOAD IN THE y-DIRECTION The boundary conditions 4.24 and 4.25 are Specialized to the finite plate shown in Figure 7(a). The infinite plate can then be considered by letting the dimensions of the plate approach infinity. As a conse- quence of 21+ m, one must set 59 60 Aln=Bln=0 for n2_2 in order to preserve the finiteness of the stress components at infinity. Then the boundary conditions become 2A 2All 11 _ — ' _ _ __ —————a(u2 "Tl-2) [(112 - “2)21 + (112 - u1)22] +—_—_a(FZ 112) [(F2 - “2)21 + (112 - u1)22] o _ C1 on x-c 23 23 11 _. 11 _. _ o _ _ - a(u -T1') (“Z-uan—afi -lJ) (Hz-ulnz - Zth+Cl on x- c 2 2 2 2 0 = Tt(x-c)+Cl on y :th 5.1 and ---—-2A11 ( -—)2+ (—- )2 +—--—-2Kll "(-— )2 a(u2-'u'2)‘ “1“2 “2 1 “2 “2 “1 2, aCfi'Z-uz) “1“2 “2 1 l i _ _ 2311 __ _ 213-11 + 112012 - 1192?] _—_a(u2 'VZ) 112012 - 111) 2 -———-—a(fi2 _ 112) 112012 -lTl)z2 on it :4 II 5.2 NO NO on y = ih Since C3 is a finite arbitrary constant, the consideration, 2 + w, applied to Eq. 5.2 results in two independent equations: A11”2(“2 ’ “1) ' Anuzmz ' “1) " B11‘1‘2(“2 ' “1) + 11“2(“2 ' “1) = 0 5.3 and 2 — 2 — _ -.-2 ., ~ —2 '. _ .-.' A11[“1(“2 .' 112) + 22(42 - 111)] - 11["1C2 .22) + 12(42 111% 72 . —- 2 11 2‘t2 ' “1’ + .11”2(“2 ' “1) t 0 5‘4 on both x = in and y ih. Then Eq. 5-1 can be written as, 1 ._ -. _ — IT _ 1: ‘5' _ ' — '1 _ “2“” — .— I .L r O ‘ -—.~. Cia(“2 1‘2’ n X z - 2 u L Coa(u u) = -Tca(iJ -'LT)+ l 2 2 0n X=-C t 2 2 2 — O _ Tt(x - C)a(112 - “2) Claw.2 - 22) + on y=ih 15 g 2 2 In order to satisfy Eq. 5.5, one must take Tta(u2 - U2) A11(“2 ‘ “1) ‘ A1152 ‘ “1) ’ B11(“2 ‘ “1) T 11‘“2 ‘ “1) = 2 5.6 and O '= Cl th 5.7 At this point, it is noted that there are only three equations available (Eqs. 5.3, 5.4 and 5.6) for the determination of four con— stants. Consequently, one of the constants must be arbitrary. One may also note that these three expressions correspond to the stress compon- ents at infinity. From Eqs. 5-3 and 5.6, one obtains T au t 2 A -B =—f—_— 5.8 11 ll - 2(p2 pl) 62 It is then convenient to express the stress components (Eqs. 4.39, 4.4” and 4.41) as follows: A - B _ 2- 2 _ 2__ 2 11 11 2 Z 21 a 2 -- Z2 Z2 3 Ox = 2R9 ( 51")- “1(“2 ”“2) Wu +“2(U? "“1) *‘ 1.1.7:?“ ‘X a“2‘2 ‘ “2_2 ‘ V72”? 21 a 42 (1 g 3' — 22 -a2 -2 2 2 +u2(u2~ul) _____;____w > 22 _. a2 2 J) A B 11 2 _ 2 _ 11 _2 2 2 2 2 5.9 A - B 2 z 1 _ 2 o 2Re 8(11 in) (112 - 112) + (112 - ul) “2 2 2 _ 2 2 _ 2 21 a 22 a .22 + (u2 - 1J1) _._... 5.10 .2 _ 2 22 a A B 2 z 11 11 — 2 o =-2Re — u(u-u) +u(U-u) xy (112112) 12 2 2-2 22 1 2_2 21 a 22 a “2‘2 “ + U (u - u ) —-—-—y 5.11 2 2 1 _2 2 22 " a; Upon applying the boundary conditions in the forms of 5.4 and 5.8, the ‘ expressions for the stress components are obtained: 63 Ttuluza L'l L42 1 xx - Re u _ u - 2 1 2 2( + ,2 2) Z2_az(z + I/2_ 2) 21 —a 21 I (.1 —a I 2 2 22 a 5.12 It [ ”221 "122 , .. o . = Re _ _ H - 3.15 y’ ”2 *‘1 2 _ 2 2 _ 2 21 a I 22 a _ ' T u p [. z z T o = Re t} 2 1 - 2 5.1:. "y 1 “2 [22 _ ,2 1’22 __ 2 1 d 2 a These expressions agree With either Savin's or Lekhnitskii's results. In order to obtain the stress field in the neighborhood of the crack tip, the following transformation is used: 2 = a + r(cos 6 + p k Sin 9) (k = l, 2) 5.15 k where r, 6 and a are defined in Figure 6. r , Considering -— << 1, one obtains a 2 2 2 2 Q g . _ _ . zk-a (zk+ zk-a)-a4/2ral\/tosb+pk51ne (La-1,2) and z .___..._15.._.__. 8‘ (k=l,2) z: _ a2 NZra Vcos 6 + ukSln U Hence, the stress components are written for the neighborhood of the crack tip as, 64 {>1 t “1“2 H1 “2 o = Re - xx. JEE' U2 _ u1 Jcos 6 + plain 9 Vbos 6 + uzsin 6 5.16 = TH]; Re 1 “2 _ "1 W .122? “2 ‘ “1 Eos e + ulsin 9 /cos a + uzsin e 5.17 O a Tr“; Re “1“2 1 _ 1 xy /§;' lul - u2 Vbos 9 + ulsin 6 /cos 9 + uzsin 6 5.18 The stress intensity factors are obviously = f- = . Kl Tt a K2 0 5 l9 i.e., the sliding mode stress intensity factor is not present in the solution. 5.3 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT SHEAR LOADING According to Figure 8, the boundary conditions, 4.24 and 4.25, with z + w, become 2A 2X l1 __ —- ' 11 _. _. ._ ._ a(u2 _11—2) [(112 - "2)21 + (:12 - ulna +a(fl'2 ”12) [(112 -u2)zl + (112 "‘11)sz - 23 _ 2B11 _ - a(1.l2 '32) (“2 - u1)22 - afiz - 112) (HZ - u1)::2 a on x==3tc y=h =-h 5.20 65 and 2All A a(u2 4:2) [“1(“2'“2)21+“2(“2"“1)22] + a(,.2 — )1: j2"""2)Z1 28 2B- _ . -— _w 11 —- .Ll ._ ._ + u--! .. --————-—— - ‘ - ._ ,. 52-, 2'_ u2(2 L11..)‘2'23 a(.12 —:2) 2(“2 1)22 ~( 2 - 2) 2( 2 2 r o _ =. CZ on X L. o = <2Tc+C on x=—c 5.21 s 2 o LTs(c-x)+C2 on y-ih The boundary conditions, 5.20 and 5.21 can further be Simplified, re- sulting in the following three equations in four unknowns 2 - 2—. .. _-— —2—_,- 72._.— A11E‘1(“2 -' “2) + 1‘2”2 “1% A11l2’1“2 "2) + “2(“2 ”13] _. 2 __ __ B11“2(”2 " “1) + B11“2“‘2 ' “1) ’ O 5 22 A11(u2 - “1) - A110:2 - “1) - B11012 - “1’ + 311(“2 - t1) = 0 5 23 A11“2(“2 ' '1) ' A11“2“‘2 ‘ “1) ' B11“2(“2 ’ “1’ + M11 2 ‘ 1) T 8(h p ) 2 2 _ _ 2 5.24 The unknown integration constant was found to be 00 = - r c 5 25 66 One may again note, as in the tension problem above, that one of the four unknown constants is arbitrary. From Eqs. 5.23 and 5.24, one obtains Tsa A - B = 11 ll 2(u2 - pl) 5.26 Upon substituting Eqs. 5.22 and 5.26 into the expressions for the stress components (Eqs. 5.9, 5.10 and 5.11), the resulting expressions become 7 2 2 2 Ta u u o = Re 3. l _ 2 xx u2 H1 22 _ a2 + (22 _ a2 22 _ a2 + 1,22 _ 2 1 z1 1 2 z2 2 a 4 5.27 T 2 z o = Re _? 1 --———3———- 5.28 yy “2 111 M22 _ 2 “22 _ a2 1 a 2 o = 5.29 T u z u 2 Re -s 1 l _ 2 2 xy 1 u2 ‘Lz _ a2 22 _ a2 1 2 With the transformation, 5.15, the crack tip stress field was found to be 3w 2 a I: Ts Re 1 “1 _ u2 xx '/2;' u2 - ul /cos 9 + ulsin 6 Vbos 6 +-pzsin 6 5.30 67 T943 1 1 1 o = :__ Re ‘ _ . [ , — yy /2r ”2 L1 Vcos 6 + u sin a /ccs 6 + p251n 6 l 5.31 T V; . u . e ' i 4 _ J " .. . . , y VZr l 2 Vcos o T t 510 a 2-: r v s sin 3 l 2 5.32 In this case, the stress inten51ty factors are K =0 K =T/E 5.33 1 2 s i.e., the opening mode Stress intenSity factor 18 not present in this solution. 5.4 CENTRAL CRACK IN AN INFINITE PLATE UNDER CONSTANT LOAD IN THE x-DIRECTION In accordance with Figure 7(b), as 2 + m, the bOundary conditions, 4.24 and 4.25, become 2Al 2A 1 M —i' -' _ 1 _____il__. —' _ ‘-' . _r- - 3012-712) [(“2 L2)"‘1+“‘2 L‘1)223‘“.:.(E2_.12)[“12 “2)z1+(”2 “fizz-J — CO on x—i’ 2B 23 1 ’ “ -'—?—-J%E—y (U2'-ul)Eé 77:7—lL-3-(UE -31)29 = 5.34 auz‘ 2 ““2‘“2 ‘ O C. on y =;th 68 2 2A11 A1 _ _ 1 _. _ _ ____..__ _ _ +___ _ a(u2-TJ'2) [111012 u2)zl+u2(u2 “1)7‘2] 8(11'2 '“2) [111012 112)zl 2B11 2B11 + 112012 -ul)22] -__a(u2 —U2) 112012 -u1)zz -———a(F2 "“2) 112012 ‘ul)22 ' + o _ i Tty C2 on x- c o 2 TLh-tC on y = h 5.35 v o — + =— Tth C2 on y h These conditions will also result in three equations having four unknowns 2 _. 2 _. .2 _2 _. 11111211012 - 112) + 112012 - 111)] - A11[“1G2 - 112) + 112012 - 111-2] ' Ta(u -T1') .2 - 2 _. _. _ t 2 2 ’ B11“2(“2 ‘ “1) + B11“2(“2 ’ “1) ‘ 2 5'36 A11(“2"”1)"A11(T"2 '“1) ’311(“2"“1)'+311(“2 "“1)‘=0 5'37 A11“2(“2 ' “1) ' A11“2(“2 ' “1) ‘ B11“2(“2 ’ “1) + 11“2(“2 ' “1) = 0 5.38 From Eqs. 5.37 and 5.38 All - B11 = 0 .5.39 Hence, upon substituting into the expressions for the stress components, one obtains I = T 5.40 o = o = 0 5.41 69 Therefore, one can conclude that the stress field is not influenced by the crack in this particular case. 6. PROBLEM FORMULATION FOR A FINITE ORTHOTROPIC RECTANGULAR PLATE CONTAINING A CENTRALLY LOCATED TRACTION-FREE CRACK UNDER CONSTANT LOAD 6.1 INTRODUCTION The motivation for considering the problem of a finite orthotropic plate with a central crack is mainly due to the structural importance of unidirectional materials in present day.design. As it will be shown, the handling of such problems is somewhat simpler than that of aniso- ' tropic problems with regard to the size of the coefficient matrix in the boundary collocation method. This simplification is a consequence of stress symmetry prevailing in a symmetrically loaded rectangular ortho- tropic plate. The consideration of stress symmetry will lead one to the conclu- sion that only the odd terms.are needed in the Laurent series expansion of the stress potentials. One also finds that the arbitrary integration constants, C? and C3, can also be determined.before.collocation. For the loads considered in the case of the anisotropic infinite plate (Figures 7(a), 7(b) and 8), the boundary conditions of the finite orthotropic plate will be.given in this.section in a form readily appli- cable to the boundary collocation method. 6 . 2 STRESS SYMMETRY . From physical considerations of the.problems.at hand (Figures 7(a), 7(b) and 8) in the case of.orthotropic materials one can.expect symmetry 70 71 in the stress components su:h that OXX(Z) = Gxxfi-Z) o (z = ” '—' 6 l W ) n" ") O (2; G 1-2) xy X] It can be shown that in order to satisfy these symmetry conditions, tne has to set ' = Q: = 0 for n ev*n gin 2n c As a consequence to this relation among the unknown :cnstanrs, the sat- isfaction of the boundary conditions is required on the y 5 0 pa:: of the boundary only. For the further consideration of special orthotropy, i.e., for o o _, , , 6 = 0 , 90 , in addition to conditions 6.1, the fellowing symmetry is also obvious xx xx 0 z = o 4“ 6,2 yy( ) yy( ) o 2) = o —z xyk Xy( ) The symmetry conditions 6-2 will lead one to the conclusion that the constants 9g, %%, gin and gin are pure imaginary quantities— Thus, the above problems can be solved with the satisfaction of the boundary conditions on the boundary in the first quadrant for specially ortho— tropic materials. 6.3 THE DETERMINATION OF THE INTEGRATION CONSTANTS First, the boundary conditions 4.27 and 4.28 are written at (c,0) and (-c,0): 72 .. 2111:3152 «1.15“ (1+ Rm Lev“ n=l u (u -u ) n 2 n 2 n 2Re % 21122-121 (Am ”313] (i) (1+ 1 - (i1) ) - (1+ 1 - (5:1) ) .-. J; Xnds] + (:3 6.4 = - / Y (18] +CO 6.5 s n l = / de] +CO 6.6 s n 2 y=0 Upon considering n odd (orthotropic materials), the addition of Eq. 6.3 to Eq. 6.5 and Eq. 6.4 to Eq. 6.6 will result in 73 C°=-l- Yd; + Yd; 1 2 n n s s 6 7 x=c x=-c y=0 y-O cg=%-de] - de] s n s n 6.8 x=c x=—1 y=0 y=0 Then for the three different cases of loading shown in Figures 7(a), 7(b) and 8, the arbitrary integration constants become Constant load in CO = T t l t 6,9 y-direction: C3 = 0 o Constant C1 = 0 o 6.10 shear loading: C. = -T c 2 s Constant load in CE = 0 o 6.11 x-direction: C2 = 0 6.4 BOUNDARY CONDITIONS In addition to the above observations, one also finds that the constant 55¢ = (112 - ul)(u2 - U2)(U2 - “1) is an imaginary quantity in the case of special orthotropy. Conse- quently, the square matrix gm? (Eq. 4.37) becomes singular due to the third column of 6' (Eq. 4.38) containing zero elements only. In order to account for the case of special orthotropy among the various orthotropic cases, it is noted that 74 yRegéRegs - yIm‘oZImég = yRe[?C?S] Considering this contraction and the determination of the integration constants, Eq. 4.38 becomes r 0 ReF y ReG 11 ll -ImF 1 2NN p Ree: Re [scars-11 * LII“ 2NN J It can be noted, here, that upon solving Eqs. 6.12 at in T 6.12 M points, the stress components are fully determinate, (see Eqs. 4.39, 4.40 and 4.41), however, the displacement expressions only allow the deletion of the first two columns of Eq. 4.38. Since this dissertation primar- ily deals with stresses and stress intensity factors, the form given in 6.12 will be used. At this point, the boundary conditions-for the three different cases of loading (Figures 7(a),.7(b), and 8).are written in a compact form readily applicable for automatic computation. 75 . . .1231: -ImF . . .ReF -ImF (Rem/C(55)} P Rel'1 '1‘“? INF 1NP 2NN 2NN 1 11 Ree? -— Imi? Re%’ ImrlNP [Re%ENN I""652111: J Ttx - Tsy -% , z in r 6.13 TLy-Tst 6.5 SOME REMARKS ABOUT THE COMPUTER.PROGRAM The consideration of Eqs. 6.13 at ‘M points on the boundary of a rectangular plate results in a matrix equation such as 4.32, where a - a = K 2Mx(NN+NP+l) (NN’+NP+i)x1 2Mxl 6.14 Then for the least squares boundary collocation method,.the problem becomes 2': 6T2? - fTK 6.15 * where ’8 is the desired solution vector. The solution of Eq. 6.15 was carried out on the IBM TSS/360 com- puter by declaring all quantities and operations in double precision 76 and utilizing the Gauss-elimination technique with complete pivoting. It should also be mentioned that the elements of the coefficient matrix, n a were scaled b the devisor Rn = c2 + hzfi . y The listing of the complete program is given in Appendix B. 7. RESULTS AND DISCUSSION or res SGLUTLONE CF THE TENSION AND saran PROBLEMS 7.1 INTRODUCTION In this section, the solutions or. he 7 U’ (P EB vr shown in Figures 7(a), 7(b) and 8, are presented and disccssed For both ihe tensiin and shear problems the effe2ts or r? he va;1ations or the dimensionless . 4. '- 'enca ion angle and the a’fi ratios (J P( *4 material ratios (3.30, 3 32), the are considered. For the tension problem, variation in the h.c ratio is also taken into account. 7.2 SPECIFICATIONS OF PARAMETERS TREATED The solutions of the three problems shown in Figures 7(a), 7(b) and 8 are based on the following specifications: 1. For all cases a: .5 inches unless otherwise specified (see Table 1). 2. The problems were non-dimensionalized by c;nsidering the can— cellation of the constant boundary Stresses frzm bcth sides of the boundary condition equations. This is effectively equiva— lent to setting the applied boundary stresses equal to unity, If other values of the applied boundary stresses are used, they are properly indicated (see Table l) 3. With the exception of the investigation of the effects of material properties on the stress intensity factors, the ratios [11 .366 {T} 78 and E11,1322 = .257 E 2 2612 12 were used throughout the entire analysis. These ratios corre- spond to the properties of fiberglass. In studying the effects of varying material properties on the stress intensity factors in both the tension and shear problems, the following parameter ranges were considered: E .05_<_E22$l 11 E E .001; 11/ 22 £00 (E11 )2 -——-— v 2G12 12 5.3.; h=1; 6 =45° c 3 c . In studying the effects of varying plate aspect ratio, h/c, and crack length to plate width ratio, a/c, on the stress intensity factors in the tension problem, the parameter ranges were chosen as follows: C) A nlm M uflna .25.: nls‘ 79 6. The effects of variation of the orientation angle, 6, were also analyzed for both the tension and shear problems considering the following parameter ranges: E. c 2 0< :3 Os 0 O) 5.90 D. c 7.3 EXAMINATION OF THE MAPPING COLLOCATION METHOD The available results which could qualify the validity of the mapping-collocation method for centrally cracked finite rectangular orthotropic plates are those of the centrally cracked infinite plates, given in section 5, and the results of Kobayashi, Isida and Sawyer (refs. 1, l9 and 20) for finite isotropic rectangular plates. In order to gain a certain degree of confidence in the mapping—collocation method as deve10ped in the preceeding sections, both the results for the in- finite orthotropic and finite isotropic plates were verified. One knows that in the infinite plate solutions A = B = 0 for ln 1n n 2.2 and furthermore only the difference A11 - B11 is needed to ex- press the stress components and the stress intensity factors whether the plate under consideration is anisotropic or orthotropic. Hence, the very first verification of the mapping-collocation method would result from the use of the constants, Re[g2eg] and 3’ d plate dimensions as compared to the crack length. However, this type of , only, with fairly large substantiation of the mapping-collocation method rests on the a priori assumption that a large finite plate is the same as an infinite plate. 80 A better method to verify the infinite plate results involves the use of the "full" form of the stress potential, ¢l(;). Upon consider“ in" a 100"x100" plate with a one inch crack, it was found that the stress intensity factors were quite insensitive to the number of terms of the Laurent series expansion of ¢l(§). For example, the use of Re[%:%;] and @% (three unknowns) resulted in Kl = .70712 JTH. and K2 = .00000 while using 83 unknowns resulted in the stress intensity factors, Kl = .70713 VIE. and K2 = .00000 with Tt = l. The corre- sponding analytical results for an infinite plate are, Kl = V75 JIH and K2 = 0. It was, however, also noted that the satisfaction of the stress boundary conditions markedly improved with the use of more and more terms in the Laurent series expansion of ¢l(§). Thus, the results given in section 5 were verified for both the tension and shear problems. The verification of Kobayashi's Isida's and Sawyer's opening mode stress intensity factors for various finite isotrOpic rectangular plates in tension was carried out by setting the ratios and £11-, 2 2G 12 12 The value of the orientation angle had no influence on the results. The cases considered are shown in Table l, where the various stress inten~ sity factors are tabulated and can readily be compared with each other. 81 As can be seen from Table 1, agreement with the known results is ex- tremely good, as a matter of fact one cannot speak of a higher degree of accuracy than the one obtained when numerical methods are used. In order to economize the mapping-collocation method with regard to the use of the number of terms in the Laurent series expansion of ¢l(C), two different truncation numbers NN and NP were assigned as limits of the sums in the various expressions. The numbers NN and NP are odd for orthotropic materials and they are actually correspond— ing to the number of the negative and positive exponent terms in the Laurent series expansion of the stress potential function, ol(c). Upon investigating various ‘NN/NP ratios for the finite plate problems, it was found that a ratio of NN/NP = 5/3 resulted in fairly fast conver— gence with rather well-satisfied boundary conditions. This NN/NP ratio was retained throughout the whole investigation of the problems given in Figures 7(a), 7(b) and 8. One of the most important indications of the degree of accuracy of the solutions obtained by the mapping-collocation method comes from the examination of the boundary stress between collocation points (ref. 2). Since the least squares method of collocation results in the satisfac— tion of the boundary conditions only in the approximate sense, the mag- nitude and frequency of oscillations of the boundary stresses serve as an indicator to the "exactness" of the solution of the stress boundary value problem. Thus, in each case, when the problem of a centrally cracked finite plate was solved, the boundary stresses were also exam— ined. For example, in the case of a square plate in tension with 82 a/c = 2/3 and 6 = 45°, the largest error in the boundary stresses was found to be about 4% and it occurred in a sudden manner at the corners of the plate. Since this dissertation deals with stresses and stress intensity factors of centrally cracked finite rectangular orthotropic plates, the pattern of convergence of the stress intensity factors, K1 and K2 was also investigated. It was observed that the insensitiveness of K1 and K to a change in the number of terms of the Laurent series ex- 2 pansion of the stress potential, ¢l(t), is an excellent indicator to how well the boundary conditions are satisfied by the least squares approx- imation. Referring to Figures 9, 10 and Table 2, one observes two types of convergences: 1. Using approximately twice the number of equations as unknowns or larger, i.e., 2M.z.2(NN + NP + 1), will result in nonoscilla— tory, steady values of K and K . These values of K and 1 2 1 K2 are still subject.to.variations as.the.number of.unknowns changes. 2. The second type of convergence is that when K1 and K2 become insensitive to an increase in the terms of the Laurent series expansion of the stress potential, ¢l(t), i.e., when a further increase in the number.of unknowns does not.change the values of K1 and K2. Thus, it is noted that for a square.plate in tension with .a/c = 1/2 and 6 - 45°, the number of unknowns NN + NP + l 8‘41 is already suf- ficient to result in converged values of K1 and K2. However, one may also note that the boundary conditions are somewhat better satisfied 83 when NN + NP + l = 83, despite the fact that the stress intensity fac— tors have already converged. Hence, one can draw the conclusion that for increasing plate dimensions less and less terms of the Laurent series expansion are needed in order to obtain converged values of the stress intensity factors. This fact was also observed rather radically in the case of the 100"x100" plate. Upon considering various materials, plate sizes, a/c ratios and orientation angles under the loading conditions shown in Figures 7(a), 7(b) and 8, two conclusions became obvious regarding the convergence of the stress intensity factors: 1. Convergence of the stress intensity factors is definitely affected by the material parameters. It was found, for example, 0 that for a square plate in tension with a/c = 2/3 and 6 = 45 the use of the dimensionless material conStant Ell/E22 (E11 _ v )2 2612 12 = .0001 E with any EZZ- ratio resulted in oscillatory values of K1 and 11 K2. In this case, the number of unknowns was taken as 83 which gave excellent results for higher values of the dimensionless material constant. 2. Convergence of the stress intensity factors is also affected by the a/c ratios. For a square plate in tension with E /E E E 11 22 2 = .257 E33-= .366 ( 11 __v ) 11 2G12 12 84 a/c = .9, 6 = 45° with the use of 83 unknowns, it was found that the boundary conditions were rather badly satisfied, the. method probably resulting in unconverged stress intensity fac- tors. For smaller a/c ratios well-converged values of K1 and K2 were obtained. 7.4 SOLUTION OF THE TENSION PROBLEMS OF AN ORTHOTROPIC RECTANGULAR PLATE WITH A CENTRAL CRACK As shown in Figures 7(a) and 7(b), two distinct tension problems were considered in this study. The solution of the problem for loading applied in the x-direction (Figure 7(b)) for various h/c and a/c ratios, orientation angles and material constants resulted in zero values for the stress intensity factors. This fact is of significance because it shows that just as in the isotropic case, the stress field in the plate is not disturbed by the presence of the crack.7 It should furthermore be noted that in the case of Figure 7(b) the crack is aligned with the direction of loading, and for this particular config- uration the plate behaves as a finite orthotropic rectangular plate loaded by a constant boundary stress in the x-direction without a crack. The consideration of the other tension problem, i.e., when the constant boundary stress is applied in the y-direction (Figure 7(a)), resulted in various parametric studies involving material properties, a/c ratios, h/c ratios and orientation angles. 7.41 Effects of Material Properties on the Stress Intensity Factors According to expressions, 3.30 and 3.32, an orthotropic material can be characterized by two dimensionless ratios constructed from the engineering material constants, C 2, and v . As 3 E11’ E22’ 1 12 85 consequence to this method of characterization, a parametric study was undertaken in order to determine the effects of these two dimensionless ratios on the stress intensity factors. Upon considering a square plate with a/c = 2/3, 6 = 45 and Tt = l, the dimensionless ratios were varied and Table 3 was obtained. From Table 3, the family of curves shown in Figures 11 and 12 were con- structed. Figure 11 demonstrates the effects of the dimensionless material ratios on the opening mode stress intensity factor while Fig— ure 12 shows the sliding mode stress intensity factor as a function of these two ratios. Considering Figure 11, one arrives at the following conclusions: 1. The complete family of curves is asymptotic to the E11/E22 _ 0 (E11 _ V )2 2G12 12 E /E ordinate, i.e., as 11 22 + 0, K + w. E 2 l < 11 __\, ) 2G12 12 2. The curves are asymptotic to various Kl'-va1ues as E11/E22 E 2 + cc; and they are asymptOtic from below the 11 ) -———-- v (.12 . asymptotes. 3. As the E22/E11 ratio decreases, the opening mode.stress in- tensity factor, K increases in value. 19 86 4. For each constant E22/E11 ratio, the stress intensity factor, Kl, will attain a minimum value which is not the K -asymptote. 1 5. Upon examining the parametric curve for E22/Ell a 1, it seems that the smallest possible minimum value which can be attained by K1 is K1 = 1.131 VIE: occurring at E11/E22 _ 2 (3111“,) 2612 12 This value of K1 is very close to the value of the opening mode stress intensity factor, K for an isotropic material 1, (X1 = 1.134 753.). The consideration of Figure 12 results in conclusions somewhat different from those for the Opening mode stress intensity factor. First of all, one may note that the presence of the sliding mode stress intensity factor, K2, is strictly due to the orthotropy of the material and to the finiteness of the h/c and the a/c ratios. It may be noted that the maximum magnitude.of the sliding mode stress intensity factor shown in Figure 12 is about 14% of the value.of.the opening mode stress intensity factor for the same E22/E11. ratio. The closer examination of Figure 12 results in the following conclusions: 1. The complete family of curves is asymptotic to the E11/E22 ('Eii _ v >2 2012 12 = O 87 ordinate, i.e., as E11/E22 (Ell _ v )2 2G12 12 2. The curves are asymptotic to various K -values as +0, K+-m° 2 E11/E22 m. E11 _ V >2 2612 12 3. As the E /E ratio decreases, the absolute value of the 22 ll sliding mode stress intensity factor, K increases. 2’ 4. For each constant EZZ/E ratio, the sliding mode stress 11 intensity factor attains its minimum absolute value at E11/E22 = m. (E11 _v )2 2012 12 5. The smallest possible minimum absolute value of the sliding mode stress intensity factor, K2, is zero as one would expect, i.e., for an isotropic material K2 = O. 7.42 Effects of Plate Width. and Plate Size on the Stress Intensity Factors In order to investigate the effects of various plate widths and plate sizes on the stress intensity factors, K and K the dimension- 1 2’ less material ratios were set at E E33-= .366 11 and 88 E /E 11 22 2 = .257 2G12 12 with an orientation angle, 6 = 45°. These ratios correspond to fiber— glass properties. The variation of the h/c and a/c ratios resulted in Table 4. In accordance with Table 4 graphical representations were constructed as shown in Figures 13 and 14. The family of curves for the opening mode stress intensity factor, K1, (Figure 13) indicates the following: 1. As a/c + 2/3 the value of the opening mode stress intensity factor becomes larger and larger for each constant h/c ratio. 2. As a/c + O the value of the opening mode stress intensity factor approaches K1 = JTE'JIEZ, the value of K1 for an in— finitely large plate, also for each constant h/c ratio. 3. As the h/c ratio decreases, the opening mode stress intensity factor, Kl’ increases in value. 4. As the h/c ratio becomes very large, the opening mode stress intensity factor approaches K1 = /:§ JEEL for the range of values 0 < a/c 5.2/3. This seems to indicate that for an in- finite strip, only very high a/c ratios will effect the stress intensity factor, Kl, such that it will deviate from the K1- value of the infinite plate. The consideration of Figure 14 leads one to.the following conclusions: 89 H . As a/c + 2/3 the absolute value of the sliding mode stress intensity factor, K2, becomes larger and larger for each con- stant h/c ratio. 2. As a/c + O the sliding mode stress intensity factor approaches the value, K2 = 0, i.e., the value of K for an infinitely 2 large plate. This also occurs for each constant h/c ratio. 3. As the h/c ratio decreases, the sliding mode stress intensity factor, K2, increases in absolute value. 4. As the h/c ratio becomes very large, the sliding mode stress intensity factor approaches the value, K = 0, i.e., the value 2 of K2 for an infinite plate for the range of values 0 < a/c 5.2/3. The presence of significant values of K are due to the value of the 2 orientation angle (6 = 45°), however, the variations in K2 are defin- itely induced by the variations in both the a/c and h/c ratios. Thus, one can conclude that the finiteness of the plate with a/c > 0 accounts for the increase in the values of both stress intensity fac- tors as compared to the case of the infinite plate. 7.43 Effects of Orientation Angle and Plate Width on the Stress Intensity Factors For the investigation of the effects of the orientation angle and the ‘plate width on the stress intensity factors, the h/c ratio was taken as one (square plate) and the dimensionless material ratios were held constant at the values specified above. The variations of the orientation angle, 6, and the a/c ratio re— sulted in Table 5. The graphs constructed from.Table 5 are presented in Figures 90 15 and 16, where Figure 15 shows the effects of various orien- tation angles, 6, on the opening mode stress intensity factor, K1, for constant a/c ratios, while Figure 16 shows the sliding mode stress intensity factor, K2, as affected by changes in 6 and a/c. An examination of Figure 15 results in the following observations: 1. As the a/c ratio increases, the opening mode stress intensity factor, Kl, also increases. As the orientation angle, 6, increases, the opening mode stress intensity factor, Kl, becomes larger and larger for each con- stant a/c ratio. The curves indicate zero slopes at 6 = 00 and 6 = 900. For each constant a/c ratio, the maximum value of K occurs 1 at 6 = 900 and the minimum value at 6 - 0°. The smallest possible minimum value of the opening mode stress intensity factor is /75 JIEI which is the Kl-value for an in- finite plate. The conclusions concerning Figure 16 are as follows: 1. 2. The existence of the sliding mode stress intensity factor, K2, is strictly due to the values of the orientation angle, 6, dif- ferent from 00 or 90°. As the a/c ratio increases, the sliding mode stress intensity factor, K2, also increases in absolute value. For each constant a/c ratio, the absolute value of the slid- ing mode stress intensity factor reaches a maximum at different values of 6. For example, the maximum absolute value of K2 for ale 8 2/3 occurs at about 6 = 50°, with K2 being -.065 JR. 91 4. The smallest possible absolute value which is attained by K2 is zero as expected for an infinite plate. 7.44 Stress Distribution on the y = 0, x > a Line As an illustration of the typical stress distribution on the line y = 0 for the range 1 < x/a 5.3/2, the values of the stress components were tabulated in Table 6 for a plate with h/c = l, a/c = 2/3, 6 = 450 and plotted in Figures 17 and 18. Figure 17 shows that both normal stresses, Gxx and Oyy’ become asymptotic to the x/a = 1 line and oxx becomes zero at x/a = 3/2 which is a point on the boundary. This is expected since Oxx was specified to be zero on the x/a = 3/2 boundary. The other stress component, Oyy’ has the value of Oyy = .68 at y a 0, x/a = 3/2. The intersting behavior of the shear stress is.shown.in Figure 18, where in the.neighborhood of the crack tip Oxy < 0 and Oxy + - w as x/a + 1. This behavior of ny in the neighborhood of the crack tip could be expected from knowing that K2 was found to be negative for the same case (see Figure 16). Then as x/a increases Oxy will reach a maximum value (Oxy = .153) at about x/a = 1.35 which is 15.3% of the applied stress.. It is noted that the presence of oxy is due to the general orthotropy of the material. It should also be observed that oxy becomes effectively zero on the boundary, which it should be for proper satisfaction of the boundary conditions. It should be mentioned at this point that for each .K1- and K2- value obtained,.the complete boundary value problem of.an.orthotropic rectangular plate containing a central crack had to be solved. In each case, in addition to the computation of the stress components on the 92 boundary, the stress components on the y = 0 line were also recorded. It gives one confidence in the method to know that for isotropic and specially orthotropic materials, the shear stress component, Oxy’ was effectively zero on the y = 0 line. It should also be noted that for problems with EZZ/Ell = l, the shear stress component, sky, on the y a 0 line was found to be practically zero for any value of E11/E22 (E11 _ V )2 2012 12 7.5 SOLUTION OF THE SHEAR PROBLEM OF AN ORTHOTROPIC SQUARE PLATE WITH A CENTRAL CRACK In addition to the tension problem discussed in the preceding sec- tion, the problem of an orthotropic square plate containing a central crack and loaded by unit shear stress (Figure 8) was also considered. As in the case of the tension problem, the shear problem was also solved for various dimensionless material ratios, a/c ratios and orientation angles. The h/c ratio was held constant (h/c = 1) throughout the en- tire analysis of the shear problem. 7.51 Effects of Material Properties on the Stress Intensity Factors The results of the parametric.study involving.variations of the dimensionless material ratios while considering a/c - 2/3 and 6 = 450 are summerized in Table 7 and Figures 19 and 20. An examination of Figure 19 reveals.the following effects of the dimensionless material ratios on the opening mode stress intensity factor: 93 l. The complete family of curves is asymptotic to the E11/E22 E 2 = 0 (11_v) 2G12 12 ordinate, i.e., as E11/E22 E 2+0, K1+—oo. 2G12 12 2. The curves are asymptotic to various values of K1 as Ell/E22 (E11 _v )2 2G12 12 3. As the E22/E11 ratio decreases, the absolute value of K1 increases. 4. For each constant E22/E11 ratio, the opening mode stress intensity factor, Kl, attains its minimum absolute value at E11/E22 = ( E11 - v )2 2612 12 5. The smallest possible minimum absolute value of the opening mode stress intensity factor, K1, is zero as one would expect i.e., for an isotropic material Kl = 0. The effects of the variations of the dimensionless material ratios on the sliding mode stress intensity factor, K2, are shown in Figure 20. According to Figure 20, one arrives at the following conclusions: . As the . It seems that the smallest possible minimum value of K 94 . The complete family of curves is asymptotic to the E11/E22 E 2 = 0 2G12 12 ordinate, i.e., as E /E 11 22 +0 K +00. E 2 ’ 2 (41... v 2012 12 . The curves are asymptotic to various K2-values as E11/E22 2 2612 12 E22/E11 ratio decreases, the value of the sliding mode stress intensity factor, K increases. 2’ . For each constant E22/Ell ratio, the sliding mode stress intensity factor attains its minimum value at E11/E22 = (E11 _v )2 2012 12 2 IS .982 JTEZ and it occurs at E11/E22 = w. (E11 _v )2 2912 12 For an isotropic material K2 = 1.022 Jihl. 95 7.52 Effects of Orientation Angle and Plate Width on the Stress Intensity Factors The variation of the orientation angle, 6, and the a/c ratio re- sulted in Table 8 and Figures 21 and 22. Considering Figure 21, one arrives at the following conclusions: 1. The existence of the opening mode stress intensity factor, K1, is strictly due to the values of the orientation angle, 6, different from O0 or 90°. 2. As the a/c ratio increases, the absolute value of the opening mode stress intensity factor, Kl, also increases. 3. For each constant a/c ratio, the absolute value of the opening mode stress intensity factor, K , reaches a maximum at different 1 values of 6. For example, the maximum absolute value of K1 for a/c = 2/3 occurs at about 6 = 46°, with K being 1 -.073 J17... 4. The smallest possible absolute value of K is zero as expected 1 in the case of an infinite plate. The effects of the orientation angle and the a/c ratio on the sliding mode stress intensity factor are shown in Figure 22. According to Figure 22, these effects can be summarized as follows: 1. As the a/c‘ ratio increases, the sliding mode stress intensity factor, K2, also increases. 2. For each constant a/c ratio, the sliding mode stress intensity factor, K2, reaches a maximum and this maximum occurs at various 0 values of 6. The minimum of K2 .always occurs at 6 = 0 . For 96 example, for the curve, a/c = 2/3, the maximum value of K2 is 1.072 (EH. occurring at 6 = 42° while the minimum value of K2 is .994 7551 at 5 = 0°. 3. The curves indicate zero slopes at 6 = 00 and 6 = 90°. 4. The smallest possible minimum value of K2 is V75 JEEL, which is the K -value for an infinite plate. 2 There is an important observation which concerns the effects of the var- iation of the orientation angle in the shear problem. As far as 6- dependence is concerned, the orthotropic plate behaves differently in tension and in shear. In the tension problem it was found that the max- imum K1 values always occurred at 6 = 90°, (specially orthotropic material) while in the shear problem it is obvious that the maximum K 2 is 6-dependent. However, it is also observed that in both the tension and shear problems the minimum values of the significant stress inten- sity factors occurred at 6 a 00 which designates the other type of special orthotropy. 7.53 Stress Distribution on the y = 0, x > a Line As an illustration of the typical stress distribution on the line y = 0, 1 < x/a.g 3/2, the values of the.stress components were recorded in Table 9 and depicted in Figures 23 and 24 for a plate with a/c = 2/3 and 5 - 45°. Figure 23 shows the normal stress components, Oxx and oyy, as they vary for l < x/a.$ 3/2.7 It is noted that Oxx + - w and o ‘+ - w as x/a + 1, but oyy becomes positive and attains the value YY of qyy - .72 at x/a - 3/2. In agreement with the specified boundary condition, Oxx becomes zero at y = 0, x/a = 3/2. 97 The distribution of the shear stress, Oxy’ is shown in Figure 24, where oxy + w as x/a + l and it has the specified value, ny = l, at x/a = 3/2, y = 0 on the boundary. It should also be mentioned that both Oxx and 0 were found to be zero on the line y = 0, l < x/a 3,3/2 with EZZ/Ell = l and E11/E22 for any value of E 2' The special cases of isotropy and 11 ) --—-— v (2912 12 special orthotropy resulted also in zero values for the normal stresses, oxx and ny’ on the y = 0, l < x/a 5.3/2 line. Tt = in; O‘O‘O‘O‘O‘O‘ 98 Table 1. Comparison of opening mode stress intensity factors for rectangular isotropic plates with central cracks 1000 psi Kobayashi Sawyer Isida Author Ref. 1 Ref. 20 Ref. 19 c a K1 K1 K1 K1 in. . 'in. psivEEI psi/EEK ' psivffil * psi/EEK 3 .25 502.4 504.7 502 502.1 3 .50 720.6 725.0 718 719.2 3 .75 904.8 910.1 900 900.3 3 1.00 1083.0 1073.7 1075 1073.3 3 1.50 1496.0 1401.5 1436 1454.6 3 2.00 2074.0 1952.2 2000 2002.8 Table 2. Typical pattern of convergence of the stress integsity factors for tensile loading (h/c = 1, a/c = 1/2, 6 = 45 , T =‘1) a = .5 in Number of equations 16 24 32 40 48 64 96 128 176 t . Number of Number of Number of Number of unknowns unknowns unknowns unknowns 15 25 ' 41 83 K1 K2 K1 K2 K1 K2 K1 K2 J13. JIE. JZH. ‘ JTH. VIE. J35. /15. J35. .9580 -.0267 .9522 -.0279 .9514 -.0281 .9568 -.0326 .9593 -.0336 .9511 -.0283 .9560 -.0325 .9590 -.0334 .9511 -.0284 .9560 -.0326 .9590 -.0334 .9591 -.0334 99 up. come 55:35 32% «SE @558 m5 .8 859m>c8 8 53:3 .833 - .a 82m: .2 he .03 "a N: "as ._ £23582 25:2 22:58 6 32:52 as a: as a: as as a 8 8 a o . _ _ _ so _ _ _ _ _ _ _ 2 n 2265:: _o .8532 .IIAY LCHI Immo. 13a. mm u 235:: Lo .2522 1'0 Av/O 1 8o. immo. O c C )0, mm u 235:: _o 8:822 :V 1 2:65;: _o .8252 .1 8o. a» u! (1x) Jone; Ansuaw! ssaus apow fiuguado 100 .2 . : .0? .0 N: . a: ._ . as 9:82 222 as .568 £225 3g... 25:. 2:5... 2: 3 859228 8 €sz _833 - .2 2:5 22:58 8 285:2 as a: s: a: a: 82 8 s s a sac. _ _ _ _ _ W a _ _ _ .1 ”8.- m 2 n 2322: _o 355:2 II? M. 5 I as. m. a m. I. 08.- .15.. m. w L E..- w. m... D. I04 1 mm..- m mm . 2255.2: 3 285:2 IO OIIIIIOI m. I B... My Illol Jo 0' mm 1 2255:: ac 285:2 S n 2265:: _o 288:2 vmo .- 101 ooo. ooH.H Nao.- ooo.: NNo.- NmH.H ooo.- oaH.H ooo.- NNN.H HoN.- ooo.: e ooo. omH.H Nao.- ooH.H “No.1 NoH.H ooo.- moH.H ooo.- HN~.H mo~.- ono.~ oo.ooH ooo. MMH.H Nao.- oma.fl amo.- HoH.H ooo.- omH.H Hoo.- ~H~.H ooo.- o~o.~ oo.o: ooo. mma.H Nao.u NMH.H oNo.u ‘oma.H omo.a omH.H ooo.- oo~.a ooo.- oao.a oo.o ooo. HmH.H Nao.u NMH.H omo.u oma.H Hmo.u mma.a ooo.- oom.a ~H~.- ~oo.: oo.~ ooo. oma.a «Ho.- nos.a oNo.u HoH.H ooo.- Ana.a ooo.- ~o~.a NH~.1 oom.a oo.H ooo. omH.H Nao.u oma.a omo.u osH.H omo.n oma.a ooo.- ooo.: ooo.- Hom.H ma. ooo. Hoa.a «Ho.u moH.H ooo.- ooH.H ooo.- ooH.H Noa.u oH~.H mam.u ooo.: on. ooo. mmH.H oao.- hmH.H ~oo.u ~oa.a ooo.- oaH.H ooa.u o-.a om~.- oom.H AN. ooo. m-.H “Ho.u o-.H ooo.- -~.H o~o.- Hm~.a omH.n o-.H o-.u mHo.H mo. mmo.- Hoo.a ooo.- ooo.: Nao.u ooo.: ona.n omo.H oo~.1 aoo.a ooo.- oma.H Hoo. .mmx .mm\ .mm\. .mmx .mm\ Jme .mmx .mmx .mmx .mmx .mm\ .mmx NH; -_mmmm me He me He me :2 me fix me He Ne HM Na ans: Him H u wmm. o. u.wwm o. n www. o. a WWW N. a WWW no. ".mwm .aa m. u a m m m m m m u A: u a .ono u o .m\~ u 0\m .H u u\:o wowvmoH wafimcmu you mucuomM huamcmuaw amouum Gnu co mowuumnoun Hmfiumuwa mo muoowwm .m manma 2_- \J'I 2“ 4:. Opening mode stress «intensity factor (K1) in in CI 102 Asymptotes: K1 " 1. 227 Kl ' 1.144 E —2—2-- o 05 E e 11 l E22 - 0.2 E11 # l. 2 "" . E rlsotroplc 22: l 0 I material E ' ll 1,1L____ I L ‘ l l .01 .1 10 100 l Ell/E22 £11-, 2 2612 12 Figure 11. - Effects of material properties on the opening mode stress intensity factor for tensile loading (NC = 1, alc = 2/3, 6 = 45°, Tt = 1). Sliding mode stressftensity factor (K2) in in :4— 105 Asymptotes: K2 = '0. 201 K2 = -0.090 K2 = '0. 012 .l l'-—' mlrfi t—IN t—iN II O N Weakly orthotropic E22 _ . — - 0.8 maternal“ E I 11 #1 0 l l l .01 . 1 1 10 100 E11’522 2 E11 _ V 2012 12 Figure 12. - Effects of material properties on the sliding mode stress intensity factor for tensile loading (NC = 1, alc = 2/3, 6 = 45°, Tt = 1). 104 o emu. qmo.1 Nmo.H wmo.l wmo.H ooo.: omH.H qu.| ONH.N qu.| Noe.m .W 0 «ms. Nmo.| wmm. Nmo.| Hum. mmo.| mmm. NOH.I qu.H mwN.| omq.N .M o mHN. mHo.I CNN. mHo.I mum. mHo.I me. ooo.: mNo.H H¢H.I wee.H .% o HE. ooo.- m2. ooo.- m2. ooo.- 92. do... o3. ooo.- moo; m C Now. 0 men. 0 mom. 0 NON. 0 Now. 0 Now. 0 are. No. me We .6. we. .2 We. .2 as. me. we. a m. NM HM NM HM NM HM NM HM NM HM NM HM MMHIU «lo NIU HIU NIU QIU 1: In I: 1: Hug Hug. .53."... AH n uH . mq ".mv wanmoH oHHmamu pom muouomm huHmcoucH wmouum co oNHm oumHe was apUHB opmHm mo muommmm .o «23 Opening mode stress intensity factor (K1) in y/I—n 105 2'2"— h/c=0.5 20*— 1.8— 1.6— 1.42- 1.2F' h/C=1 1.0— hlc= 2 .8— h/c=133 .6 l 1 l 0 1/6 113 1l2 2/3 Crack length to plate width ratio, a/c Figure 13. - Effects of plate width and plate size on the opening mode stress intensity factor for tensile loading (6 = 45°, Tt = 1). 106 .2 n J .omo . 9 268. 222 Loo 388 £235 385 $2: 9:23 of co one. ago new £23 3% .o moot“ - .3 2:5 ea 5:2 £2; 2% 2 £32 :3 EN N: mm u 2: MuoE mduoE 2. 3.- mo .- NH. 2!} u! (3)1).101321 linsuaiu! ssaus apow fiugpus 107 ooo. Nmm.fi ooo. mmo.H ooo. cow. ooo. “mu. 0 “oh. om ooo.- Nmm.H goo.u mmo.a Noo.u 0mm. Hoo.s Rm“. 0 “ON. mm mHo.u omm.H oao.u mmo.a ooo.: mwm. ~oo.s emu. o No“. mm mqo.u mo~.H «mo.- mmo.H mHo.n Hum. moo.n own. o nob. o“ o~o.u New. ooo.: mas. 0 men. we mmo.u «mo.H H~o.u om». ooo.- sq“. o No“. qo mmo.u mmo.H o No“. No mmo.u NHo.H o~o.u omm. ooo.: «ch. 0 no“. 00 «00.- qm~.H mmo.u Hoo.H @Ho.u Hem. moo.u Hen. o No“. on “mo.u mmm. wao.u mmm. o “oh. Nm moo.n mo~.H o No“. on «00.: owH.H mmo.I mmm. mHo.I mam. voo.l mmn. 0 new. we Noo.u mmH.H o “on. oe mmo.: mma.~ o “on. on mmo.n Noo.H oao.n «mm. moo.n ooh. Hoo.u mm“. o “oh. om oao.n ~mo.~ Qoo.u mfim. Hoo.u «m5. ooo. “N“. o No“. m «00.: omo.H ~oo.u ohm. ooo. mwh. ooo. “N“. o No“. N ooo. cmo.H ooo. ohm. ooo. man. ooo. NNN. o No“. o .mm\ .mm\ .mm\ .mm\ .mm\ .mm\ .mm\ .mmx .mm\ .mm\ mwn NM HM Nx fix my . HM Nx ax ax fix o M U N U M U 0 U U Mum mum mum mum cum 5“?» u flfl a a .H u u\nv wcfivmoa mawmcmu you mucuomm muwmcwucw $33 05 no 3me :Oaumucwwuo vcm 56H: 3me mo 303mm .m «dame Opening mode stress intensity factor (K1) in {In 1.4 ms a/c - 2/3 alc-1I3 '8_ ale-116 '7 wc-O .6 l l J I 0 20 40 60 80 Orientation angle (6) in degrees Figure 15. - Effects of orientation angle and plate width on the opening mode stress intensity factor for tensile loading (h/C ‘ 1, Tt ' 1). 109 -.07 — a/c = 28 Li ‘ 06 r- .9 AN if -.05 —- E E. _ 04 __ alc = 1/2 § 33 g; -.03 - ‘” alc = 1/3 E -.02 D5 .E .2 7. -.01 0 20 40 60 80 Orientation angle (6) in degrees Figure 16. - Effects of orientation angle and plate width on the sliding mode stress intensity factor for tensile loading (h/C =1, Tt=1)' 110 Table 6. Stress distribution on the y = O, x > a line for tensile loading (h/c = 1, a/c = 2/3, a = 45°, Tt = 1) a = .5 in. X 3' Oxx 07y Oxy 1.016 7.33 9.59 - 482 1.046 3.54 5.53 -.227 1.076 2.37 4.27 -.131 1.106 1.74 3.58 -.069 1.136 1.33 3.12 -.021 1.166 1.04 2.77 .020 1.196 .81 2.50 .055 1.228 .63 2.26 .085 1.258 .49 2.05 .111 1.288 .37 1.86 .131 1.318 .26 1.67 .146 1.348 .18 1.50 .153 1.374 .11 1.33 .150 1.410 .06 1.17 .136 1.440 .03 1.01 .108 1.470 .00 .84 .063 1.500 .00 .68 .000 Normal stress 111 10— 1.0 1.1 1.2 1.3 1.4 1.5 x/a Figure 17. - Normal stress distribution on y = 0, x > a line for tensile loading (hlc = 1, alc = 2/3, 6 = 45°, Tt = 1). Shear stress 112 -5 l l l l l 1.0 1.1 1.2 1.3 1.4 1.5 xla Figure 18. - Shear stress distribution on y = 0, x > a line for tensile loading (hlc = 1, alc = 2/3, 6 = 45°, Tt = 1). 113 N00. 000. HHO.H NNO.H 5N0.H mm0.H mmo.H mHH.H n0q.H Hub. N 000. 000. 000. 000. 000. 000. 000. 000. NH0.1 .mmx mwm. nmm. NHO.H mN0.H 0N0.H nm0.H 000.H 0NH.H mm0.H MN. N n00.1 000.1 HH0.1 ma0.1 ma0.1 0H0.1 0H0.1 H~0.1 0m0.1 mm... n00. H00. 0H0.H 0N0.H ~m0.H H00.H 000.H MNH.H N00.H 0H0.1 0H0.1 0N0.1 0N0.1 0m0.1 mmo.1 ~m0.1 m00.1 mn0.1 mam. M00.H 0N0.H nm0.H N00.H H00.H 000.H mNH.H 00¢.H 0m0.1 mm0.1 500.1 000.1 000.1 000.1 000.1 m00.1 cmH.I mm0.H 0M0.H mm0.H 000.H «~0.H N00.H 000.H NOH.H eme.H 500.1 N00.1 m00.1 000.1 H0H.1 n0H.1 0NH.1 mma.1 mmm.1 00H.H HOH.H m0H.H 00H.H 00H.H NON.H 0HN.H NNN.H mmn.a 00H.1 NmH.1 H0H.1 00N.1 MHN.1 MNN.1 ch.1 00m.1 000.1 NV. f N H . u wanna. Nx NM NV. Nu 00H 0.N 0.H mm. on. mm. 00. ~00. wcwvmoa ummnm ~00 mucuumw zuwmcmuafi ammuum wnu so mmauumaoua Hmwumuma mo muummwm 'ty factor (K1) in Opening mode stress intens1 {1? 114 Asymptotes: K1 - -0. 148 K1 - -0. 057 K1 ' '0. 007 1;. j .l w -. 1 Weakly orthotropic 22 __ materiaI-\ ~5— - 0.8 .01 . 1 l 10 100 E11"‘322 £11-, 2 2912 12 Figure 19. - Effects of material properties on the opening mode stress intensity factor for shear loading (hlc - 1, a/c - 2/3, 6 - 45°, TS - 1). 115 1. 4 r— Asymptotes: T —=0.05 .H N :—a ._| Sliding mode stress intensity factor (K2) in ‘/1-n 1.0 Isotropic material-1’ E 0 E2? = 1.0 11 .9 l I l J .01 .1 1 10 100 E11’522 2 £11-, 2012 12 Figure 20. - Effects of material properties on the sliding mode stress intensity factor for shear loading lhlc = 1, a/c = 2/3, 6 = 45°, T.S = 1). Effects of plate width .5 1n. deg 20 34 38 42 44 45 46 48 50 54 58 62 64 66 68 70 85 88 90 K1 HE. A? OOOOOOOOOOOOOOOOOOOOCO K2 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 .707 116 Table 8. and orientation angle on the stress intensity factors for shear loading (h/c - 1, Ts - I) aim K1 .000 .000 .000 .000 -0001 -.004 -.004 -.001 .000 .000 K2 .720 .720 .720 .722 .725 .727 .727 .727 .727 .727 nlm K1 .000 .000 .000 -.002 -.005 -.007 -.007 -.008 -.009 -.010 -.011 -.012 -.012 -.011 -.004 -.002 .000 nah-4 K2 ER. .760 .760 .761 .770 .779 .781 .781 .781 .782 .782 .782 .782 .781 .781 .779 .779 .779 nlm .000 -.002 -.011 -.022 -0024 -.024 —.025 -.025 -.024 -.019 -.006 -.002 .000 .e C K2 K1 AT. 63 .839 .000 .839 -.003 .841 -.008 .861 -.037 -.064 .881 -.069 -.072 .883 -.073 .883 -.073 -.073 .882 .882 -.o72 .879 .869 -.041 .861 -.010 .861 -.004 .861 .000 hflh) h‘ F- he P0 P‘ P‘ he s: he F4 h‘ .477. .994 .994 .997 .037 .068 .072 .072 .072 .071 .070 .066 .027 .005 .005 .005 Opening mode stress intensity factor (K1) in {E l O 00 I O \I 1 O C» 1 O U1 1 C) b 1 O u: 1 O N I O ._.1 117 alc = 2/3 a/c =' 1/2 alc = 1/6 alc=0 20 40 60 80 Orientation angle (6) in degrees Figure 21. - Effects of orientation angle and plate width on the opening mode stress intensity factor for shear loading (h/c= 1, TS = 1). 118 1.1- alc = ZIB 5...: O I J T ‘f Sliding mode stress intensity factor (K2) in Vin alc = 1I3 a/c = 1/6 1 1 1 alc = 0 . . 7 _I l 1 I 0 20 40 60 80 Orientation angle (6) in degrees Figure 22. - Effects of orientation angle and plate width on the sliding mode stress intensity factor for shear loading (hlc = 1, TS - 1). 119 Table 9. Stress distribution on the y = O, x > ao line for shear loading (h/c = 1, a/c = 2/3, 6 = 45 , TS = l) a = .5 in. x 3. Oxx Oyy Oxy 1.016 -3.01 -O.60 8.77 1.046 -1.55 -0.36 5.14 1.076 -1.09 -0.28 4.02 1.106 -.085 -0.23 3.43 1.136 -0.68 -0.19 3.04 1.166 -0.56 -0.15 2.76 1.196 -0.46 -0.11 2.53 1.228 -0.37 -0.07 2.34 1.258 -0.29 -0.03 2.17 1.288 -0.23 0.03 2.01 1.318 —0.18 0.09 1.87 1.348 -0.12 0.17 1.72 1.374 -0.08 0.25 1.58 1.410 —0.04 0.35 1.44 1.440 -0.02 0.46 1.30 1.470 0.00 0.58 1.15 1.500 0.00 0.72 1.00 Normal stress 120 l l I 1.0 1.1 l. 2 l. 3 1. 4 1. 5 xla Figure 23. - Normal stress distribution on y = 0, x > a line for shear loading (hlc - l, alc - 2/3, 6 = 45°, T - 1) S o 121 Shear stress 1 l I 1 1 1.0 1.1 1.2 1.3 1.4 1.5 xla Figure 24. - Shear stress distribution on y = 0, x > a line for shear loading (hlc = 1, a/c = 2/3, 6 = 45°, T = 1). s 8. SUMMARY AND RECOMMENDATIONS The problem of multiply connected regions in plane linear aniso- tropic elasticity was discussed by both.Lekhnitskii.and Savin (refs. 4, 6), however there was no practical method of solution developed from their discussion of the problem. The motivation of the development of the mapping-collocation method.for finite.anisotropic.regions with cen— trally located traction-free cracks is due to Bowie's successful results (ref. 2) for finite isotropic regions with centrally located traction- free cracks. The method developed herein closely parallels Bowie's method for the isotropic problem, i.e.: 1. 2. The crack is mapped into the unit circle. The boundary conditions on the crack are satisfied exactly by using Muskhelishvili's.continuation principle. . A form of representation is.assumed for the unknown stress potential function. . The boundary conditions on the outer boundary are satisfied by the method of least squares boundary collocation. During the development of this method of solution for finite aniso- tropic regions, the following intermediate results were obtained: 1. The possibility of.a.disp1acement function formulation was ,.clearly demonstrated. It was shown that Marguerre's displace- ment function is directly derivable from a displacement func- tion in anisotropic elasticity. 122 123 2. For orthotropic materials two dimensionless material ratios were defined which were constructed from the engineering material constants, E11’ E22, G12 and v12. These ratios are of convenient forms when a parametric study is contemplated involving stresses and stress intensity factors. 3. Upon applying Muskhilishvili's continuation principle across the unit circle, one of the complex stress potentials were ex- pressed in the form of the other such that the zero traction conditions on the crack were satisfied. This result was di- rectly applicable in the mapping-collocation method. 4. It was shown that the infinite region results of Lekhnitskii and Savin can be obtained from considering a finite rectangular region and extending its outer boundaries to infinity. As final result, the complete formulation of the problem solution of a finite anisotropic region with a centrally located traction-free crack was obtained. The mapping-collocation method as developed herein for anisotropic regions with centrally located traction-free cracks was then applied to a large number of rectangular orthotropic plate problems. The effects of varying material properties, orientation angles, crack length to plate width and plate height to plate width ratios on the stress inten- sity factors were thoroughly studied for both tensile and shear loadings. The parameter ranges considered in the solutions of problems of rectan- gular orthotropic plates with centrally located traction-free cracks are summarized as follows: 124 1. Variation of Material Properties for both the Tension and Shear Problems; E .05$E22£l 11 E11/E22 .001 __ __ (E11 )2 -——-- v 2012 12 2-= 31 3-= 1, 6 = 45° c 3 c 2. Variation of Crack Length to Plate Width and Crack Height to Plate Width Ratios for the Tension Problem; a 2 Ofiw 00 amo maoummm wmumuou can 0cm Hmcwwfiuo ocu a“ mommocmwwum oflummao map 0cm mooamfiaaaoo afiummam msu wawumaou mcowumavo Hm mnu umcu mcfiuoz .mcowumsvm cowumLCOMmcmuu onu xwflfiaENm umLBmEom HHNB cowuwuou m 5020 .mcoum>m woumuom oLu 0cm HchmHuo mam mo mommlax mnu coozuon mawcm mcu mcwmp 0 C 0 :Hm ll E? 0 moo ommxs 132 ma 0H ma nfill lllmlll! flmmmvmmoj We saw ch Amemvmqo cat Na 1 NE c8 @3quva Fun 851 NEL M0mmmv0mo fl:.1 Ne GEN GENH AmvamNo :6. Na Na .3661 .. a .. 3 m0mmmv0mmy was 1 ma ma 1 CNS CNEN «GEN fissmmvsso ems - ma Nae - ma NaaN- aNeN Amwmvmmo GNEI NGEI ma CNS Aqmqumo Nae cNEI GNEI ma Amamvmao awe was «:8 ma .mqaqua&l F masl GNB mCI Nca 6 EN1 GEN GSNI a 8N1 133 (- 1 1 F l" 7" C34(S34) m "“ C34(S34) 1 1 = 16 c (s ) n m c (s ) L35 354 _ J k 35 354 I I C33(S33) = C33(S33) 17 It is noted here that the transformation equations are arranged in six groups. The components of each group interact with one another but are completely uncoupled from the other groups. This information is useful in studying elastic symmetries. AEOLOTROPIC MATERIALS Aeolotropic or monoclinic elastic symmetry is defined by the invar- iance of the elastic stiffnesses and compliances under the transforma- tion given by -r q l O O t = O 1 0 18 O O -1 (ref. 14) c .J This transformation specifies a material which has the xl-O—x2 plane as its plane of elastic symmetry, i.e., the properties of the material at +x3 are equal to those at -x3. In order to derive the elastic stiffness and compliance matrices for an aeolotropic material, consider the stress and the strain tensors under the above transformation. The stress tensor obeys the transformation ij = tiktjtokz 19 C 134 The strain tensor transforms in a similar manner I Eij = tiktjzeki 2° Carrying out these operations, one obtains I I I I 011 = O11 023 = ‘023 811 = E11 823 = ’523 I I I I 022 = O22 031 = ‘031 622 = 822 831 = "531 21 I I I I °33 = °33 012 = G12 633 = 833 612 = 612 Upon changing to the appropriate indices in contracted notation, Hooke's I Law is written in both.coordinate systems. For example, 0 and o l 1 become I I I I I I I 01 = 81161 + 31262 + 31383 + 31464 + 81555 + 81686 22 01 = 81151 + 8128 + S e + S s + S e + S e 23 2 13 3 l4 4 15 5 16 6 Substituting the relations 21 into 22 and comparing 22 with 23, one must take for an aeolotropic material S14 = S15 = 0' By considering all the equations of Hooke's Law, it can be shown that for an aeolotrOpic material the elastic.stiffness.and compliance matri- ces must be of the form 135 , S11 S12 S13 0 ° 516} S12 S22 S23 0 0 S26 S13 S23 S33 0 0 S36 8 = 24 o o 0 s44 845 o 0 o 0 $45 855 0 E16 S26 S36 0 0 S69 and 1 C11 C12 C13 0 0 C161 C12 C22 C23 0 0 C26 C13 C23 C33 0 0 C36 C = o 0 0 C44 C45 0 o o 0 045 055 o 1316 C26 C36 0 0 C661 Hence, it is concluded that an aeolotropic material is.characterized by 13 independent constants. It may also be noted that for the state of plane stress or plane strain, aeolotropy is the first assumption. ANISOTROPIC PLANE STRAIN The state of generalized plane strain is obtained.by assuming u1 = “1(x1’ X2) 112 = u2(xl, x2) 26 U =UO 135 , S11 S12 S13 0 ° 316' S12 S22 S23 0 0 S26 S13 S23 S33 0 0 S36 S: 24 o o o 344 845 o o 0 o 845 555 0 E16 S26 S36 0 ° 369 and 1 C11 C12 C13 0 0 C161 C12 C22 C23 0 0 C26 C13 C23 C33 0 0 C36 C: o o 0 cm4 C45 0 o 0 o 045 055 o 1316 C26 C36 0 0 C66, Hence, it is concluded that an aeolotropic material is.characterized by 13 independent constants.. It may also be noted that for the state of plane stress or plane strain, aeolotropy is the first assumption. ANISOTROPIC PLANE STRAIN The state of generalized plane strain is obtained.by assuming 1: I 1 ' “1(x1’ x2> u = u2(x1, x2) 26 U =11 O 136 where no is a constant. Upon using the infinitesimal strain- displacement relations in conjunction with the generalized plane strain assumption, one obtains 633 = s = e = O 27 When Hooke's Law for aeolotrOpic materials is applied,.one finds that the normal stress in the x3-direction can be expressed in terms of the other stress components. Consequently, it can be eliminated from Hooke's Law, and the resulting stress-strain relations can be written for the state of anisotropic plane strain as, e. = a,,o 28 1 13 j where C13° 3 a = c . —-———¥1— (i,j = 1, 2, 6) ij 1] C33 The constants, a , are called.the reduced compliance coefficients. If 13 the elastic stiffnesses are used,.Hooke's Law becomes 0 = S .s,. (i,j = 1, 2, 6) 29 i e., the elements of the elastic stiffness matrix remain unchanged. Considering both forms of the stress-strain relations, the rela- tionships among the elastic constants can immediately be obtained. .Taking = S o 30 1 = 51351 ijajkok and using Kronecker's delta, results in Sijajk = 61k 31 137 Hence, for the state of anisotropic plane strain, the expressions which relate the various elements of the reduced elastic compliance and the elastic stiffness matrices are given as follows: Elastic stiffnesses in terms of reduced elastic compliances: 2 2 _ “22“66 ' “26 “11“66 ' “16 s — s = 11 [GI 22 [a] _ “16“26 ' “66“12 _ “16“12 ' “11“26 s — s — 32 12 161 26 [01] C1 01 - 0. " 2 S = 12 26 22“16 S = “11“22 “12 16 [a] 66 [a] Reduced elastic compliances in terms of elastic stiffnesses: 82 s2 a g S22366 ‘ 26 a g S11366 ' 16 11 Is] 22 [s] a = 316326 ' 312866 a = S16512 '-S11826 12 Is] 26 Isl 33 s s - s s s s - 82 a = 12 26 22 16 a = 11 22 12 13 Is] 66 |s| In these expressions Ia] and [SI designate the determinants of the reduced elastic compliance and the elastic stiffness matrices, respectively. For a stable solid, the strain energy density function must always be positive definite for any value of the strain.(stress).components, except when all the strain (stress) components are zero. From matrix theory, the necessary and sufficient condition for the quadratic form i.e., the strain energy density function, 138 1 W = 2 aijoioj 34 or w-is es 3 2 “11 to be positive definite is that all the leading principal minors of the matrices, a and 8 .must be positive. Hence,.the.materia1.constants must satisfy the following inequalities: a a “11 “12. “16 11 12 > all 0, > 0, 012 022 “26 > 0 36 “12 “22 “16 “26 “66 and s 3 S11 S12 S16 . 11 12 $11 > O, S, S .. > O, 312 S22, 526 > O 37 12 22 S16. S26. S66 These relations in conjunction.with.relations 32 and 33 1mp1y that the diagona1.elements are all positive. ANISOTROPIC PLANE STRESS .The state of anisotrOpic plane.stress is defined.by the.following assumption: 0 = o = o = 0 38 Upon the use of these assumptions in Hooke's Law for aeolotropic.mate- rials, it is found that the normal strain in the x3-direction can be eliminated from the strain relations. The resulting.stressvstrain rela- tion for the state of anisotropic plane stress becomes 0 = B .s. 39 where S S ..=S..'—iij'3— (isj=1929 6) 13 13 S33 In this expression, the constants, B , are called the reduced elastic ij stiffnesses. When writing Hooke's Law in terms of the elastic compli- ances, one obtains ei_= Cijoj (1,3 = 1, 2, 6) 40 Hence, in the case of anisotropic plane stress, the elastic compliances are not reduced. The relationships among.the reduced elastic stiffnesses and the elastic compliances are given as follows: Reduced elastic stiffnesses in terms of elastic compliances: 2 2 8 = “22“66 ' “26 B = “11“55 ‘ “16 11 |CT 22 |c[ e = “16“26 ‘ “12°66 B = “16“12 ’ “11“26 41 12 [cl 26 ICT c c - c c c c - 02 B = 12 26 22 16 B = 11 22 12 16 [cl 66 [cf Elastic compliances in terms of reduced elastic stiffnesses: 140 2 2 C = “22“66 ’ “26 C = “11“66 ‘ “16 11 lel 22 [BI c = “16“26 ' “12“66 c = “16“12 ' “11“26 42 12 [BI 26 [8T 8 B - B B 8 B - 82 C = 12 26 22 16 c = 11 22 12 v ,, 16 [BI 66 [8] The consideration of the positive definiteness of the strain energy density function results in inequalities similar to those obtained for the state of planestrain. ' B B “11 “12 “16 ll 12 511 > 0» B B > 0. 812 822 826 > 0 43 12 22 “16 “26 “66 and C C “11 “12 “16 ll 12 C11 > 0, C C > 0, C12 C22 C26 > O 44 12 22 “16 “26 “66 In this case, it may also be concluded that the diagonal elements of the reduced elastic stiffness and the elastic compliance matrices are positive. ENGINEERING CONSTANTS From simple tensile and torsion tests, one can obtain constants which can readily be related to the elastic stiffnesses and compliances. These constants, customarily.ca11ed engineering constants,.are given in ref. 4, and they are reproduced here for reasons of completeness. 141 For an aeolotropic material, the elastic compliances are expressed in terms of engineering constants as follows: C 1.1. C =_"_12_=_21 C =_1_ 11 E11 12 E11 E22 44 023 V V “22 = El_' “23 = ' Egg_= ' 872' “55 ='El" “5 22 22 33 13 c ._1_ C “13.11-32 C ._1_ 33 E33 13 E33 E11 66 012 “31 12 “23 31 n12 1 n1 12 “45= G ‘ G C16=E ‘0 23 13 11 12 T) I) n n C 6 = 12,3 = 3,12 026 3 12,2 = 2,12 3 E33 G12 E22 G12 where: I . E11, E22 and E33 are Young s moduli of elastic1ty with respect to the coordinate directions. 623, G13 and C12 are the shear moduli.for planes which are par- allel to the coordinate planes,:mf0-x3, xl-O-x3 and xl-O-xz. 021, v32, 013, v12, v23 and v31 are the Poisson ratios which characterize contraction in the direction given by the second sub- script due to extension in the direction given by the first subscript. The constants, u31,12 and u23,31, are called the coefficients of Chentsov. These constants characterize shear in the coordinate planes. For example, characterizes shear in the planes “31,12 142 parallel to the x1 and x2 coordinates inducing shear stresses in planes parallel to the x3 and x1 coordinates. 023 1, etc., are called the coefficients of mutual influence of the first kind and n 3, etc., are called the coefficients of 1,2 mutual influence of the second kind. The mutual influence coef- ficients of the first kind characterize stretching in the direc- tions parallel to the x coordinte, etc., induced by shear 1 stresses in the x2-0—x3 plane, etc. The mutual influence coef- ficients of the second kind characterize the shears in the planes parallel to the x2 and x3 coordintes, etc., which are due to normal stresses along the x1 coordinte, etc.. Using these engineering constants, the elastic compliance and the re- duced elastic compliance matrices for the anisotropic plane stress and plante strain states become 1 ' ”13V31 _ v21 + v31"23 n12,1 + v13”12,3 E11 E22 “11 a = _ V12 + V13V32 1 ' v32°23 n12,2 + v23“12,3 E11 E22 E22 n12,1'+"13"12,3 n12,2‘+"23”12,3 1 ’ n3,12”12,3 __ E11 E22 G12 46 143 and F’ _. 1 _ V12 n12,1 11 E11 “11 V n C: _E2_1 _1_ 7152.2 ,7 22 22 22 n n ' 1,12 2,12 1 “12 “12 “12 The corresponding reduced elastic stiffness and elastic stiffness matrices can readily be obtained from relations 32 and 41. INVARIANTS For a generally anisotropic material there are a number of invar- iants associated with the compliance and.stiffness matrices with respect to a rotation about the x3-axis (ref. 14). The proofs of these invari- ants follow directly from the transformation equations 12, 13, 14, 15, 16 and 17. 'These invariants can, of course, be specialiZed for the various cases discussed by taking the appropriate terms as zero. 11 = C11 + C22 + 2C12 I2 = “66 ' 4“12 I3 = “44 + “55 48 I4 = “13 + “23 2 2 I5 ‘ “34 + “35 I6 = “33 144 and J1 = 811 + $22 + 2812 “2 = “66 ’ “12 J = S + S 3 44 55 49 J4 = “23 + “13 2 2 “5 = “34 + “35 “6 = “33 It may also be noted that the compressibility of an anisotropic material given by 1 + I6 + 21 5“ is also invariant. ORTHOTROPIC MATERIALS An anisotropic materia1.possessing two.perpendicu1ar.planes of elastic symmetry, say, x3 = 0 and. x = 0, wi11.automatically have the 1 x2 = 0 plane also as its plane of elastic symmetry. Such a material is called orthotropic;.characterized by 9 independent material constants. The compliance and stiffness matrices.formed by these constants must be of the form 145 F’ _. 011 012 013 0 0 0 “12 “22 “23 “ “ “ c c c 0 0 0 c = 13 23 33 51 0 0 0 044 0 0 o o 0 0 055 0 0 o 0 0 0 C66 and r— __ 811 812 $13 0 0 0 812 322 823 0 o o 813 823 s33 0 0 0 s = 52 0 0 0 844 o 0 0 0 0 o 355 0 0 0 0 o o s _ 69, ORTHOTROPIC PLANE STRAIN The consideration of an orthotropic material in the state of plane strain results in “16 = “26 = “36 = “45 = “16 = “26 = “36 = “45 = “ “3 consequently a = “26 = 0 54 Then the reduced elastic compliance matrix, a, and the elastic stiffness matrix, S, will be of the form 146 r“ 2 C _ “13 C _ “13“23 0‘“ 11 C33 12 C33 2 C C C “ = “12 ‘ 63 23 “22 "E“1 “ 5“ 33 33 L. O 0 C“§ and r- 1 “11 “12 “ S = $12 $22 0 56 0 0 S66 In terms of engineering constants, these two matrices become F 1 ‘ “13“31 “21 + “31“23 E E 0 11 22 _ “12 + “13“32 1 ' “32“23 a - E E 0 57 ll 22 o o G1 1 12. and KVE11(1 ‘ “23“32) Kv“11(“21 + “23“31) “ “ = KvE22(“12 + “13“32) KvE22(1 ‘ “31“13) ° 5“ O 0 G e 1% where K = (1 — v v - v v - v v ~ 2v v v )-1 v 12 21 23 32 31 13 12 23 31 ~147 ORTHOTROPIC PLANE STRESS For the state of plane stress, the reduced elastic stiffness and the elastic compliance matrices satisfy the following forms: P 82 S S q S _ 13 S _ 13 23 0 11 S33 12 833 2 S S S e = s12 - —l§—32- s22 — gfié- o 59 33 33 O 0 S L 66- and ’- 1 C11 C12 0 C = 012 C22 0 6O .9 0 C663 These matrices can be given in terms of engineering constants as, “' “11 . “21“11 o 1 1 ’ “12“21 1 ' “12“21 “12“22 “22 B = 1 - v v 1 - v v 0 61 12 21 12 21 L O 0 G12- and 148 P ‘1 1 “12 — _ _E_._. O 11 11 \) c = - E21. E%—- o 62 22 22 o o '6“— L 124 GENERAL ORTHOTROPY IN PLANE PROBLEMS In orthotropic problems, the consideration of different material and reference axes forces one to use the equations transforming the material constants from the material axis system to the reference co- ordinate system. The material constants are customarily defined with respect to the material axes; hence, the elements of the elastic stiff- ness and compliance matrices will vary in accordance with the prescribed transformations between the two coordinate systems. In the case of plane orthotropy, one finds that the elastic stiff— ness and compliance matrices will become full matrices upon transforma- tion; however, the number of independent material constants remains un- changed. In the literature, one speaks of general orthotropy when the reference axes do not coincide with the axes of orthotropy. The transformation equations which apply.for plane orthotropy are given in group 12 with the specialization “16 = “26 = “16 = “26 = “ Consequently, the elements of the elastic compliance matrix can easily be obtained for plane orthotropy as, 149 ' C0845 sin46 sin226 (1 2\)12) C = + + G --——-— “1 “11 “22 4 12 “11 c' g _ “21 + sin226 <“ + “12 + 1 + “21 _ 1 > 1“ “22 4 E11 “22 “12 C) II + + - ' C0845 s1n48 sin226 ( 1 2“21) c 22 E22 E11 4 12 E22 63 ' 20 2 C66 = 311,225 (EL + El— + E21) + co: 26 11 22 22 12 F 2 2 ' _ sin 6 cos 6 l v21 Cl6 - sin 26 E - E + (“G - E “)cos 26 22 11 12 22 b P . , 2 2 v 026 = sin 28 £§§—§-- 5%3-§-- §é¥—-- E2§)cos 28 22 11 12 2 b where 6 is the angle between the material axes of orthotropy and the reference coordinate system. The elastic stiffness, the reduced elas- tic compliance and the reduced elastic stiffness matrices can then be constructed in accordance with the appropriate relations given before. B. COMPUTER PROGRAM The computer program whose listing is reproduced on the following pages, will perform the required operations in the following order: 1. The coordinates of the collocation points are calculated. 2. The elements of the load vector, K, are computed for each collocation point. 3. Lekhnitskii's complex material parameters are obtained. 4. The elements of the coefficient matrix, at are calculated and the multiplication, 6T6, is accomplished. 5. The modified load vector,<§TK, is obtained. 6. By calling the subroutine, DGELG, the least squares solution vector, V1“, is found. The DGELG subroutine utilizes Gauss' method with complete pivoting. 7. The stress intensity factors, K1 and K2 are obtained. 8. Finally, the stress components are calculated on the boundary and on the y=0 line. The input variables are designated as follows: TX, SY one half of the constant specified boundary tensile , stress; TX on the y=h, -c < x < c and SY on the x =.i.c, 0.< y < h part of the boundary. SX, TY minus one half of the constant specified boundary shear stress; SX on the y=h, -c < x < c and TY on the x = i:c, 0 < y < h part of the boundary. 150 151 NYM, NXM, number of collocation points on the following sec~ NXP, NYP tions of the boundary: NYM for x = -c, O $.y $.h; NXM for y=h, -c < x 5.0; NXP for y=h, O < x 5_c; NYP for x=c, O < y < h. A half-length of crack C half-width of plate H halfuheight of plate DELTA orientation angle T“? r ‘1 V12 Poisson's ratio 3 ' I “11’ E22 Young s moduli of elasticity . _ 612 shear modulus la NN, NP truncation numbers of infinite sums associated with the negative and positive exponent terms of a Laurent series, respectively The results recorded in the output are the coordinates of the collo- cation points, the elements of the least squares solution vector, the stress intensity factors, K 1 and K2, the stress components on the boundary and on the y=0 line at discrete points. 0000100 0000200 0000300 0000400 0000500 0000600 0000700 0000750 0000800 0000900 0001000 0001100 0001200 0001300 0001400 0000100 0000200 0000300 0000400 0000500 0000600 0000700 0000800 0000900 0001000 0001100 0001200 0001400 0001500 0001600 0001700 0001800 0001900 0002000 0002100 0002200 0002300 0002400 0002450 0002500 0002600 0002700 0002800 0002900 0003000 0003100 0003200 0003300 0003400 0003500 0003600 0003700 0003800 0003900 0004000 0004100 0004200 0004300 0004400 0004500 0004600 0004700 0004800 0005000 0005100 0005200 0005300 0005400 0005500 0005600 0005700 0005800 0005900 0006000 0006100 0006200 0006300 0006400 0006500 0006525 600 IJSZ COMMON/DATA/YIZOO)’XIZOOI9ASIGXXIZOOI9ASIGYYIZOOIQASIGXYIZOOIOSYQSX9TY’TX9C9H9AQDELTA’ 4V129612pE229E11.HFAD(7),NYP,NYM,NXP,NXM,NN,NP REAL*8 X.Y,SY,SX,TXyTch9H9A90ELT49V1266129E22,E1l,ASIGXX.ASIGYY.ASIGXY DIMENSION ETEM(100.lOO) REAL*8 ETEM NAMELIST/OUTCON/NP,NN,NYM,NXH.NXP,NYP,- i“TXvTszxosY9C1H,A;DELTA,V12,GIZ,1322,E11 READISyDUTCON) N4=NP+NN+1 WRITE169600) FORMATIIHI) HRITEI7,0UTCON) CALL CALCN(ETEM9N4) STOP END SUBROUTINE CALCNIETFMvNQXI DIMENSION ETEMH(1009100)9XKMI20079AUXI100) DIMENSION ETEMINQXoNQX) DIHFNSION EMIZOOQIOOIOCMIIOOIQKMIZOOIQZCONI67’ZC0NN1679 *ETKMIIOOI9ZCMIIOOIQZCMZIIOOI COMMON/DATA/YIZOOI9XI200)’ASIGXXIZOOIQASIGYYIZOOI'ASIGXY‘ZOO)OSYOSXOTYQTXQC,HQA'DELTAQ 5 VIZQGIZVEZZOEII'HEADI?)9NNYP’NNYMQNNXP'NNXM’NN'NP REAL*8 F11.522.Gl2.VIZ.DELTA.A.A2,H.C.X.Y.TX.TY.SX,SY.ETEMH, EM,XKM,CM,KM,AUX,ETKMyTEMloTFM29C229C119C669C129N0AN9R905LYM, ”FLXMgoELXPgDELYPoDELX9DELTARQDR,TEMCyTEMS'RKoEPSQCMKoCMKpl9 CMKPZ9CMKP3QDELY9XXQYYQSIGXX9SIGXY95IGYY90REAL90AIMAGpETEMgSUM *9ALPHAZ,BEIAZQKngZyADR'RECCCSQASIGXX,ASIGYYgASIGXY COMPLEX*16 U10.U20,Ul,U29U18.UZB.UIDB,U2089IMAG.F1N,F2N.GIN.GZN.ZI- *922921892289UC0NI,UC0N29UCON39UC0N4QZCONIyZCDNyZTEMyZCUNN92CM'ZCH2- #ySIGXXS’SIGYYS,SIGXYSQZTEMBQTERMI,TERM29UC0N59UISOOUZSOQZSTRS' #QUZRMUIVUIMU29U2MU239HIlNyHI2N9XJ1N9XJ2N9UCON69UC0N79UCON8 *9SQTZI9SQTZZQSOTZZBQZCONNIozCUNNZyZCUNN3oZCUNN49ZCUNN592C0NN69CSUH- *yFllocll9CC9ZC09KSUMNP9KSUMNN92109220922809SXXNNQSYYNN9$XYNNQSXXNP- *QSYYNP'SXYprTEMICQTEMZC DATA OR /.0174532925 / DATA IMAG/(OovloI/ NAMELIST/RUGO/N4gNPvNNvMoHZvaH,AQRQDELTA9V1206129FZZQEIIv- *TX.TY,SX.SY.C11.C12,C226C66,U10,UIDB,UZD.U2DB.U1.UlB.U2.UZB,- *UIMUZyUZMUZByUZBMUI9UCON29UCON3yUCON49UC0N50UCON69UC0N79CC NAMELIST/RUGI/ZIvZZyZTEMyZCON NAMELIST/BUcleIN9F2N961N9GZN92CUNN ITRIG=O IUTRIG=O N4=N4X N:NN NYM=NNYM NXP=NNXP NXM=NNXM NYP=NNYP AZ = A¢A N15 NYM +1 N25 NIS+NXM N35 = NZS+NXP DELYM ' H/FLOATINYM) DELXM C/FLOATINXMI DELXP C/FLOATINXP) DELYP = H/FLOAT(NYP I XIII = ‘C VIII = 0.0+0 M=NYM 00 1 1:29" XII): 'C YIII = VII-17+DELYM YINIS) =H XINISI = ‘C NSP1=NIS+I M = M+NXM DO 2 I=NSP19M VIII: H XIII = XII'II+DELXM YINZS) = H XINZS) 3 0.0+O NSPI = NZS+I M = M+NXP 00 3 I=NSP19M YIII=H XII): XII-II YIN3$I 3 H XIN3S) = C IFIIDELThoGEo #560 + DELXP 1.).AND.(DELTA .LE.89.)) GO TO 350 0006550 0006575 0006600 0006700 0006800 0006900 0007000 0007100 0007200 0007300 0007400 0007500 0007600 0007700 0007800 0007900 0008000 0008100 0008200 0008300 0008400 0008500 0008600 0008700 0008800 0008900 0009100 0009200 0009300 0009400 0009500 0009600 0009700 0009800 0009900 0010000 0010100 0010200 0010300 0010400 0010500 0010600 0010700 0010800 0010900 0011000 0011100 0011200 0011400 0011500 0011600 0011700 0011800 0011900 0012000 0012100 0012200 0012300 0012400 0012500 0012600 0012700 0012800 0012900 0013000 0013100 0013200 0013300 0013400 0013500 0013600 0013700 0013800 0013900 0013950 0014000 0014100 0014200 0014300 0014400 0014500 0014600 0014700 0014800 0014900 0015000 350 15 153 0ELYP=H/FLDAT(NYP+1) Y(N3$)=H—0ELYP NSP1=M38+1 M+NYP 00 4 I=NSP1.M Y(I)= Y(1—1)-DELYP X(I1 = C CONTINUE J=1 DU 5 1:19M KM(J)=TX*X(I)+TY*Y(I) J=J+1 KMIJ)=SX*X(I)+SY*Y(I) J=J+1 M2 = 2*M IF(ITRIG.E0.1) C11= 1./E11 C12=-V12/F11 C22= 1./F?2 C66 = 1./012 R=DSORT(C*$2+H$$2) H = 60 T0 200 TEMIC = (C66+2.*C121/(2.*C11) TEMZC =CDSHRT(TEM1C**2-C22/C11) U10: IMAG* CDSORT(TEM1C+TEM2C) H108: -U10 U20:IMAG$CDSQRTITFM1C-TEM2C1 UZDRz -U20 DELTAR= 0FLTA*0R TEVC = DCUSIDELTAR) TENS = DSIN(0ELTAR) U1 = (UlnfiTFHC - TEMS1/(TEMC+UID*TEMS) U150: UI*U1 U2 = (U20*TFMC - TEMS1/(TEMC+UZD*TEMS) “250: U2*U2 UlB=0CONJG(Ul) 028=DCUNJG(HZ) 02M028=U2-U28 UCONI = (028-011/(U2MU281 U28MUI=U28-Ul U1M02=Ul-UZ UC002=U1*02MU28 UCON3=UZ$UZRHUI UCON4=U28*U1MU2 UC0N5=UC0N2*U1 UC006=UCON3*02 UCON7=UCON4*028 HCON8= (UPB-U181/(02MU28) CC=UCON5+UCON6+UCON7 HRITE1698000) J=1 00 25 1:1,” Z1 = X(I) + U1*Y(11 22= X11) + U2*Y(I) ZlB=DCONJG(Zl) ZZR=0CONJG(ZZ) SQTZ1= CDSQRT(Zl—A)*COSORT(ZI+A1 ZCON(1)= (21+SQTZI) ZCON(4)= (Zl-SOTZI) SOT22= CDSORT(ZZ-A)*CDSORT(ZZ+A) ZCON(2)= (Z2+SOT22) ZCON(5)= (Z2-50TZ?) SQTZ28= DCONJG( SOTZZ) ZCON(3)= (ZZB-SOTZZB) ZCON(6)= (228+SOT228) K=4 JP] = J+l F11=(U2MH?B*SOTZI+U2RMUI*SOTZ2-U1MU2*SOT228)/R Gl1=(UCON2¢SOTZ1+UCON3*SOTZZ-UCON4*SOTZ28)/R FMIJ,1)=0.D+0 EM(JP1.11=Y(1)/R FM1J921=0REALIF111 EM(JP1.2)=DREAL(GII) EM(J.3)=-0AIMAG(F11) EM(JP1931=—DAIMAG(GII) IFINP.EO.I) 00 TO 17 00 15 KK=39N992 00 16 JJ=1.3 ZCDNN(JJ)=(ZCON(JJ1/R)$tKK FIN = UZMHZR*ZCONN(1) +028MUI*ZCONN(2) GlN=UC0N29ZCONN(1) + UCON3*ZCONN(2) EMIJ0K1= DRFALIFIN) EMIJP19K1= DREALIGIN) K:K+1 FWIJ9K1 = ‘0AIMAGIF1N1 EMIJP1.K)= ‘DAIMAGIGIN1 K=K+1 +UIMUZ #ZCDNN( 31 + UcnmgmzCHNN(3) 0015050 0015100 0015200 0015300 0015400 0015500 0015600 0015700 0015800 0015900 0016000 0016100 0016200 0016400 0016450 0016475 0016500 0016510 0016520 0016550 0016570 0016575 0016600 0016700 0016800 0016900 0017000 0017100 0017150 0017200 0017300 0017400 0017500 0017600 0017700 0017800 0017900 0018000 0018100 0018200 0018300 0018400 0018500 0018600 0018700 0018800 0018900 0019000 0019100 0019200 0019275 0019300 0019500 0019600 0019700 0019800 0019900 0020000 0020100 0020200 0020300 0020400 0020500 0020600 0020700 0020800 0020900 0021000 0021100 0021200 0021300 0021400 0021500 0021600 0021700 0021725 0021750 0021800 0021900 0022000 0022100 0022200 0022300 0022400 0022500 0022525 0022550 1016 1015 6501 6497 6498 6001 26 265 27 275 6502 390 6503 999 6512 60 640 997 9975 6508 998 9985 6509 6511 30 15M: IFINN.EO.1)GO T0 25 DO 1015 KK=3,NN,2 DO 1016 JJ=4v6 ZCONN(JJ)=(ZCONIJJ1/R)**KK FZN = U2MU28¢ZCONNI41 +U28MU1¢ZCONNI51 +U1MU2 I“ZCUNNI61 GZN= UCON2*ZCONN(4)+UCON3*ZCONN(51 + UCON4*ZCONN(6) EMIJyK) = DREAL(F2N) EMIJPloK) = DREALIGZN) K=K+1 EMlJvK1 = v DAIMAGIF2N1 EMIJp19K1 8 -DA1MAG(GZN) K=K+1 J=J+2 FORMATIIHL. 8HK VECTOR IlHK) HRITE(6,64971 FORMATIIHLv'X ARRAY' llHK) NRITEI6,6001) (X(1101'19M1 HRITE1696498) FORMATIIHLo'Y ARRAY'llHK) HRITEIbvéOOl) (Y11101319H1 HRITEt6.65011 , HRITE(6.6001) (KM‘11’I'19M2) F0RMAT11X18615.71 DO 265 J319N4 DO 265 13J9N4 SUM 3 00 DO 26 K=19M2 SUM = SUM + EMIK911 * EMIKyJ1 ETEMIJ,I)=SUM ETEM(I,J)=SUM DO 275 I=qu4 SUM=00 DO 27 J=loM2 SUM: SUM+ EH1J11)¢KM(J1 ETKM11135UM IFIIOTRIGoE0.0) GO TO 6512 HRITEIbvbSOZ) FORMAT (1H1927H E TRANSPOSE ‘ E MATRIX // ) DO 890 I=lgN4 HRITE(6960001 (ETEM119J10J'19N41 WRITE1696503) FORMAT (1HL027H E TRANSPOSE * K VECTOR // 1 HRITEI6,6000)(ETKM(I1.1-1,N4) DO 999 IBIgN4 DO 999 J=1vN4 ETEMH119J1=ETEM110J1 DO 60 I=qu4 CM(1)=ETKM(1) DATA EPS/1.E-15 / CALL DGELGICMoETEM 1N4ylv59591ER) FORMATIIHLo'IER= '915) HRITE(696401 IER IFIIOTRIG.E0.0) GO TO 6511 DO 9975 l'lvMZ SUM=O. DO 997 J=11N4 SUM=SUH+EM(I,J)#CM(J) XKM111=SUM HRITEI696508) FORMAT(1HL 22H E t C VECTOR ll) HRITE(6,60001(XKM(1191=19M21 DO 9985 1319N4 SUM=OI DO 998 J'19N4 SUM=SUM+CH(J)*ETEMH(IcJ1 AUXII1=SUM HRITE16.6509) WRITE169600011AUXII191-19N4) FORMATI1H1,12HETE*C MATRIX 1 K=3 ADR=AIR KSUMNP=Oon+O RECCCS-CMII) ZCD=CM(2)+CM(3)*IMAG NRITEI69619) RECCCSvlCD HRITEI7.6191 RECCCSoZCD IF!NP.EQ.1) GO TO 32 HRITE16962O) K=4 DO 30 J'3oNP.Z KP1=K+1 ZCM(J)=CM(K1¢1MAG*CMIKP1) HRITEIbobZl) JvZCH1J) KSUMNP=KSUMNP+AOR*¢J‘FLOAT1J)*ZCM1J1 K=K+2 IFINN.EO.1) GO TO 33 “Elfin-h 0022700 0022800 0022900 0023000 0023100 0023200 0023300 0023400 0023500 0023600 0023700 0023800 0023900 0024000 0024100 0024200 0024300 0024400 0024500 0024550 0024600 0024700 0024800 0024900 0025000 0025100 0025200 0026300 0028400 0025500 0026600 0025700 0025000 0025900 0026000 0026100 0026200 0026300 0026400 0026500 0026600 0026700 0026000 0026900 0027000 0027100 0027200 0027300 0027400 0027600 0027600 0027700 0027R00 0027950 0027900 0020000 0028100 0028200 0020300 0029400 OUZRROO 0020900 0(1?()0(10 0029100 0029200 0029300 0029400 0029450 0029800 0029600 0029700 0030100 0030200 0030300 0030400 0030500 0030600 0030700 0030800 0030900 0031000 0031100 0031200 0031300 0031400 0031500 1.55 HRITE(6.622) 00 31 J=3,NN,2 KP1=K+1 ZCM?(J)=CM(K)+IMAG*CM(KPI1 HRITE(6.621) J.ZCM2(J) KSHMNN=KSUHNN+AOR#*J*FLOAT(J)*ZCM2(J) 31 K=K+2 620 FORMATtlHL/6x.'N'.l6X,'Cl') 621 FORMAT11X.12,7X.2615.7) 622 FORMAT(1HL/ 6X.'N',16X,'C2') 619 FORMAT(1H1/1HLy4X.'REfiL(CC*CS1=',GlS.7,4X,'CD='.2015.71 33 ALPHA2=0RFAL(U21 BETA2=0AIMAO(U2) TEM1=2./050RT(A) ZTEM=U1MH2*H2MUZR*(ADR*ZC0+KSUNNP-KSUHNN) TEMS=DAIMA01ZTF01 Kl=-TEM1/RETA2*TFMS K2=-TEMI*DKFAL(ZTEM)+TEMI*ALPHA2/RETA2*TFMS HRITE16,644)K1.K2 HRITE(7,644)K1.K2 644 FORMAT(1HL.'K1='.015.7,'K2='9615.7/(1X.IS.2GIS.7)1 DATA NYMS,NXMS,NXPS,NYPS/50,50,50,50/ NYH=NYMS NXM=MXMS NXP=NXPS NYP=NYPS ITRIG =1 00 T0 100 200 wRITE(6,67S) 625 FORMAT(1H1/1HL,19X.'Z'918X,'SIGMA XX',8X,'SICMA YY',RX,'SIRMA XY'/- *l) XX=X(I) YY=Y(I) 21 = XX+ 01*YY 22 = XX+UZ*YY Z18=0CO0JG(ZI) ZZR=0CONJ01221 SOTZI= C050RT(Zl-A)*CDSORT(ZI+A) ZC00(1)= (21+80T71) ZC0014)= (Zl-SOTZI) SOT72: CDSORT(22-A1*C050RT(22+A) 2C00121= (22+SOTZ2) 2CON151= (ZP-SflTZZ) SOT7PH= nCUNJG( SOTZ2) ZCON(3)= (Z20-SOT2201 ZC00(6)= (220+SOTZZB) SXXNN=0.0+O SXYNN=0.0+0 SYYMN=0.0+O SXXNP=0.0+0 SXYNP=0.0+0 SYYNP=0.0+0 1F(HP.E0.1) 00 T0 45 DO ‘91 J=39Np72 DU 40 JJ=193 40 ZCONH(JJ)=(ZCON(JJ)/R)** ZCONN(1)=ZCONN(1)/SOT21 ZCONN(2)=ZCONN(2)/SQTZZ ZCONH(3)=ZCONN(3)/30TZZB 2::SP = SXXNP+XJ*ZCM(J)*(UC0N5*ZC0NN(1)+UC0N6*ZC0NM(2)-UC0N7*ZC0“N- m 53;;P = SYYNP+XJ*ZCM(J1*(U2MU28*ZCONN(11+U28MU1*ZCONN(71-U1MU2¥7CO- *4 41 22;33)= SXYNP+XJ*ZCM(J)¢(UC0N2*ZC0NN(11+UC0N3$ZC0NN(2)-UC0N4¢7C00V- 45 “12:0N.E0.1) 60 T0 46 00 43 J=3,NN,2 00 42 JJ=4,6 42 ZCOHN(JJ)=(ZCOH(JJ1/R)**J ZCONM¢41=ZCONN(4)/30TZI zc00N(5)=zanN(S)/SQTZZ ZCONN(6)=ZCONN(61/SOT228 XJ=J SXXNN = SXXNN+XJ*ZCM2(J)*(UCON5*ZCONN(41+UCON6*ZCONN(51-HCON7*ZC00- *N(6)1 SYYNN = SYYNN+XJ$ZCM2(J)*(02M028*ZCO0N(41+HZRMUI*ZCONN(51-01002t7c- *0NN(6)) 43 SXYNN = SXYNN+XJ*ZCM2(J)*(UCON2*ZCDNN(4)+UCON3*ZCONN(51-HCOH4*7CON‘ *N(6)) 46 210:21/SOTZI Z20=ZZISOT22 2280:228/50Tl28 ZTEM=ZCD*(UCON5$ZlD+UCON6¢ZZD~UCUN7*ZZRD)IR SIGXX=2.*(RFCCCS/R+0REAL(ZTEM+SXXNP-SXXNN)) " I 1.. 0031600 0031700 0031800 0031900 0032000 0032016 0032032 0032040 0032100 0032200 0032300 0032400 0032500 0032600 0032700 0032800 0032900 0033000 0033100 0033200 0033300 0033400 0000100 0000200 0000300 0000400 0000500 0000400 0000700 0000800 0000900 0001000 0001100 0001200 0001300 0001400 0001500 0001600 0001700 0001800 0001900 0002000 0002100 0002700 0002300 0002400 0002500 0002000 0002700 0002000 0002900 0003000 0003100 0003200 0003300 0003400 0003500 0003600 0003700 0003000 0003900 0004000 0004100 0004200 0004300 0004400 0004500 0004600 0004700 0004300 0004900 0005000 0005100 0005200 0005300 0005400 0005500 0005400 0005700 0005800 0005900 50 630 77 6000 (In rsfigwtfintflrahcwtfirxncwf\ntfiryhcwtinCWrrficerntfirmncwrintflr)n;fif1rafi<fir10c1r\n:wrintficwfi 156 ZTEM=ZCD*IU2MU23*Z10+UZBHU1*ZZD'U1MU2*ZZBDI/R SIGYY=2.*DREALIZTEN+SYYNP~SYYNNI ZTEM=ZCD*(UCON2*ZID+UCON3*ZZD'UCON4*ZZBDI/R SIGXY=-2.*DREALIZTEH+SXYNp-SXYNN) ASIGXXIII=SIGXX ASIGYYIII=SIGYY ASIGXYIII=SIGXY HRITEI6963OIII0XII)vYIIIoASIGXXIIIyASIGYVII)oASIGXVIII'I'IyMI FORMATIIX9I49551507I DATA NCP/IOO/ IFIITRIGoEQoZIRETURN MzNCP DEL=2.*C/FLOATINCP'II Y(1)=0. X(1)=-C 00 77 I=2.NCP YIII3OO XII)=X(I'1)+DEL ITRIG=2 GO T0 200 FORMATI 8E15o71 END OOOOOOOOOOOOOOO...00.0.0.9...OOOOOOOOOOOOOCOOOOOOOOOOOOOOOOOOOOOOODGEL SUBROUTINE DGELG PURPOSE TO SOLVE A GENERAL SYSTEM OF SIMULTANEOUS LINEAR EQUATIONS. USAGE CALL DGELG(R.A9M.N.EPS.IER) DESCRIPTION OF PARAMETERS R - DOUBLE PRECISION M BY N RIGHT HAND SIOE MATRIX (DESTROYED). ON RETURN R CONTAINS THE SOLUTIONS OF THE EQUATIONS. A - DOUBLE PRECISION M BY N COEFFICIENT MATRIX (DESTROYED). M - THE NUMBER OF EQUATIONS IN THE SYSTEM. N - THE NUMBER OF RIGHT HAND SIDE VECTORS. EPS - AN INPUT CONSTANT WHICH IS USED AS RELATIVE TOLERANCE FOR TEST ON LOSS OF SIGNIFICANCE. IER - RESULTING ERROR PARAMETER CODED AS FOLLOWS 1ER=0 - N0 ERROR. IERt-l - N0 RESULT BECAUSE OF M LESS THAN 1 OR PIVOT ELEMENT AT ANY ELIMINATION STEP EQUAL TO 00 IER=K - WARNING DUE To POSSIBLE LOSS OF SIGNIFI- CANCE INDICATED AT ELIMINATION STEP n+1, WHERE PIVOT ELEMENT WAS LESS THAN OR EQUAL TO THE INTERNAL TOLERANCE EPS TIMES ABSOLUTELY GREATEST ELEMENT OF MATRIX A. REMARKS INPUT MATRICES R AND A ARE ASSUMED TO BE STORED COLUMNWISE IN M*N RESP. M*M SUCCESSIVE STORAGE LOCATIONS. ON RETURN SOLUTION MATRIX R IS STORED COLUMNWISE TOO. THE PROCEDURE GIVES RESULTS IF THE NUMBER OF EQUATIONS M 15 GREATER THAN 0 AND PIVOT ELEMENTS AT ALL ELIMINATION STEPS ARE DIFFERENT FROM 0. HOWEVER WARNING IER-K ~ IF GIVEN - INDICATES POSSIBLE LOSS OF SIGNIFICANCE. IN CASE OF A HELL SCALED MATRIX A AND APPROPRIATE TOLERANCE EPSv IER-K MAY BE INTERPRETED THAT MATRIX A HAS THE RANK K. ND WARNING IS GIVEN IN CASE M=l. SURROUTINES AND FUNCTION SUBPROGRAMS REQUIRED NONE METHOD SOLUTION IS DONE BY MEANS OF GAUSS-ELIMINATION WITH COMPLETE PIVOTING. DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL 0....OOOOOOOOOOOOO0.0.00.9...OOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOOGEL SUBROUTINE DGELGIRcAonNvEPSoIER) DIMENSION AIIIyRIII DOUBLE PRECISION R.A.PIV.T89TOLQPIVI IFIMIZ3923TI DGEL DGEL DGEL DGEL "C’Otl~l°\fl&dflho~ pl..- 57 0006000 0006100 0006200 0006300 0006400 0006500 0006600 0006700 0006800 0006900 0007000 0007100 0007200 0007300 0007400 0007500 0007600 0007700 0007800 0007900 0008000 0008100 0008200 0008300 0008400 0008500 0008600 0008700 0000800 0008900 0009000 0009100 0009200 0009300 0009400 0009500 0009600 0009700 0009800 0009900 0010000 0010100 0010200 0010300 0010400 0010500 0010600 0010700 0010800 0010900 0011000 0011100 0011200 0011300 0011400 0011500 0011600 0011700 0011800 0011900 0012000 0012100 0012200 0012300 0012400 0012500 0012600 0012700 0012800 0012900 0013000 0013100 0013200 0013300 0013400 0013500 0013600 0013700 0013800 0013900 0014000 0014100 0014200 0014300 0014400 0014500 noon On On an (TOE-To ‘IO‘U'II‘ 14 15 16 17 18 I9 157 SEARCH FOR GREATEST ELEMENT IN IER=0 PIV = M M : I»1::=M NM=N*M 00 3 L=I.MM T8 = DABSIAILII lFITH-PIV)3.3.2 PIV=TB I=L CONTINUE TOL=EPS*PIV AII) IS PIVOT ELEMENT. MATRIX A 0.00 PIV CONTAINS THE ABSOLUTE VALUE OF A(II. START ELIMINATION LOOP LST=1 TEST ON SINGULARITY IFIPIV12392394 IFIIER)7.5.7 IFIPIV-TOL)6.6.7 IER=K-l PIVI = l.DO/A(I) J=(I-1)/M I=I-J*M-K J=J+I-K I+K IS RON-INDEX, J+K COLUMN-INDEX OF PIVOT ELEMENT PIVOT ROW REDUCTION 00 8 L=K9NM.H LL=L+I T8=PIVI*R(LL) RILLI=RIL1 RILI=TB AND ROW INTERCHANGE IN RIGHT HAND SIDE R IS ELIMINATION IFIK-M)9.18918 TERMINATED COLUMN INTERCHANGE IN MATRIX A LEND=LST+M-K IF(J)12.12;10 II=J$M OO 11 L=LSToLENO TH=AILI LL=L+II AIL)=A(LL) AILLI=TB R00 INTERCHANGE AND PIVOT RON REDUCTION 00 13 L=LST,MN,M LL=L+I TH=PIVI¢A(LL) A(LL)=A(L) A(L)=TB IN MATRIX A SAVE COLUMN INTERCHANGE INFORMATION AILSTI=J ELEMENT REDUCTION AND NEXT PIVOT SEARCH PIV = 0.00 LST=LST+I J=O DO 16 II=LST9LEND PIVI=-AIII) IST=II+M J=J+I DO 15 L=IST9NM9M LL=L-J AILI=AILI+PIVI*A(LLI T8 = OARSIAILI) IFITB-PIVIISTISTI4 PIV=TO I=L CONTINUE DO 16 L=K9NN9M LL=L+J RILLI=RILLI+PIVI*RILI LST=LST+M END OF ELIMINATION LOOP BACK SUBSTITUTION AND BACK INTERCHANGE IFIN-1173922919 IST=HH+M DGEL DGEL DGEL OCEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGEL DGFL DGEL DGEL DGEL DGEL DGEL DCFL DGEL DGEL OCFL DGEL 00FL SP. 70 71 72 73 76 77 86 R7 88 94 95 97 98 107 108 114 115 117 118 138 139 140 141 0014600 0014700 0014800 0014900 0015000 0015100 0015200 0015300 0015400 0015500 0015600 0015700 0015800 0015900 0016000 0016100 0016200 0016300 0016400 0016500 0016600 0016700 0000100 0000200 0000300 0000400 0000500 0000600 0000100 0000200 0000300 0000400 0000500 0000600 20 21 22 23 158 LST=M+1 DO 21 I=29M II=LST~I IST=IST-LST L=IST-M L = AIL) + .500 DO 21 J=II1NMgM TB=RIJI LL=J DO 20 K=IST.MM.M LL=LL+1 TB=TB-A(KI#R(LL) K=J+L RIJI=RIKI R(K1=TB RETURN ERROR RETURN IER=-l RETURN END FUNCTION DREALIZ) DIMENSION 2(2) DOUBLE PRECISION ZoDREAL DREAL=ZI11 RETURN END FUNCTION DAIMAGIZ) DIMENSION 2(2) DOUHLE PRECISION ZoDAIMAG DAIMAG=ZIZI RETURN END DGEL 160 DGEL 161 DGEL 162