BTERATWE‘SOUJTIONS 0F PLANE EMSYOSTATK PROBLEMS Thesis fax the Degree of Ph. D. MIG-1.3%!“ STATE UNWERSITY Chesfier L“ Davis 1965 {MR m3. THESIS 0-169 MIC” / 7‘" I, / I‘IIIIII‘ 'I-‘I III/I’M I II 1. II; 3 129 1 .. L ”II [III ‘- III) I / I ,' I I“ [III/A758 I ,/ I‘I/I/I/I/ .5869 This is to certify that the thesis entitled ITERATIVE SOLUTIONS OF PLANE ELASTOSTATIC PROBLEMS presented by Chester L. Davis has been accepted towards fulfillment of the requirements for __P_h_.I_)_,__ degree in Miss ”7 n (7 / ’ 'l , {ézigccoaewmxik_ Ma'or rofessor g E. Malvern V Lawr no Date November 18, 1965 an... _- ._A—. &A LIBRARY Michigan Sm: University ABSTRACT ITERATIVE SOLUTIONS OF PLANE ELASTOSTATIC PROBLEMS by Chester L. Davis Three objectives of this thesis are: to compare the efficiency of three iterative methods of solving biharmonic-finite—difference equations, to report on ISOPEP, a system of computer programs designed for the numerical solution of biharmonic plane elastostatic problems, and to demonstrate the utility of this system of programs and the advantages and limitations of numerical solutions by examples. Three matrix iterative methods considered are point successive overrelaxation, the alternating direction implicit method, and the cyclic Chebyshev semi-iterative method. These are compared in terms of computer storage and time required for the solution of a model problem. Numerical results indicate successive overrelaxation is best unless the mesh is refined so the number of points exceeds 350. Then the alternating direction implicit method is superior. ISOPEP, a system of FORTRAN II subprograms for the iterative sol- ution of plane elastostatic problems, is explained. Documentation of ISOPEP, including listings of source decks, specifications for input and the output from a sample problem, is provided. Discrete values of the stress functions and stress components for six ISOPEP problem solutions are provided. These problems include three notched tensile specimens, an infinite plate with a square hole and a semi-infinite plate with a uniformly distributed load along a Chester L. Davis portion of the edge. The numerical solutions of the six example problems indicate that a high speed digital computer with a large main memory provides an effective and economical means for the analysis of plane elastostatic problems. Good agreement as shown in the com- parison of the numerical and eXplicit stress solutions for some of these problems. Use of numerical solutions for the investigation of stress concentrations is shown in several of the examples. ITERATIVE SOLUTIONS OF PLANE ELASTOSTATIC PROBLEMS By Chester L. Davis A THESIS Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Department of Metallurgy, Mechanics and Material Science 1965 ACKNOWLEDGMENTS I wish to express my appreciation to Professor Lawrence E. Malvern for his suggestion of the problem.and his significant con- tributions to the study. Also, I wish to thank Professors William A . Bradley, George E. Mass, Charles S. Duris and Robert H. Wesserman for their service on my graduate committee. Financial support was provided by the National Science Found- ation, the Consumers Power Company, the Ford Foundation and a Michigan State University Graduats_Research Assistantship. The facilities of the Michigan State University Computer Laboratory and the University of Toledo Computation Center were used. The substantial support of these contributors is acknowledged with gratitude. I wish to thank Mr. John W} Hoffman, Director of the Engineering Research Division, for his assistance and support. To my wife I wish to express my appreciation for her aid in typing and counsel throughout this endeavor. ii TABLE OF CONTENTS Page I. INTRODUCTION . . . . . . . . . . . . . ...... . . . 1 II. THE DIFFERENTIAL EQUATIONS . . . . . .’. . . . . . . . .h A III. DERIVATION OF THE DIFFERENCE EQUATIONS . . . . . . . . . 9 IV. ITERATIVE METHODS . . . . . . . . . . . . . . . . . . . 20 A. POINT ITERATIVE METHODS B. METHODS OF BLOCK ITERATION v. COMPARISON OF THREE ITERATIVE METHODS . . . . . . .h. . 52 VI. SOLUTIONS OF PLANE ELASTOSTATIC PROBLEMS . . . . . . . . 59 VII. SUMMARY AND CONCLUSIONS . . . . . . . . . . . . . . . . 93 BIBLIOGRAPIIIY O O O O O O O O O O O O O O O O O 0 0 O O O 98 AWDICES O 0 O O O O O O O O O O O O O 0 O O O O O 0 O 101 iii TABLE 1. 3. h. 5. 7. 8. 9. 10. 11. 12. 13. 1A. LIST OF TABLES Comparison of iterative subroutine characteristics . Comparison of machine time and number of iterations required for convergence . . . . . . . . . . . . . . Summary of problems solved with ISOPEP program . . . Comparison of the ISOPEP program for the IBM 1620 and CDC 3600 computers . . . . . . . . . . . . . . . Stress function at points adjacent to the semiu Circula'r mesh 0 O O O O O O O O O O 0 O O O O O 0 O O sttress component at points adjacent to the semi: Circular nOtCh O 0 O O O O O O O O O O O O 0 O O O 0 Y—stress component at points adjacent to the semi- CirCUlar nOtCh e e e e e e e e e e e e e e e e e e e Xlwshear stress at points adjacent to the semi- circular nOtCh e e e e e e e e e e e e e e e e e e e Stress function at points adjacent to the V-notch . . X-stress component at points adjacent to the V-notch Y-stress component at points adjacent to the V-notch XY-shear stress component at points adjacent to th. V’DOtCh e e e e e e e e e e e e e e e e e e e e o Stress concentrations factors for the V-notch tCHBil. BPCCimQD e e e e e e e e e e e e e e e e e e Comparison of boundary values of C75 obtained frmm the iterative solution with those from Savin's solution . . . iv Page 52 5h 56 58 7k 7k 75 75 76 77 77 78 80 91 LIST OF FIGURES FIGURE Page 2.1 Boundary forces . . . . . . . . . . . . . . . . . . . . . 6 3.1 Rectangular mesh . . . . . . . . . . . . . . . . . . . . . 9 3.2 A regular point . . . . . . . . . . . . . . . . . . . . . 10 3.3 (a) Unequal mesh spacings (b) Subregion l . . . . . . . . 13 3.11. Dualmeshes....................... 17 3.5 Cell . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 3.6 Twelve cells for P5 . . . . . . . . . . . . . . . . . . . 19 h.1 The model problem . . . . . . . . . . . . . . . . . . . . 20 h.2 Mesh for the model problem.. . . . . . . . . . . . . . . . 20 h.3 Successive overrelaxation - Namber of iterations vs.63 .-. 30 h.h alternating direction-implicit method - Number of 1t.r‘t1°n. '391' e e s e e e e e e e s s s e e e s s e e 3“ ‘05 Two-m. block. 0 O O O O O O O O O O O O O O O O O O O O M h.6 Cyclic Chebyshev semi-iterative method-Number of itCr‘tian. 7's 9 e e e e e s e e e e e e e e s s e e e s ‘8 h.7 Comparison of iterative methods (a) Iterations vs. number of interior points . . . . . . . 50 (b) lichins time vs. number of interior points . . . . . . 51 6.1 Boundaries of the model problem.. . . . . . . . . . . . . 60 6.2 Distributions of the stress function and 6} for hobl“ 1 O O O O O O O O O O O O O O O O O O O O O O O O 62 6.3 Problem.2 - Semi-infinite plate . . . . . . . . . . . . . 63 6.4 Comparison of numerical and exact solutions for hablm z 0 O O O O O O O O O O O O O O O O O O O O O O O 6A LIST OF FIGURES FIGURE Page 6.5 Distributions of the stress function and a; for ”Oblom 2 O O C O O O O O O O O O O O O O O C O O O O O O 66 6.6 Solution domain for Problems 3, l. and 5 . . . . . . . . . 67 6.7 Distributions of the stress function and d} for PTObl'fllBeeseeeeeeeeeeseeeeeeesee 68 6.8 Distributions of the stress function and a} for Pr°blm 1‘ O O O O O O 0 9 O O O O O O O O O O O O O O O O 69 6.9 Distributions of the stress function and a": for HObICMSeseeeeeeseeeeeeesseeeees 70 6.10 Stress concentration in a semi-circular notch . . . . . . 72 6.11 Stress concentration in a Vinotch . . . . . . . . . . . . 73 6.12 Stress concentration in a rectangular notch . . . . . . . 79 6.13 Stress concentration factor It for a notched flat bar . . 81 6.1h Boundary points for the semi-circular notch . . . . . . . 83 6.15 Treatment of boundary points . . . . . . . . . . . . . . 83 6.16 Interior points near a boundary . . . . . . . . . . . . . 8h 6.1? Problem 6 - Infinite plate with a square hole . . . . . . 85 6.18 Approximation of the square hole . . . . . . . . . . . . 87 6.19 Domain of solution for Problem.6 . . . . . . . . . . . . 88 6.20 Stress function and 0; distributions for Problem 6 . . 90 vi LIST OF APPENDICES APPENDIX Page A. CONVERGENCE OF ITERATIVE METHODS . o . . . . o . . . 101 B. ISOPEP _ A FORTRAN PROGRAM FOR THE ITERATIVE SOLUTION OF PLANE ELASTOSTATIC PROBLEMS . o . . . . . 115 vii I. INTRODUCTION Determination of the stresses in a plate under conditions of either plane strain or plane stress is a fundamental problem for the structural engineer. The literature of the classical theory of elasticity includes exact solutions of numerous problems. For an account of the mathematical theory see Timoshenko and Goodier (1951) and Mbskhelishvili (1953). However, exact solutions are available only for those problems which have rather simple geometric shapes and boundary constraints. Finite-difference equations have been used as an alternative for the analysis of more complex practical problems. This dissertation compares several iterative finiteodifference solution methods with respect to their efficiency in treating plane elastostatic problems. The finiteudifference solution of a boundary-value problem is a two-step procedure. First the partial differential equation and the associated boundary conditions are replaced by difference equations which relate the discrete values of an approximating function at a finite number of points. 1 regular mesh of lines is superimposed on the domain of the given boundary-value problem. A finite system of linear equations is formed by writing a difference equation for each node of the mesh. The second step is the solution of this system.of simultaneous equations. It is often necessary to solve for values at a thousand or more points. The relaxation technique was employed by Southwell (1946) for the solution of the biharmonic difference equations. Though Southwell and his colleagues solved a number of complex engineering problems, 1 2 the relaxation method has not. been readily adapted for digital computer solutions. Iterative methods have been used extensively for computer solutions of systems of equations. These methods make repeated use of simple algorithms which at each application provide an imoved approximate solution at one or more of the mesh points. The exact solution is the limit of the sequence of the adjusted point values. Though there are many iterative methods, the three which seemed to offer most promise for the solution of the biharmonic finite-difference equations were: the technique of point successive overrelaxation introduced by Frankel (1950) and Young (1951.), the alternating-direction implicit method of Cents and Dames (1958), and the cyclic Chebyshev semi-iterative method proposed by (h-iffin and Yarga (1963). The comparison of these iterative schemes is the first objective of this dissertation. Using the computer time required to solve a given problem as the measure of the efficiency of an iterative method, the numerical results obtained indicate that for the biharmonic equation there is a critical mesh spacing 11* such that the point successive overrelaxation iterative method 1. the best of the three methods tested for h 2. 11* mm. for h< h“ the alternating-direction implicit method is best. lach of these methods uses a parameter for accelerating the con- vergence. Formulas for the bounds of optimum parameters are included for rectangular regions. For the more general case of irregular boundaries there are no convenient relationships for estimating the optimum para- meters at the start of the iterative solution. The procedure for deter- mining the accelerating parameter is different for each of the three meth- ods. There existed the possibility that the choice of the most efficient method might be more dependent on the method used for determining the 3 acceleration parameter than on the performance of the iterative method. To investigate this possibility the optimum parameter for the solution of a model problem.was determined on an empirical basis for each method. The results were consistent with the other comparisons of the three methods. However, this study did reveal that the machine time required could be reduced from 18% to 42% by starting the solution with an optimum parameter. Another objective is the preparation of a set of computer programs for solution of plane elasticity problems. Six FORTRAN-II routines have been written for the Control Data Corporation 3600,and slightly modified versions have been run on an IBM 1620. The solution of a particular prob— 1em.requires the preparation of a pair of routines for the given boundary conditions. The main program which provides linkage of these subroutines is called ISOPEP. A full description is provided in Appendix B. The third objective of this dissertation is the demonstration of the utility of the ISOPEP program. Solutions in terms of stress funew tions and stress components are provided for six problems: (1) A square plate with uniformly distributed loads on portions of two edges, (2) A semi-infinite plate with a uniformly distributed load applied on a saga ment of one edge, (3) A flat-plate tensile specimen with two semim circular notches, (A) A flat-plate tensile specimen with two V-notches. (5) A flat-plate tensile specimen with two rectangular notches, (6) An infinite plate with a square hole. Several of the problems were select- ed as examples of the numerical calculation of stress concentrations. Two of the problems are included for comparison of the exact and numerw ical solutions. The numerical solutions obtained with sufficiently small mesh intervals provide stress components which are in good agreement with values from.the exact solutions. II. THE DMD]. EQUATIONS Plane stress is the state of stress approximated in a thin plate which has loads applied only on the boundary and parallel to the plane of the plate. In a three-dimensional Cartesian coordinate system a state of planes stress exists if the stress components 03,7;2 ) 7}: are sero at every point. Consideration of static force equilibrium under conditions of plane stress leads to the equilibrim equations. 3o? } ‘53—+—7-’;1+X- (1) L? + 11v +Y: . where I and I m the components of body force per unit volume. For plane stress the Hooke's law relationship between stress and strain is, E G): l’ d} '- VO' +- BK 5’ f (2) 56,8 0} *1/07 4- Eofi‘ E 33.9 I 2 ( H1073“; . where 3 is foung's modulus, Vie Poisson's ration, §(x,y) is the differ- ence bet-sen the local current temperature T, and the original tempera- ture To, and (is the coefficient of thermal expansion. it a point in the plate the strains are defined: EX=§§r, 67:15,};y7-‘fi1-fi3f‘ (3) where u and v, the components of the displacement, are continuous functions. Since the three strain components are expressed in terms of the two functions u and -'9 they cannot be taken arbitrarily. By differentiating and combining Equations (3) it is possible to obtain the compatibility equation, z .55.... 751+ my = a"! (h) At each point in the plate Eq.'s (l) and ()4) must be satisfied. 0. B. 1117 derived a single differential equation in terms of a stress function, ¢(x,y), which will satisfy both the compatibility conditions and equilibrium equations. The stresses are determined by the follow- ing: a u. =‘g—3d-V, W=§X£+V9 Wiv=-§—xg4‘; (5) where Y(x,y) is the potential of the bow forces. The Hoohe's law relationships (2) are substituted into the compati- bility equation (h): 231%;- Vo-y + Each 3—1[o--;:-1/c; +E«€]=2(1+;/)§;y (6) It is advutageous to introduce the equilibrium equations by forming the st- of the derivative of the first with respect to x and the derivative of the second with respect to y. 2—? ”Ere +43%- £95.23 +£75.51.” 251+ =0 °’ va=_[a ax +95%%+%¥1T In.“ (S) and (6) can be combined so the term containingT :yis elimi- (7) mated. 51??” + E «£1 + 5715*; was.) - «leg—3 gag—+3 Th1. IWMiesto («$14. +§?)[dx+d’+fqrd+(1+1/)[ [3%(—+§ QY =0, In terms of the two-dimensional Laplacian operatorv2 this can be written V‘[r+o- o—+o(E€]+(I+‘/)[%§- +93%, = Substituting the expressions given in Equations 5 for the stresses and 6 assuming body forces can beucxpressed in terms of a potential function “1.1) the equation can be written V’4>+Eacv1(‘+(z—V) V2V= o (a) If the body force potential 1.: harmonic, the equation becomes W + Ear V‘é‘ . o (9) and if‘i'iti‘additifdn"iteniperature changes are negligible, it becomes V‘tb = o (10) where in Cartesian coordinates, 9 +47 ’46 a Vl’d’ = $79 +2999y=~+ " a y ‘I- The solution . of a plane stress problem, when the body forces can be expressed in terms of apetential, thus consists of finding the stress function fl which satisfies the appropriate one of lqfls (8), (9) or (10) and the prescribed boundary conditions. _ As noted by Bakelnikeff (1956) , even in a thin plate the stresses vary somewhat through the thickness . The two dimensional problem formulated here, strictly speaking, applies to the averages through the thickness and this is. often called a state of generalised plane stress. For a state of plane strain instead of plane stress the governing equatiens have a similar formulation . The appropriate equations are obtained by replacing a! by V/(l-ll), s by l/(l-Vz) and 0( by (I (1+1!) 1. has (2). (5). (7). (1) and (9). Hi and? are the components of external loads per unit area acting on the boundary and (I and 5 are the direction angles which the normal nabs withthe x and y-axes respectively, then neglecting bow forces ‘7 X! Figure 2.1 Boundu-y forces X: 1032 + M Ky “1‘72 mm} + )7 73.3 (11) where = cos 0:. 3 731 = coo-{3 In terms of the stress function the stresses on the boundary are given by Xflfl ,5. ”ii—3%,}, Y=W~°§$ri§$£ (12) Introducing coordinate axes s and n, tangent and normal respectively to the boundary, the boundary conditions can be written X= %:$1‘3§"* xlafifé z 3'J(§§) 7:.- egg? 37%;}: a -33“; 5%) These can be integrated along the boundary é—f =- - I You i (if). 3435- - W + at. ... sat- hes-e . integrating along the boundary yields 4> Jasmin” who I Yests + (i?) )(7- in (fines) + 49. (15) where (our -g;3, 51M: - substituting into the normal derivative of ¢, #= §$+va g9; :-— Cosxf Y4; +51» or] X45 +( @f’kosmfigg) 33nd (i6) Iqfls (15) and (16) determine ¢ ande ~9- at every point of the boundary in (13) (1h) terms of the boundary stresses 1, I and the constants of integration (%%)a,(%3 )Oand ¢°. These constants of integration m be chosen arbi- trarily since they do not appear in the expressions for the stresses. If the constants of integration can be chosen so ¢ is symmetric with respect to a line of physical synetry for the plate, then Eggs) 2 o , a, (v‘dmgw : o (17) for points (x,y) on the line of symmetry, where n is the normal to the line of symmetry. The problems herein considered will have boundary conditions as given by Equation (17) on a line of symmetry of the body and of the form u” ’ = 5am), 1%” = flow) (18“) which follows from Equation (114) or WW) 313 (Km) 9 1 #3212: F (m) _ (1810) which follows from lqfls (15) and (16) on an outer bomdary. The func» tions f1, f2, f3, and fh are valid only on the boundary and could be expressed in terms of the single parameters 3 n is the outward normal. For problems with displacements specified on the boundaries the Iavier Equations should be used rather than the biharmonic equation. This pair of coupled second order partial differential equations for plane stress conditions, obtained by replacing the stresses in Equation (1) with expressions in terms of strains and then substituting deriva— tives of the displacements for the strains as given in lquation (3). can be written: 171E37I7[Z%+(’ V)‘a_1 ay +"*");§‘5‘9’]+X=0 W[’.§u +I.. 377 Uo%{+Iex 91).er + In «tafwgogj: u,.= U. - (U006 (Uy)L+[qu:)+2(UI>+ (:uy )]2, (+39: 3—,)U 3:7 ”37‘ +3?)+Uo‘éf; “)5U051+(39X + Dy)‘ U0 Io Q! U“? U’O‘ZIUO = [(Ux‘)‘, +2 (0+ij) (Uyl’ a)i]%.hk+ [(U I‘ll”. 4(UX331)¢ +‘(UX (1)9 . he +‘I(nys)o+(U 90] if + gig—'7 M, where Mhsmaxhlx ”MO, 40 n + m = 6. Similarly, Us* Um ZUD‘: [(Ux1)o+2(ny)e+( Uy‘lfi fiHHUx") flux?) 0‘01'6(sz7) “WWW (oUII’IiJZ-fivf 557” M5 where 145. < maxlanym'812, n + m = 6. Adding we find Up U8 + U,O+ U,2- ‘IU :2[( me) )022h! + (UxI)2—-h:+( +(2)Uy 57+ “HUME 6 +(Uxay‘)ofi+ + 2'2"? —"'"M6 where M6 < mex|ux.y,.|, n + m = 6 along all lines through the mesh pomt. 12 Substituting from above (Uny‘L' ”I'I'IIUvUe +Um+Uu+ IUa~2ILI+IlgIIJj+L§g ~_- léfi‘M‘ I (22) Hence '1‘ I 9. -, ‘7’?) = eafi‘tzfgefifiro ia approxinated by (0x90 +2 ”1%),“ U190 9" M UI+UI+U9+ Lin "' 21"”? Us" qe“ Um) , . ~8(U, IUZW5 +U.,) +on,J+R,,—-—o (23) Ihere R7, LS 'fijé-th‘ . Kantorovioh and mm (1961:) give a, _<_ $9- 112)! where l in obtained by replacing all the derivatives of higherarder in lode (20) 9 (21) and (22) by their maximum absolute values on the meeh lines Joining point 0 and the surrounding pointee , ' The remainder 95 providee ‘a bound for the di, ecretication error, that ie, the difference between fixer) and the truncated eeriee approximation” U(n,y) at the nodal pointee Thin indication of the order of the approxi- nation in a major advantage of the Tyler" e eeriee nethod for deriving the finite difference equation. Inflation Techigue _ The use of integrals for the derivation of finite difference equa- tions has not been need an extensively as the other methods. Verge (1962) providee a general introduction to the method and indicate: it has been need in nuclear reactor design computer. codes for several yeereo Use of this method for the biharmonic equation is made by Griffin and Verge (1963). The aatiefaction of the partial differential equation (8) at every point of a region R is equivalent to the satisfaction of the integral equation IVV1(V2¢+BV*C€)J‘4> :- 0 (21:) 'R 13 where B = (1mV), C -CKE, for every arbitrarily chosen subregion Ai' Hence, for any arbitrarily chosen subregion Ai bounded by 31, Green's theorem gives 2 )de 3953- ‘ I Cfl'fi-o f/v(v2~

iii= if; iigu, Uy-M ML] ifm- Mg] (29) If the mesh spacing is constant, Equation (29) reduces to the form of Equation (23). ‘The evaluation of the right side of Equation (23) requires 12 addi- tions and 3 multiplications. (Ihen the mesh spacing is uniform.the multiplication by h"2 is not-necessary.) The evaluation of the right side of Equation (29) requires 29 additions and 25 multiplications if the ratios giéa. are computed once and stored. For the IBM 1620 the time required for 10 floating point additions is equivalent to that required (for one multiplication. Thus on the 1620 Equation (23) could be evaluated for approximately 7 mesh'points in the same time required to evaluate the right W use of Equation '(29) for a single mesh point. Any iterative method or solution of the finite-difference equations requires the evaluation of one of these expressions at every interior node at least once during each iteration. Several hundred iterations may be required. The rate of convergence of an iterative scheme which uses Equation (29) would have to be 7 times the rate of convergence of a scheme which uses Equation (23) for the same number of mesh points before Equation (29) would be preferred. In some problems it is possible to reduce the number of mesh points significantly through the use of arbitrary mesh spacing. Regardless of whether it is possible to establish the superiority of Equation (29) over Equation (23) for use at 15 every point, Equation (29) provides an excellent method of handling irregular boundaries and the changing of mesh size from one subregion to another within R. Variational Formulation Application of the variational method for deriving finite-differ- ence equations can be found in Courant and Hilbert (1953) and Forsythe and ‘lasow (1960). Ehgeli, Ginsburg, Butishauser and Steifel (1959) show its use in deriving biharmonic difference equations. The variational formulation of difference equations for plane elastostatic prdblems is especially convenient for prdblens given in terms of two displacement functions. Griffin (1965) has shown the advantages of this approach. The basis for this derivation is the Principle of Stationary Potential Energy which states: Fran the set of continuously differentiable displacement distributions which satisfy the given‘boundsry conditions of an elastic body, the displacement distribution which actually occurs is the one which makes the potential energy stationary. n For a plate subjected to plane stress, take 1,! as the body forces per unit volume and if the surface forces per unit area. Then the strain energy V, per unit volume is V: fiféf+ £37 + 21’Exe'y + 11:12! ‘6}; The strain energy per unit volume can be expressed in terms of displace~ ments (u,v) if these are continuously differentiable functions of (x,y) by substituting 5x2??? €y=§a§ra Ey=(%1+35% v= m%a[ei->*+er+~2~—%i + "flee—Ir] ‘3” ‘1 change or variation of total strain energy will occur'with any arbi- "16 trary variatiu of the displacements 8n, 8v, “(7 = Sffdexdydz = f/ SVdXd)’ where for convenient?” the s-dineneio: is taken as one and it is assmd that v does not vary with s. The virtual work 6"“;" done by the external forces under the virtual displaoauts Sn, {v is given by 5 We“ affix“ +Y5v)Jny+ {(5582: +YJrHJ the potential energy 1. defined ‘ Gaffvdxdy - [fix-u +Yr) )clxdy- -f(Xu+YI})JS (31) Applying theR principle of stationary potential energy and considering bow forees ad surface forces constant SQ = 0 A m =1)” svelxdy - [{sz +Y§r)Jny ‘JKXéatY-Eflda’ (m This states that the-Zion‘s in strain energywill be «mil to the mark dementhebetbyentenalforoee for an arbitraryvirtnal dieplnoe- meets (n.8v,.fh1swillbetrneonlyif w arethe actualelastie displse-snte predued by the external forces. . the usual procedure in the calculus of variations is to consider the eonditions imposed on the integrand of Donation (3!) and thee derive a lineu- partial differential equation of the form W} = 0 __ (33) there “(a,v). Ilse additional bondary conditions arise which are called natural boundary conditions. Instead of finding the differug and mum of the fern of lquation (33) and then act-emu. emu differences the variational form of location (32) will be solved snorioany. There are two major advantages of this approach. first the nestles v which satisfies Iquation when: automatically «um thenatu'alboudaryoeaditione sethat fin-thereonsueratieneffie 17 natural boundary conditions is not necessary. Second, the discretization of the problem leads to a linear system of‘equations that has a coeffi- cient matrix which is symetric and positive definite. These relation» ships are important because they are part of the criteria for the convergence of iterative methods for solving linear systems of equations. (See Appendix A.) For a plate in the region R of the xy plane and unit thickness the total potential energy is Q - 2,, ,4, —5—~[/M33‘ ”fa-351+“ V){(§%)+(%19, }+ " 11W; #5?an - [/(querMny - 2’3 (23, affine (3“) The stationary value of potential energy is a Smininum for stable equilibrium under specified boundary conditions. a rectangular mesh and its dual are imposed on the region B. The y dual mesh lines (indicated nith ‘ i If dash lines) are' parallel to and b—[w — T __ t . spacedhalfwaybetween the lines —.%--.._§-- .:L ,__4:__ of the primary mesh. .L J. : : _-l__._1_..‘.... '_...J_._+..... The variational. problem : i : : Equation (32) is given in tons L—E—.._-.E}.-i .é. “1r“ of a continuous displacement ' l ‘ 1 4-” distribution u(x,y) and v(x,y). Figure 3.14 Dual meshes This is replaced by a problem in which the displacements have discrete values “1 and vi at the mesh points. The potential energy of each mesh point is approximated ,and the sum taken over all mesh points represents the total potential energy 0(u1, v1) . The problem is reduced to finding the unknown displacements “i and 7i at each node which will make Q(ui,vi) 18 stationary. This requires ‘1’: =0, 3’-”233"-“~—” “‘ . (35) g 80’ Ja‘,2,3)..~——77) '3 where n is the number of nodes at which “i is unknown and m the number where vi is when. Equations (35) represents a system of n+m linear difference equations. , hgeli, Stiefel et al (1959) show that a variety of quadratic functions may be used to approximate the potential energy. The polygons formed by mesh lines and dual lines will be called cells. A primary mesh point will occur at one vertex of a cell. A regular interior point will be the? canon point of four adjacent cells as shown in Fig. 3.5. The potential energy integral (31.) is approximated for each cell under the assumption that the functions a and v and their 20 h. derivatives are uniform over the f ----- . -2- "I l 3. “all 5 1 hol L cell. The functions take the ‘ i i l w l discrete values at thenode r0. L__-______j . . . h... The derivatives are approximated a; J» for cellB by Figure 3.5 Cell 6 . 9X 30’ 5 y - hoz ‘ Th0 approximation for the potential energy is Q-OZI gflEW[V(‘kJ’*%)+ o . (1.16) This scheme is known as the Richardson iterative method (or method of simultaneous displacements, point Jacobi, or point total—step method.) See Varga (1960). It requires all elements of the vector iterate U“) for the computation of 9-01-11) . Other iterative methods introduce the new values of the eluents of the vector as they are determined, and these are used in the computation of successive elements. Methods of the latter type use only half the storage required for Richardson’s method. Iterative methods are characterised by the repeated application of a computational scheme which yields an approximation to the exact answer as a limit of the sequence of successive vector iterates. 1 basic question which must be answered affirmatively for any useful method is 3 does the sequence of vector iterates converge? To answer this an error for each vector iterate is defined Eh) = Eh) _ P. n_>_ o Subtracting Equation (I45) from (I46) we find EIm-tl) . (1:. .. 5) 2(3) Repeated application of this relationship, starting with the initial vector 1.1“) gives . ems-(14)" 155°) n> . (m considering a single element 310') of 30-), if the H; 010!) . U,- it is necessary that $111“) and 135031“) exist and 31.313031“) as 0. This condition will hold for all elements only if 26 $134; - i)‘ of” = 9 74 a _. (us) for any arbitrary vector E (c). Q is the null vector of 1: elements. Iquation (ha) is valid than only if gig-we - i)‘ - .9. (1:9) where 9 is a square null matrix. . ' hilne (1953) shows that Equation (19) is the necessary and sufficient condition for convergence of an iterative method and this is assured if all the eigenvalues of the matrix 2:, e _I_ - a are less than one in absolute value. _ .. Iindsor (1957) has shown by derivation of the eigenvalues of 21 for a rectangular plate that Riehardson' s method is not convergent for the biharmonic equation. In the remainder of this section iterative methods will be reviewed and those’which are especially useful for solving the biharmonic equation will be identified. a. Log-:31 - Iterative lethods . The general form of Equation (1:3) for the biharmonic difference equation on a rectangula- region with n colums and p rows of interior mesh points is of the form ii-:_ , , so where l is a m: sparse, non-singular matrix, Q and g are column vectors of 1: components and k '. up. The derivation of the matrix form of some iterative methods is simplified if we tales A a Q 'I' 2 + E . (51) where g is a strictly lower triangular matrix with the elements ‘13 (i) 3) below the diagonal of _a_ and all other elements zero. ]_J_ is formed of the diagonal elemmts ‘ii of _1_ and g is a strictly upper .27 triangular matrix with the elements ‘13 (i < 3), Q, 2 and g are lock matrices. Richardson' 3 Method U(m+l) - - - (52) where 2‘1 is the inverse of Q, a diagonal matrix with elements l/an. For the model problem, Dal, n25, ps3 the point value is Inn)” on) m) UM) ( '27) U33- =53 MU}? 3+3? 3+U +£33.31 .UI .+-U ”"~ +U-,333 m3} m H +Ui+:}n]- .05[U”3 +Ubj;z+U,(;2j«/-IJ Mafia], (53) 1525.7), Isl—P. Though the introduction of exterior mesh points may not be the best. way to account for the normal derivative boundary condition one advantage of this approach is that the term F13 will be automatically included when the calm subscript of an term is not within the range 1 to n and when the raw subscript of any term is not within the range 1 to p. Hence the 1'13 can be dropped from Equation (53). Gauss-Seidel Ile___t___hod ‘ Varga (1962) identifies this also as the Liebmann method, point single min method and the method of successive displacements. The main difference between this and the Richardson method is the i-ediate intro- duction of the adjusted point values into a single vector iterate. The matrix form illustrates how this is accomplished. (2+9) 2-21-11 9. 2“”) = (earl i—(aurl s i“) (51‘) For the example problem, the point values are computed with Us”: 3 ru-+u33:;’+u3‘::,’3 11:33-33:23: (1.31:7: 33-3333; 33 can-5.33:3" 333" L332;- U313; ‘5" 28” Here too, the 1'13 will be automatically included if point. values have been provided on the boundary and on the first exterior mesh line. ‘ Starting atUil and sweeping along columns successivelyche required (n+1) iterate values will be available as specified by Equation (55). Boundary values computed initially remain fixed. Exterior point values can be recomputed at the end of each iteration. This method is convergent for the biharmonic difference equation but the rate of convergence is so slow that it has limited usefulness. for a problem with 66 interior mesh points 1098 Gauss-Seidel iterations satisfied the same convergence criteria as 22h iterations of the successive overrelaxation method. Successive Overrelaxation This is also known as the extrapolated Liebmann method, Parter (1959) . systematic overrelaxation and the extrapolated Gauss-Seidcl. Young (1951:) proposed an acceleration of the convergence of the Gauss— Seidel method based on an examination of the changes introduced in a given vector iterate and then introducing a multiple of the change at each point. The point value obtained by the Gauss-Seidel method will be designated 31““ ),1 5 i 5 1:. Then Uimu) ___ Uilez'. ”{Uilan Ufm} : ‘ I _ 00) ugh») + w Uitm-a) (36) The quantity 6) is the relaxation factor. For overrelaxation, l< “<2. ‘hen e a l the method is Gauss-Siedel. The matrix representation must account for the prior application of Equation (56) at all preceding points. Hence, D Ulm0+ §LZ""”')= f _Hylm) (57) U’”*” = UW-t ML] ""21! W) (58) 29 Eliminating Em” between qu's ,(57) and (58) we find the matrix tom of the successive overrelation aethod (Q+w§)g‘”’*" =[u-w)p- tutu Um-t cuff (59) The component U13 of the vector iterate at a point where 15 13 n,l g j g p is given by “a, b”) mu “fl its“) 3+1) (m H), (m) Um (7n) (7») My) he) FUH’jH: -bUl'flfltl) .Q 5(U2-2 J +U51-2+Ui+z,,+Uy+2UiJ (60) The rate of convergence of Equation (59) depends on the value of a). The optima value of w is given by the formula 2 w = —-————— b [+1/ I .. 92. (61) where f3 is the spectral radius of the point Jacobi. method. a," =. IE‘"’II/||§“”"’ll 62 m «° = w A»: N. _ _‘ ’ Any nor. of the error vector 2 could be used. One readily cmpnted 1. “Em H= 2:: hi ml a , i (see Appendix A) a canon proc;dure for the deter-nation of the relaxation factor is to set it initially to one; and, after a number of iterations, say 100, use Dee's (62) and (61) for calculating a not w. Subsequently m can be recoaputed every ten or weety iterations. Forsythe and Issue (1961) pp. 368-372 describe two alternative approaches and indicate that the determination of a good estimate of Us early in the computation is an investigation in which there is continuing interest e The relationship between the estimated spectral radius and the nmer of iterations required for convergence of the model problen of Fig. h.1, p.20 is given in Fig. ho3. 30 U1 .8 g : Iodel Problem 5 :‘s eOOOOI 6 x 13 flesh The "' e05 SEC/IT. § NO. OF ITERATIONS l-' w .8 ? :9? 598 .. 597 :96 [95 I91; ,_ ' Estimate of Spectral Radius 9 Figure h.3 Successive overrelaxation Number of iterations vs. P B. lethods of w Iteration The iterative methods considered thus far used an explicit formula for the calculation of each component of the vector iterate. Is it possible to use direct methods to find s block of components of the vector iterate? Consider the model problem in terms of submatrices, Equation (13). Taking arbitrary values for the column submatrices Eh and g3, the components of 32 can be determined by solving E 92 (“1): .122 - .3. Hi ('0 - a 23‘" (63> share 92 a. {U21, U22, 1123, Ugh, U25}. Since a is an (nxn) matrix, substantially smaller than _1_ which is (Ink) the direct solution of Equation (63) will not require unreasonable blocks of computer memory and is an acceptable procedure. Using this method the point values of the vector iterate are not determined explicitly one at a time; instead 11 components are determined simultaneously. Hence, this procedure is called the simultaneous displacement method and is classified as implicit. . . v .IV. '1 4 ." 11?». tlll‘lvr' vi. Ill-Ills!!! I! I‘ll 31 The matrix form of certain block methods assures faster average rates of convergence. See Appendix A, Section 5. lrms, Gates and.Zondek (1956) and Keller (1958) have investigated block methods using the components of the solution vector on a line as the basis for partitioning the matriximo For the biharmonic difference equation, "two line“ schemes and the alternating—direction implicit method have been studied and appear to have advantages over other block methods. See Parter (1961A) The AlternatingéDirection Implicit Method Peaceman and.Rachford (1955) found that the rate of convergence of a "line" method could.be substantially improved if after sweeping all rows using the simultaneous displacement method and a relaxation factor, the next sweep of all the mesh points‘was made by columns. Coats and Dance (1958) derived a convergent, alternatingndirection iterative method for solving the biharmonic equation. This method is similar to the alternating-direction method for solving Laplace's equation proposed‘by Douglas and Rachford (1956). The derivatives in the biharmonic equations are replaced using central difference approximam tions 9X? ~ ij = 8:Ul'j =Uz'J'-2 RUUZJ-I +6U _LfUz J-H +UIJJ+l 93V~zLflU1U= 8;sz—Ui-)2J‘ PUU-IJ+6U:J-4U?+l) J+UZ+11J .2142. u w ”1 2’3 L~IEU 2851 U "'2 L+U" 2[Uz’-I,J +UX+IJ+ 1,.J'l 1.1 X 9‘1 1'“: X J 2J g (J + UHIJ.“ + UH J'H +Ul+')JI”1-UI'U"J‘] '4 WW) 2-. 5x U23 + 2 5x151; Uz'J' + 3‘} Uz'j Is shown for the model prdblem, Equation (h2), introduction of the difference equations reduces the problem.to the solution of a linear 32 system in k unknowns— 53> Id ll l’zj which can be written (5 + 9. + 2) u E where l}: flij’ [g gjij and L gjij respectively represent the components of g g, 911 and £11 at the mesh point (i,j). For the alternating direction implicit method, Equation (1:3) is replaced by a pair of matrix equations '11=<£~ rs) E+r§ where r is any positive scalar. The iterative scheme in the form pro- posed by Peaceuan and Rachford would appear 3 +1) 1.1““) = (_I_ - rll g) 31“)+ rlll 2 (6h) (rm-H H ‘H (rmi'l 9- + rn+1 .P- + D Elba) '-'-' (I " rn+1 fl) Ehfié) + r.+1 E. 565) 9- - rmH The first sweep of the mesh is by rows and only Equation (614) is solved for the vector iterate 9542'). This is an implicit method since all the components of Eh'w) in one row are determined by_ (614). It is necessary to retain all the components of pf.) while solving for 6‘99. During the second sweep Equation (65) is solved one column at a time for the components of 2“”). The solution along one row or one column is obtained by direct elimination. Conte and Dames use a factorization techniquemhich is well adapted for the case of Fig.4.2, where there are, at most, five unknowns along a mesh line. The Douglas-Rachford method is a variant obtained by changing (65) to a form which does not contain E. Substituting for g _t_I_(‘+'/2) from (61;) we obtain (rm-l Q + I'mi-l 2- + -I-) E-(mfl) 1’ 2 ._ Haw/2) "° (1 " I.1114-1 9 " rm+1 2) 2(a) 33 Coats and Dames use a simplified form not containing 3 use) = “(as ._ rm [9 Hum) .. Q Hm] . (6'6) The use of Equation (66) for the column solution provides a. Simpler computational procedure than the one associated with Equation (65). The determination of optimum acceleration parameters requires consideration of another form of the system of linear equations. This is found by combining Eq. (6h) and (66) into a single equation. 2W1) = are; 11‘“) + an r; (67) where gr c (_I_+ rg)'1[(r§+;_)’1'(;-rg- r2) + r51] r = (l + marl '(ra + 9‘1 . A The difference between the nth vector iterate 2(3) and the solution vector 2 is the error vector 1!“). It can be shown (see Appendix A) that gun) =_ 5,.“ 50-) (55) If g“) is the initial, arbitrarily selected vector iterate then 3(0) 3 2(0) .. 1.] and §(l) . firm 0 firm-l eoeeeefirz oar]. 2(0) For a convergent iterative scheme Elm) 92 as m 960 , and convergenceof this method can be accelerated by the choice of rm for each iteration. However, the relationship of this iterative method to others considered is less complicated for the Specialease when a single value is assigned to the scalar r. Then i g _ 15‘“) 4353‘?) f _ (69) and the eigenvalues of the cents-Dames matrix gr provide the basis for detemining the optimum acceleration parameters. See Fig. tub vhich shows the relationship between the acceleration parameter r and the amber of iterations required for convergence of the model problem. Fig. 4.1, p. 20. 500. Model Problem U“ 5 a .00001 2’ 4p 6 x 13 Mesh o boom Tine . .098 SEC/IT. 22 0’ 300- DJ F. -200- LL (3 .100n (D Z I I t I I l l 0 b0 60 _. 80 10° 120 ”A .Ualues of Acceleration Parameter r Figure heh Alternating direction implicit method Number of iterations vs. r Choice of optinnn values of rm for a square plate has been considered by Conte and Dames (1958) and Faimatherand Mitchell (1961;). For this particular geometry the eigenfunctions of Equation (68) can be expanded in the form E3) = a sin (pfiih) sin (qfljh), (p. q .. 19 2. °°°°. 104) Consider an amplification facecr ) = if” Ft ELI” If this is substituted into Equation (68) we find (1 a 16r +3 Se 33)? (16 rml))\p,qe1+16rm1r(t§;h,g;) + 256 ,mlg :3? (70) where Sp : sin P3212 , 5Q 2: sin i192, h is the mesh interval. The error associated with the initial vector iterate can be 35 expanded in the for. u— (O) h-’ (0) ' n . ' E" 2%-97’1?‘ SMOOTH/05m (3/106), , _ th and the error vector at the n iteration can be written (a) d (1») , , . ‘ E5] zg’C”, 5/”(f7NIIISIM (fly ’1’, where no = (a) It follows fro- !quation (70) that Of (16r)! A3351 for all p and q if r! is positive. Hence, after 11 iterations each compontlt of the error decreases by a factor m . 1U: (167')! A?) 3/ if r is positive. The minimum value of (16r) A found from I X P"! Equation (70), is zero and occurs when I lér = S’s-5;. This indicates that an appropriate choice of the 1:? can be nade so all components of the error vector will vanish and the exact solution can be obtained. Rather than attempt to find this optimum I} , a less conplicated procedure consists of choosing the r]? which optimise the rate of conver- gence of the method. _ Consider a slightly different expression for an amplification factor - -__-_- (I—lérglf’sle 2- ’67'Afz "’ (Hursfiiflz Since 36 '6'; +33? 2. 25,2592. we find from Equation (70) fl I16r~7| 2 pl IérM for all p and q. The analysis of the problem of finding an. optimum set of acceleration parameters is treated by Doug Lea and Rachford (1956) . The factor by which the error is decreased after I iterations is 2. a 2 2,0936?) ’67:?) = I:— jyy/MIQ) = l: [2%f:+ 5:?)1] It is necessary to find the set oer( (=1, 2 ,-—-—-t) giving a maximum value of 2t which is as 8.311 as possible for (p,q :3 l, 2, 3 ---k-l) . The Douglas-Rachford solution treats this as aChebyshev minimal: problnpwhilo- Conte. &1DCI08 reconend a set of acceleration parameters of the fora 1611, “((14), X: l, 2, -f-, t where 0< or < 1. This permits the determination of an upper bound on 2t l-( -TL Z. (5' 3,, I671) < Titan): [W e 1 (71) The formal procedure, outlined by Conte and Danes, starts with a choice of Pt ((1) which permits the determination of the number of cycles of t double sweeps of the mesh required for the selected reduction factor. However, the choice of Pt ( or) is subject to the empirical observation that best results are obtained for or< 0.2. The number of iterations, t, per cycle is conputed from t _>_ I + Wee/5:341) (72) 2; = T s X= 52,3,-~,t (73) 37 gives the value of rk for the gthiteration in the, cycle. Is an example of the method consider a 20 x 20 grid on which it is desirable to reduce the initial error by a fector 10-6. ”Selecting or a 0.2 makes Pt” 0.01. The ntmber of cycles required is determined from (Pt)n == 10"6 _ . Thus, 11 g 3 is the number of cycles. 'The number of iterations per cycle is > [05521919 ~ {I‘— l+ 087.2) ~ 7'35 ortg8and Factcrisation Technilue . Il’he difference equation (6h) can be written in the form Mb“ (11+Vr)(*'+Vz) (we) I’M-Va): , Pg.) ‘6.)- 4U]- J6] +(6 +— 7M) U5). L/Lj).§b/. + Ui+2)J Fé') (71‘) where Uh”) m) ('m)_‘_ (7» ’33' '(?.,l,,—, "WI 3-W( I+IJ+U5-I,J')+ )8(UI', 1,13: +U1IJ:‘)) m (on) (M (7" 3") ’” - 2 ( Ll“, J+I +U -I,J'H+ UH l-l “(limp/'1) U1) Iii-2- ~ U61"; similarly, difference equation ( 66) can be written in the fore (n+I)_ Um“) (”11‘") Uhs-II) I‘m-H) ‘ U,,_z +l6 6+%—IU,,, LILU,,,,.+U,,,,,= F2, (75) where l’” 4%) Urm) (1n)+ (ML (7»; re = —,— ,, . ":‘wa.’ +6 Us 4-4113"; Use. Usins (7h) at all mesh points in a row leads to a 'quidiagonal" system of linear equations of the form 38 " cl -h 1 "U11 l ”F111 7 ~15 02 4| 1 U12 P32 1 --h 03 -h 1 U3 r33 1 -h 0,, at 1 111,, rm, 1 =21 6n -l ‘1‘ Ui,nol FBn-l l -h on Uipn rfln D .J l- where the C's all have the value (6 + l ). 1 system of linear equa- rn+1 tions of the same general form is obtained when (75) is applied at all mesh points in a single column. The factorization technique is an effi- cient direct method of solving a quidiagonal system of equations. Take '1 = 01, 'q 8 Cq " gq_2 _. quQ"1 BO :09 Blg-h/Cls qu-(h-tdq gq_1)hq, Bn=o~ ‘0 I 0’ 81 a: l/Cl, . {q : lhq , (n a ‘11-]. I O t - “b - Bq-z h0 =3. 0, h1 *Li'm/QI’ hq = [FR - no:2 - dq hq_1 1/ .WQ’ 2g q 5n Starting with the specified initial values the n values of h are computed. Then the values of U are computed by back substitution Un=hn q = ha " Bq ”w ~' " 8.9 ”q+2 ’ q = “‘1’ “‘2’ The alternating-direction method was initially applied to the ‘1 U ,1 solution of parabolic and elliptic partial differential equations of the second order. For Laplace“ 8 equation with sufficiently small mesh spac- ing it provides a significant increase in the rate of convergence per iteration over the point successive overrelaxation method. This advant- 39 age must be weighed against the requirement for the double computation in a two sweep scheme and the need for storage of the (n+2?) vector iterate. Compared with recently developed block successive overrelaxau tion methods, superiority of the rate of convergence is not rigorously established. Varga (1962) notes that the convergence of the Peaceman- Rachford and similar alternating direction methods has been established only for rectangular regions. Though there has been some success in applying alternating direction implicit methods to more general regions, the convemence in the general case is yet to be Justified. There is no general theory for the determination of optimum ‘ convergence para- meters except for rectangular regions. These same limitations apply to the solution of the biharmonic equation and in addition Keller (1961) demonstrated for several block methods that corresponding biharmonic schemes converge more slowly than Laplace schemes. Semi-iterative Methods A system of linear equations 5 y. " E: can be solved by an iterative method of the form 31“”) g g E“) + E m a o (76) it a .-.-.- _I_ - g is a positive-definite nan matrix. See Appendix A. as m as, 11(3) converges to the unique solution of the system of equations. A semi-iterative method uses an algebraic combination of solution vector iterates EU!) as a means of increasing the rate of convergence. Starting with an initial estimate 9(0), the error vector associated with the vector iterate 2(11) is given by 30“) . U0“) - U . m 1“ E(°) )40 A linear algebraic combination of the vector iterates Eh“) is introduced 1““) :3 f: pjtml 2(3) 1:120 (77) J=0 where the coefficients pj(m) are selected so that each [(13) is a weighted average of the 2(3) and a better approximation to the solution vector y than 2011). For the special case 2(0) 2:, g it is necessary that X | »" f filmbl (78) i=0 Then 10“) g g for all m _>_ 0. The requirement of Equation (78) will be . imposed on the constants p.) (m) for any arbitrary choice of 11(0). The error vector associated with 3:0“) is denoted iii/(m) 2"!“ g :6“) - y, = fpjon) 11(3) - g (79) _ .. . i=0. Since the Constants pj(m) must satisfy Equation (78) ’ I... '1 5"" — z: srwfl” ~ ( i F} WW I . -0 EW= z: firmIflw—g] E/(m’= Z 8’ (M)_E_:(J) ‘ 0m): j (a) ( gm M )E (80> 0. H 0 If we introduce a polynomial in a component of 115, defined by Pm(u) s: 22: p3(m) u:3 m a O (81) J=o then we can write (80) in the form 1.35m) ___ Pm (E) 2(0) where Pmogg) is a polynonial in the matrix 2°. The condition imposed by Equation (78) on :ijm) requires Pm(l) g 1. Using the definitions of matrix and vector norms given in Appendix in A, we can write IIE"’"’H =IIBfM)E°’ll i N am: HU’H ' m For “I309” < l the average rate of ccnvergence of the semi~iterative IV :3 method for m iterations is defined as ._ .— E MN 11 [ 13(3)]: _XnJJFLi If the polynomial is selected so Pn(u) e um, the vector Em) is identically 11f.) Ind the average rate of convergence can be readily simplified. ghnqfl. ~121LEmlfllfl. W ._-_._ RPM”? 7n 3' --'m v" 4 and, as shown in Appendix I, it is equivalent to the convergence rate of the basic iterative formulation of the problem. The average rate of convergence of I“) will be optimised by finding the minimum of II Plum/l under the restriction 241) = 1’. See Varga (1962) Chap. 5. Chebyshev Semi-iterative lethods Golub and Varga (1961) identified a polynomial which satisfies the requirements for optimising the average rate of convergence. They used [ 2ub- (b+a) Pm(u)-:m[2-(b+a)'| 5 11112.0 m H J Cos (m Basal 5), =1 5 3 5 1: m .. 0 where on“) : Cosh (m Cosh-’1 s), ‘ > 1, are the Chebyshev polynomials The Chebyshev semi-iterative method for the problem specified in Equation (76) has the form (n+1) (n) (maul) (zeal) tn 1 =wm1[!.1 +z-z 1+1 (6:) 1:2 where (See Varga (1962) p.138) 2 Om l ‘ w,.1, “n+1- 1 p , m_>o and p is the spectral radius of 1. lots should be madeof the absence of the vector iterate E“). It is not usedim the computation. The vector iterate [(#1) is formed directly from the preceding, vector iterates 35') and 10"“. It is necessary to store both I“) and 10"“. Spectral radius is defined in Appendix A, Sec.2. E2129. £22: 2 man 11:22 Consider a slightly different approach to the solution of lquation (70). Let Id II [N [’1 I + , t or equivalently IN IS Id '1". !-!!+9 an where ‘-" [fl 1'" [1! 2J' , 9' E] If _I_ is an an convergent matrix, Equation (83) has a unique solution IO II E and the subvectors E and E are equal to the solution of Equation (76). Since the solution of. (83) would require“ the manipulation of.(2n12n) matrices, it is not reco-mnded as a practical method. However ,' it does provide the basis for an improvement of the suit-iterative method. Issue that E in Equation (76) is m nxn, Hermitian, convergent matrix with the special fen 4:. it] where the sub-strides are square and gl- is the conjugate transpose of g. 1:3 A matrix [which can be expressed in this form is said to be weaklycyclic of index 2. l weakly cyclic matrix 5 of index k :is such that 51' has real“ non-negative, eigenvalues. For the assumed form of g the successive over- relaxation method for Equation (76) can be written (Varga, (1962) p, 11,9) 21““) = w [a 22“) + :1- 91 m] + 21““). gab-pl): “@(fll)+§_fi2(m)]+ge(m), mgo _ where g and g; are partioned- into 111, He and :1, 32 respectively. 1'_h_g m Chehz‘ “ shev Sent-iterative M The Chebyshev semi-iterative method corresponding to this special matrix ! can be written 1 (Verge (1962) p.150) 9.1"”) = am [3 22“) + £1 - 21““1’] + QM)» 22‘”1’=w.+ [s2“’+za-U“’9J +95“)? 1:1 Since the iteration parameter a) “1 is a function of the number of iterations n, the cyclic characteristic of _l_ per-its the following fore 5‘2”” :- n+2 “1L- am” ”1;- 42‘2”] +952” --:° IPhis is hows as-the cyclic Chebyshev semigiterative method. bus to skip- ping half the vector iterates the rate of convergence otthis method is twice tint of tin Chevyshev semi-iterative method. Basically though, this method is Just a variation of the successive overrelaxation method. Varga (1562) shows that the average rate of convergence of the cyclic Chebyshev semi-iterative aethod is better than the average rate of convergence of the successive overrelaxation method. The cyclic Chebyst semi-itera- tive nethod will be used to solve the model problen and the convergence rates of successive overrelaxation, the alternating direction implicit method and the cyclic Chebyst seniaiterative methods will be compared. hh M Chebyshev §eni-iterative Hethod for the Biharmonic Elation Consider the solution of the biharmonic equation for a rectangular region with a grid of n interior mesh points in the x direction and p interior mesh points in the y» direction. The difference equation can be written in the form f J §H:£ . ‘ (8h) ' where N a up, _A_ isan Hid! +~=w l coefficient matrix, .3111. the solution column vector and {is a colunn vector N\° A l I IP""T—- which accounts for the specified boundary conditions. The for-ulation'of the Figure in; m line blocks _ .. _ __ model problem in terms of sub-utricesfiquation (113), can be extended for the general case .. .. ‘- - r' - I l emtgggg, 9. 21 Eli I l 21.22.29. 9. 92 :2 ---r-1--I-fi-e~d p- -- . I ,. OI'BsusI! 0 U = r ::'..‘.."..:::;..-:. :5 :9 L. : : 4 99. ----- _I. 2.! 2,, lip L. _ _ _J _ J The l, g, 1; are nxn matrices and El: and El: are n conponent coluln sub- natrices corresponding to the solution vector components along one row. hfi Griffin and Verge (1963) select a two res block of the couponcnts of g and show that the assoc1utod matrig Q is of the form required for the cyclic Chebyshev scheme. The number of rows of interior mosh points; p, must be cvono This permits the selwciifin ot t e E two row column suh~ matrices for the p rtioning of y, {“ " , 7 E 5 Q 0 2; #i T '1 - .1: 2S, .9 2‘2 Lé T 0 K L o T e F‘ - - - ., —3 - -3 (85) 9. 2.2? I: T F' L _. intJ _ht . ‘where 1!. 2 SF. 9. E. B. ,. / E = D E = 9 ET 2 9 Q :3 (1‘2" B n B I o I q , / _,._ .. _ —. - q 2 The difference equations (8h) can be expressed in the form T L - ' - T - 6 ~29’Eq E-q-l Elm-1 (8) or since §.is a 2nx2n symmetric positivexdefinite matrix, which assures the existence of the inverse Lfl, it can be written: T g gfl F' - L“l KT T - Lfl g T "q "q "" "‘ 'q-l ”gel ‘which establishes the fern required for use of the cyclic Chebyshev semi-iterative method. Assuming the solution is known at all mesh points except those in the two-row'hlock q, Enuation (86) provides the basis for the determine: tion of the remaining unknowns. The forn.of the matrix §_is conveniently simplified if the components of Tq are selected alternately fron.the two TOWBe 1:6 Take I? = {Uh};Ua,5+|2Ua,j1U2,J-Hs “““ ”‘1 UmjiUmJ‘“; ’ than 20 -8 ~3 2 1 O .3 20 2 m8 0 1 0 -3 2 20 ~8 =8 2 1 0 1 O l L,: 0\\ l\ __o o -8 l 2 O O 1 O 0 2 o -8- 1 2 5 = o O 0 O l 0 o 0 2 0 ”'5‘ 1; 2 \“ \\~-~ \ \\‘\ L v ‘\\.\ " ‘ e \~-. - \ \\ ‘ C. 2 f*o The matrix E is symmetric with nonmzero elements appearing only in _ (2nx’2n) “\as ‘~1 (2nx2n) the main diagonal and the eight adjacent bands. Hence 9 the linear system of equations (86) can be solved by direct methods. The square—root method can be used towadvantage. See p. 1.8. After obtaining a method for the solution of (86) for_one block, the full cyclic Chebyshev semiaiterative method is introduced h? ' I\(1»III)_ (1233)] (am) 0 5n; g; “($1 - TW ri’... T 51:21 I, a z <8“ TQM“): (2)11-1)+ *l2M+')_I1-I-(2m I) > ‘wlw I2R+l+ w"”+l[-I 2m: -- 2H! J’ 777 ‘ O and “(2314-2) - mm“) I (1)144!) — > LI“ =2.f:’k —KT _._—2", gig/3+, 3 ”-0, (88) am I i) um‘f' H2?" 42) (23;) ¢ -t .J’s / , “1:22 :1; wzmn[- Tzk Tzn ]’ ’ \ A The Chebyshev method consists of solving the difference equations (87) over the first two-row block, and all subsequent odd-numbered blocks, then difference equations (88) are solved over all even-numbered blocks. The iteration parameters (02 m+1’ (0sz are computed recursively using (0‘1“:[I— Eli-us] ‘For 5?. 2 share f’ is an approximation of the spectral radius of A. If _1_ is the coefficient matrix for a rectangular region, _ an approximation to the maximum eigenvalue, given by Griffin and Varga is inwéfyi where a?“ . 76: (it)? 5.1%“ + '5") + 3 will (89) h is the uniform mesh spacing in both the x and y directions and a and b are the dimensions of the plate. The relationship between the number of iterations and the value selected for f is shown in Fig. h.6, from the model problem of Fig. n.1, p. 20. he lodel Problem 2501, 8 : 000001 6 x 13 Mesh m. Time = .16 SEC/IT. Z 200, 9 E o: 150. Lu ._ #4,— — lmw La. 0 5 O' 0.. Z .99. .96 . .9“; - .95 Estimate of Spectral Radius {9 Figure 13.6 Cyclic Chebyshev semi-iterative method 3 _ Number of iterations vs. ’0 zaageaeama l direct-method algorithm for solving a system of linear equations, the square-root method, which is suggested by Paddeeva (1959) as one of the most efficient, is applicable only to symmetric systems. In Equa- tion (36) which can be written ag=a I the matrix 9 is nine-diagonal symetric and positive definite. Hence, 1;. can be expressed .as the product of two triangular matrices, one of which is the transpose of the other a=¥s —9 ' "I 811 312 313 81h 815 0 - - - _ o 823 ”an 325 826 - - - - 822 \\ \ O \ \ \\\ \\ S : 0 \ \ I. ’2n-h,2n O\\:\~.. 92n-3 ,2n “ \‘x‘ t’2n-2, 2n 0\ \\: 2n-l,2n s2n,?.n $— _ “cording to the rules for matrix multiplication the following relation- ships hold 213-51; 8132* ‘21 ’23 + " ' ' ‘ ; " +811 313' 1‘3 (9o) 111"]2.i+§21""'f"“11’ 1‘3 The I“ are ooppnted by recursive use of equations (90). Since L3,: . a4“ a s4”: . a",J . o, '1131/111“, ”133%: ZSJSh 11 31/711 “:53 ' 1) 1' '1J=MS"S“ , '1+h.>_dzi '13 .-_- 0, DJ _ The solution of _I_.. 11: g g is obtained by solving two systems fast . . éz-z Il’he components of the vector 1! are computed with the recurrence formulas Number of Iterations 50 5-: p1 p1 .— 2:_ 821 m2 . “1...... ms 9%: ’ Similarly the components of‘l are detegmined with m2!) mi - :— Sig-T2- T2n 82n,2n , T1 811 , i 2n 5000- 30004 Successive Overrelaxation X ___ Cyclic Chebyshev Semi-iterative Method n 3000“ """""' _______ Alternating Direction Implicit Metgod o 1000] 800‘ 600‘ 300? 200‘ L4 L 4‘ V fi" 0 100 200 300 LOO 500 Number of Interior points (a) Iterations ve. number of interior points Figure A.7 Comparison of Iterative methods Execution Time in Minutes 10.04 8.0~ 6.0< 5.01 Loo‘ 300‘ 2.0+ 51 r f” Successive Overrelaxation X / _____ Cyclic Chebyshev Semi-iterative Method a I Alternating Direction Implicit Method 0 * x . L L4 . A 1 4 0 1 00 200 300 A00 500 Number of Interior points (b) Machine time vs. number of interior points Figure h.7 Comparison of Iterative methods V. COIPIRISON OF THREE ITERATIVE METHODS A prime objective of this study is the identification of the rela- tive merits of three iterative methods for the solution of the biharmonic difference equations. Commenting on the problem.of selecting the best iterative method for solving engineering probless Varga (1962) p.2h5 makes several observations. (1) For the general case there are no theoretical arguments which rigorously establish the superiority of any one method. on the numerical experiments of many investigators. (2) The present evaluation of these methods has been based (3) The numerical results indicate that for each two-dimensional second-order partial differential equation boundary-value problem there is a critical mesh spacing he such that the two line cyclic Chebyshev semi-iterative method is superior for all mesh spacings h 2 he, while for h < he a.multip1e- acceleration-parameter PeacemanéRachford method is better. Table 1. Comparison of iterative subroutine characteristics I _.Y_ Number of Storage of 1 Time ‘sec.[itgratéon} subroutine FORTRAN arrays for Name Statements hSO points Mesh_pts. Mesh pts. 80R 32 one .050 .223 am 89 1221. .098 .585 CHEB 162 1126 .167 .612 80R - Successive overrelaxation method ADI - llterna'ting direction implicit method 6333- Two line cyclic Chebyshev semi-iterative method The evaluation of these iterative methods for biharmonic difference equations is more complicated.thsn for Laplace's difference equations. For Laplace's equations, Golub and‘Varga (1961) indicate that the cyclic Chebyshev semi-iterative method requires effectively no more additional 52 53 arithmetic operations or vector storage than other iterative methods. This is not true for the biharmonic difference equations as shown by the comparison in Table l of the three iterative subroutines included in the ISOPEP program of Appendix B. The overriding consideration for many is the actual machine time required for a solution which satisfies a specified convergence criterion. Convggence criteria V _ _ The norm H E“) ll II of the error vector E“) is often used for terminating the iterations; the norm is defined by H Em) “II I ._ ‘Eémfl A value is assigned to a parameter 52,- and when . . Ila“) Hus a the ‘solution vector iterate gun) is accepted as the solution. Griffin (1963) indicates that for the two-line cyclic Chebyshev semi—iterative method, if the nuber of two-line blocks is large the average difference between corresponding values of U1 in successive iterations will be approximately 2 8/1, where k is the number of mesh points at which the stress function is unknown. The criterion used in the 1801!? program of Appendix B , which includes the subroutines 8G, ADI and cm, is I‘ (I) la: 5 5 with 8. 105. for the model problem '01 |m s 3.6. The three iterative methods were used to solve the model problem of Section VI - 1. “try conditions were used for the stress function at the plate centerline. rm. is a sunny of the results in Table 2 and comparison of the methods in terms-of number of iterations and machine time required for different mesh intervals in Fig. h.7. Sh Table 2. Comparison of machine time and number of iterations required for convergence P—r Number NEST-er of E 55 lesh " of ' Iterations _ .7 Execution time, ”It“ :3: 13663:? son I 0031 %m%‘ig’fii“ 1/10 5:11 50 138 93 91 .070 .108 .116 1/12 6:13 72 19h 9s 99 .161 .156 .283 as 1/12 60:13 72 336 128 123 .280 .205 .31.2 1/20 10:20 190 636 317 - s- 1.08 1.35 - is 1/20 10:21 200 500 33h 357 .87 1.1.7 2.08 1/30 15:30 105 1981; 712 - s 7.25 6.73 - a 1/30 15:31 1.50 21.116 7M 1053 9.10 7.29 10.75 e There must be' an even number of rows for the SCSI subroutine. el- Optimum relaxation parameter used 32;: of numerical: results . Numerical experimentation can provide insight for_the appraisal of iterative methods and may provide a basis for theoretical investigations. The results obtained are based on solutions of the model problem only and conclusions should be qualified accordingly. 1 number of observa- tions are presented for consideration. For the model problem the cyclic Chebyst semi-iterative method is iteratively faster than point successive overrelaxation for all mesh sises considered. As shown in Fig. h.7(a) it is also iteratively faster than the alternating direction implicit method for mesh spacing h>l/16 9 or for less than 125 points. In terms of machine time required, suc- cessive overrelaxation is best for less than 350 mesh points or h'31/26. If h <1/26 the alternating direction implicit method is best. The rélaxation factor or acceleration parameter was selected on the basis of the" treatment given in the discussion of each of the three methods in Section IVmB. For successive overrelaxation an initial value “)0 g 1.5 was selected, and this was changed after every 10 iterations. As shown by the second and third problems in Table 2 this will not necessarily assure a good approximation of the optimum value ofw-b. When the cyclic change was used for the 6x13 mesh size problem approxia mately 70% more iterations were required then when the optimum re] exam tion factor was used. Another indication of the variation which can be expected when using this procedure is shown by the scatter of points in Figures ho7(a) and (b) for the 50R subroutine. This contrasts with the curves fitted to the points for the other two methods. It should be noted that a random selection of values for a) will produce greater charge in the number of iterations required for convergence for the successive cverrelaxation method than for either of the other methods. Comparison of Figures 11.3, 11.11 and 11.5 shows that the alternatingmdirection implicit method and the cyclic Chebyshev semiaiterative method do not impose as severe a penalty on overestimation of the appropriate parameter as does point successive overrelaxation. For a considerable range of values above the optimum, the number of iterations increases only slightly above the minimum for these two methods , while the minimum in Fig. 11.3 for point successive overrelaxation is much sharper. Another important consideration is the computer storage required by each of the three subroutines. The data in Table 1 clearly identifies the successive overrelaxation method as the one which requires least 'storage for both the arrays of data and the sequence of instructions. Because of the relative simplicity and minimal storage requirements the successive overrelaxation method was incorporated as the main iterative method in the ISOPEP program. 56 Since the number of iterations required for convergence appears to increase significantly with the number of mesh points, attention was directed toward the possibility of improving the rate of convergence by first solving the problem with s relatively coarse mesh, then interpolating between thee vslues to generate a better initial solution vector for s finer mesh. some of the results are given in Table 3. mu 3. ' 01'. .. elv t ISOPEP ‘0'. l subroutine CHANG was written for this interpolation and ' Estimate Prob. Converg. Ilesh Number of Time Time No. criterion sise Iterations ‘75-‘15?- LHours) of f 101 10'"5 5x10 118 3.31 .9571 102 10'5 10x20 521 51.9 .9801; 103 10"5 15x30 15h8 5 38.3 .99112 201 10‘5 3x6 39 .85 .719h 200 10-5 6x12 811 2.62 .9556 203 10-5 12x21. 86 6.31 .9906 205 10"5 7x18 276 16.9 ‘1 .9737 206 10-6 11x28 1551 5»- 29.7 ‘7 .9900 303 10"5 6x12 225 2.33 .9506 306 10-6 21:18 1038 h 1.9 .9793 801 10-5 6x12 179 2 .30 .9670 I102 10-5 12x21. 126 5.67 .96118 has 10"6 2141118 166 1 3h.7 : .8862 501 10-5 6x12 135 I 1.80 .9h76 502 10"5 121211 295 12.59 .9781 506 10"6 2111118 1157 5 117 .6 .9858 601 1:15 15x15 176 1 20 .9802 601 10'7 name 7072 57 _ 11.2 .9998 57 Problems 101, 102 and 103 in Table 3 were solved using the subroutine CHANG. The final solution of Problem 103 required a total time oi” 6 minutes, 311 seconds and s. total of 2187 iterations for the three problems. The same problem solved with initial values of the vector iterate set to .1, required 7 minutes 15 seconds and 19118 iterations as listed in Table 2. The other sets of problems listed in Table 3, such as 1101, 1102, hoe, were solved using the iterpolation procedure as the mesh was refined. For Problem 1108. only 166 iterations were required after the last refinem- ment of the ash while for Problem 506, which has the same number of interior points, 1157 iterstions were needed. Both of these problems have more than twice the number of mesh points used in Problem 103 yet the solution time is approximately the same for 506 and substantially less for 1108. mobilities g; the Computers used - Two computers have been used in‘this investigation. One, the Control nets Corporation, 3000 at llichigsn _State University Computing Center, is very fast and has a large main memory of 32,000 words of 12 decimal digits eech. The floating. point multiplication of two numbers with lO-digit mentissas requires less than six microseconds. The other computer is the 18! 1620. at the University of, Toledo Computation Center. The 1620 performs e. floating-point multiplication of two numbers with eight-digit nentisses in approximately 12 milliseconds. The main wry provides for storsge of 20,000 decimal digits which is about 1/19 of the 3600 memory. However. the min memory is supplemented with en 13! 1311 Disk Pile which-provides 2,000,000 decimal digits of secondary storage. Several problems were solved on» both computers end e. comparison is provided in Table 11. The slight difference in the number of iterations 58 required for convergence is possibly a consequence of the difference in the mantissas of the floating point numbers used. The difference in the speed of the two computers contrasts sharply. In fact for problems with a hundred or more mesh points the cm 3600 is approximately 1200 times faster than the 13]! 1620. Lt rates of 3375.00 and 830.00 per hour for the 3600 and 1620 respectively, the cost of solving problem 1.102 would be $8.82 on the 3600 and $572.10 on the 1620. m qualifications of the comparison should be noted. First, the time required for compiling a representative combination of the ISOPEP W subroutines on the 3600 uses between 20 and 30 seconds. Compilation time on the 1620 ranges between 10 and 15 minutes. Forthis purpose the 3600 is only about 30 times faster than the 1620. A second consideration is the time saving possible on the 1620 through use of a machine language successive overrelaxation subroutine. The introduction of such a sub- routine for the iterative solution of Laplace' s equstion resulted in a 50$ reduction of machine time required. mu. 1;. Comparison of the 130m program for the IE! 1620 and one 3600 Computers Number of CDC 5635 Tfie IBM 1355 '— Prob. flesh , Iterations Total beaution Time No. Size 3600 1620 Min. Sec. Kin. Sec. Hours lethod 1.101 51:10 118 116 25 3.35 I .97 SOR 2 .101 5x10 90 9O 29 6.0 l . 78 ADI 1.102 10x20 636 622 l 25 1 5.0 19.07 son 2.102 10x20 317 311 2 2 1 21. 27.20 ADI 1.103 15:30 1981; 19714 7 36 7 15. 158.00 son VI. SOLUTIONS OF PLANE ELASTOSTATIC PROBLEMS The second of three objectives of this dissertation is a general computer program for the iterative solution of plane elastostatic problems. The program ISOPEP includes six subroutines which contribute to this objective. The program is not completely general, however, since sub~= routines _for the boundary conditions must be added for any specific problem. 1 flow diagram, listings of FORTRAN source decks, description of input preparation, and the output of a sample problem are provided in Appendix B. The analysis of physical phenomena with mathematical models, the analytical or numerical solutions of the models and the subsequent test» ing by comparison of the solutions with experimental observations, consti- tute major pursuits of scientists and engineers. The creation of the ISOPEP program contributes to the analysis of plane elastostatic problems and this computer program is considered a major part of the dissertation. Duonstration of the utility of iterative methods for the solution of several problems is the third objective of this dissertation. Problems number two and six are included because analytical solutions are available for comparison with the numerical solutions. Problems three, four and five are representative of the practical problems for which analytical solutions are not available. The treatment of boundary conditions for irregular regions is included in the discussion of problem three. 1. T_h__e 29331 Eoblem . The square plate and edge loads as shown in Fig. 6.1 are synetrical with respect to the x—axis. By proper choice of constants of integration the stress function will also be symmetrical with respect to the x-axis. Starting at the origin and proceding counter-clockwise the boundary conditions are obtained from Equations (11;) and the requirements for 37.91370 59 60 Ali F e HP - ‘T ______________ a V ’a 3 *— ‘——« fl-d ”—1 k—4 .44 H1 5 c 1 a: a 7 (x Figure 6.1 Boundaries of the model problem The constants of integration 1, B, C, D. E, and F in the following equations are determined so as to produce the desired symmetry of 0 and consistent values of ¢ at each of the points b, c, d, e, f, g in the functional expressions for 95 along neighboring line segments. Two of the constants are arbitrarily chosen to make ¢ 3 0 and if: 0 at b. Along ‘53 a“? ' a _ nay-'0' 3945-0; along-Ed. ‘ §§=—ry+(3$)c=-m, gamma 4) = 'P92/2 ”7 : ‘Pyz/a ‘IzIPcz1 ,‘ alongd—e %§=(}$)J= “'44P," ‘3‘?"0 . _._ 45: "9QP9+B 2 ‘— 1+pr + o/anP 2' along ef M= (M)e=—,‘+a¢>5 =0 59 as ) <1) I: C 2 (436') : -o.oaa‘p; 61 was?! §$=*‘+PV+D $2.0, 4: = ~2ry‘+Dy+E = ~27°y1+126dry‘-32“a?; maxi? .. ad a .. ... ‘37: (5)3“0’ 53?" ’ ¢=F.= 00 Values of the difference-fequation approximation of the stress function 013 and nodal point values of thenormal stress component 0:11 are given in Fig. 6.2. The point values of the stress function are approximated 2 3...... — where C -'-’ 32 . The constants p and a have been assigned values of unity. 2. gemi-infinite 23333 33333 uniform 393.3 93 33:3 156mm" ' ‘ The exact solution of the problem of a semi-infinite plate with a uniformly distributed load of intensity 'Po acting on an interval of the bemdary -e 5 x 5 c, is given by Timoshenlao and 00on (1951). The stress function has the form 9 m l (r20 - r191) lumerical solutions of this probl- have been reported by Veyo and Hormbeck (1961:), who used point successive overrelaiation of the biharmonic finite-difference equations, and Pisacane end lalvern (1963), who used _a_ numerical mapping technique for application of the luskhelishvili complex variable method. For comparison with the other numerical solutions a . square region, 2 units by 2 units, of the semi-infinte plate will be considered. The intensity of the distributed load P0 is chosen as {one and the value of c isfi. For the assigned values of c and Po the stress {motion can be written 62 FIGURE 6.2 Dlstributlone of the Stress Function and 0. for Problem 1 COLUMN 0 0 ' 1 0 7 C 1 10 11 13 11 :‘00000 101110 1.0300 3.1910 20.000 3.0710 109.00 1.1990 07200 000000 07300 11. «210'u"1106'"'1100'"°1165"332566°.3210."Hudunflod""1210"00000.“:«100 6: 31.1101 1.1111 1.1910 1.1120 2.0111 2.0129 1.9100 1.1111 .1120 .0021 :1200 -5,“ 000:,“ 000:}.‘,00: £5461...:!;6 ...:5“‘..-:,z‘,...:}’i’..::‘16‘..: 28“’OOOE°1.’ 1.1119 1.0111 1:1102 1:1110 2:0001 2.0012 1.9101 1.1011 .0911 -.0001 .1200 ‘.32tt133:3221233:32221.33 32112333: 12113313221133z3211133: 31.2533z3111233: 3.1113333..:. 31.0202 1.1020 1:2121 1:0101 2:1290 2.1112 1.0199 1:1021 .0101 -.0002 .1200 '-32:1233:31:1133:31:2133:: 212233: 31:1133zz212133: 3221;33:3112133: 331211333321113333.... 1.1100 1: 2120 1:1000 2:9111 2:0001 2.2010 1:0000 1:2010 :0121 - :0120 .1200 3-3111 3333211133:321:133zz211133zzzss 3333121 33:3111133: 3112 33::112133::1121333:.7:s 131010 1:1010 1:0211 220210 251010 2:1190 1:1110 1:2100 30100 -:0190 .1200 ’ -:122233=3122133::211133::112133::222133: 3J12a233 321113333121133::112133z3112133331... 1:0101 1:0000 250019 250992 250119 2.0912 1:0110 1:1100 :0000 0.0219 .1200 '-3221233=3121133:3112133z3111133::212233::212233::22213333222133zzxtst333311113333:... ’2.0100 2.0010 2.1190 2.1000 2.1102 1.9919 1.0020 1.1119 .1002 -.0122 .1200 3J1115" 3111233:3111132 3221533331.: 33:3121 33:3222133z3221233:3122 333311113333:1ez 2:1012 2:0120 2.1111 2:0109 2:1911 1:0919 1:1290 1.0010 .1010 -.0101 :1200 3° -3122133: 3112233:3111133: 3111;33:3121233zz2o2s32 3212133::212 33:32:513232221333:21.. 2.1121 2:0009 2:0011 2:2001 2:0111 1.1090 1.0100 130110 .1102 —.0011 S1100 “.:121133::112133::115132H131233 3112533z3111232 3212133=:112133=:121133::1151333:25.: 231000 252000 2:2090 2:0001 1:9090 1:0100 1:1019 :9000 :1109 -:0110 :1200 3.3111133::111133:311.133z314212311113 :31121 3:3112133z3112133:32.1:33:32as23333...1 2.0110 2.0101 1:990? 1:0902 1§1011 1:1192 1:2010 :9109 :0002 -.0021 .1200 " 3121 "131116"131111"131516"1:ii. 3313111133:3211232 “21523313212233 322133333.... 1.0110 1.0102 1:1091 1:0010 1:1011 1:1912 1:1191 :0091 20010 -.0121 £1200 3-:221233zzat1133::222233z3111132 3111133::125133::212z33::121233::111133331112333z.... 1:1009 1:1100 151219 1:0011 1:1010 1:2110 130119 E1090 E0002 ~50000 :1200 " 3222233::222133::211133::.21133:3111133: 31112333:22.133°3122133::122133::21223"3.... 1:2011 1:2001 1:2101 1:2111 1:1121 1:0120 :9010 E0190 .1101 -:0902 :1200 3‘ 3211233=3222233z3211133:32222'333222233 :3112133:3121233:3111132 3121233:::21 3333.... :9990 :9902 :9001 :9010 :9221 :0110 :1091 :1111 .1010 -:1111 .1200 "-:222233=:211133::212233:322123 :3221233:3112 33:3111133:3211 33:3112 33:312213333.... .1100 .1110 51101 .1010 :0001 .0091 :1011 .0190 .2110 -.1110 :1200 *' -32112323222133::211133:311:233z32.1133:.3211233::111233::212233::211233::11123333.... .000. ...90 3.002 ..soo ....7 ..301 ..072 .131. .1..a -3191: 21200 *’ :66!!°":6611"::6315'"36611"136111"=36111":31119"531119'231111°'1§3‘11'":.122 .2221 .2210 :2219 :2100 :2100 .2119 g2122 .1991 .0910 -.1102 :1200 1°. :6“‘000066‘,000: .6“ 00006650000066. 000:6" 00006,; 000:,6‘ 000:“: 000:,‘§ 000:’,°, .0011 :0011 :0010 $0011 :0110 €0101 :0010 :0109 .0101 -:1091 51200 3‘ 3.1213333221133332.113333221233332.2233:. 3251133zza12133zz112133:::11233::22213313.... 22 :0000 0.0000 0.0000 0.0000 0.0000 0.0000 0:0000 0:0000 0.0000 -.1000 .1200 38666'°6:6666"836666°'616666"6:6666° 6.6566"5.6666"6.6556'°"1565°°1'6666" toooo 2 -== #2[{(x-2)‘+y3}ta2" iii/3*) 3‘ {02011232714 1%)“g *5"? Values of the stress .. - _- ,00‘0 function on the bound» (DP arise of the square region which are inside of the semi-infinite plate can be determined from the exact solution . The values of 91 .2 points one nesh interval outside the square region can be similarly computed or T__’_""' ‘1 Figure 603 derivatives of the exact 3_ Problem 2 - Semi-infinite plate obtained from the normal solution at the boundary. Due to symmetry it is possible to determine the solution by considering one half the square region. Along the y-axis which 1. the .111. of sy-etry, 211. nor-.1 derivative, g} g o, is used to determine the values of 0 at points exterior to the solution 0min. Alon; 3310 use Equations (11:) and find 2. .- at %—$="R“" 373%:‘0: _. ... 9 _. '33?- X+B° 793412.39 (15 = -£X2+Bx+ C . Sy-etry of ¢ with respect to the y-axis requires 3.0. Two of the constants of integration may be arbitrarily selected. Let A u. c n 0. 6192 Thonnlongo—n- §2Q=-X, c$13$=09 ¢=-PLX2' é¥=(%¥)a=-’4o §%=(§$)Q=O2 ‘ ¢=—j¢x+D=-hx+V32 Thovolno o£D=1/32 numwmo row-mt gag-U32 intho 1103135 moo-1m m- ¢ um. "53 moi-5'. Iota tho hat to tumor tho mu Motion - é- : ond 1/3! on nooouory on n oonooquonoo of tho arbitrary oohotion of I and c. _ Tho vnlnon of tho nor-:1 oonponont of mono 0")" nnd tho ntrooo {motion nro 31m in 11:. 5.5 for 9 M8 nub. rho oonvorxonoo , oritorion nood m 10", nnd on mum «unto at tho solution notor tor tho 111128 noon m obtunod fron intorpolotion bottom tho oolntion moor M for o 7:11: noon. Tho total nubor of notation roquirod; V’ PH .1-» 02‘ ' y - 1/19, h - t/n. Exact Solution V Numerical Hopping 0 Succoooivo Ovorrolmtion of Finite-differenco Equations m. 6.1: cup-risen of mnoriool uni mot 50le for Problon t 65 for both problou too 1830. A flat distribution of 0.1 ot an interior pointoyoo {nod for tho initiol ootinoto for the 7x11; nooh. Voyo and Bamboo]: roportod opproxinotoly 3000 itorotiono woro roqnixod for tho 111128 looh with on initinl volno of -.0625 oooignod 'ot oll intorior point- ond o oonvoo-gonoo oritorion of .5210'7. though point voluo of tho oolntion of tho problon oolvod by Voyo ond Hombook voro not ovoilnblo for oonporioon, it no not poooiblo to diocorn ow oignifioont littora- onoo in tho rooulto moon in Fix. 6.14 Ind oinilor grnphiool rooulto in thoir "poflo _ _ > . rho «notion of tho itorotivo solution at tho mot. «liftoronoo oquotim iron tho mot oolntion in greatest in thooo. rogiono whoro tho otrooo groiionto oro motoot. Tho otrooo nodionto oro lone ot tho tint row of intorior pointo_ bola! tho boundu-y onbjootod to tho diatributod load. This dotribntion of 0'? for tho onct ind finite -— ution-moo solution! at y» l/lh io ohoon in 11¢. 6.14. Tho nonoriool uppu. oolntiono 0,! mooono no mm no o1oo ohooh. It is opporoht thot tho um upping toohniqno .providoo only ought]; bottor roonlto in tho rogion oirootly mdor tho oppuoo zooo thou-o rtrooo' notionto *oro highoot ond noh pooror roonlto .noor tho bonnduy x g 1. Both nuoriool nothodo provido bottor ”promotion to tho onot oolntiono ot ran of nooh pointo oorrooponung to y a ah. Problonob Egg-Iota“ tonsilom. rho problolo of prootiool intoroot to tho onginoor froqnontly hovo more!» bond-tion Ind ottontion in ofton oontorod on Aotrono oonoontron- tion. on innotuotion o: otrooo oonoontrotiono hoo hooo oohomtoo by oolving oooh notohod tonailo opooi-on problon for ditioront nooh mom Throo diftoront tnoo of notohoo homo boon oonoidorodx 66 FIGURE 6.5 Distributions of the Stress Function and cry for Problem 2 COLUMN Row 3 6 5 5 7 5 9 1o 11 12 13 16 15 16 17 .00000 -.ooz55 -.o1ozo - .02295 - .06017 -.05503 -.o7589 - .09575 -.11156 - .1296o .16732 -.1o517 -.15355 -.2ooo9 -.21e75 ~116666°1i16666'1116666 11165666 11666661116666"116666'°116666 1116666°'116666°°'16666 '1" .0666 11:666611116666"11oooo . C O : z . . : C I z . : C C O z : :600251 -600506 -601256 -602682 -606086 -605821 - 607593 -609376 -611161 -612966 61o132 -616518 -610306 -1zoo9o -62101? .6:6666' :66666'::66666'::66666'::16 616'::66666'::6166 '::66666'::66661'::66666' ':66666'::66666'::66666":66666"66.... o o o o o o o : 6 o 6 o o o o o o O O o 0 -100522 -1o1ooo -1o1755 ~15259o -1otato -105977 -1o7o9o -1o9ooo -111zoo -1129eo 116759 -:1?5:o ~115521 -1zo1o7 o121591 1166666°1166166'1166666°1166666'1166565'111696 '1165665'1166666'1166666'1166666"166666"-16661 1166666'11666661'160259 ~1o1555 -1o1505 -1ozoos -1o3555 -1oo773 -1ooza. -1o7910 -1o9ooo -111327 -115075 116556 -11oooz -115575 -1zo151 -121925 1H66665' '.6666 '11 6666'1166666'116666i°1166fi '111666 '1165656'116666 '1166666"161666'1166666'1166666'1166666"156011 O 0 -102672 -102657 -1o:235 -1o~1oo -10533o -1oo725 -105255 -1o9oo7 -111535 -1152oo 116972 -115717 -115671 -126253 -121999 1H56666':W66666°':66666'1166656'1166566'1166566'11i56i6°1116665'1166666'1166665"166666'1166666'1166665'1166665"-1oo~5o -1o3650 -1o:oo3 -1oo116 ~1o5951 -1ooooo -1o7253 ~1057oe -:10239 -111562 -115696 115152 -116593 -115521 -126352 -122112 116666 '°'6o666°"66666"'66566‘1166566°"‘6656'116056 '1116666°116655 '116656 "16656 '1166651'116666 ’1161566'1161155 3105655 -1o¢o9o -1osooo -1o5755 -1oo7oo -1o793o -1o9252 -11o7oo ~11223o -113525 115655 -117155 -115556 -1ao5oo -122271 0 oo o oo o o Ooo I o o o ooo O O -166666'1W666'1E66666 116566611135676'-.66166'1166766'1:15566'1E16666'11o6566°°:66665 -.66616 11666661-166666'1oozozs O 0 : : : : I : : : : : : 0 -165571 -1o5o11 ~105025 -1oooo9 -1o7551 -105562 -1o9599 -111259 -112713 -115260 115522 -117oo5 -1191o1 -1ao751 ~122651 ¥6656 '116666 '1166566'3 6666°1165fi '116666 '1163 6'116666 '1166666°11i656 "16666 '1166656'116666 '1166566°1Eoz969 -1oo529 -105555 -1o7o3o -1o7537 -105653 -1o9655 ~11ooo5 -111556 -113256 ~11~725 115265 -117521 -119o55 -121o7o -12276: 111666662 H66666'166666'1161166'1166566'116616 '1166H66'116666'1H665i5'166666'1166566'1165666 1166666'1166666'1153906 “607606 -60771, '600060 “600616 -60936I -61029§ -611370 '612573 ~613002 '615276 6167~0 -618260 -619825 -6Il¢20 '623056 11166666216661121616662166666 1166661'1166666'1166656'1166665'1116666‘111666 °°166666'1166656'1165166'1166656'1166561 C 0 O z 0 O O O O C : : : : : z z : z : : O O C '00.6’0 -.00196 'o09609 -o09620 '610315 'o16176 -.1218~ 'ol3317 -olfi557 'o65..7 567296 '010797 -oZOZTQ ':360’2 '623525 1166666 111665 6 11166666111666661166666 1116666 '11 656611166656'1H6566111666661116661611166665111666661116666611605729 'o097.’ .009.33 -o10172 ‘:106‘5 -066290 'o62095 'o13030 'olhlo.’ 'o15233 '51‘556 66.7097 'o69309 '02077. ':':390 ':230‘3 11166666116 6566611 W66 '116355611666665 1165 6666516666611166“ 11666656116666661166266911666662" ~169166 16655661166555? ~110555 —110975 -11126o -111ooo -112257 -115039 -115927 -616936 -11oo53 -117252 115552 -119911 -121329 -122797 -126369 1116666621 665 ’1166666'1166566'1166656'11 6666'116666 '11 66662116666'1H1656 "116666’111166 '116666 '1166666’1107255 ‘611995 0612079 '612329 -6127§O -613303 -61b008 -61§Ib§ -615800 -616.60 -610015 619251 -620550 -621921 -6ll)50 -626020 1,:ooo‘oooooooooooo ...5...’;;I...:5.6......... o-o-ziooooooo oooo ooo ooooooooo Ooooo ooooooo - 325 6 -.32505 -.6165 3 -. 7931 -.2557 3261 -.6656o -:65506'1 166661 .13965 1116166 -.1o555 -.o9I66°11o793o O U C C . C C : : : : : : O : : : I O O : O 1:115105 -113155 -113626 -113505 -11o333 -11o997 -115757 -11oo92 -1177oz -115566 119959 -121265 —1225?? -123965 -12537o oo o o o o o o o lo -1665662 :6666 '11 6666'1166663'1166666°1166666 1166666°1.66695'1.ii666'-1i66 '6 '166655'1166655 -111666 1166666" -1oo~95 . z : : : : O O z I O O 6 I : 1;61¢226 -61¢297 -616517 -61001I -615376 -616003 ~616750 -617610 -618571 -619625 620161 -621971 -623267 -636502 -625960 o o ooo o oo o ooo o 0 -’66666'- .6666! 1165836‘-.66566".26§£6-:25706'-.11666'- .16663'1215617".16167"116.76 1116662 -111666' 1116161’1165976 O 6 : : : 6 6 6 I z : : : : 1:61936! -615662 -615610 -615999 -616629 ~617022 '617732 -6105b9 -619661 -620h15 621565 -622730 ~623961 -63525) -626590 00 000000 000 00000000... ooooooooooooooo oooooo ooooooo oozi ooooooo oOooooooooooobooo o 000 1165165 -:26616 -:666 95 -:2555o -.65155 -.22755 -:61666 1116521 -:i7506 -:66666 io5oo -.631oo -.11756 -.1o6oI -.o93oo o o o o o o z 6 6 o o o o o o o o o o o '6IOQOQ ‘66‘529 -61672Q '6170Q6 -617¢91 -610056 -66|729 -61950. '62038b -621350 622397 -6Z3§18 -6ZQ707 '635950 '621263 1116666 “116666 '1166666'116666 '116665 '116165 ’116661 '1116666°1165666'11 6666"116666'1116666'1111666'1116561°11o9557 : : : : : : z z : : O : z : : 23117557 -117559 -117o33 -115139 -115551 -119o97 -119759 -126655 -121521 -122257 125252 -12~556 ~125652 -1aoo95 -127959 0 o. oo oo o o ooo o o no. oo 0 I o o 1666662166666'1166666'1166963 116666611.2ios6'1.l6fi '1.655352 165676'1116666°'. 15906 -166666 1112666 116066611109963 0 C . : : : : z z : : : z : : : O . -11o711 -11577o -1159oo -119255 -119635 ~620168 ~120711 ~121o72 -122275 ~123155 125132 -12517o ~125256 ~127655 -125556 1116666 111666661116 26661116666511166666111 6666'116616 11166665°11666i6°1166665'°166616°11162551116666611166666'1110139 ~619l37 -619093 -620060 -6203!7 -620720 ~621207 ~121793 -122675 -623255 -626090 625036 -626036 -621109 -63|2§4 -629536 ‘3666666': :61661'::66 666'::66666'::66666'::66666'::16666'::66666'::16666': :66666":16666'::66666'::66666'::16616'::1.... ~62096¢ -621010 -621111 -621§41 -621001 -622273 -622036 -623687 -625226 -625050 625941 -626910 -627955 -619055 -630212 '1666666': 66666'::16666"'1o616"'16.66"'116.6' 1666 "'16 66':: 66:11' 16666"'66666':: 6666'::66666'::66666'::..m -622092 -62214! -622295 -622960 -622090 -623365 -623082 -62b509 -625221 -626012 626879 ~621|17 -625822 -629OIO -631012 22:66616':: 666'" 6666'::16166'::66666'::6666'::16661': 16666"'1.666'::6666"'1616""16.16'::66111' :11111'::..... '623221 -62326’ '623QIS '623657 -623993 -62#¢21 -626930 -6225Q1 -626226 -626969 627026 -6ZI713 -629706 -63o160 -631036 z116611621666662166616' 11i6666"—116666 1115656 1166666" -1i6695 '1116665 '1116696 '166666'1116 M6 6'11 666511166666'1116571 O C C : : z : : O O . C O 6 : : : : z-62Q’SO -62§3’6 '62b536 -62Q70. -615091 -625502 ~626000 -626500 -627261 -627977 62.7.6 ~62966¢ -630607 -631660 -632671 0.0.0....0.0.0.0....OOOIOOOOOOOIOOOOOOOOI......OOOOOIOOOOOOOOOOOOOOOOOOOCOOOOO0.0000............OCOIOOOOCIOOCOC H. U 51y 67 Problen 3 _- Seu1circular notches, Problem 11 - V-notches and Problem 5 =— Rectsnguler notches. Sketches of the specimens are included in Figures 6.10, 6.11, 6.120 E 1-66---” - --- ----g-o—H 23--------—---—-—-—7n H Figure 6.6 Solution domain for Problem 3, h and 5 Boundary conditions for the three problems are identical. Intro- ducing e Cartesian coordinste system as shown in Fig. 6.6 integration constants will be chosen so the stress function is symetricel with respect to the x and y axes end it is necessary to solve: each problem for only one quarter of the whole plate. It theend '53 e uniformly dis.- tributed tensile loed of Po units intensity is applied. Using Equations (11;) we find “>461 94>... '::FI, E“Ey+(%$l,’ But (6'35 3 0, if d is sy-etric with respect to the x-nds; hence 43= B‘f‘é + B = Elf/a. The constants L and a are arbitrarily selected so sero. Then along 671’ 93—)? = 0 5 ’3'? : 24 7; I ¢ Z ¢C : 2 a2 P0 ‘ Along the ho axes of sy-etry the nonel derivative conditions are used to compute values of the stress function in two mesh ran exterior to 27 .1250 .1253 25 .1091 .1061 23 .0833 .0031. i£&§'°‘255ii"°2655i"'2 21 .0626 .0632 19 .0“26 .0830 17 .0012 .0239 .0258 COLUMN I 20002 :0010 :0070 I009. I 0.0000 0.0000 13 -.0182 ‘I0179 —.0170 ‘I0156 ‘.007‘ . 'I0003 11 -.0321 315"! 3 I I . . ‘.0205 I I -I0159 -2019. 3&1§"i I I -:0317 -Ioaos I ‘.0250 : -:oazo -Io.ov I -Iozes I I -.0305 -I0025 0.1. 0355 3357"i . '.0316 I I . I I I 'I0056 . I ’.0992 ‘.0497 -10355 -10399 ”1:57.56”125356"i “.0539 : -I°53“ I ‘.0520 I I I ’.0690 “.0553 ‘I0560 '.053I -.0512 '.0609 I I -0036. I I ROH FIGURE 6.7 Distributions of the Stress Function and o, for Problem 3 0 6 0 b. 0 1. 0 1. 0 5 0 fly 0 6 0 0 0 3 0 7 0 2. 0 6 0 7 0 by 0 9 0 6 6. 3 Ra 0 E4 2 K. 0 ‘4 1 Eu 0 5. 6 5' 6 G4 6 H. 6 S. 0 fl. 1 ‘5 6 K: 6 ‘3 2 I: 3 2 2 .2 I 2 7 .2 9 2 1. 2 1 2 9 2 6 2 2 2 7 2 1. 2 5 .2 7 2 9 .2 0 .2 0 1 0 .1 1 .1 2 1 I. 1 6. 1 6 1 6 .1 7 1 0 1 0 1 9 1 9 .1 9 1 9 .1 0 .1 0 I II III III III II III II II II III III III III II III I. III III II III III III II II III II III III III II III.II III II III III I I I I I I I I I I I I I I I 61 I1 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 2.9 3.0 ‘.6 SI“ .6I1. 7.6 7.6. 0.1; 0.9 900v 9Ifiv 9IR‘ OIO- OI? Olfiv 006 6.0 6.3 6..“ 6.0 6.6 6.0 6.7 6.5 6.6 6.7 6.2 6.0 5.3 5.9 5.1 5.7 0.9 0.nv 0.0. 0.9. 0.0 0.6 0.2 0.0 0.3 0.7 0.1. 0.6. 0I6 0I7 0.9. 0.9 1.1 1.3 1.6. 1.6 1.5 1.6 1.7 1.7 1.0 1.0 1.9 1.9 1.9 1.9 1.9 1.9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII I I I I I I I I I I I I I I I I I I I I I . . I I . I I . I I . I I I I I I I I I I I I I I I I 0.3 1.7 5.7 9.2 2.6 5.1 0.9 0.9 2.5 3.9 5..» 6.6 6.3 7.6 7.... 0.3 3.0 6.6 6.5 6.6 5.5 5.0 5.0 6.1 6.9 6.1 6.9 6.2 6.1 6.6 6.0 6.6 0.7 0.6 0.6 0.1 0.7 0.3 0.7 0.2 0.5 0.9 0.1 0.6 0.6 0.7 0.0 0.9 0.3 0.6 0.5 0.6 0.6 0.7 0.7 0.0 0.0 0.0 0.9 0.9 0.9 0.9 0.9 0.9 IIIIIIIIIIIIIIIIIIIII‘IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII.IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII I I I I I I I I I I I I I I I I I I I I I I I I I I . I I I I I I I I I I I I I I I I I I I I I 0.6 0.9 6.3 6.7 02 6.3 2.3 6.3 0.6 3.7 62 0.... 0.6 1.0 2.2 3.7 6.0 6.6 5.5 6.7 7.5 7.9 8.0 8.3 9.3 9.5 9.5 9.2 0.6 0.9 0.9 0.6 60 6.6 6.9 6.3 67 6.0 6.5 6.6 6.9 6.1 6.3 6.5 7.6 7.7 7.0 7.9 05 0.6 0.6 0.7 07 0.0 0.0 0.0 0.0 0.9 0.9 0.9 0.9 0.9 09 0.9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 209 6.2 9.0 1.1 2.9 2.6 0.6 7.6 6.6 9.6 6.0 7.2 0.7 3.3 67 5.6 5.0 6.0 7.6 9.2 0.6 1.1 2.5 2.9 3.2 3.6 6.5 6.6 5.5 5.6 5.1 5.7 I :1 6 2’ 6 I“ I I6 5 0! $5.9 81.0 5 :1 5 :5 5 7“ 5 33 filob ‘4I7 K'Ia 5 XV ‘4I9 0.0 0.8 0.0 0.0 0.6 0.8 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII.IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 0 35 0 :1 6 30 11.0 9 36 91.7 6 :9 E’I6 6..“ 5.01 0 :2 14.7 “6.6 1.06 3 l5 hvo6 7 .3 9.3 1.9 3.0 6.5 6.2 7 .1 0.2 9.6 0.7 0.0 1.3 1.7 2.1 2.5 2.0 .2.2 .2Io .3.0 15.0 3.7 3.7, 3.7 3.7 3.7 5.7 “.0 “.0 .0.6 .9.9 6 A, 6.9 0.0 0.0 0.9 0.9 0.9 0.9 0 .9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 0.9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII I1 I1 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 3 I“ 0 23 1..5 1..3 14.9 flvoo 6 96 9 :1 I.I0 I.I9 9 23 6 :5 5.I7 6 :1 0 22 1..7 22 6.7. 7.2 9 .5 1.3 3.7 B .3 5.3 7.5 0.0 0..» 9.1 0.9 0.9 1.9 1.9 1 .9 .....u. 1.1 1 .0 2.6 2.8 2.4. 2.2 2.1 2.0 2.0 2.0 3.9 3.9 39 3.9 0.1 0.1 0.1 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.9 0.9 09 0.9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIOIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII .1. Il IL I1 .1 I7. I1 I11 I1 III III. I1 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 2.... 7.6 5.1. 1.0 5.3 6.0 5.5 1.9 5.5 0.9 0.7 7.0 6.0 9.0 3.0 6.0 1.2 1.6 “.0 7 .1 9.9 1.1 3.8 5.8 6.2 7.8 0.7 9.0 0.2 0.6 1 .3 1.0 001 0.5 0.1 0I7 003 1.1 1.8 1.6 1I5 1.3 1I2 .101 «(III 2.0 2 .0 2.0 0.3 0.2 0.2 0.1 0.1 0.1. 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII . I1 I1 I1 .1 71 I1 .1 .1 I1 .1 I1 I1 I1 21 '1 I1 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 6 .0 1.6 9.5 9 .9 2.5 1.7 3.... 2.6 0.6 2.7 5.0 5.0 3.6 0.7 I .1 7.7 2.0 9.2 5.2 2.7 0.8 2.3 6.3 6.6 7.3 9.3 0.7 1.3 2.1 3.2 3.6 3.1 IIIB o :5 o :0 o :5 nUIg nUI6 o :5 o :U o :0 AVIb 7.I~ 1.I3 IIIZ \IIl 7II0 1.I° 0.3 0.3 0.2 0.2 0.1 0.1 .0 .1 0.1 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII .01 _.1 ..1 _.1 . ‘1 .1 .1 I1 .1 .1 .1 I1 .1 I1 ,1 I1 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 3.... 6.9 1.2 9.6 0.5 3.2 9.7 0.0 0.0 6.6 9.6 1.6 0.6 0.6 3.1. 6.1 1.6 7.1 6.2 0.5 3.2 5.3 2.7 0.6 1.0 2.3 3.2 5.6 6.9 67 7.0 7.2 2.3 1.0 1.3 1.0 0.6 0.0 0.6 .025 0.0 0.0 0.6 0.4 0.2 0.1 0.0 0.0 0.8 0.3 0.3 0.2 0.2 0.2 0.1 0.1 0.1 0.0 0.0 0.0 0.0 0.0 0.0 0.0 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII ..1 ..1 _.1 ..1 .01 ..1 _.1 _I1 I1 .1 .1 .1 I1 91 '1 I1 I I . . . . . I I I I I I I I I I I I I I I I I I I I I I I I I 6.1 0.0 1.5 7.0 5.5 7.9 1.6 0.6 0.3 1.6 6.6 5.5 5.0 3.6 9.6 3.9 7.0 3.1 0.3 6.7 3.2 0.0 0.2 5.7 3.6 2.0 0.6 0.3 1.6 2.1 2.0 3.2 2.5 2.1 2I6 1 ..1 1.7 1.3 0I9 OIS OIZ °I9 0I7 OI.) 0.3 002 001 OIO 0.6 0.4 0.3 0.3 0.2 0.). 0.1 0.1 0.... 0.0 0.0 0.0 0.0 0.0 0.0 0.0 II....II.....I.....II..I......I...............II.._.III..........IIII.....I..IIIIII...IIIIIIIII. .IIIIIII ..1 _.1 ..1 _.1 ..1 ..1 ..1 ..1 ..1 ..1 ..1 .1 .1 I1 91 .1 I I I . I I I I I I 6 I I I I I I I . I I I I I I I . . I I I I 5.9 5.6 7.0 2.1 9.1 9.2 2.2 9.5 0.3 0.0 5.0 2.2 1.0 2.2 3.1 7.0 .1.0 .IIS .3.0 .0.“ 6.9 .3.6 .1.6 0.0. .6.6 .9.6 .3.0 .2.0 .109 0.1. 0.7. 0.3 3.6 2.2 2.0 2.3 1.0 1..» 1.0 0.6 0.3 0.0 0.0 0.5 0.3 0.2 0.1 0.0 0.6 0.4 0.3 0.3 0.2 0.2 0.2 0.1 0.1 0.1 0.0 0.0 0.0 0.0 0.0 0.0 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII .01 ..1 _.1 _.1 .21 ..1 ..1 _I1 .01 ..1 ..1 ..1 _.1 .I1 21 .1 I I I I I I I I I I I I I I I I I I I I I I I . I I I I I I I I 705 7.6 9.0 3.2 02 0I1 37 9I6 0I9 °I0 6.3 1.03 °I6 III 56 1o]. 2.... 0.9 6.6 1.9 0.6 5.... 2.0 9.3 7.9 6.9 6.3 3.0 2.0 1.6 0.1 0.3 3.7 2.2 2.0 2.3 1.9 1.5 1.1 0.7 0.3 0.0 0.0 0.6 0.6 0.2 0.1 0.0 0.6 0.6 0.3 0.3 0.2 0.2 0.2 0.1 0.1 0.1 0.0 0.0 0.0 0.0 0.0 0.0 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII . 1 . 1 . 1 _ 1 . 1 _ 1 . 1 _ 1 - 1 . 1 . 1 o 1 . 1 - 1 . 1 — 1 9 1 3 5 7 9 1 3 5 7 9 1 3 5 7 9 1 2 2 2 2 2 3 3 3 3 3 6 6 I 6 I I .1250 1.0000 .1050 1.0000 .0868 IIIIIIIIIIIIIIIIIIIIIIIII I 1.0000 00703 .0000 "i .0555 .0625 .0312 .0217 .0130 .0070 I0036 .0000 .0000 69 FIGURE 6.8 Distributions of the Stress Function and 0‘; for Problem 4 COLUMN 3 3 7 9 11 13 13 11 19 21 23 23 27 Rg'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII IIIIIIII. 1.3023 1.3203 1.3166 1.6963 1.7391 2.3066 3.6332 2 2 2 2 2 2 2 I I I I I . I -00‘5I ‘OOIII -0060. '603I’ -0026’ -001.6 60°10 0°20. IIIIIIIIIIIIIIIIOIIIIIIIIIIIIIIIIIIIIII 0 6060000666... 6 1. 166 1.3317 1.3963 1.3316 1.7936 2.2361 2.6331 .6169 2 2 2 2 2 2 2 2 -20639 —20627 -20392 -20331 -20263 -20126 20029 20213 .0616 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII IIIIIIIIIIIIIIIIIIIIIIIII 1.3361 1.3776 1.6606 1.3676 1.7729 5.9929 1.6333 1.0333 .3076 2 2 2 2 2 2 2 2 2 I I I I I -2o613 -.06o3 -.o366 -20306 -20213 -.oo96 20031 .0226 .0620 .0623 IIIIIIIIIIIIIIIII'IIIIIIIIIIIIIIIIII II IIIIIIIIIIIIIIIIIIIIII II:I ”I I 1. 006 1.6166 1.6699 1.3626 1.6135 5.7266 1.3619 1.0919 65552 .3371 2 2 2 2 2 2 2 2 2 2 I -20366 -2o373 -20336 -20213 -20163 -20066 20076 20263 20626 .0621 .0633 “1:11.11"1:1111°°1:1111'°1:1111'°1:1111'13111"1:1111"1:1111°°°:1111°'°:1111"221.. I : I I I I I I I I I I I I I I I I I I I I -.0332 -.o339 -.0302 -.0239 -.0130 -.0036 .0101 .0262 .0660 :06)! .0636 .1061 ’ IIII HIIIIIII HIIIIIIIIIII'IIIIIIIIIIIIIIII III II “I IIIII IIIII ”I HIIIII ”I I 1.6396 1.6635 1.6393 1.6176 1.6131 1.6136 522556 1.0635 :7555 25552 :3055 .1103 I I I I I I I I I I I I I I I I I I I I I I I I -.0316 -.0306 -.0266 -.0206 -.o116 -.0001 .0126 .0261 .0633 .0660 .0637 .1062 .1230 1 IIIII IIII HIIIIIIIIIIIIIOIIIIIIIIII “III IIIIIIIIIIIIIIIIIIII:IIIIIIIIII IIIIIIIIII'IIIIIIIII 1. 6292 1.6255 1.6307 1.6239 1.3919 1.3263 1.1937 1.0216 6117 :6552 .6062 .2060 .0101 2 . 2 2 2 2 2 2 2 2 2 2 2 I I I I I I I I I I I n-20200 -20266 -.0231 -.0111 -.0066 .0021 .0130 .0300 .0661 .0666 .0661 .1063 .1230 I I I I I Q I I II I I I 355 2252655522525565 2523555 25.5552225:2 25 2522525552 222522 002522225555222265252222559322225155 22.1260 I 2 2 2 2 2 2 . 2 2 2 2 . I I I I I I I I I I -20263 -20233 -.0191 -.0136 -.0036 .0066 20173 .0316 .0660 .0636 2°222 .1066 :1250 , IIIIIIIIIII IIIIIIIIII IIIIII. .IIIIII I IIIII HIIIIIII IIIIIIIIIIIIIII HIIIIIII HIIIIIIIIII 1. 3692 .3655 1.3360 5.3296 5.2630 5.2139 52 5561 .9692 .6679 .7026 .3319 .6110 22111 I 2 . 2 2 2 2 . . 2 . . 2 2 1-20211-20199-20165 -20101 -20026 20012 20196 20336 20692 20666 20669 21063 21230 H. .556IIIi:555,II!:ii£IIIi: .iai...izizii..!:i; ‘IIIlasaaIIIzasaiIIIIa‘ésIIII,5‘IIIIIII‘iIIII; IIIz‘o39 -20119 -2o166 -20133 -20019 -20002 20093 20213 20369 20302 20611 20632 21066 21230 a’x:i“i"i:iafla..1.”I...i: 5;;,IIi:’ai’II‘:‘::,IIiIIIi‘IIII550IIIIIIiigIIII,7iIIII:;;;;III:;;‘; II:,l96 -20130 -20160II-201o7 -20033 20020 20113 20230 20362 20312 20677 20333 21061 21230 5 666666 ”I 666666666666666666666666666666 6666666666666 66666666666 HIIIIIOIIOIIIO 600.0 102I71 162555 .162320 1.20., 1.170. 16115.. .10016, 6 “I7 .IIII 6.033 67230 0663‘ :6179 2 . 2 2 2 2 2 2 2 2 2 2 2 1-2012t -2o116 ~20062 -20030 20061 20136 20263 20316 20320 20662 20636 21061 21230 122555225225552252 55522252 .5555225.' 5555225 2555 22525525222255552'22555222225525222255552222525522227006 I I I I 6 6 6 6 I 6 : 6 2 : 6 6 6 o I 6 6 o 6 6 6 o 6 6 6 6 6 6 6 “60101 "60091 ”60061 ”60°50 00060 60550 6025. 603II 60527 00“. 00.60 010‘. 61250 29 6666666 6666666666666666666666666666666666666 666 6666666666666 6666666666.66666666666666666666 53175I 1.1720 165655 161~36 161152 ‘007‘l i6°221 69651 69091 6.592 6.57. 67562 67702 6 6 6 6 6 6 6 6 6 6 6 6 I 6 6 6 o o 6 6 I o 6 I 6 -20061-2oo71-20061 20001 20016 20163 20269 20393 20336 20690 20662 21066 21230 1 IIIIIIIiIII I.IIIIII.IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII.II.I 1: 1662 :1555" 5255552252555522525555" 525366 1.0160 .9660 .9206 .6629 .6367 .6330 226292 I I I I I I I I I I I I I I I I -20066 -20036 40023 20022 20069 20113 20279 20601 20339 20693 20663 21069 21230 ” OOOOIOiOUO 0'0 000060 00. 6600000 00 000000’ 000666066066 600060666600. 66666600660600... 15.1161611 52 55 55 956 5.0766 520655 1.005 .9613 ”5 .9031 .6662 .6163 .6793 I I I I I I I z I I I I I 33 40066 -20039 -20011 20033 20101 20163 20261 20601 20363 20696 20663 21069 21230 125555225255572 5255252252555522525555225255552252552222 22555522225552222252552222555 2222555522229213 I I I I I I I I I I I : O o O o 0 6 o 0 I I o I I -20036 -.0026 020000 20061 20111 20196 20296 20612 20361 20696 20666 21069 21230 ’11 III IIIIIIII II II .IIII'IiIIII III IIIIIIIIIIIIIIIIIII. IIIIIII IIIIIIIIIIIIIIIIIIIIIIIIIIII 6561.0665 525 625 5:03731.0655 5.0232 .9979 .9705 .969 .9376 .9331 .9399 .9332 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -.0023 -.oo16 .0010 .0036 .0119 .0201 .0300 .0617 0330 .0100 .0666 .1030 .1230 9 IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII II IIIIIIIIIIII: IIIII IIIIIIIIIIII IIIIIIIIIIIII 1.0669 1.0632 1.0669 1.0639 1.0330 1.0113 25555 .9132 9311 .9319 :5556 .9623 .9606 I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I -6001? '60008 6001. 60063 60526 60206 0030‘ 60610 60552 60701 60067 61050 61250 IIIIIIIIIII .IIII .II. IIIIIIIIIIIIIIIIIII'IIIIIIIIIIIIIIIIIIIIIIIIIIIIII. III II I 1. 0326 1.0555 2.55552 252 5526 1252552 25.0112 .9929 .9766 .9663 .96 3 .9102 .9793 .9976 I . . . 2 2 2 2 2 2 2 2 2 2 -20010 -20001 20026 20069 20131 20210 20306 20622 20336 20702 20661 21030 21230 , IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII.II IIIIIII IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII. II 1.0191 1.0212 1.0227 1.0223 1.01752 2525555 .9921 .9616 .9767 .9137 .9623 .99552 25.0061 I I I I I I I I I . I I I I I I I I I I I I I I I I I I I I I I I I I . I I I -.0003 .0002 .0029 .0013 .0136 .0213 .0310 .0626 .0333 .0103 .0666 .1030 .1230 IIIIIIII .0 IIIII II.III.IIII.II ”II III II II... IIIIII IIIIIIIIIIIIIIIIIII III I. 1. 0102 1. 0555 2555522 1.0160 1.0105' 525555 29960 .9655 .9529 .9636 .9913 .95552 25.0060 2 2 2 2 2 2 2 2 2 2 2 . 2 ~20002 20006 20032 i20076 20131 20213 20311 20626 20333 20703 20666 21030. 21230 7 III II. III. ..II.IIII.IIIIIIIIIIIII HIIIIIII II. III. IIIIIIIIIIII ”III II I. 120065 125552 1255552 20011 1.0036 1.0010 25556 .9925 25501.9931 .5570 .55552 25.0037 2 . 2 . 2 2 2 2 2 2 2 . 2 .9 .0000 20006 20036 20071 20136 20216 20312 20623 20333 20703 20666 21030 21230 1255552 221255522 25255552252: 552522 2555522 2555 222255552222555522225555 '2’2555522225555225255552 2520020 I . . 2 . . 2 2 2 . 2 . . 2 I I I I I I I I I I I I I .0000 .0006 .0036 .0076 .0136 .0211 .0312 .0623 .0333 .0703 .0666 .1030 .1230 51 III. ...:...6...:...6..i.656 IIiIII.IIIiIIIIIIIIIIIOIII 1.000 000 000 .0000 .0000 1.0000 5255552252555522525555225255552252555522520000 u 0‘: 70 FIGURE 6.9 Distributions of the Stress Function and er, for Problem 5 COLUMN Row 3 5 7 9 11 13 15 17 19 21 23 25 27 -.0020 -.o0o3 -.0552 -.0007 -.0307 -.0191 0.0000 1. '5555'°i35555"535555"535555°' 35555"535555"539503 . . 3 3 3 . 3 s-30010 - 30002 -30551 -30000 -3o300 -.0191 030000 00 00000.0... 0.... 00.000.00.000... 0 .0 1.9257 1. 9027 1.9795 2.0227 2.05733 5. 0335 539052 3 3 3 3 3 . 3 -30013 .0590 -30507 ~3o003 -3o305 -30190 030000 1000000 .00....C0000000I0.00.0000000000..0..00.0 135510 1. 9099 1.9520 2.0115 2.0730 2.0902 2.0100 . 3 3 3 3 3 3 0 0 . 0 0 0 0 9-:°0°3 ”00501 '00,}. -00057 -9030} '00109 000000 0000.. .00000000000.000.00.00000000 00.00.00. 0 1:0307 1. 0500 1.9007 1.9029 2.0905' 2.1709 5.1219 . 3 3 3 3 3 3 -30507 —30572 -30525 -30007 -30335 -30107 030000 l100000000....000D0.0000.00.00.....000000000000000 157600 5.707; 10.307 109225 290709 202973 203626 3 3 3 3 3 3 3 -30505 -30551 -30500 ~30031 -303 0 -3o193 030000 oooooooo:ooo 900000900 000000000... 00 900059000 1. .0993 7155 1.7530 5.0272 1.9005 5.3020 3.0090 . . 3 3 3 3 3 -30537 -30523 .-30000 ~30007 -3o305 -30170 030000 .0200 .0010. .0025 .0033 .1001 .1250 5 0000000000000 ..0000000000.000.000.00000000000000000 00 .00 000 00. 00000000000000... 1.0020 1.0553 1: 0001 1.7302 1.0300 2.0027 0.2013 03550000 0.0000 5.5555 .0000 0.0000 0.0000 3 3 3 3 3 3 3 3 3 3 3 3 3 -30502 -30000 -30000 -30370 -30270 -30105 30010 30212 30010 30020 30032 31001 31250 H.000... H00 ”9 0000 9009909009909000000900000.900990 concoooooooooooooonzo 0.0.0000 H0000. ”OI. 130019 1.0155 1. 0020 1.0900 1.7009 1.9310 1.0000 .5500 .1530 .m5 .0250 -.0037 -30503 . . 3 3 3 3 3 3 3 . 3 3 . . . . 0 0 0 0 0 . 0 0 0 -.0002 -30000 -.o007 -.0330 -.0200 -.0112 .0003 .0225 .0021 .0025 .0032 .1001 .1250 19 o ..i.... 09.0.. ...... o. ... 00. 0.0.0.000...000000000000:0. 0000: 0000000000 H00... ”0 I 1:5555 005 '1.0125 1.0525 5.7050 5.7000 1.0509 .0320 055 000 .0900 .0135 -.1130 0 . 0 0 0 0 0 0 0 : 0 0 0 0 0 0 : z : 0 0 I I 0 0 0 0 0 0 0 0 0 -.0010 -.0005 -.o300 -.0290 -.0201 -.0o70 .0071 .0203 .0031 30029 .0033 .1001 .1250 00.00.00.000 0000009 000000000000... 00 09.99099009905090009090900999090.90900050 000...... .0000 1 5005 1.5505 1.5755 1.5950 1.5975 3.5200 1.2907 .9210 .5000 .3510 .5000 .0053 -.1117 o o o o o O 0 0 0 0 0 0 0 0 0 0 : z : 0 0 0 0 0 0 0 0 0 0 0 0 -.0372 -.0359 -.o320 -.0250 -.0101 -.0o03 .0099 .0203 .0003 30035 .0035 .1001 .1250 135555 'i35iii"i35555"i35555 '533535"! I33555°°i': 51 35"‘355 55"‘3355 "°33555"'3 535‘ "35333" 330001 0 0 . . 3 3 3 . . 3 3 . 3 3 . 0 0 0 0 0 0 0 0 0 -0032. -0031, ‘00275 ‘00211 ’00122 '90009 00127 00203 :0... 0000’ 00.39 :5002 :1250 099000.000.0.900.000.00090900.00.00.05090000000955.0000...0000000509 00.0000 0.0000 00.00 I .00. 1.0039 1.0713 1.0720 1.0577 1.0117 1.3100 1.1027 .9000 .7505 055 .0030 .5055 .0001 0 . 0 0 0 . 0 0 0 0 0 0 0 O O 0 0 . O 0 0 : 0 0 . -30201-30209 -30232 -30171 -30005 30023 30153 30303 30070 30051 30003 31003 31250 00000000. to. 00.o.0000000000000.090000000-00.0000000000000000000000000.000.900000000000000000000 1. 0120 133300 1.5100 1.3910 1.3015 1.2539 1.1200 .9090 .0009 .0370 .0092 .3023 .1952 3 . 3 3 3 3 3 3 3 3 3 3 3 0 0 0 0 0 0 0 0 0 0 0 0 0 29'90239 '00227 '90192 “00133 ”00°51 00°53 9017. 90322 000.3 90059 90.07 030.0 03250 135553"i35353°'535535°'I35553"535555'°1355533'535555"'35555"'35555"'33535"335553"'33355°"33200 0 0 0 0 0 0 0 0 0 z : 0 . 0 0 0 0 z : ,1-30199 -30100 -3o150 -30097 -30019 30000 30201 30300 30090 30007 30050 31005 31250 135555"i35555°'i35555' 5355553'535555"535335"i35555"'35555"'35553"'3'355'°'33355"°35555"°3051o 0 0 I 0 0 . : z : 0 : : 0 0 0 . 0 . 0 0 0 0 0 0 0 0 0 0 0 0 -.o103 -.0152 -.o120 -.oo05 .0009 .0100 .0222 .0350 .0507 .0070 .0050 .1000 .1250 33 0000.0 .0 .000 .00. 000000.000 0.0 .00.... 0.000000000000000... 0000000 000.00.000.00..0..00 135555 132503 5.2095 1.2270 531552 5.1330 1.0591 .9700 .0055 .7035 .7053 .0255 .5002 3 3 3 3 . . 3 3 3 3 3 3 3 0 0 0 0 0 0 0 0 0 0 . 0 -30131 —.0120 -.0009 -.0030 .0035 .0120 .0200 .0370 .0517 .0000 .0057 .1007 .1250 351.000.0600... 00 Hooooooooooooooooooooooo H00... 000.000....OOO0000.......OOQOI0000.00.00.00000000 :2022 1.2077 1.2021 1.1039 1.1513 131009 1.0501 .970 .9000 .0209 .7037 .7002 .0000 . 3 3 3 3 3 3 3 3 3 3 3 3 O c o o o o o O o -.0102 -.oo92 -.0001 -.oo11 30050 .0100 .0250 .0303 30520 .0005 30059 31000 .1250 H0... 0.. Moo-0900000909990090.9000. H.000... ono.0000.o0000000000000000.000000.900000000000000 131555 1.1030 1.1590 1.1000 1.1102 1:00001.0317 .9777 .9101 .0591 .0150 .7737 .7525 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -.oo70 -30000 -30037 30011 30079 30105 30271 30390 30530 30090 30002 31000 31250 000000000 000 005.900a90.000cocoo.u00090-00099000099099.00099.000000.000000.00900000000000.0000... 1.1190 131205 1.1215 1.1101 1.0093 1.0001 1.0215 .9001 .9302 .0095 .0000 .0330 .0279 0 0 0 0 0 0 0 0 0 0 0 z . 0 0 0 0 0 0 0 —30050 -30005 -3oo10 30030 30090 30100 30203 30003 30500 30090 30000 31009 31250 ‘1 0.00 00000.0.0000....0000.0.00.00.00.000000000000.0.000000....00000.0000000000000...00000000000000 1. 0052 1. 0099 1.0001 1.0790 1.0002 1.0027 1.0130 .9032 .9009 .9100 .0992 .0007 .0910 . 3 3 3 3 3 3 3 3 3 3 3 3 0 0 0 0 0 0 0 0 0 0 0 0 h) ':0030 '90020 .°:°°°° 00000 00110 90393 00293 90011 00500 :0097 00.05 .1009 .1250 00.000 00.000 0.000000000000000000000 000.000.000.000. 00-00 .0 00.0 .00 0000000000 0 1. :5575 135555 13 0593 31.0539 1.0033 1.0200 31.0070 .9007 '.9029' 35003 ' .5510 .9205 '339013 3 . 3 3 3 3 3 3 3 3 3 3 3 0 o o 9 o 0 o o O O o 0 o -.0021-:oo12 .0010 .0059 .0122 .0203 .0301 .0017 .0550 .0700 .0000 .1050 .1250 , 0000 0. 0.0.00 0.0000 0.0.000000.000000.00.000.00000000000000000.00000000. 0.0000. 000 0 135330 135555 1.0355 1. 5320 1.0259 1.0107 1.0037 .9900 .9750 .9011 .9557 :9555 .9770 3 3 3 3 ‘ 3 3 3 3 3 3 3 3 ‘1-30009 -30001 30025 30009 30131 30210 30307 30021 30553 30701 30007 31050 31250 13555 "i35i55"535555'°535i55"i': Hi53'°535553'° 35555°"35553"°35555"'3555 55"‘3555 "'35 53" 39997 0 z 0 0 0 : 0 0 : 0 0 0 0 O o o o o O o o I 0 0 O O -.0002 .0000 .0032 .0075 .0130 .0215 .0311 .0020 .0555 .0702 .0007 .1050 .1250 00.000.000.000... 00. 00000000000.00.0.00.000000000000000.000000000000000 000000... 0000000 00 1. 000 3 1.0003 135537 1.0007 1.0030 1.0020 9 .9909 .9901 .9550 .9939 .9932 1.0052 0 0 0 . . 3 3 3 . 3 3 3 3 3 3 3 0 0 0 0 0 0 0 0 0 0 0 0 .0000 .0000 .0030 .0070 .0130 .0217 30312 .0025 .0555 .0703 .0000 .1050 .1250 1 000.000.000.00... .00 ”O 00 0.0000000. 09.900.00.00... 00 o. 0000. 0000000000 1.0000 1.0000 130055 1.0555 1.0000 5.0000 1.0000 13 5555 535555 535555 1.5555 1.0000 1.0000 "9|: 71 the solution domain. Along the xsaxis, .3? :3 0, while along the yaaxis g? = 0. Along the notch edge if or To: and 3? 53%va %§= (%:,4’)4=2a73 (b = ZQBy +C = 2aE(y—a) The value of c is determined by the value of ¢d in the squetion used along 3'5. The lost expression for g! is valid for :11 3 problems but for Problem 5_ it simp1ifies to g s 0 along 3?. . ' The oolutions given in Figures 6.7, 6.8 end 6.9 were obtained with boundsry oonoitions for which Po 9.! 1 end It “3;. Though a 12:21; mesh is used for listing the solution in all three figures, the results are token fro: the solution obtained for e. 2MB Iesh. Stress concentretions occur at the base of the notch st y a e. For the V snd seni-oirculsr notches the concentration occurs elong the y—sxis. For the rectangular notch the concentretion occurs st 1 u s, y g s. The oonsrisons of stress ooncentrstions for different choices of nesh sise ere given in Figures 6.10, 6.11, 6.12. B. I. Peterson (1953) gives sn snot stress oonoentretion of 3.08 st the bees of thee-i- oirsulsr netshes in e tensile specinsn infinitely long. The stress eon- osntrstion vould be slightly higher in e specimen of finite length. Bouthsell (1956) determined s stress oonoentretion footer of 3.0 using s 'relsxsticn' solution for e neoh of 116 points. .Solutions' obtsined vith the ISOPEP progre- produce stress concentretions of 2.92 and 3.211 for neshes of- 65 end 1038 points respectively. It the bees of the V-notoh or the corner of the reotsngulsr notch even slight strsins induce stresses of high-ugnitude end the- stresses detersined from s solution of the bihsrsonic difference squeticns would herdly represent sn sctusl ptvsiosl stste of stress. However, stress B i? L in A ..l 3, « I: * +—1 .__l +4 3651’ ‘ Emu-eta :4 O U3 U3 :2 2. ‘ E—u‘ LC 4 A ififflrw 1 5” --«--———- ISOPEP - 1038 Mesh PtSe J m~-~-- SOUTW‘IBLL - 116 Mesh PtS. --------- ISOPEP - 65 Mesh P‘tS. 1. " w _ g s ; A :1 - ; . . . 1*; y 0 1/6a 1/3a 1/2a 2/3a 5/6a a ORDINATE AT 2: = 0 Figure 6.10 Stress Concentrationin a semi-circular notch distributions obtained under the assumption of an ideally elastic body are useful in the analysis of the plasticity problen. Exact solutions of these problems would provide V infinite stress concentrations for these two problems. The numerical solution provides a set of discrete values of the stress function. The value at each nesh point is an average value associated with an area surrounding the point. The area in general is proportional to hz. As the uniform nesh space h is decreased, the stress function solution at a point of high stress concentration would provide a better approximation of the high stress. Figures 6.11 and 6.12 illus- trate this phenomenon for the V and rectangular notches. STRESS ox H g ‘3 9." L \ f a_ . Fr 4—4 H (—-u ....) 8.1 3,. 1' L4; 4—4 r--> ‘_1 -—> *“ ,///l :: § ”’ 7.’ * pl 3,. ‘ oil—.4 ISOPEP 6 {_ — 107A Mesh Pts. 0 ... _ __ -i 267 Mesh Pts. ________ 66 Mesh Pts. r'8.65 1/6a 1/3a 1/2. 2/3a 5/6a ’ ORDINATE AT x.= 0 Figure 6.11 Stress concentration in the V-notch a 78 TABLE 5 STRESS FUNCTHON A16801N15 ADJACENT TO THE SEME-CflRCULAR NOTCH ne-e—e—e-sjg NO‘WJ‘WN. ———u———— eQQLUMN WP]; 13 111 15 16 17 18 3 -302569 -301823 ~300971 03 00000 4 ~302560 ~301816 ~300966 83 88080 5 -302533 -301795 ~300953 . 300002 6 -302h88 ~301757 -380928 300009 a -3OZ#23 -301703 ~388890 380025 ' -302380 ~301631 ~300833 308056 301842 9 -302239 ~3015h1 -308759 300187 381857 ' 10 -302123 -301h35 -388678 300173 301093 302085 ’ 11 -301992 ~3©1316 -3@0566 300256 3011h9 302107 303125 12 ~301851 -301186 -308&50 308353 301223 302155 3031h0 13 -301700 -3010#6 -3®032h 300863 301313 302221 303180 18 -301Sh1 -300898 -300189 300582 301115 302302 303239 15 -301378 -3007&5 ~300889 300788 301523 302392 303309 16 -301212 -300590 388098 300837 301637 302888 303386 17 -.01o&§J-.0013§ 00239 .ooagafl .01153 .02588 .03168 19 go 21 V_22 23 21 25"" 30h166 30h180 305209 ' ' . 30h216 30522h 386253 307291 ' ‘ 308266 305255 306269 387297 308338 309375 310% 308328 305296 306298, 307310 308339 309376 310 .0h389 .OSth .06326 .07338 .083h9 .09379 .101 6 TABLE 6 X-STRESS COMPONENT AT POHNTS ADJACENT TO THE SEMB-CURCULAR NOTCH M h _ 6 CQLJJMNB 2 1 1 1 ‘ 1 17 1 19 20 21 2 3 23th 23785 33238 8 23831 23701 23952 5 23369 23629 23576 6 23275 23506 23111 3 23158 23345 23061 23031 23130 23193 ' 9 13917 13955 13911 13386 ' 18 13805 13801 137h5 13668 13817 11 13708 13678 13611 13586 13367 13186 ' 12 13628 13586 13516 13188 13213 3955 3695 13 13563 13518 13118 13338 13178 39th 3653 3839 ' 18 13508 13h55 13380 13278 13131 3980 3711 3860 3238 ' 15 13860 13h©5 13331 13229 13999 3936 3786 3543 3351 3612 16 13h17 13362 13288 13193 13078 3931 3769 3598 3131 3286 17 1.379 1.323 1.251 1.162 1.853 .927 .786 .638 .891 .358 TABLE 7 Y-STRESS COMPONENT AT POBNTS ADJACENT TO THE SEMI-CERCULAR NOTCH CQLHMN 13 10 15 16 17 18 19 20 21 22 ROW ' ‘~ . . .H H .3 3 3329 3209 0. 000 0 3328 3199 3000 5 3361.263 3110 6 3395 .326 3220 7 3025 .012 .305 8 .000 3002 .000 9 3350 .303 .357 .087 10 .312 3326 .372 .051 ‘ 11 32713295 .308 .038 .580 ' ' 12 3222 3200 .286 3350 3003 3590 3323 13 3169 3185 3210 .257 .318 3005 3500 10 3116 3126 3105 .175 .216 3268 3321 .352 ‘ 15 3068 3070 3088 .108 .136 .170 3205 3226 3219 .180 16 .025 .030 .000 .055 .076 .102 .127 .107 .152 .135 TABLE 8 XY-SHEAR STRESS AT POINTS ADJACENT TO THE SEMI-CHRCULAR NOTCH 0 13 10 1S 1% 17 1§_‘ 19 29 21 22 3 03000 03000 0 -3100 -3153 6 -3266 -3391 g -;300 -0‘95 -3307 -306 ' 9 -3312 - .051 -3605 10 -3312 - .030 -3578 11 -3302 -.399 -3S11 - .633 '12 -3287 -3365 -3008 -3 35 - .628 - .578 ‘3 -;275 -0337 -01‘01‘ ..o 70 -;530 ’:562 - 10 -3266 -3319 -3370 - .026 -3070 -3091 -3069 15 -3261 -3307 -3353 -3396 -3031 -3007 -3035 -.385 -.303 ' 16 -.258 .298 .338 -.375 -.003 -.017 - .009 -.376 -.317 -.202 76 TABLE 9 STRESS FUNCTHON AT POENTS ADJACENT TO THE THE V-NOTCH QQLUMN _ 12 HB 111 LJS lg H7 IQ 73 C) i ‘d‘dddda -;02fl87 -.02fl62 -;02095 -208997 -;on876 -;0fl7h@ -00“ 595 “50310152 -;03286 -00‘ ‘ 28 -;00969 -;00811 -;00656 -000353 ‘.011570) -;0fl5h0 -;oflh66 -;©U362 -301] 2h© -;Oflfl06 -300965 -;008fl8 -300669 -;005fl8 -000367 -;00298 -000070 ;00072 .002fl3 -900853 -;008U6 ‘200737 -0 ©0631} -0005fl 8 -;00392 --000260 -;00023 500035 .00355 .00295 .00h3h ;0057O 9007©h .00834 0.00000 g00035 .003®h .00fl9fl .0029fl 'ZOOAOfl .005fl7 .00638 .00763 .00890 OOUOUW .0flflh3 .0fl268 .0fl390 \ #9 ’0fl@60 '.0fl266 50357fl 003 O .022 -__?IT_— .OUOhfl .Oflflflfl .0flfl8fl .0fl36fl .OflQGh .Ofl682 0fl79h :0fl906 .02®fl7 .0292? 201005 ;02®83 502097 .02fi36 202193 .0226h .023h5 .02h33 .02526 :02621 .02738 ;028fi5 502911 .03125 ;03336 ;03167 .OBZUA .03273 .033hl 303hn5 '.03h9h .03575 '03557 .03739 :03832 agddfi .0h966 .Ohfl75 .0h201 :0h239 .0h288 ;0h3h5 30hh07 .Ohh72 .0h539 .0h607 ease Nmmrwnaammugl‘Nmmkwnaemmw¢mrw .0h67h .05208 ;05205 .05236 .05268 .05398 .OSSSh ;05%®# .05#57 .QSSflU .05565 .06250 .06255 .06272 .06298 .06330 .©6366 .06h®6 .06hh7 .06h88 .0729fl g07296 .07309 .07329 .07353 0073811 .07h00 .07hh0 303333 308336 .083h6 008360 .08378 008397 .08hfl7 ;fl0hfl6 ;30h37 ;30421 gfl0h26 .IOhBI TABLE 89 XMSFRESS COMPONENT A8 POTNTS ADJACENT TO THE VuNOTCH _ _ _« COLUMN ‘__ ._ 83 8h 85 86 87 88 89 20 28 22 Row 3 2.306 3.858 8.653 A 2.350 2.952 3.557 8.628 5 2.288 2.608 2.633 8.539 .874 6 2.809 2.269 2.839 8.55C .958 .642 7 8.992 2.088 8.853 8.573 8.05 .708 .567 8 8.898 8.888 8.670 8.900 8.086 .798 .559 A82 9 8.725 8.672 8.582 8.335 8.098 852 638 .956 .337 80 8.622 8.562 8.886 8.280 8.085 .880 .699 .525 .375 .273 88 8.539 8. 975 8.372 8.235 8.072 .902 .736 .588 .938 .308 12 8.070 8. #07 8.385 8.896 8.059 .983 .767 .625 .492 .367 13 8.483 8. 352 8.268 8.868 8.085 .989 .789 .668 .538 .02 Th 8.366 8.306 8.230 8.837 8.032 .928 .805 .690 .577 .h67 85 8.326 8. 269 8.898 8.885 8.028 .928 .887 .783 .680 .509 16 8.292 8.238 8.872 8.096 8.088 .928 .826 .732 .637 595 87 8.263 8.288 8.850 8.079 8.003 .920 .839 748! .668 .576 TABLE 88 Y-STRESS COMPONENT AT POTNTs ADJACENT TO THE V-NOTCH COLUMN 83. 8h 85 86 87 88 89 20 28 22 Row 3 8.376 8.705 8.622 A 8.033 .983 .787 .874 5.665 .529 .882 .725 .682 6 .h28 .388 .289 .868 .573 507 7 .273 .282 .287 .325 .h26 .969 .482 8 .879 .887 .863 .233 .383 .365 .380 .337 9 .885 .098 .887 .867 .229 .278 .305 .388 .273 _ 80 .066 .059 .076 .885 .86h .206 .236 .250 .258 .286 88 028 .026 .088 .072 .880 .886 .875 .892 .899 .896 82 -.002 - .003 080 .035 .066 .097 .823 .848 .809 .850 83 -.029 - 027 -.086 .003 .028 .955 .078 .099 .808 .807 In -.058 -.049 -.O39 ~.023 -.002 .089 .080 .056 .067 .078 85 -.O69 -.066 -.058 -.Ohh -.027 - .008 .008 .028 .035 .088 I6 -.O8h -.080 -.073 -.O68 -.047 - .038 -.085 -.008 .0098 .08; 78 TABLE 82 XY- SHEAR STRESS COMPONENT AT POINTS ADJACENT TO THE V-NOTCH The finest mesh used for the notched tensile specimen problems "a We gram.includes a subroutine STRESS for the calculation of any or all of the stress components a}, a; and The corresponding mesh interval h is 1/L8. T;y. C08.UMN ROW '3 8“ 852 _86 87 flak 89 2O 2 22 3 0.000 0.000 0.000 h 0‘“, 0.000 -0672 5 - 0°96 -ofl2h _06S© ‘09®2 6 O. 000 -;228 -0563 ‘6773 “.698 - 7 "' .080 -0268 ‘0h98 -0650 -0653 -055h 8 ‘:‘]Jh -e286 "ekhe -056“ -0586 ’.5'00 -3160 9 -o“69 -0289 -0'O09 “.1895 -9525 -950“ ‘5‘th -0369 ‘0 “.891 -0286 -0378 'ekhs -0u75 ‘0“67 -9h30 “037“ -0301 ' In -020“ -028. -0352 -0“06 "oh3h -e"33 -e“07 -e36h “.380 -02“, '2 “.28” ”027‘. -9330 -0375 -0399 -el§@2 -6385 -6352 -0307 -025“ ‘3 -0236 -0267 -e3‘flk -0350 -9378 -0376 -036“ -03,8 -030“ ".255 I“ -.289 “.260 “e299 -9329 -03“? -0353 “63NS -0325 -029h -025“ ‘5 '.22® -.255 '.287 -.301 -.328 -.333 “6328 -.381 -.286 -.251 ‘6 -0120 -0250 -e277 -e298 -0301} “630$ ‘03“: "'e299 -0277 -02~7 The ISPOEP pro- Values of all three stress coupon- ents and the stress function in a region surrounding the point of stress concentration are provided in Tables 5 to 12 for the semi-circular and V- notches. See Appendix B for the finite-difference equations used. These equations are valid at interior points only. Calculation of stresses at points on an irregular boundary are not included in the subroutine ‘J ._ a 7: 6.28 6. , gr .— +—-> «J .9 4—4 —-> w t4 -t X 5° ;_ f —+ ‘— f H 4059 ' (— 4 *—-> ‘L _i , 4° ‘ Li. STRESS ox 996 Mesh Ptso - ‘_ 2A6 Mesh Pts. 1; _'3"’_- ...._.__m 60 Mesh Pts. O 1/6a 1/3a 1/2a 2/3a 5/6a a ORDINATE AT x.= 8 Figure 6.12 Stress concentration in a rectangular notch Though a solution obtained with a coarse mesh does not provide a good.meesure of the stress at a point of stress concentration, the results shown in Figures 6°10, 6011 and 6012 indicate a decrease in mesh interval 11 produces only slight changes in the stresses computed at all other pointso Use of an average'value of the stress function over an area surround? in; a mesh point depends more on the magnitude of the area than on the geometry'of the region. For example at the'base of the Yanotoh the area surrounding the mesh point is taken as (bout half the area for an interior mesh pointo “Does this not imply that the solution obtained at this point is a better approximation for s.tensile specimen'with a r30 rounded rather than a sharp corner at the base of the notch? Some evim dence for this conjecture is provided by comparison of the stress con» centration factors fovmd from the numerical solations of Problems hOl‘9 1102 and [£8 with the stress concentration factors provided by Peterson (1953) for notched flat bars in tensiono Peterson provides values of the stress concentration factor Kt for tensile specimens with deep notches with parallel sides and a semicircular base as sheen in Figure 6°13o The stress concentration is defined 653(3): norm where 0:0“ is the average stress across the minimum cross section of Kt: the specimen" The curve plotted in Figure 6013 shows the relationship of xi to r/h where r isthe notch radius and D is the'width of the bars This curve is based on the Neuberatheory solutions tabulated by Peterson (pp 2647) for d/D s 05 where d is the minimum distance across the bar at the base of the notch. The stress concentration factors for three different numerical solutions of Problem )4 are listed in Table 133 Table 13 Stress concentrations factors for the V-notch tensile specimen # w Problem :32 h 0;: Ofnbrn Kt r/D hol 6x12 1/12 h0362 2 2.181 1/12 1102 121211 I/Zh 60083 2 3 901:2 l/Zh hoe 21am 1A3 80653 2 _ b.326 1/h8 Assuming the numerical solution approximates the solution for a bar with circular notch such that r g 119 the stress concentration factors Kt for problems 401 934029 108 are plotted as three points labeled in Figure 60139 and the Neuber theory solutionso There is good agreement between the numerical solutions Peterson reports that for notches with inclined sides having an included angle (X and a circular are at 31. the base, the notch angle has very little effect if 0° < (If 90°. 0 ”mm """“““'<' -.; ,r 6" “it a” I “'T (1 D P 4. i :F 5o .. / ’X~_-i X“ // U L‘\ r S V '\ r Mino radius is) c: L” .- .3 R; \ H o.) C'. G.) is.) 8 3° m w 3 u U) Y 29 ‘ Neuber theory solutions 7 Numerical Solutions 2;? 0 .01 .02 .03 .04 .05 .06 .07 .08 .09 .{0 r/ D Figure 6.13 Stress Concentration factor Kt for a notched flat bar Treatment g§_1rrg‘g;ar boundaries The uni-circular notch ie typical of the irregular boundaries often encountered in practical problem. Griffin and Verge (1962) chow how different mesh intervals can be used so that each horizontal. grid line crossing the irregulnr boundary intersects a vertical grid line on the boundary. There are several adventegee to this approach, but it doee increase the number of arithmetic operatione needed for the calculation of each vector iterate. For this reason e uniform neeh 82 interval is usedin ISOPEP, andthe irregularhanndaryieamroxinated with straight line segments, introduced arbitrarily, as shown in - Fig. 6.1L, connecting nodes of the primarpmeah. , “(gmhe dual-mesh is formed by lines halfway between lines [of theapttimary mesh.) The points labeled B are on the boundary and valueeof "31.1. are conputeduein; Equation (28) of Chapter III. The only complication is the evaluation of VZUi for points on the boundary. A few nameless till illustrate the technique. - ..L 3 a VZUO “ Ho 9g #45 ‘oeH ‘4 d 1e" 17" ' / L: 16» C x ‘ 158—3.! /l/ Th fie A J; 13 f l ‘ J l I I. l -- _. 1+ 5 6 7 8 9 10 1.1 12 J-—y figure 6.11: Boundary points for the semi-circular notch It an interior point where no . h? the line integral is approxilated (see Fig. .393 and Equation (27)) so a Wfiiw we: +‘Uz-Uo)'fff + wa- uni—3 +(u4-Uo)%—z] For the notch bounduy 9.15- a _ ay"/a, 33—0. 83 at a boundary point such as P153 (I a 15, J . 3), Figures 6.1!; and 6015<‘)a [1st —-=—-¢cose *(§-§L$M]3%, = g , =I35° 2- ‘ (b)'P18,11 . VZUO=F[UI+Uq-2U,+%Jo ,, A/C At the point ’15 6’ Films 6.11; and 6.15“); o I ‘ , c 'b L}; :15 z (‘f-U Us )9 f5 Sid!) z'(U'.-U°) . a flack/f may)? W] a y 4 cl c{ L-gcl‘szpii'é‘yUa-Uo) ., ; 3%®~M! [_UuU- 2U 52(in m] 02' + MW vzuo= ”#054917" 1+ Figure 6.15 Treatment of boundary points . ll . IillI i“ (in! “111.31! ‘1! . 8h The evaluation of V2U(15,6) is used in equation (28) for the finite- difference equations at the point Pm’ée It is not possible to write finite-difference equations for point ’15,,6 using Equation (28). .Hence the evaluation of “(15,6) and the discrete values of the stress function at other points labeled with a c in Fig. 6.11; are obtained by. interpola- tion. Fox (1950) treats this problem with the Gregory-{Newton forward interpolation forsula. A point 0 exterior to the boundary is introduced as shosn in Pu. 6.16. Iith the value of 03 and III; . (fig—)3 specified on the boundary, it is M g |~—h .4 possible to eliminate U o and O 1 2 3 +1; V x determine 01 at the first interior point in terns of the boundary conditions and one or more interior Figure. 6.16 Interior points points as follows near a boundary - (x-xa QQJX-XJX-x A3U.(x-x. x- - A’LL - - - U‘Uo'i-AUO h + 2! +4 +T—W4-WWW4' """‘ . Let s - (XE - xo)/h then U8 2 Ua + (U,-U,)S+(U2-2U,+U.)‘i%;!fl + (U3- 3U1+3U,"U.)"" {SEEM-2) + (Ut- ‘+U,+eu2-Lw,+u.)fl‘il}+§im +-— —-- , U2:(U'_U°)+‘UZ-2UI+U0)(S‘*)+(U3”3U3+3U,-Uo)'1§'2;3§!£:—2 + (Ut-‘+U3+suz-I+u,+u.)iL-fir’—-L€"' ‘:’-=1 - +---- 9 or US = mu, +D2U, +1302 + D,U,+D5u., +H, AU; = EnUO + EZUI + LUZ + Etua + E: U? + H2. where the D1 and 81 are fourth degree polynomials in I obtained by retaining differences up to fourth order while H1 and 32 involve fifth 8:3 and higher differences of U as well as higher degree terms in s. Neglect:- ing H1 and Hz the external point value can be eliminated and the expres« sion of U1 is u.=[%z— as + (gs — %)Ua + (%-§%)U,+(-Ef-Bf>ut]/(-§é- 15-2) . 0 Though it is preferable to retain differences of at least the fourth ' order, similar formulas can be written for higher or lower order differ» ences. Hoblem ._6_ - grectfiular hole _in an infinite plate The solution for an infinite plate subjected to uniform tensile stress Po at x 2.00 and x a ~00 is provided by Bevin (1951) for several different approximations of the boundaries of a rectangular hole. The stresses on the boundary of a square hole obtained from one of Savin's solutions will be used for comparison with the stresses obtained from a numerical solution. The problem solved numerically is for a 2 unit by 2 unit region surrounding the hole. The hole dimensions are 2/3 by 2/3. Values of the stress 11 ' JL function and normal derivatives on the outer boundary of the Pa 35 E *- . selected region are determined from Savin's solution. The nunerical solution of the biharmonic equation provides a basis for the determination of £9. lo" the stresses on the boundary of the hole,and these are compared with stresses from Savin's solution for a curvilinear Figure 6.3.? Problem 6 Infinite plate with a square hole 86 approximation to the square hole. Savin's solution is obtained using luskhelivshili's method. It is known that the stress function 95 (x,y) can be written ¢ (m) - u. [‘z‘mz) + 5(a)] where Re means the real part of the complex expression, (X (z) and 5 (a) are analytic functions of the complex variable a = x + 1y. '2— a: x-iy. Savin's solution is obtained by conformal mapping of the region exterior to the hole onto the interior of the unit circle in the complex ¢-plane by a mapping function 2 a w (g‘ ). The mapping function, approximated by three terms of an infinite series, and the complex stress functions a (fi‘) 5 (I: (1“? fl and the derivative of g (3‘) 5 @[OKS‘fl are given by Savin (1961) pp 51-53 as com) = RH? -— ZL¢3 +32% 3'71} (fig) = BREW + o.426§'+0.046(’3+o.008$15+ 0.00‘I' V719 ..~_.§_c/ CS")..- _:_ 0.55%?-0.657i'3—o.026<’5-o.0299'7l LP(S')_, d? "' BRLQS, + /+ 0.55“, —-O,/25878 ’ where 735:4- i)? .-_- (9e19, (f: ,7? ), ( E, 6) are rectangular and polar coordinates respectively in the complex q plane. Points on the circum- ference of the circle correspond to points on the boundary of the square hole, but when the series for Ca) ( $’) is truncated, the correspondence is not exact. The Cartesian coordinates of any point in the plate in terms of the polar coordinates in the complex plme are given by x = R (+C059-g-3C0539 +—§.—>g7-C0,S' 76) , , = .3 (fysine + game -— sumo. The boundary of the hole defined by these equations when p g 1 is not an exact square. The edges are slightly bowed and the corners are approximated with circular arcs of radius r e 0.02hsa, where a is the distance along the x-axis from one edge of the curvilinear square to 87 the other. H is the length of one side of the square hole. The value of _R a .39321508 is used in the problem solved numerically. This assures the transformation of the hole “Y —————--—~—~— boundary at a distance of 1/3 from the origin onto the bounds- arycfthe circle (Delinthe l 1 g; I complex plane. Thus a g 2/3, a l ’ l and the corner radius of curva- E I ture for Savin's solution is : ' 1 l r .— o.0163. 0 v" . ' Figure 6.18 Approxima- The equations given by tion of. the square hole Savin do not permit the direct evaluation of the stress function along the boundaries 1 g l and y ..-. 1. Derivatives can be found from ‘7'— %¥2+‘2—§§ == «(3’) + tagger) + W?) ‘ See, for example, Savin (1961) p. 6 or Huskhelishvili (1953) Po 183. Values of (3 which correspond to points on x a l and y a l, are deter:- mined for 0° 5 9 5 90° at intervals A9 a 3°. For these values of and 91% and .3; are calculated.“ The availability of a basic set of complex variable subroutines and statements in 3600 ' roam simplify these calculations. . The values of, the derivatives of ¢ thus determined relate to points unequally spaced along x e l and y a 1. A five point Lagrange interpolation formula is used to find values of the derivatives at equally spaced points. Along the x—axis, quadrature of g—g provides the point values of at up to an additive constant. The g; values specify a normal derivative. 1. similar procedure is used to obtain values of ¢ and normal derivatives along the boundary y g l. 88 The general problem of a multiplyeconnected region is discussed by Griffin and Varga (1962). The constants of integration can be arbitrarily chosen on an exterior closed boundary. However, for each interior closed bomdaly it is necesealy to determine “wig-21“” one point P]: on ax’ the boundary in such a w that the components:6 of displacement u.v and the rotation a), . J. (.._ 3' - 351;.) will be single-valued. The three additional mlmowns for each interior boundary are related to the point values of the stress function in the region surrounding the hole in a manner which permits use of iterative methods. Since the purpose of the “Y . ”sent nmerical solution was | d b B only to test to. efficacy of the -' finite-difference method in re- 33' producing details of the rapidly- " varying stress near a stress- &.L—F 3 concentration point by comparison ’ with Savin's solution, the proce- a, L, X dure of Griffin md Verge for a '6 26 l ' ” Figure 6.19 nomain of solution multiply-connected region was not 4 ., for Problem 6 followed. Instead the results of Gavin's solution were used to choose the integration constants on the boundary x a lin such a m that {1 -,.° at the point ( l/3,o ), and the assumption that d is sy-etric with respect to the x and y axes will satisfy the required conditions of single-valued- ness. The symmetry conditions then took care of the other constants on the hole boundary, and it was possible to solve the problem Iith one quarter of the. plate as. a simply-connected region, see Fig. 6.19. The x-axis is a line ‘of. symmetry. Hence gg— - 0 along the x-axis. The 89 ywaxis is also a line of symmetry. Hence éég— = 0 along the y-axis. Along Z? 2 _.JLQ ._ _ 312... 07""3y3-_o’ 7701“" QXAy—Os 22_ .21 - a _A _ ay ' (99):.‘0: x“(ax q— ,q ¢rB=o The value of B is zero since it has been assumed that ¢a = 0. Along 3? y- 3X‘—0 5):“ jangyo -' 09 §$= (Hie-‘0’ 3.2-“ 53%): ¢=C=¢F=O Hence i-O. Along R and Ed boundary values of the stress function and its normal derivatives are obtained from Savin's solution. Iterative solutions were obtained for mesh intervals of 1/15, 1/30, and 1/h8. Values of the stress function and a; at the inter- section of every fourth grid row and grid column of the M8 mesh are given in Fig. 6.20. A comparison of the stress concentrations found in the numerical solution and those given on p. 53 of Savin (1961 ) is provided in Table 11.. Values of 6'; given by Savin are the values of 6} in a rectang- ular coordinate system which has the origin at the point under investi- gation and the y-axis tangent at the given point to the curve 6 - 1. See Savin p. 8. The table includes comparable stresses, of, for the iterative solution with a convergence criterion of ‘IO'6 and a second set 65 for the solution obtained with a convergence criterion of 10-7. Both are included because this appears to be an exception to the general case in which stress concentrations increase as the convergence 90 FIGURE 6.20 Stress Function and 6‘: Distributions for Problem 6 COLUMN no. 3 1 11 13 19 23 21 31 33 39 63 61 31 .00000 -.00266 -.00886 —.01661 -.02600 —.03201 0.06057 -.o6119 c.05652 g ’ ...-IIIIIIIIIIIIIIIIIIOIIIIIIIIIIII IIIIIIIIII ..IIIIII ......IIII... .0000 .0110 .0009 .2013 .3216 .6300 .3230 .3911 is I I I I I I I I I I I I I : I I I I I I I I I I I I .00000 -.00260 -.00033 -.01309 -.02366 -.03136 -.o3016 -.o6311 -.03223 7 eeaeeeeeeeeeeeeeeeee eeee 000000;. eeeOee eeeeeleeeeeete 0.. as .0000 .0266 .1165 :2366 .366 69692 .3666 .613 36 I I I I I I I I I I I I I I I 200000 -200233 -200131 -201369 -202002 ~202663 -203311 -203933 -206320 11 eeeeeeeeeeeeeeeeeeeeeseieeaeeeeea ee 00 eeeeee eeeeeeeeoe eeeeeOeaneI .0000 .0660 .2110 .3332 .6300 :3393 .6036 .6396 Q I I I I I I I o e e e e I e e e 200000 -200139 -200662 -200031 -201316 -201613 -202326 -202060 -203363 1, eee-eeee aesaeeeeeaeeeeOeeeaseeeae eeaeeeeeeeeoee eoeeeee 00006 20000 . 123 .6011.3106 mi .6616 :6933 .9293 I I I I I I I I I I z I : z I .00000 0.00000 0. 00000 0. 00000 0200000 200163 200130 200033 -200191 -200303 -200031 -201239 -2o1633 19 e eeeeeee 00600.0 .6eeeeeeeeeeeeegoeeeeoeeeeeeeeeeeeeeeeeeeeeeeeeeeeee eeeeeeeeeeee 1.66662 266992 2121136 0363 6 6330 1.1013 .9119 .0660 .0119 .1999 .1992 .0032 . 2 2 2 2 2 2 2 2 2 2 2 0 a a .00313 200302 200609 200610 200061 .01161 .01606 201301 201603 201361 201066 0200900 0200660 II... III... IIIIII IIIIIIIIIIIIIIIIOIIIIIIIIIII I IIIIIIIIIIIII IIIIIIIIIIOI 126116 1262962 21. 6609 1.6390 21. 3063 1.3063 1.1916 1.0336 .9936 .9263 :6969 .0022 I I I I I I I I I I I z I I I I \ I I\ I I I I I I 202211 202296 202310J 202323 202190 203139 203610 203100 203031 203063 203011 203692 203323 21 e e eeeaeeeeeeeeeeeeee eeeeeeeeeeeeeeee eeaeeeeee eee eeeeeee so es 000 eeeee 000006000000... 1.6120 1.3006 1.6010 1.6210 1.3391 1.3220 1.966 1: 1326 126936 1.0196 :9130 .9600 : : : : : \ : : . : : . : : I I 203010 203033 203166 203369 203661 .06030 .06393 206106 206936 201013 201132 201121 201036 2212661222123929221236662212322922222661221226722212232622121996221212222212676I22120296222299952222 200139 200101 200912 209136 209.39 209001 210116 210322 210000 211026 211161 211262 2112.0 22 1291692212 365622122666221229952212 26662 222952212266622222996222226962212699622126639221269562222 e I e e e e e e e e e e 6 e e e e e e e e e e e e e e a .13370 :13615 213551 213715 .16076 .16626 .16797 .15155 15673 .1573? .15939 .16001 .1616? 9 eeeee ee eee eepeeeeeeeeeeeeeeeeeeeee..eeeee6006000600eeeeeeeeeeeeeeeeeeeee eee aeeeee'e ea 1. 2623 12 2321 1. :2273 1.2139 1.2053 1‘C953 1.1660 1.1653 1.1395 1.11 .0826 1.056 0 e e 0 e e e e e e e e a e e e e e e e \\ 0 e e e e e e '0 e e e e e e e e .10.66 .10910 .1906! .19259 .19565 .19880 .20260 .20597 .20929 .21221 .2166! .21653 .21792 eeeeeeee eee eeeeeeeee eeeeeeee eeeee.eeeeeeceeeeeeeeeeeeeeeeeeeeeeeeeeeeeeeeeeee see-eeeeeeeeeeeee 1.1936 1:1959 1.1051 1.1707 1.1731 1.1683 1.1612 1.1696 1.1326 1.1121 1.0900 1.0690 0 e e e e 0 e e e e e e e e e 9|. 0 e \ . e e \6 e e I e e e e e e e e e e e 225190 .25233 .2535 .25562 .25031 .2616? .26690 .2683. .27172 .27677 .27766 .2796? .20166 .1 ae eeeee eeeeeee 0‘ eeeeeeeee eeee eeeeeeeeee eeeeeeeeeeeeeeeeeeeee eaeeeaaeeealeee9000600600000 1.1599 w151.153‘1.1506 1.11616 1.1655 1.1616 1.1360 1.1229 1.1002 1.0916 1.0763 0 e e e e e e e g . e e e e e e e e e e \\ e e e e s C e O 6 a .32320 .32359 232671 .32665 232916 .33209 .33532 .33866 .36196 236503 .36702 235026. .235232 ‘ IIIIIIIIIIIIOIIOQIIIOIIII.IOIIOOQQIOOIICIOQII ......IIIIIIIIOIIII.I.IIIIIIIIIIIIIIIIIOIIII... \ | COLUMN //’Ex “0% 13 16 11 10 19 20 21 22 I /’0y \ .ggooo "861 010 -:03052 -. .colgo .03 34/ / .1, I, \ :,z§50ezazgg’aeee:iazaeeeaéeaeeeeel 3/ X/ x y .0000 .0620 .1612 .2231 .2 9—” 1 3 200000 -200013 -200061 -. 230916 -22 I U 0’ (DO (I: ~M— N l .00 0660 10 II eeeeeeeeeeeeeeeeeeee anOeeealese ..620 .3091 -.5029 .1693 -.1 :0000 0100 0 192 ’:3000 e 00 -. 000 -2 22 .0 030 - 66 0 22916 .3631W6 203 6 l7 eeseeeeeaeeeaee 000:00000006 -.0916 -.g162 -.1112 -.16 6 0 e 725 e 172 0397 0‘3 9 a e e O 00 006 .00 13 . 0021 00 26 I 0 23900 .0326 .8600 1206 2 1. eeeeee eeeeeeeeeeeeo 0000000 9 -.1321 -.o306 -.1292 -.1121 l .0000 .3113 .6639 .3031 .30 2 . 8881° 2- 88° 598881 528823 53%? .‘181? 9 aeeeeeae: 6e e :0 see 2 a: es s. ’ -8888 - 8’"81::"iiiéi"'-292§'°-‘sz1 I I I I I I 200066 2 00 200133 20021. 200260 200313 201.9901 2. 0 1.1333 1.3603 1.6031 1.2166 IIIIIIIII IIIIIIIIIIIIIIIIIIIIIIIIIIIII Ms - 7 11:2 use we: I000 I 7 I a I I1 I‘96 I I I I I I .00116 . 02 .00366 . 0613 0 619 o 3.0 1.0003 1. . 1.6611 1.2607 126291 123 90 21 eaaiaeeee 62222166922222612.212162522 e1988 23916 2 0 23663 29330 23231 23163 290309 2 31 280605 2 0619 203133 23002 2260 60‘ 1e I 10 565 60 132 le‘ .9 1e .72 eeee eeee eOee 0666600 00000 eeee .3336 . 5 .;.12 .661; :3996 -.581 2 2 . 3 22 3 001 9 2o 610 2 61 2 933 1 11 201 91 2 1161 I I I 23.6690 12 3 1.2865 122631 1.6316 1. 063 229962222 52222199622 2262992222 .§599221.1266 eZIII e 7 e299. .2920 e .7. a 325 91 Table 16. Comparison of boundary values of cr' obtained from.the iterative solution with those from Savin's solution. Savin's Solution ISOPEP Solutions 9' x; y' Cfis C51 (TE 1; y 0° .3367 o «.936 -1.005 -.902 .3333 w o -.991 -.895 .3333 .0833 -.973 -.883 .3333 .1250 -.966 -.861 .3333 .1667 —.897 -.828 .3333 '.2083 -.801 -.786 .3333 .2500 -.152 -.118 .3333 .2917 35° .3361‘ .2952 -.566 .170 .185 .3333 .3125 60° .3352 .3166 .605 65° .3296 .3296 6.368 6.78 6.52 .3333 .3333 50° .3166 .3352 6.660 .3333 3.602 3.281 .3125 .3333 55° .2952 .3361 2.888 2.710 2.665 .2915 .3333 2.176 2.057 .2500 .3333 1.979 1.869 .2083 .3333 1.880 1.775 .1667 .3333 1.819 1.722 .1250 .3333 1.780 1.690 .0833 .3333 90° 0 .3367 1.760 1.718 1.668 0 .3333 The stresses OE and 0" are compare 10 to 05 for iterative solutions with the convergence cr teron of 10' and 10-7 reapectively. 92 criterion is decreased. The values of cf; and ng along x = .3333 for y‘<.3333 are values of (j; computed from the point values of the stress function. Along y = .3333 for x<.3333 the values of UK and 073 are the values of 5; computed from the stress function solution. At x.= .3333, y = .3333 Which corresponds to 6 = 45°, the values of CyELgiven by Savin would be along a line which makes an angle of 45° ‘with the xpaxis. Values of’cy;} C5; and ’1gy computed at the corner of the square were used to compute values of (TR and C75 in the same direction as 0'3. The values, 0‘6 - 4.368, a; = 4.78 and 53 = b.52 are the largest stress concentrations occuring at corresponding points in the solutions considered. The fact that slightly different values are found for the numerical solutions might be expected from the different representation of the boundary at the corner in the numerical solution and Savin's solution. VII SUMMARY AND CONCLUSIONS The objective of this study has been three-fold; (1) The identi- fication of efficient iterative methods for the solution of plane elastostatic problems; (2) The preparation of a system of computer programs for solving this class of problems; and (3) A demonstration of the use of iterative methods. An investigation of the numerical solution of elliptic differential equations resulted in the selection of three matrix iterative methods as the alternatives which should be considered for inclusion in a digital computer program for the solution of plane elastostatic problems. Computer programs have been written for the point successive overrelax— ation method, the alternating-direction implicit method and the cyclic Chebyshev semi-iterative method. Solutions of a model problem for various mesh intervals and convergence parameters are used for compara ing the methods. The model problem.is a square plate with uniformly distributed loads on portions of two edges. See Section VI-1. The superiority of one iterative method over another may be judged by comparing the number of iterations required to satisfy a given con~ vergence criterion. The results given in Fig. 4.7 (a) show the cyclic Chebyshev semi-iterative method is iteratively faster than the other methods when the number of mesh points is lessthan approximately 125. This number of mesh points corresponds to a mesh interval h = 1/16. For h o for an 9,... Hell - 0. M Hall = lClll all where C is any ...m, (c) ”2+1”: ”g H + “I“, the triangular inequality. There are three norms of interest. For the vector §g= {X1, X2, X3 ---- In; , we define Io Ill HI = max lxil, IL HEHII = |x1| + lle 1' -... +lxnl: III. ”gunI =1/lx1l 2 +--- + |1rn|2 , Euclidean norm Since norms are often used as bounds it should be noted that urnI .<. ”sun .<. n “anI 103 Fadeeva (1959) states that a necessary and sufficient condition for a sequence ;(m) of vectors to converge to a vector §_is that “_X-(m) =§H-—> O as m——>OO for any norm satisfying conditions (a), (b) and (c); he gives a proof for the three norms I, II, and III. The norm of an n x.n matrix §,is a nononegative number [lg H which satisfies (a) IIAII> 0 1r Meandnai: =0. (b) Hcall = ICIHaH (c) Hysll Still! +1133.“ (a) Heal! .<. Hall Hail The necessary and sufficient condition for convergence of a sequence of nxnmatrices 5"“) is that Ham) - £||~+O as m—eOO. Any norm which satisfies (a), (b), (c) and (d) can be used to establish (m) all—*0 then I|g(m)|l->l1e!lo convergence. If ['5 As in the case of vectors it is possible to introduce matrix norms in a variety of ways. But for the purpose at hand it is convenient to introduce a matrix norm which satisfies some special requirements. One of these is the requirement for the matrix norm to be compatible with a given vector norm. This condition holds if MAX.” 5 Hall H ..X.l! . Compatimutx condiiiaa It is possible to determine a matrix norm which satisfies (a), (b), (c), (d) and which is compatible with a specified vector norm. For example, as a matrix norm compatible with the Euclidean vector norm X we define the s ectral norm 3A” of a matrix A as the H" “111’ p I’ "III ... least upper bound of the norm of the vector 5,; as ; runs over all vectors of norm unityg 11111111111 1.11:1 101» llélllll = max ”5. 39,111 With H E’HIII =1 2. Eigenvectors and eigenvalues An eigenvector (also called proper vector, characteristic vector or latent vector) of a linear transformation A is any non-zero vector x such that as= A; where A may be a complex number and is known as an eigenvalue (or proper number, latent root or characteristic number). If a vector x is an eigenvector it must satisify aux: +a12X2+ ‘ " ’ ’ “" amX”=;\x' QZle +QZZXL+ - - -- - - -- aznxn=XX1 amx. + anaxa+ - - ‘ ' - - ' ‘1an ’2’“? a system of homogeneous equations, which can have a non-trivial solution only if an ‘1 dig a,” am an") 0211 0.7” an; dnn‘)‘ The determinant is equivalent to an nt'h degree polynomial in A which is called the characteristic polyomial of the matrix A. For a real symmetric matrix all the eigenvalues are real. Definition: If _A_ is an n x :1 complex matrix with eigenvalues )1, 1 g i 5 n, then 6 (£1) = max | A1 I is the sEctral radius of _A_. 105 3. Properties of Matrices Theorem.A.1: If g'is any arbitrary matrix then Hill 2. ed). [Varga (1962) p.10] Theorem.A.2: If §.is an n x.n matrix then M = we?) where A? is the conjugate tranSpose of g, [Varga (1962) p.11] Definition: An arbitrary n x.n complex.nonsingular matrix 5 and an n- dimensional vector §_can be combined as g? g_; to form.a homogeneous second degree expression in terms of the components of §_which is called a quadratic form. If the value of the quadratic form is positive for all non-trivial g then the form and the matrix glare said to be pgsitixg definite. Definition: If §,is a n x.n complex.matrix such that any two elements situated symmetrically with respect to the principal diagonal are complex conjugates, aji ==§ij, then §_is called a Hermitian matggx. Notable characteristics are: 1. Diagonal elements must be real. 2. Any symmetric matrix with real elements is a Hermitian matrix. 3. The eigenvalues of a Hermitian matrix are real. 4. The eigenvalues of a positive-definite Hermitian matrix are positive. Definition: A square matrix.g_is diagonally dominant if “15?) 3 filau‘l J=l, Jae: If the strict inequality holdsfor all i, §,is said to be strictly diagonally dominant. Theorem.A.3: If §.is a Hermitian n x.n strictly diagonally dominant matrix with positive real diagonal entries, then §_is positive definite. 106 [Varga (1962) p. 23 1 Theorem Ll}: If A is an n x n Hermitian matrix, then || 5” = pg) [Varga (1962) p.11] Theorem Lg: If A is an n x n complex matrix, then A is convergent if and only if (3(5) <1 [Verge (1962) p. 13] Thus, we see that if the matrix A in Eq. (A-1) is symetric and positive definite, then 5 is convergent if the spectral radius of A is less than 1. In the next subsection the convergence rates of several iterative methods will be considered. A. Convergence Rates of Iterative Methods The system of linear equations (A-1) can be written in the form (1 - 1.4.) H. = E. (AL-3) where g is an n x n complex matrix and I_[_ is the identity matrix. If (:_L - g) is non-singular a unique solution vector 11 exists. Consider an iterative method for the solution of Eq. (4-3) 2(”1)-!U(“)+£ Ian-1.2.3 For any vector iterate 9f“) the difference Um) (A-h) - Q - go“) is a measure of the deviation of the vector iterate from the solution vector. The (hit) error vector _E_ can be written in terms of the preceding error vector §(m) g g E-(m-1) This is obtained by subtracting (A-B) from (A-h) and it can be readily shown E(In) a E” §(0) From Theorem A.5 we conclude that the error vectors g0“) will tend to zero for an arbitrary E“) if and only if g has a spectral radius €(_Ig)< 1. Theorem A.6: If A and g are n x n matrices, then (1) l|£||>0 or i=9; 107 (2) if k is a scalar, “kg“ = [k l [15,”; <3) Hue” .<. Hell + Hall3 m Mean s “an Hill: (5) Hi!“ _<_, ”A“ ll _qll for all vectors y_ of n components. Also, there is one nonzero vector _v in the n-dimensional vector space such that “A!” - “A“ Hill [Varsa (1962) m9] Since the solution vector iterates g0“) and the associated error vector iterates g“) are n-dimensional, Theorem 11.6 can be used as justification for New” - ”ram“ 0, ° ||f||<1 then R(.!m)' ~41! [llgmllk’ ]- Jbgbfi-IL- (A-S) is defined as the average 22 _of convergence £9; :_n_ iterations of the matrix 3. [Verge (1962) p.62.J The convergence of g is said to be iteratively faster than the convergence of 5 when Mg") > R (f). The significance of Rm") as a measure of the average rate of convergence may be seen by consideration of the average error’ reduction factor per iteration 0: HEM)” k" a” 113°)” J 108 For the case when “gull <1 0' s H!“ 11"" = ““3” (.-.) where e is the base of the natural logarithm. Thus this) is the exponential decay rate for the upper bound of the average error re- duction 0‘, per iteration in an m-step iterative process. Let NIn - 1/R(_I~f‘). Substituting into (A-6) we find O’Nm .<. 6 from which we conclude NIn is a measure of the number of iterations required to reduce the norm of Q“) by a factor e. Theorem A.7: R(Am), the average rate of convergence for m iterations of an n x n convergent matrix A has a limiting value of -In E(A) as m increases without bound. [Varga (1962) p.67 ] Definition: The asmtotic 22L! 2;; convergence Rx (A) is Rag (A) = -1n 9(5) Corollgz: Let A be a convergent n x n matrix, then Roe (A) .>_ ME“) for any positive integer m. For Hermitian matrices A and g, the spectral radii may be used for comparison of rates of convergence, since by Theorem [1.6, 9(a) <9(1.3)<1 implies ”an” My?) , so that of two convergent Hermitian matrices, the one with the smaller spectral radius will have a faster average rate of convergence for any 111. It should be noted that though H g“ ||—-> 0 as m—boo for two iterative schemes with matrices A and _B_, it is possible in general that for a selected value m1, matrix A may be iteratively faster than g but for a second value m2, 3 may be iteratively faster than A. .109 Identification of the convergent iterative methods for biharmonic difference equations is aided by consideration of the following three theorems. The form of matrix A'in Eq. (A-1) is the determining factor. Theorem A.8: If A_is a strictly diagonally dominant n x n complex matrix, then the associated Richardson‘s and Gauss~Seidel matrices are convergent and the corresponding methods are convergent for an arbitrary initial vector 11(0). [ Varga (1962) p. 73. ] The matrix ADassociated with the biharmonic difference equations is not strictly diagonally dominant and Windsor (1957) has shown that Richardson's method is not convergent for the biharmonic equation. Theorem A.9: If A_= Q +-g;+-gf is an n.x n Hermitian matrix where Q is Hermitian and positive definite and (Q +'§) is nonsingular, then the Gauss-Seidel iterative method is convergent if and only if A,is positive definite. (g? is the conjugate transpose of'g). [Varga p.78.] Theorem A.10: If i - p + _q + (2* is an n x n Hermitian matrix and g is Hermitian and positive definite, then the successive overrelaxation method is convergent for any arbitrary'gflo) if and only if‘A is positive definite and (Q +609.) is nonsingular for 0< (.0 < 2. [Varga (1962) p. 80.] 5. Optimum relaxation factors Definition: An n x n matrix g'which has zeroes and ones for elements and only one non-zero element in each row and each column is called a permutation matrix. Definition: Given A_an n x.n complex matrix with n >.1, if there exists an n x,n permutation matrix 2 such that [*0 lb '11 ll L$2 £22 j where A11 and A22 are square matrices of order k and (n-k) respectively, then A is called reducible. Otherwise A is called irreducible. Theorem A.11: If A is an irreducible n x n matrix, then: (1) A has a positive real simple eigenvalue equal to €(A) the spectral radius; (2) Increasing the value of any element of A will increase Q (A); (3) Corresponding to the eigenvalue (MA) there is an eigenvector with all its elements positive. [Theorem of Perron and Frobenius, Varga (1962) p.30] Definition: If A is an n x n irreducible matrix with non-negative elements which has a single eigenvalue of modulus €(A), it is said to be primitive. If A has I: eigenvalues with modulus of e (A) then A is gygl_:I._c 2f. _igggx A, k _>. 2. Each eigenvalue of modulus €(A) is a simple eigenvalue. [Verge (1962) p.35 ] Definition: A square matrix A of order n is said to be M EM 2A index A (k > 1) if an n x n permutation matrix A exists such that 2 9 2 5., A 0 O 0 '1' -2)1" ‘” " 252’s a .9 9 2 l 31' ' | \\1 | I 0 O A 0 L-' -' '-&ar- J where the null diagonal submatrices are square. A may or may not be reducible. A matrix can be simultaneously weekly cyclic of different indices. [Varga (1962) p. 39 ] Again consider the matrix equation 111 i932". where A is an n x n complex matrix which can be partitioned A“. 51,2 —’ "’ "" ALN 52,1 52,9. - - " £2,” .4. a (A‘7) _ AM! Am - - — Arm. where the square diagonal submatrices are nonsingular. Let Q be formed of the submatrices Ai’i, 5.1" Q 9. 9 i2), 9. Q- o o A ‘ L" " ”M; then A is also nonsingular, and the matrix equation can be written (i-mi+22-£ or _1 (Au-8) u-a2+2 F 1 where A - _I_ - _B_- A is the iteration matrix for the block Jacobi iterative method 2(m-1-1) ... 2 £011) + D-1 .E The matrix A is called the block Jacobi matrix of A. Definition: If the block Jacobi matrix A is weakly cyclic of index p, then A is Ecyclic in the partitioned form (A-7). [Varga (1962) p.99 J Definition: The p-cyclic matrix A is consistently ordered if all the eigenvalues of the matrix 112 _B_(k) as kfl + R413") .3. are independent of k for k 7‘ 0, where A is the block Jacobi matrix with sero diagonal elements, A and A are respectively strictly lower and upper triangular matrices such that A - A + A. [Varga (1962) p. 101] If the partitioned form of A is block tri-diagonal, it is con- sistently ordered and 2-cyclic. The optimum relaxation factor 60b which maximises the asymptotic rate of convergence of the block successive overrelaxation matrix for p = 2 is given by (4)6: -__3___. I + 1/1 - (”(5) where A is the block Jacobi matrix. [Varga (1962) P. 110 1 For the cyclic Chebyshev semi-iterative method the acceleration parameters are given 2 C .1 l/( ‘ (1),, = #12 , ‘ A . e (5) Cm(/é(5)) where A is the block Jacobi matrix. [Varga (1962) p. 138] The alternating-direction implicit method for the case of a fixed acceleration parameter r is a slight variation of Eq. (67) p.33 2(n+1) _ B:— U(n)+ 4!. If there is a block Jacobi matrix AR associated with Ar, the asymptotic rate of convergence will be a function of 62(h). Thus we find that the optimum relaxation factors for successive overrelaxation and the cyclic Chebyshev semi-iterative method as well as the acceleration parameter for the alternating direction implicit method are functions of the spectral norm of an associated block Jacobi matrix. In fact Verge (1962) shows that the three iterative methods considered have the same asymptotic rates of convergence for Laplace's equation solved for a square. 113 There appears to be a dearth of good approximations of (3 for irregularly shaped regions. It is common practice to use approximations based on the numerical results obtained in the iterative solution. Consider point successive overrelaxation. Forsythe and'Wasow (1960) p. 250 Show that the eigenvalues }d_of the matrix of the simplest iterative scheme, the point Jacobi method, are related to the eigenu values 721 of the matrix of the successive overrelaxation method by 71. A: . This is applicable to the eigenvalue )(1 of maximum modulus and hence to the spectral radius of the successive over- relaxation matrix. The Optimum relaxation factor as derived by Forsythe and wasow p. 253 is 2 2 L05 8 m : [fl (Aw-9) The error vectors from one iteration to the next are related Hamil 5. He II ll§(m")ll where A,is the appropriate matrix for the iterative scheme being considered. The dominant eigenvalue Th is the limit of II Am)" / [[§(m"1)“ as m increases without bound. Any vector norm of A,may be usede One computational procedure for estimating Cub consists of starting the fprcblem solution with w- 1, then after a number of iterations approxi— mate ‘fi; 7(1 ” ”Em ll / “Em-1)“ (A-10) and solve (A-9) for an approximate value of 00b. Another approach to approximating cub consists of selecting various values of Ce), running through several iterations for each and then by comparison of results, select the best. The procedure used in the ISOPEP code consists of computing a set 11b of values for 6.) obtained by applying Equations (A—9) and (1-10) every 10 or 20 iterations. This has the advantage of being an automatic procedure, and reasonably good results have been obtained for a number of problems. Approximating (db numerically as the solution proceeds is more computing art than science... Better methods for determining optimum relaxation factors would contribute much, to the usefulness of iterative ‘thOd' e APPENDIX B ISOPEPBA FORTRAN PROGRAM FOR THE ITERATIVE SOLUTION OF PLANE ELASTOSTATIC PROBLEMS The system of computer programs named ISOPEP was written as a general system for the analysis of plane elastostatic problems. It includes a set of FORTRANwII subroutine subprograms and a main program which provides linkage of a selected subset of these six subroutines and additional boundary-value subroutines which must be provided for each problem. Thislappendix includes a brief description of the general sub- routines, six examples of boundary-value subroutines, instructions for preparing input data for the program, a sample of the output, a flow diagram.and listings of the FORTRAN Source decks. I lasts This name is used for the system.of programs and also for the main program which provides linkage between the subroutines. The system was designed for use on an IBM 1620 with a‘20 K'main memory and a CDC 3600. The number of subroutines linked by ISOPEP may be reduced to increase the storage available for arrays. Only the Call statements in ISOPEP and all Dimension statements need to be changed. Only one of the three iterative method subroutines SORLX, ADI or CHEB is normally used. 9.92m OUTIN is the input and output subroutine. It provides for the initialisation of arrays. When problems using a relatively coarse mesh 115 116 spacing have converged, the point values are stored for later use in calculating an initial stress-function distribution for a finer mesh. 212119. The initial stress function distribution for the second or sub- sequent mesh refinements of a problem.can be computed by interpolating between the values obtained from an earlier coarseemesh solution. The problems may be solved consecutively, with the earlier solution saved in memory, or the preceding solution may be read from.punched cards. Execution of this option is controlled by the input of an appropriate value of the control number MESH. mass The stress components 0;, 0'; and TV will be computed and punched out on cards if specified by one of the options determined by the input of a control number NSTRS. The difference equations used to calculate the stress components at each interior point are 0;” —“'-'-’- [ U2;J'+1"2Ui,l +Ui,J-11/;.f- ' 0}“. =5 [ U}+hj"'2Ui,J' + Ui¢bJ1/hl 3 ”Km 5“ "'[-U2‘+1,J‘+T U2'+1,J-T Ui-1,J'+I+Ui‘bj"/Ih‘ ' £031.! AQAAA is the subroutine for the point successive overrelaxation iterative method. This has been adopted for general problems with irregular boundaries. As written, it is limited to simplybconnected regions for which all mesh lines,parallel to one of the Cartesian coordinate axes, are continuous segments connecting two boundary points. ADI The alternating-direction implicit method has been written as two 117 subroutines ADI1 and ADIZ. The second is executed immediately after the first. The specification of two subroutines makes more storage available for arrays when the programs are run on the disk-oriented IBM 1620. The solution domain is limited to rectangular regions. gggg The cyclic-Chebyshev semi-iterative method subroutine has been run only on the CDC 3600. It is written in FORTRAN II and could be divided into two or more subroutines for a computer with a limited main memory. Only rectangular regions can be treated with this program. It is necessary to specify an even number of interior mesh rows for this subroutine. PB N BD The values of the stress function at points on the boundary must be computed only at the beginning of the problem.aolution. This is part of the initialization of the solution array. A subroutine of this type must be provided for each problem. Examples are included for 1 5 N 56. PB N EX Derivative boundary conditions must also be provided for each problem. These may be treated by the introduction of exterior points or with special finite-difference equations for each point on the boundary. Both approaches are illustrated in the examples included for 1 g N _f 6. This subroutine must be executed during each iteration. lgpgt’Preparation The input to the problem can be provided on a single 80 column card unless point values of the stress function U13 are provided. The deckeof stress function values would be placed immediately after the 118 ITERATIVE SOLUTION OF PLANE ELASTOSTATIC PROBLEMS FLOW DIAGRAM 0F ISOPEP START ;;. I INITIALIZATION l KPE-1 INOT-1 ,fe CUTIE , Read g CALL OUTIN . Input Data ‘ L * PRINT Title of Solution Method Nbv \ PBNBD Calculation 0f Boundary CALL PBNBD CHANG Values U3 ::::::]< , Interpolation Routine for Initial 013 PBNEX ‘Calculation ,of Exterior Values UExt STRESS Calculation of Stress Components CHEB Routine for _, Cyclic Chebyshev ‘ Semi-iterativelhth. SORLI Successive Overrelaxation 'Iterative Routips ADI, Alternating - I Direction Implicit INOT , 2 _lterative Routine ‘ Ies:__ CNSTRS 0 3L 4 Y3; A fl OUTIN J Print (or Punch) m CALI. OUTIN M Cznerlgledgn Output } jNo ‘1 ll. ( (I'll, all... 119 control card which includes identification, dimensions and control numbers for input and output options. The stress function values must be provided in agreement with the FORMAT (I5,7F10.7/(5X,7F10.7) ). Output is punched in this form in order to simplify restart procedures. The FORMAT specification (F10.3,F10.8,2F10.4,101A), is used for the first card. The following list gives the use and symbolic name of each of the 1A numeric entries on the card. Column Symbolic Numbers Name Function 1-10 PRNO Problem.number for users identification. 11,20 DE Convergence criterion, usually in the range 1Om5=10°7. 21-30 REA Relaxation factor 314.0 SPYIJ H _E ”II Norm of the error vector in the preceding iteration when restarting. Set to 1.0 initially. h1-Lh HZ Number of mesh intervals along the xeaxis. AS-AB MY Number of mesh intervals along the yaaxis. 49-52 IF (2) If y = 0 is an exterior boundary (3) If y‘= 0 is a line of symmetry 53-56 JP (2) If x.= 0 is an exterior boundary (3) If x - 0 is a line of symmetry Note: Two mesh lines exterior to the domain are reserved for derivative conditions if the boundary is a line of symmetry otherwise only one line is reserved. 57-60 N N-so: Only one card of input is required. The initial value of the stress function at all interior points is set to 0.1. N>0: Count of number of iterations completed. A deck of Column Numbers Symbolic Names 120 Function 61-61. 65-68 69-72 73-76 N<02 NOUT ND ND point values of the stress function must be provided immediately following the first data card. This is the control number which terminates the processing of a sequence of data sets. A final data card should always be provided with this entry. Hhximmm.number of iterations permitted. This is the choice of the user and provides for punched output for restarting. Output will be printed after NOUT iterations and every subsequent set of NOUT iterations. 3,NT: The initial value of EPA will be used for all iterations. . ND'< NTz‘ The relaxation factor will be computed and DESK (1) (2) (3) changed at the end of every ND iterations. No mesh refinement calculations Read in values of the stress function from.a prior problem and interpolate to find an initial stress function distribution for the current problem. Use the stress function stored in.memery from the preceding problem.as the basis for interpolating te find an initial stress function distribution for the current problem. The user must be sure there is at least one problem.apecified on input cards pre- ceding this one. 121 Column Symbolic Numbers, Names Function A 'Am. i'.‘ -.,L-,_...i.- 14ae.-‘.uu .u‘.‘ “‘9‘ lL-~.—“n:—‘|’_u ‘Ammn 4.. “.-.-.. v I“.— Lr.W _ _. ‘_.._" --.-‘_..-‘__.. ..— 77—80 NSTRS '< 0 Punch Specified stress components INSTRsl = 1 Print. (and punch) 0.3., INSTRSI = 2 Print (and punch)O-Sf and o; INSTRSI .. 3 Print (and punch) Txy9 o; and 0; The second control card is used only if MESH = 2 on the first card. The four integers read in with a (L15) FORMAT Specification pr0w vide the dimensions of the problem solution prOVided on the cards following this control card and the number cf mesh spaces along the yaaxis of the problem to be solved. Column Symbolic Number Name Function ___ 1-5 IP Number of mesh columns of input data 6~1O JP Number of mesh rows of input data 11-15 M1 MY for the input solution 16-20 M2 MY of the problem to be solved. Example gf_Input Data First half of the first card Card columns O 1 2 3 A 123456789012345678901234292829 M2 67§20 1.190 .00001 1.5 1.0 Second half of the first card Card columns 5 6* 7 8 12345678291234567890123h56j8©01234567830 5 1O 3 2 1 200 200 10 1 -3 122 IBM 1620 OUTPUT OF SAMPLE PROBLEM ¢¢JOB 5 ##FORXSB 1 0078“ CORES USED SUCCESSIVE OVERRELAXATION 0 SOLUTION OF THE BIHARMONIC EQUATION C.L. DAVIS PROBLEM CONVERG. MESH SPACES MESH NUMBER CRITERION MX MY IF JF STRESS 1.190 .0000100 5 10 3 2 1 -3 RELAX. SUM OF ITERATIONS POINTS NOT FACTOR ERRORS COUNT MAX.N0. CONVERGED 1.5000 1.0000 0 200 50 1.6352 1.7275 10 200 “5 1.6325 1.h1ko 20 zoo #5 1.5618 “.7068 30 200 “S 1.u7o7 .1865 #0 zoo #5 1.8761 .0372 50 200 #1 1.558 .0136 60 zoo ht 1.598 .0082 70 zoo #3 1.5832 .0085 80 200 #0 1.5579 .0019 90 zoo so 1.5258 .0006 100 200 28 1.8878 .0002 110 200 .9 1.190 .0000100 1.h878 .0001 5 10 3 2 116 200 200 10 1 -3.. - . .. STRESS FUNCTION 2 3.6000000 3.0200000 2.8800000 1.9800000 .7200000 -.7200000 3 3.5072552 3.3395735 2.7961993 1.9099653 .696h0hh -.7200000 5 3.302709“ 3.1353795 2.6327760 1.798069“ .6555203 -.7200000 5 3.0285702 2.8790977 2. #277819 1.6699096 .6106627 -.7200000 6 2.6981876 2.5725967 2.1885523 1.5266087 .5626013- -.7200000 . 2.3065270 2.2096326 1.9063169 1.3587928 .5068087 -.7200000 1.8h37695 1.7788512 1. 5669582 1.1528863 .h365079-.7zooooo 9 1.3121758 1.2792279 1.1621542 .8963012 .3hk3551 -.7200000 10 .7h826h2 ..7h100h7 .705771h .5860380 .2259561 -.7zooooo 11 .2501268 .25h2h91 .2625820 .2h90h07 .0918619 -.7200000 12 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000 -.7200000 ISOPEP 15207 NEXT COMMON C L DAVIS '1'] I ll! f,.ll ._ - NVgCDUDGFQChU1aflflhJ 19=5CMO 53,82,82 82 CALL STRESS GO TO 50 END INSERT SUBROUTINE OUTIN - CHECK ARRAY DIMENSIONS INSERT SUBROUTINE CHANGE - CHECK ARRAY DIMENSIONS INSERT SUBROUTINE STRESS - CHECK ARRAY DIMENSIONS INSERT SUBROUTINE SORLX - CHECK ARRAY DIMENSIONS SUBROUTINE PB6BD BOUNDARY VALUES FOR REGION OF INFINITE PLATE WITH A SoUARE HOLE TYPE COMPLEX PHI PSI,PDU OMEC,DOMG DOMGB,DPHI OPHIB PSIB,zET DIMENSION U(64,6A),IK(6H ,lL(6h),DI5).P(31,BII,UI(6A) DIMENSION UXD(6A),UYD(6A , DIMENSION DUX(62),DUY(62),X(6z),Y(62),UBx(3I),UBY(3I) COMMON U,IK,D,MX,IF,IA,M|,MY,JF,JA,MJ,N,NOUT ND,NT,NDL,RFA,RF1, I SPYIJ,PRNO,DE,DI,KPE,INOT,MESH,NSTRS,MINC,IP,JP,AMI,IL,P,UI COMMON UXD,UYD LAGRANGE INTERPOLATION FORMULA FOR UNEOUAL INTERVALS GRANF (XA) =UA*( 5(fi§;¥2)*(XA-X3)*(XA-Xh)*(XA-X5)I/((X1-X2)*(XI-X3)*(XI-Xh)*(XI-X5)) + iifié;XI)*(XA-x3)*(XA—xu)*IXA-x5))/((X2-XI)*(x2-x3)*(x2-XA)*(xz-xs)) 2(XA;¥I)*(XA-XZ)*(XA~X4)*(XA-XS))/((X3-X1)*(X3-X2)*(X3-Xh)*(X3-X5)) +UD , ' g(Xé;XI)*(XA—X2)*IXA-X3I*IXA-xs))/((x4-XI)*(xu-x2)*(XA~X3)*(x4-x5)) +U _ 9IXA-XI)*(XA-Xz)*(XA~X3)*(XA-X4))/((xs-XI)*(x5-x2)*(X5-x3)*(xs-XR)) chXF(RH,TH)A(COSF(TH)/RH-RH**3*COSF(3.*TH)/6.+RH**7*COSF(7.*TH)/ I5 . *R chYFERHSTH)a(SINF(TH)/RH+RH**3*SINF(3.*TH)/6.-RH**7*SINF(7.*TH)/ 15 . * -R 22 FORMAT (zOIA) - 50 FORMAT (ABHO INFINITE PLATE WITH A SQUARE HOLE C L OAVIS//) 6O FORMAT (IIHO THETA = ,FIO.5,6HRHO a ,F8.2,IsH NO BOUNORY PT./) PRINT 5O 2 A=MX DIaI./A RaI.I796h523/3. PEEBOBLI'E59265 THETA=2.*PI LIM=37 I33 DANGaLIM-I DANGauPI/(2.*DANGI XMAXaI.0 YMAXfiII o O DLROa-.05 DO 6hU lafloLIM RHOaI.0 IF(THETA-7. *PI/h I 62D, 6IO, 6IO 6IO RIDX:XMAX-ABSF(CORXF(RHO, THETA)) IF(RIDX) 6I6, 63D, 6I2 6I2 RHO=RHO+DLRO IF(.IO—RHO) 6IO ,6IO,6IA 6Ih PRINT 6O,THETA, RHO . GO TO 630 6I6 Do 6I8 K=I 50 DRO=((CORXF(RHO, THETA)/XMAX-I )*RHO)/(R*(RHO**7*CO$F(7. *THETA)/7 - Iz. *RHO**3*COSF(3. *THETA)/3. )IXMAx-I. ) IF(ABSF(DRO)-. OOOOOOI) 63D, 630,6I8 6I8 RHO=RHO- DRO . - VA GO TO 6Ih 62D RIDYaYMAX—ABSF(CORYF(RHO, THETA)) IF(RIDY) 624 ,630,622 622 RHO .RHO+DLRO _ IF(O. I-RHo) 620, 62D, 6Ih 62h DO 626 K-I ,50 DRO=((CORYF(RHO, THETA)/YMAx-I )*RHO)/(R*(RHO**7*SINF(7. *THETA)/7. - I2. *RHO**3*SINF(3 *THETA)/3. )/YMAX-I. ) IF(ABSF(DRO)-.OOOOOOI) 630,630,626 626 RHOsRHO-DRO . GO TO 6Ih 630x I)=CORXF(RHO, THETA) Yl)=CORYF(RHO, THETA) ZET=CMPLX(RHO*COSF(THETAL RHO*SINF(THETA)) OMEGuR* I° /ZET-ZET**3/6 +ZET**7/56. ) DOMGaR* .I25*ZET**6- 5*zET**2-I IzET**2) DOMGBa CMRLX(REAL(DOMG),'-AI MAGIDOMG)). PH1=R*(. 25/2ET+ 426*2ET+.Oh6*zET**3+. OO8*2ET**5+. 00h*ZET**7) DPHI=R*(. 028*ZET**6+. 0h*ZET**h+ U38*ZET**2+. h26-. zs/zET**2) OPHIB: CMPLX(REAL(DPHI), -AI MAG( DPHI)) . PSla-R*( 5/2ET+(.5&8*2ET- A57*ZET**3-. 026*ZET**5-. 029*ZET**7)/ I (I.+. 5*2ET**A-. IZS*ZET**8)) PSlBa CMPLX(REAL( PSI), -AI MAG(PSI )) PDUaPHI+(0MEG*DPHfl B)/DOMGB+PSIB DUX(I)aREAL(PDU) DUY(I)=AIMAG(PDU) 640 THETA aTHETA + DANG C BOUNDARY POINTS REGION OF AN INFINITE PLATE READ 22, (IL(J), J33, JA) UBYEI ;:—.Osh5259 UBX .323200z LUMaLlM/Z A, DO 650 Rafi, LUM KaLIM-I UBX(I+I)=UBX(I)-. S*(DUX(K)+DUX(K+I))*(X(K+I)- X(K)) USA 65D UBY(B+I)2UBY(E)+.5*(DUY(I)+DUY(fl+I))*(Y(I+I)-Y(l)) ' BOUNDARY POINTS ON THE SQUARE K=IL(3) DO 652 l=3,K UEK,I;=0.0 652 U I,K =0.0 YVAR=0.0 K=3 DO 662 J=JF JA 654 IF(YVAR-Y(K5) 66D,66O,656 656 lF(K-LUM+2) 658,658,660 658 K=K+fi GO TO 65A 660 XI=Y(K-2) x2=Y(K-I) X3=Y(K ) Xh=Y(K+I) XS=Y(K+2) UA=UBY(K-2) UB=UBY(K-l) UC=UBY(K ) UD=UBYEK+I§ UEBUBY K+2 U(lA,J)-GRANF(YVAR) UA-DUX(K-2) UBaDUX(K-l UC=DUX(K ‘ UD=DUX(K+I) UEHDUX . A IIh X-X+DI RETURN END ALONG Y=O ...-g O0 0“” ALONG X-I ALONG Y-2 C SUBROUTINE Iasza BOUNDARY CONDITIC-N: I DIMENSION U(28 52), II DIMENSION ARI7I COMMON U, IK, D Mx, 1 SPYIJ, PRNo, DI, DII EQUIVALENCE '(AR(II, I) I2 22 FORMAT (2014) A-MX DI-I./A lB-lA-l READ 22, I PRINT 22, JBuIK(3)+I DO ZIO J-JB, JA U(IA J)-. 125 Y-O. 6 DO 211 1-3, IB U(I, JA)-. S*Y*Y fi2)-Y 2IO KrE, U(I EQUIVALENCE(D$3,U( I59) EQUIVALENCE (DHI, UI I6 IhZ IHE NOTCHED PLATE I:EZI.D(5) ILISZ) II MI ,MY, JF, JA MJ, N NOUT N INOI MESH, 'NSIRS, NINC, IP, 29 )5. (D51, ,U(Ié7)) IIDS 3”(054, ,U( I60) (055 ,U( IK J), J23, JA) IK J), J=3, JA) THE VALUE OF U ON THE CIRCULAR ARC U(|03)-05*(Y-0 25) Y=Y+DI I U22, I#)-U 22 ,3 ) U 26,15)-U 26, 3 ) x- O. 6 JB-JB-I DO 216 J-3, JB THE VALUE OF Y CM TH U(Ml, J)=, 211 E CIRCULAR ARC 5-SQRTF(. O625-X*X) THE VALUE OF U ON THE CIRCULAR ARC U(IA, J)=. S*(U(MI, J)-. 25) 0: JPM 2. I PROB. NO. 3 NT, NDL ,RFA,RFI, , m1 IL P,UI U( 61 (156 ) )), (DS6,U( 162)) THE VALUE OF DELY/H FROM THE MESH POINT TO THE CIRCULAR ARC 212, 212 .25)/D1 IF(J- 3) 216 DEL=(U(MI J IF(DEL- 1. I DELBDEL’ i0 GO TO 213 U(IA-I, J)-DEL X=X+DI “4 ,h UI6: 8 78 UI8,II =U I8; 3 U 19.12)=U 19.3 DO 2I7 K=2,11 KK-IA-K . U(KK ,h)=U(lA-I, K+3 DSI-OS*UE'A-I'5) I— 2I2 213 21h 215 216 ) ) ) ) 2I7 DSZ- 053.2.*U(1A-I,7)-I.. Dsu-. S*U(IA-I.9) I Dss.I. 5*U(IA~I.9)~.S ZIS. ZIA, ZIh ) lhB DS6I.S+DSh DHI=5.*DSS/3. AR l)=;5+DSl , AR u =.25* 053+656+I.O) AR S)-.5*DSS*DHI A AR 6 .;5+u IA-I,9)-AR(5) AR 7 -.5+O IA-I,IO) . RETURN - ‘ n H- END. , ‘ ' SUBROUTINE PBBEX { C EXTERIOR POINTS FOR NOTCNEO PLATE DIMENSION U(28,52),IK(52),O(5).IL(52) COMMON UIK,D NX,IF IA M|,MY JF,JA MJ,N,NOUT NO NT,NDL,RFA,RFI, 1lgPIAJiPfiNO,Dé,DI,K‘E,‘NOT,M§SH,NS§RS,NINC,I5,J5,AMI,IL,P,U1 220 no 223 l-3 IB ' U(I-MJ)-U(i,JA-I)«~ IF(i-IK(3)—I 222,222,223 222 u I,I;-uiI,5 ,. .,. .. . u I,2 -u I,A 223 CONTINUE c INTERIOR POINT ADJACENT TO CIRCULAR ARC ON HORIZONTAL LINE - 7 JQI‘O 60 T0 22A 218 IQZO -.- 2| Jul-I if - ' 22 .CA-I.-u(l.h) . ... CB-CA* CA-I. I2. CC—CB* CA-2. I3. CD-CC* CA-3. In. . BB-CA-25 . , BC-i 3.*CA-6.)*CA+2.)/6. . DAl-l.-CA+CBeCC+CD .. . .. DAB-CB-3.*CC+6.*CD . DAh-CC-h.*CD DBI-BB+BD-BC-l; ' oaz-I.-2.*OB+3.*OC-A.*BO O33-aB-3.*BC+6.*BO . OBN-aC-4.*OO ‘ a .1./(DA2/DA1-DBZ/DBI) u I J)-BE*(U(1 3)/DAI+(DBBIIBI-DA3/DAI)* U(I,J+I) IT ?é?/DBl-DAh/6Al)*U(l,J+2)+(BD/DBI-CD/DA1)*U(I,J+3)) a + H . H . H H IF (J-I ) 2I8,2I9 221 , 221 IF(l-2h 225.22A,226H 225 [-23 _. .. .. .... INA J=I5 GO TO 222 226 IF(I-ZS) 221+, 2213, 227 227 OD 229 J=1,JA . ”.‘z'é ikfii'fi'i IF(fi—IN(3)5 229, 229, 228 228 u(MI J)-U(IA-1, 3I+DI 229 CONTINUE C INTERIOR POINT ADJACENT TO CIRCULAR ARC ON VERTICAL LINE DUB. 5 l-IS CC-CB* CA-2. ls; CD‘BCC* CA-Bo ll}. DB-CA-.s BC.{§3. .*CA-6. g*CA+2. )/6. BD- (h.*CA-12 .)*CA+22. )*CA-6. )/2u. DA1-1.-CA+CB-CC+CD 1. DAz-CA-z. *CB+3, *CC-A. *CD DA3-CB-3.*CC+6.*CD DAA—CC-h.*CD , DBI-BB+DO-BC-I, DBZ-I.-2.*BB+3.*BC~#.*BD 033-33-3.*BC+6.*BD . DDAaDC-A. *BD BE-I./(DA2/DA1-DBZ/DBI) U I ,J)-BE*(U(IA, J)/DAI+DUD*DI/DBI+(D33/DBI-DA3/DAI)*U(I-I, J) 13 gé?/DBl-DAh/DAI)*U(I- -2 ,J)+(BD/DDI-CD/DAI)*U(I-3 ) II -|- . IF(J-7) 230, 230, 232 2321 I.-I6 - IF(J-9) 230,230, 233 233 Ui'S'Iw -U 15.5 .-. fig 2AA RETUfiN . . END_ SUBROUTINE PBABD — c BOUNDARY CONDITIONS FOR THE v-NDTCNED PLATE PROD. NO. 3 DIMENSION U(28, 52), IK(52), 0(5), P(25, 13) ,UI(25) IL(52) COMMON U, IK, D ,Mx, IF, IA ,MI ,MY JF, JA MJ, N 'NOUT N6 ,NT,NDL,RPA,RPI, I SPYIJ, PRNo' DE, DI, KPE, INOT,MESH, NS*RS, NINC, Ifi, JP ,AMI,IL,P,UI 22 FORMAT '(2OI&)' H . A-Mx .- Dl-1./A lB-lA-l READ I 22, ,(IKEJ), ,J-3, JA) PRINT 22,(IN J), JuB, JA) JB.IK(3)+1 DO 210 JaJD JA 210 U(IA,J)-.12§ J-s. 230 CA-l.-U(lA-‘ J) CB-CA*§CA-1.;/2. C IRS Y= ~O. 0 DO 215 |=3, IB U(I JA)=. 5*Y*Y IF(I-IK(3)) 215,215,212 BOUNDARY VALUES 0N NOTCH EDGE 212 KBI-M/N' U(IoK)=os*(Y‘025) 215 Y=Y+DI RETURN END SUBROUTINE PBAEX EXTERIOR POINTS FOR THE V-NOTCHED PLATE DIMENSION U(28, 52), IK(52), 0(5), P(25, 13) ,UI(25), IL(52) COMMON U, IK, 0, MX, IF, IA, MI :MY, JF, JA, MJ, N ,NOUT, N0 NT, NDL, RFA, RFI, PROB. NO. # I SPYIJ ,PRNO,DE,0I, KPE, INOT, MESH, NSTRS, NINC, IP, JP, AMI, 220 OD 223' I=3,IA . U(I MJ)=U(I,JA-I) IF(I-IK(3)-I) 222,222, 223 222U I, WI)=UEI I ,2)=U I 22A CONTINUE DO 228 J= 1,JA U(1,,J;=U 5,J) U(Z, J 3U “J J) IF(J-IK(3)I 228,228,226 226 U(MI J)-U(IA-I,J)+01... 228 CONTINUE . . RETURN END SUBROUTINE PBSBD BOUNDARY CONDITIONS FOR A PLATE WITH A SQUARE NOTCH PROB. NO. DIMENSION U(28, 52), IK(52), 0(5), IL(52) COMMON U IK, D MX, IF, IA, MI ,MY JF, JA, MJ, N ,NOUT N0 NT ,NDL, RFA ,RFI, I SPYIJ, PPNO' DP, DI, KPE, INOT, MESH, NSTRS, 'NINC, IP, JP,m 22 FORMAT '(2OIL)' IB-IA-I .- A-Mx , DI-I./A READ JB=|K(3)+I DO 205 Ja3, JD 205 U(JB, J)-O. O 00 210 J-JB, JA 210 U(IA, J)-.I25 Y: O. O DO 215 I=3, IB U(I JA)=. 5*Y*Y IF(I- JB) 215,215,212 212 U(I, JB)-. 5*(Y-. 25) 215 Y=Y+DI . . RETURN END 2. (|K(J), J=3, JA) PRINT 2, (IK(J), J=3, 'JA) IL, P, U1 ,IL,P,UI C d Ih6 SUBROUTINE PBsEx C EXTERIOR POINTS FOR A PLATE NITH A SQUARE NOTCH DIMENSION U(28,52)§IK(52), 0(5), IL(52) COMMON U IK, D Mx, IA MI 'MY JP, JA NJ, N, NOUT ND NT, NDL, RFA ,RFI, I SPYIJ PRNO, DE ,DI,KPE, INOI ,MPSH, NsIRs' ,NINC,IP, JP ,AMI, IL, P ,UI JB-IK (§)+l IB.IA- 220 DO 225 _, U(I‘NJ IF( 2220 UI' (I, JA-I ) K )-I 222,222,223 ,1 -U I,5 2 -U I ,h GO ID 225 223 U(I JB—I)-U(| ,JB+I) 225 CONIINUE . DO 228 J-I, JA UEI ,Jg-Usfi J; 2 Jw IF(J -JB) 22 226,226 226 U(MI J)-U(IA-I ,J)+DI. 228 CONTINUE H RETURN END BOUNDRY CONDITIONS CAN BE INTRODUCED USING EQUATION (28) AS INDICATED IN SECTION 6 - TREATMENT OF IRREGULAR BOUNDARIES. THE ADDITIONAL INSTRUCTIONS REQUIRED FOR SORLx FOR PROBLEMS 3 A AND 5 ARE LISTED UNDER THE HEADINGS SOR3 SORh AND SOR5. RESPECTIVELY, THE CARDS ARE INSPRTED.BETNEEN THE SUBROUTINE SORLx STATEMENT NUMBERS 133 AND Inc. 1 IL?” (3 00 0000000 . SOR3 IAI'IgITION TO SORLx FOR SEMI-CIRCULAR NOTCH II 4. IF(J-lK(3)-2) 230,230,130 230 OD 23I LI-I 5 . A -LI/3-I 5* *(-I)**Li . . L-J- (Ll-l) /3)*(-I)**( LI I- 12 231 2(5I.2.(U(K+I,L)+U( K-I, ,L).+U K ,L+I)+U(K ,L-l)-h».*U(K,L)) GO TO (233,235,236, 237,238 ,239.2AO,2AB, 25O,2sO,252,252,256, 260), K 233 Ifichoéo 23A D(2)-(TEMP+U(I,J)+.5*(U(|+I,J+I)+U(I+I,J-I)+DI)-2.*U(I+1,J))*ANV GO TO 232 235 ANV-l../AR§12 TEMP-DSI* U I+I,J+I)-U(I+I,J)) GO TO 233 H H H 236 ANY-I. /AR(2) TEMPI-DSZ*(U(I+I,J+1)-U(I+I,J)) TEMPu-TEMP-I-‘I'EMPI. .. GO TO 233 237 ANV-I./AR(3) I47 EEMIE-2§2Pn+' .S*(U(I+I, J+I)-U(I+I, J)) 238 0(2g- U(I J)+U(I+I, J-I)+U(I+I, J+I)-3. *U(I+I. J)+(U(IA. J)-U(I+1, J)) 0530(2) GO TO 2A2 239 022);(U(I J)+. 8*U(I I+I J+I)-I. 8*U(I+I,J)+. 5*DI)/AR(4) 35-.(2§ GO TO 2&2 ZhO D(2)-U(I, J)+U(U+I, J+fl)+U(9+fl, J—fl)-3. *U(I+I, J)+(U(IA, J)-U(I+I, J)) 1 /U(lA-I. J) DS-D(2) GO TO 2h2 E52)?§§(I+i, J+I)+U(B, J)- 2. *U(I+I, J)+DI /2. )/AR(7) 63.0 T0 2h2 DéZ):2. *(U(I+I, J+I)+U(I, J)- 2 .*U(I+I, J))+DI os-D(2) IF(J-13)21I2A ZSI ZSI 25I 0(2)-. *oslfla 05. I M 05.0(2 GO TO 2&2 252 D(2)-(U(I+I J+I)+U(I, J +056*U(I+2 J) -(2. +056)*U(l+1 Jg+ lg?g;(U(I+I, ,§)-U(I+I,J) /U(I+I, h)+zI.-DSG)*DI* 5)/AR{6 GO TO 221:225 256 IF(I2 257 D E) l 2h8D 250 (U I 2J +(U I—I J-I)+U(I+I J-I)-2. *U(l, J-I))*DSé+(U(I. 3)- S/U(I,u u(I,J-I))/Dsé. I'J-I/I Io 2&2 258 D 5;-(U U(-l, J)+D$6*U(l-U LJ-fl)-(l.+056)*u(| J-fi)+. *DI*DS6)/AR(h) l/Uglw 'IJ’*"("*“ J+“’*"(“*2 J)‘3 *U<'+1 J)+( U |+I.3)-U(I+1 J)) 05.0(25-U(I, J) KT .0 . 60 T0 252 260 IF(I-zh) 262 .263. 265 262 o(s)-os+u(I-I, J-I ) . A KT-I . GO 10 2A2 263 TEMPI-(U(I+I, J-I)-U(I J-I))*DSZ TEMP-TEMPI+§U(I- I,J- -I5-u(I,J-I))*. 5+. 5*( .5-052)*DI so T0268 - ' 26h IF(l-26) 265,266, 252 265 TEMP-TEMPI . . TEMPI-(U( I+IJ J-I)-u(I J-II)*DSI IBIP.-IEMP+IEMPI.+ 5*(652—0511 )*DI ANV-I./AR(2) so T0268 ( 1K0 G0 2 266 268 2&2 24h 2h3 245 230 231 232 234 235 £22 1&8 KT-O TEMP--TEMPI+. 5*DSI*DI ANV-1. lAR(1) . D(5)-(TEMP+U(I J)-2.*U(l J-I)+.5*(U(I-I,J-I)+U(I+I, J-I)))*ANV YlJ-(h. *o(3 -o(I)-o(2)-o(h)-o(s))*.05*RFA . Y2-ABSF(Y|J . . SYIJ-SYIJ+Y2 U(I ,J)=U(I, J)+Y|J IF(YZ-DE) 2h3, 244, 2K4 KPEBKPE+1 IF (KT) 1h0,1h0, 2h5 I-|+1 _ GO TO 230 SORh********(SEE COMENTS on PAGE lhé) Ao?érlou TO SORLX FOR A v NOTCH H I +1 IF(J-IK(3)- I) 230 230,140 no 231 LI-I n.. l-éL1/3-15*(-1)**Li L-J- (LI 1)/3)*(-l)**(L1-l D(Ll)-(U(K+1,L)+U(K-1,L)+U(K,L+1)+U(K,L-I)-h.*U(K,L)) IF(J-h) 232 234, 234 . g§2%?2)*(u(‘+1'J+1)+U(I'J)+U(|+I ,J-I)-3.*U(I+I.J)+Dl*.5)/3. 50 T0 2A2 0(5 )-05 IF(J-IK(3)-I) 235. 236 236 D(2)-2.*(U(| +1,J+1)+U(|, J)-2.*U(I+I,J))+DI DS-D(2) 32128322. w. )5.“ I, W O +1+O O YIJ-(h. *D(3§- -D(I)- 0(2 )- D(h #5 -D(5);* 05*RFA Y2-ABSF(YIJ . '. SYIJ.SYIJ+Y2 2% 230 231 232 2k2 U(I ,J)-U(I,J)+YIJ KPE-KPE+1 u.. . SOR5********(SEE COMENTS ON PAGE 1&6) ?D?éT'°N T0 SORLX FOR A SQUARE NOTCH .. I +1 IF(J-IK(3)-I) 230.230.1uo no 231 Ll-155 2.- K-l—-éL1/3-15*(-l)**Li L-J- (L1-l)/3)*(-l)**(Ll- -Iz D(L1)-(U(K+l L)+U(K-I, 'leu K, L+l)+U(K. L-l)-h. *U(K, L)) lF((J-lKE3)-15 232 2A2, 22 D(2)-(U I J *2 .+U(l+l J-I-I +U(l+l J- I)-lI. *U(l+l ,J))+DI YIJ-(h. *o(3 -o(I)-o(25-o( )-o(5)$*. 05*RFA . Y2-ABSF(YIJ . . . SY'J-SYIJ+Y2 IF(YZ-DE) Iuo Zhh, zuu u(I, J)-U(I, J)+YIJ' . Zhh KPE-KPE+I Ih9 TABLE OF SYMBOLS PRNO, DE, RFA, SPYIJ, Mx, MY, IF, JP, N, NT, NOUT, NO, MESH, NSTRS, IP, JP, Mu. M2. SEE INPUT PREPARATION FOR DEFINITION. U(I,J) I J MI IA - MJ JA P(J.I) UI(J) D1 KPE INOT IK(J) IL(J) DISCRETE VALUES 0F STRESS FUNCTIONS. Row INDEX. ' COLUMN INDEX. MAXIMUM l a MX+|F+I ' Ml-I'- LAST BOUNDARY RON INDEX 2 MAXIMUM J - MY+JE+I MJ-I -LAST BOUNDARY COLUMN INDEX STRESS FUNCTION DISTRIBUTION SAVED EOR GENERATION OF INITIAL ESTIMATE OE U(I.J) IN NEXT PROBLEM (SEE CHANG P. 129). A COLUMN ARRAY USED IN CHANG EOR TEMPORARY STORAGE. MESH INTERVAL H COMPUTED IN PBNBD; D1 DEPENDS ON PHYSICAL DIMENSIONS OF THE SOLUTION DOMAIN. ' COUNT OF POINTS AT WHICH CONVERGENCE CRITERION WAS NOT SATISFIED IN THE LAST ITERATION. CONTROL SWITCH EOR SUBROUTINE OUTIN. INOT - I SIGNALS INPUT OF NEXT DATA SET. INOT - 2 SIGNALS USE OF THE OUTPUT PORTION OF THE SUBROUTINE. LAST INTERIOR POINT ROW INDEX IN COLUMN J. EIRST INTERIOR POINT Row INDEX IN COLUMN J.‘ _ IK(J) AND IL(J) ARE INTRODUCED IN PBNBD FOR SPECI- FICATION 0F IRREGULAR BOUNDARIES. [’1“!Ill 11 .1 150 SPECIFICATION 0F COLUMN TERMINAL POINTS FOR SORLX A RECTANGULAR MESH WITH THE SPECIFIED DIMENSIONS OF THE ARRAY U ENCOMPASSES THE PROBLEM SOLUTION DOMAIN. FOR A REC- TANGULAR REGION THE ITERATIVE METHOD SUBROUTINE SORLX SNEEPS THE MESH BY COLUMNS STARTING NITH J - 3, USING 35;I