Introduction In this chapter, we will describe some modeling guidelines, including generally recommended mesh size, natural subdivisions modeling around concentrated loads, and more on use of symmetry and associated boundary conditions. This is followed by dis-cussion of equilibrium, compatibility, and convergence of solution. We will then consider interpretation of stress results. Next; we introduce the concept of static condensation, which enables us to apply the concept of the basic constant-strain triangle stiffness matrix to a quadrilateral element. Thus, both three-sided and four-sided two-dimensional elements can be used in the finite element models of actual bodies. We then show some computer program results. A computer program facilitates the solution of complex, large-number-of-degrees-of-freedom plane stress/plane strain problems that generally cannot be solved longhand because of the larger number of equations involved. Also, problems for which longhand solutions do not exist (such as those involving complex geometries and complex loads or where unrealistic, often gross, assumptions were previously made to simplify the problem to allow it to be described via a classical differential equation approach) can now be solved with a higher degree of confidence in the results by using the finite element approach {with its resulting system of algebraic equations}. .A. 7.1 Finite Element Modeling We will now discuss various concepts that should be considered when modeling any problem for solution by the finite element method. 7.1 Finite Element Modeling ... 351 General Considerations Finite element modeling is partly an art guided by visualizing physical interactions taking place within the body. One appears (0 acquire good modeHng techniques through experience and by working with experienced people. General-purpose programs provide some guidelines for specific types of problems f12, I5}. In subsequent parts of this section, some significant concepts that should be considered are described. In modeling, the user is first confronted with the sometimes difficult task of understanding the physical behavior taking place and understanding the physical behavior of the various elements available for use. Choosing the proper type of element or elements to match as closely as possible the physical behavior of the problem is one of the numerous decisions that must be made by the user. Understanding the boundary conditions imposed on the problem can, at times, be a difficult task. Also, it is often difficult to determine the kinds of loads that must be applied to a body and their magnitudes and locations. Again, working with more experienced users and searching the literature can belp overcome these difficulties. Aspect Ratio and Element Shapes The aspect ratio is defined as the rglio of the longest dimension to the shortest dimension of a quadrilateral element. In many cases) as the aspect ratio increases, the inaccuracy of the solution increases. To illustrate this point, Figure 7-1(a) shows five different finite element models used to analyze a beam subjected to bending. The element used here is the rectangular one described "in Section lO.2. Figure 7-1 (b) is a plot of the resulting error in the displacement at point A of the beam versus the aspect ratio. Table 7-1 reports a comparison of results for the displacements at points A and B for the five models, and the exact solution [2}. There are exceptions for which asPect ratios approaching 50 still produce satisfactory results; for example, if the stress gradient is close to zero at some location of the actual problem, then large aspect ratios at that location still produce reasonable results. In general, an element yields best results if its shape is compact and regular. Although different elements have different sensitivities to shape distortions, try to maintain (1) aspect ratios low as in Figure 7-1, cases 1 and 2, and (2) corner angles of quadrilaterals near 90 Figure 7-2 shows elements with poor shapes that tend to promote poor results. If few of these poor element shapes exist in a model, then usually only results near these elements are poor. ]n the Algor program [12J, when IX ~ 170° in Figure 7-2(c), the program automatically divides the quadrilateral into two triangles. D • Use of Symmetry The appropriate use of syrnmetry* will often expedite the modeling of a problem. Use of symmetry allows us to consider a reduced problem instead of the actual problem. '" Again, reflective s}mmetry means correspondence in size, shape, and position of loads: material propenies; and boundary conditions that are on opposite sides of a dividing line or plane. 352 ~ 7 Practical Considerations in Modeling; Interpreting Results Y,v /l---------------....l 4O,OOO-lb A T 8 in. ' A - - - -.. X. ~:: shear U l Parabolic 1·. . .- - - - 4 8 i n · - - - - - - l - I load distribution E = 30 x Ilf psi v =0.3 t 12 . In. = 1.0 in. I. x '2 In. 6 in. x I in. elements (typical) elements (typical) (5) AR=24 48' • In. n X (4) AR=6 II. '3 In.- 3 in. x 2 in . elements (typical) h (3) AR =3.6 elements (typical) (2) AR = 1.5 2.4 in. x 2~ in. elements (tYPical) (1) AR= 1.1 Figure 7-1 (a) Beam with loading; effects of the aspect ratio (AR) ilhJstrated by five cases with different aspect ratios Thus, we can use a finer subdivision of elements with less labor and computer costs. For another discussion on the use of symmetry, see Reference [3]. Figures 7-3-7-:-5 illustrate the use of symmetry in modeling (1) a soil mass subjected to foundation loading, (2) a uniaxially loaded member with a fillet, and (3) a plate with a.hole subjected to internal pressure. Note that at the plane of symmetry the displacement in the direction perpendicular to the plane must be equal to zero. This is modeled by the rollers at nodes 2-6 in Figure 7-3, where the plane of symmetry is the vertical plane passing through nodes 1-6, perpendicular to the plane of the modeL In Figures 7-4(a) and 7-5(a), there are two planes of symmetry. Thus, we need model only one-fourth of the actual members, as shown in Figure!:! 7-4(b) and 7-5(b). 7.1 Finite Element Modeling Exact solution ---------- 0 - 5 ~ (I) Q, 353 --- :Ie (2) -10 -= '0 . (3) -15 ~ ::::> -20 C It) e -25 8eo (4) ~ -30 ::0 -35 .5 g u -40 E -4$ ~ -SO l. (5) -55 2 6' 4 8 10 Asped ratio 12 14 16 18 20 22 24 lonsest dimension shortest dimension Figure 7-1 (b) Inaccuracy of solution as a function of the aspect ratio (numbers in parentheses correspond to the cases listed in Table 7-1) Table 7-1 Case 1 2 3 4 5 Comparison of results for various aspect ratios Aspect Ratio t.l 1.5 3.6 6.0 24.0 Exact solution [21 Number of Number of Elements Nodes 84 85 77 81 85 60 64 60 64 64 Vertical DispJacement, v (in.) Point A Point B Percent Error in Displacement atA -1.093 -1.078 -1.014 -0.886 -0.500 -0.346 -0.339 5.2 6.4 -0.328 -0.280 -0.158 23.0 56.0 -1.152 -0.360 11.9 Therefore, rollers are used at nodes along both the vertical·and horizontal planes of ' symmetIy. As previously indicated in Chapter 3, in vibration and buckling problems, symmetry must be used with caution since symmetry in geometry does not imply symmetry in an vibration or buckling modes. 354 " 7 Practical Considerations in Modeling; Interpreting Results '---__~I h ,b» h b (b) Approaching a triangular shape p~ a C / a»{3 large and very small comer angles (c) Very Figure 7-2 (d) Triangular quadrilateral Elements with poor shapes . -112 in. b---IO Ib/in. 6 ' 7 24 36 42 4g 54 60 66 19' 31 37 43 '49 55 '61 I.~.-~---------%iD'----------~11 ~--- Axis of sylll.lIletry Figure 7-3 Use of symmetry applied to a soil mass subjected to foundation loading (number of nodes 66, number of elements = SO) (254 cm 1 in.. 4.445 N = lib) = = Natural Subdivisions at Discontinuities Figure 7-6 illustrates various natural subdivisions for finite element discretization. For instance, nodes are required at locations of concentrated loads or discontinuity in loads, as shown in Figure 7-6(a) and (b). Nodal lines are defined by abrupt changes of plate thickness, as in Figure 7-6(c). and by abrupt changes of material properties, as in Figure 7-6(d) and (e). Other natural subdivisions occur at re~trant comers, as in Figure 7-6(0. and along holes in members, as in Figure 7-5. 7.1 Finite Element Modeling ... 355 lA=m~ 4 in. EL~L ---"'i1. . . - - 4 in. I ..I. 4 in. __---",;.I~.S~~~~~~~~~f3 1--3in.-J 200 Ib/in. (a) Plane stress uniaxially loaded member wich fille' (b) Enlarged finite element IllOdeI of the cross-batched quarter of the: member (number of nodes = 78. number of elements 60) (2.54 em '" 1 in.) Figure 7-4 Use of symmetry applied to a uniaxially loaded mer,nber with a fillet Sizing of Efements and the hand p Methods of Refinement For structural problems, to obtain displacements, rotations, stresses, 'and strains, many computer programs include two basic solution methods. (These same methods apply to nonstructural problems as well.) These are called the h method and the p method. These methods are then used to revise or refine a finite element mesh to improve the results in the next refined analysis. The goal of the analyst is to refine the mesh to obtain the necessary accuracy by using only as many degrees of freedom as final objective of this so called adaptive refinement is to obtain equal necessary. distribution of an error indicator over all elements. The discretization depends on the geometry of the structure, the loading pattern,~ . and the boundary conditions. For instance, regions of stress concentration or high stress gradient due to fillets, holes, or re-entrant corners require a finer mesh near those regions, as indicated in Figures 7-4, 7-5, and 7-6{f). We will briefly describe the hand p methods of refinement and provide references for those interested in more in-depth understanding of these methods. The h Method of Refinement In the h method of refinement, we use the particular element based on the shape functions for that element (for example, linear functions for the bar, quadratic for the beam, bilinear for the CST). We then start with a baseline mesh to provide a baseline solution for error estimation and to provide guidance for 356 ... 7 Practical Considerations in Modeling; Interpreting Results (a) Plate with hole under plane stress y _ _ _ Axis of symmetry Irlc+-+....~r---.. x (b) Finite element model of one-quart.er of the plate Figure 7-5 Problem reduction using. axes of symmetry applied to a plate with a hole subjected to tensile force mesh revision. We then add elements of the same kind to refine or make smaller elements in the model. Sometimes a uniform refinement is done where the original element size (Figure 7-7a) is perhaps divided in two in both directions as shown in Figure 7-7b. More often, the refinement is a nonuniform h refinement as shown in Figure 7-7c (perhaps even a local refinement used to capture some physical phenomenon, such as . a shock wave or a thin boundary layer in fluids) [19]. The mesh refinement is continued until the results from one mesh compare closely to those of the previously refined mesh. It is also possible that part of the mesh can be enlarged instead of refined. F or in~tance) in regions where the stresses do not change or change slowly1 larger 7.1 Finite Element Modeling (a) Concentrated load A. (b) Abrupt change of distributed load .......... Nodal line III f. 't2 Material CD Material (3) Nodal line 1'" \ Node (c) Abrupt change of pJate thickness (d) Abrupt change of material properties 'P2 P A (e) Basic model of an implant (cross-hatched) in bone, located at various depths X beneath the bony surface., using rectangular elements, (t)' Re-entrant comer, B P w (g) Structure with a distributed load Figure 7-6 Natural subdivisions at discontinuities (h) Using dements to distribute the loading and spread the concentrated load 357 358 ... 7 Practical Considerations in Modeling; Interpreting Results elements may be quite acceptable. The h-type mesh refinement strategy had its beginnings in !20-23). Many commercial computer codes, such 'as [12], are based on the h refinement. p Method of Refinement In the p method of refinement [24-28}, the polynomial pis increased from perhaps quadratic to a higher-order polynomial based on the degree of accuracy specified by the user. In the p method of refinement, the p method adjusts the order of the polynomial or the p level to better fit the conditions of the problem, such as the boundary conditions, the loading, and the geometry changes. A problem is solved at a given p level, and then the order of the polynomial is nonnaI1y increased while the element geometry remains the same and the problem is solved again. The results of the iterations are compared to some set of convergence criteria specified by the user. Higher-order polynomials nonnally yield better solutions. This iteration process is done automatically within the computer program. Therefore, the user does not p F ~ ~ ~ ~ ~ r F ~ F (a) Original mesh (b) A uniformly refined mesh F F % % % ~~ ~ F p % / (c) A possible nonuniform h refinement Figure 7-7 Examples of hand p refinement (d) A possible uniform p refinement 7.1 Finite Element Modeling Figure 7-7 A 359 (Continued) need to manually change the size of elements by creating a finer mesh, as must 1;>e done in the h method. (Tbe h refinement can be automated using a remeshing algorithm within the finite element software.) Depending on the problem, a coarSe mesh will often yield acceptable results. An extensive discussion of error indicators and esti· mates is given in the literature [l9}. The p. refinement may consist of adding degrees of freedom to existing nodes, adding nodes on existing boundaries between elements} andlor adding internal degrees of freedom. A uniform p refinement (same refinement performed on all elements) is shown in ·Figure 7-7d. One of the more common commercial computer programs, Pro/MECHANICA 129], uses the p method exclusively. A typicar discretized finite element model ofa pulley using Pro/MECHANICA is shown in Figure 7-7e. Transition Triangles Figure 7-4 illustrates the use of triangular elements for transitions from smaller quadrilaterals to larger quadrilaterals. This transition is necessary because for simple CST elements, intermediate nodes along element edges are inconsistent with the energy 360 .. 7 Practical Considerations in Modeling; Interpreting Results fonnulation of the CST equations. If intennediate nodes were used, no assurance of compatibility would be possible~ and resulting, holes could occur in the deformed model. Using higher-order elements, such as the linear-strain triangle described in Chapter 8, allows us to use intennediate nodes along element edges and maintain compatibility. Concentrated or Point Loads and Infinite Stress Concentrated or point loads can be applied to nodes of an element provided the ete· ment supports the degree of freedom associated with the load. For instance, truss elements and two- and three-dimensional elements support only translational degrees of freedom, and therefore concentrated nodal moments cannot be applied to these elements; only concentrated forces can be applied. However, we should realize that physically concentrated forces are usually an idealization' and mathematical conve· nience that represent a distributed load of high intensity acting over a small area. According to classical linear theories of elasticity for beams, plates, and solid bodies [2, 16, 17}, at a point loaded by a concentrated normal force there is finite displacement and stress in a beam, 'finite displacement but infinite stress in a plate, and both infinite displacement and stress in a two- or three-dimensional solid body. the consequences of the differing assumptions about the stress fields These results in standard linear theories of beams, plates, and solid elastic bodie~. A truly concen~ trated force would cause material under the load to yield, and linear elastic theories do not predict yielding. In a finite element' analysis, when a concentrated force is applied to a node of a finite element model, infinite displacement and stress are never computed. A concen· trated force on a plane stress or strain model has a number of equivalent distributed loadings, which would not be expected to produce infinite displacements or infinite stresses. Infinite displacements and stresses can be approached only as the mesh around the load is highly refined. The best we can hope for is that we can highly refine the mesh in the vicinity of 'the concentrated load as shown in Figure 7-6(a), with the understanding that the deformations and stresses will be approximate around the load) or that these stresses near the concentrated force are not the object of study, while stresses near another point away from the force) such as B in Figure 7-6(f), are of concern. The preceding remarks about concentrated forces apply to concentrated reactions as well. Finally, another way to model with a concentrated forte is to use additional ele'ments and a single concentrated load as shown in Figures 7-6(h). The ~hape of the distribution used to simulate a distributed load can be controlled by the relative stiffness of the elements above the loading plane to the actual structure by changing the modulus of elasticity of these elements. This method spreads the concentrated load over a number of elements of the actual structure. Infinite stress based on elasticity solutions may also exist for special geometries and loadings, such as the re...entrant comer shown in Figure 7-6(f). The stress is predicted to be infinite at the re...entranf corner. Hence, the finite element method based on linear elastic material models will never yield convergence (no matter how many times you refine the mesh) to a correct stress level at the re-entrant corner 118J. are 7.1 Finite Element Modeling ... 361 We must either change the sharp re-entrant corner to one with a radius or use a theory that accounts for plastic or yielding behavior in the material. Infinite Medium Figure 7-3 shows a typical model used to represent an infinite medium (a soil mass subjected to a foundation load). The guideline for the finite element model is that enough material must be included such that the displacements at nodes and stresses within the elements become negl~gibly small at ioeations far from the foundation load. Just how much of the mediUm should be modeled can be determined by a trialand-error procedure in which the horizontal and vertical distances from the load are varied and the resulting effects on the displacements and stresses are observed. Alternatively, the experiences of other investigators working on similar problems may prove helpful. For a homogeneous soil mass) experience has shown that the iniluence of the footing becomes insignificant if the horizontal distance of the model is taken as approximately four to six times the width of the footing and the vertical distance is taken as approximately four to ten times the width of the footing [4-6}. Also, the use of infinite elem~ts is described in Reference [13]. After choosing the ~orizontal and vertical dimensions of the model, we must idealize the boundary conditions. "Usually, the horizontal displacement becomes negligible far from the load, and we restrain the horizontal movement of all the nodal Points on that boundary (the right-side boundary in Figure 7-3): Hence, rollers are used to restrain the horizontal motion along the right side. The bottoIl! boundary can be completely fixed, as is modeled in Figure 7-3 by using pin suppo\:ts at each nodal point along the bottom edge. Alternatively, the bottom can be cbnstrained only against vertical movement. The choice depends on the soil condltions ~t the bot~ tom of the model. Usually, complete fixity is assumed if the lower boundary is taken as bedrock. In Figure 7-3, the left-side vertical boundary is taken to be directly under the center of the load because symmetry has been assumed. As-we said before when discussing symmetry) all nodal points along the line of symmetry are restrained against horizontal displacement. Finally} Reference [111 is recommended for additional discussion regarding guidelines in modeling with different element types, such as beams, plane stress/plane strain, and three-dimensional solids. Connecting (Mixing) Different Kinds of Elements Sometimes it becomes necessary in a model to mix different kinds of elements, Sl,lch as beams and plane elements, such as CSTs. The problem with mixing these elements is that they have different degrees of freedom at each node. The beam allows for transverse displacement and rotation at each -node, while the plane element only has"inplane displacements at each node. The beam can resist a concentrated moment at a node, whereas a plane element (CST) cannot. Therefore~ if a beam element is CODnected to a plane element at a single node as shown in Figure 7-8(a), the result will be a hinge connection at A. This means only a force can be transmitted through the 362 J.. 7 PraCtical Considerations in Modeling; Interpreting Results Plane elements Plane elements (a) (b) Figure 7-8 Connecting beam element to plane elements (a) No moment is transferred, (b) moment is transferr,ed node between the two kinds of elements. This also creates a mechanism, as shown by the stiffness matrix being singular. This problem.can be corrected by extending the beam into the plane element by adding one or more beam elements, .shown as AB, for one beam element in Figure 7-8(b}. Moment can now be transferred through the beam to the plane element. This extension assures that translational degrees of freedom of beam and plane element are connected at nodes A and B. Nodal rotations are associated with only the beam element, AB. The calculated stresses in the plane element will not nonnally be accurate near node A. For more examples of connecting different kinds of elements see Figures 1-5, 11-10, 12-10 and 16-31. These figures show examples of beam. ~nd plate elements connected together (Figures 1-5, 12-10, and '16-31) and solid (brick) elements connected to plates (Figure 11-10). Checking the Model The discretized finite element model should be checked carefully before results are computed. Ideally, a model should be checked by an analyst not involved in the preparation of the model, who is then more likely to be objective. Preprocessors with their detailed graphical display capabilities (Figure 7-9) now make it comparatively easy to find errors, particularly the more obvious ones involved with a misplaced node or missing element or a misplaced load or boundary support. Preprocessors include such niceties as color, shrink plots, rotated views, sectioning, exploded views, and removal of hidden lines to aid in error detection. Most commercial codes also include warnings regarding overly distorted element shapes and checking for sufficient supports. However, the user must still select the proper element types, place supports and forces in proper locations, use oonsistent units, etc., to obtain a successful analysis. Checking the Results and Typical Postprocessor Results The results should be checked for consistency by making sure that intended support nodes have zero displacement, as required. If symmetry exists) then stresses and displacements should exhibit this symmetry. Computed results from the finite element .. "'1111 !=\hould be compared with results from other available techniques, even if 7.2 Equilibrium and Compatibility of Finite Element Results '" 363 Figure 7-9 Plate of steel (20 in. long, 20 in. wide, 1 in. thick, and with a Hn.'i'adius hole) discretized using a preprocessor prowam [15] with automatic mesh generation these techniques may be cruder than the finite element results. For instance, approximate mechanics ofmaterial formulas, experimental data, and numerical analysis of simpler but similar problems may be used for comparison, particularly if you have no real idea of the magmtude of the answers. Remember to use all results with some degree of caution) as errors can crop up in such sources as textbook or handbook comparison solutions and experimental results. In the end, the analyst should probably spend as much time processing, checking, and analyzing results as is spent in data preparation. Finally. we present some typical postprocessor results for the plane stress problem of Figure 7-9 (Figures 7-10 and 7-11). Other examples with results are shown in Section 7.7. :I 7.2 Equilibrium and Compatibility of Finite Element Results An approximate solution for a stress analysis problem using the finite element method based on assumed displacement fields does not generally satisfy all the requirements for equilibrium and compatibility that an exact theory-of-elasticity solution satisfies. However. remember that relatively few exact solutions exist. Hence, the finite element method is a vr;ry practical one for obtaining reasonable, but approximate~ numerical solutions. Recall the advantages of the finite element method as described in Chapter 1 and as illustrated numerous times throughout this text. 364 ... 7 Practical Considerations in Modeling; Interpreting Results 1000 psi 20 in. Figure 7-10 Plate with a hole showing the deformed shape of a"plate superimposed over an undeformed shape. Plate is fixed on the left edge and subjected to 10DO-psi tensile stress along the right edge. Maximum horizontal displacement is 7.046 x 10-4 in. at the center of the right edge We now describe some of the approximations generally inherent in finite element solutions. ' 1. 'Equilibrium of nodal forces and moments is satisfied. This is true because the global equation E = K 51 is a nodal equilibrium equation whose solution for 51 is such that the sums of all forces and moments applied to each node are zero. Equilibrium of the whole structure is also satisfied because the structUre'reactions are included in the global forces and hence in the nodal equilibrium equations. Numerous example problems, particUlarly involving truss and frame analysis in Chapter 3 and 5, respectively, have illUstrated the equilibrium of nodes and of total structures. 2. Equilibrium within an element is not always satisfied. However, for the constant-strain bar of Chapter 3 and the constant-strain triangle of 'Chapter 6, element equilibrimn is satisfied. Also the cubic displacement function is shown to satisfy the basic beam equilibrium differential equation in Chapter 4 and hence to satisfy element force and 7.2 Equilibrium and Compatibility of finite Element Results .& 365 2_ -..:t2:t ZT'IIl.7 1 ,-- ',_1'" 1UO"'I%l 11S42.1>11 m.llOOC 6"'_ -= 1.:D!B:!..o12 Figure 7-11 Maximum principal stress contour (shrink fit plot) for a plate with hole. Largest principal stresses of 3085 psi occur at the top and bottom of the hole, which indicates a stress concentration of 3.08. Stresses were obtained by using an average of the nodal values (called smoothing) moment equilibrium. However, elements such as the linear-strain triangle of Chapter 8, the axisymmetric element of Chapter 9, and the rectangular element of Chapter 10 usually only approximately satisfy the element equilibrium equations. 3. Equilibrium is not usually satisfied between elements. A differential element including parts of two adjacent finite elements is usually not in equilibrium (Figure 7-12). For line elements, such as used for truss and frame analysis, interelement equilibrium is satisfied, as shown in example problems in Chapters 3-5. However, for two- and threedimensional elements, interelement equilibrium is not usually satisfied. For instance, the results of Example 6.2 indicate that the nonnal stress' along the diagonal edge between the two elements is different in the two elements. Also, the coarseness of the me.§h causes this lack of interelement equilibrium to be even more pronounced. The normal and shear stresses at a free edge usually are not zero even though theory predicts them to be. Again, Example 6.2 illustrates this, with 366 £. 7 Practical Considerations in Modeling; Interpreting Results ~----------------7r--_5~lb CD ,...-- lOin. t _ _ .J ~------------------~~S~lb 20 in. Ex.ample6.2 a, 0',. ::; 301 psi r----'-----::~ • or..,.. = 301 psi .!\= ax = 995 ~ 2.4 psi ~ t. r = 2.80 ~si = - 2.4 psi :.--a~ psi 0; = 'OOSpa '1:.., 440 psi =:: ~ '~PS;99Spa = 397 psi 'Ij,,, = 0'., :: Stresses on a differential element common to both finite elements, iIIus£rating violation of equilibrium - 2.4 -1.2 psi Srress along the diagonal between elemenlS, showing normal and shear suesses, , (I,. and'tnl' Note: aft and t.t are not equa.l in magnitude but are opposite in sign for the two elements, and so interelement equilibrium is not satisfied Figure 7-12 Example 6.2, illustrating violation of equilibrium of a differential element and along the diagonal edge between two elements (the coarseness of the mesh amplifies the violation of equilibrium) free-edge stresses O'y and '1:xy not equal to zero. However, .as more elements are used (refined mesh) the O'y and'txy stresses on the stressfree edges will approach zero. 4. Compatibility is satisfied within an element as lOrlg as the element displacement field is continuous. Hence, individual elements do not tear apart. 5. In the formulation of the el~ent equations, compatibility is invoked at the nodes. Henc:e) elements remain connected at their common nodes. Similarly, the structure remains connected to its support nodes because boundary conditions are invoked at these nodes. 6. Compatibility mayor may not be satisfied along interelement boundaries. For line elements such as bars and beams, interelement boundaries are merely nodes. Therefore, the prec:eding statement 5 applies for these line elements: The constant-strain triangle of Chapter 6 and the rectangular element of Ghapter 10 remain straight-sided when deformed.. Therefore, interelement compatibility existS for these elements; that is, these plane elements deform along common lines psi 7.3 Convergence of Solution A 367 without openings, overlaps, or discontinuities. Incompatible elements, those that allow gaps or overlaps between elements, ~an be acceptable and even desirable. Incompatible element formulations, in some cases, have been shown to converge more rapidly to the exact solution [I}. (For more on this special topic, consult References [7J and 18].) 1 7.3 Convergence of Solution In Section 3.2, we presented guidelines for the selection of so--called compatible and complete displacement functions as they "related to the bar element. Those four guide. lines are generally applicable, and satisfaction of them has been shown to ensure monotonic convergence of the solution of a particular problem [9J. Furthermore, it has been shown {I 0] that these compatible and complete displacement functions used in the displacement formulation of the finite element method yield an upper bound on the true stiffness, and hence a lower bound on the displacement the problem, as shown in Figure 1-13. Hence, as the mesh size is reduced-that is, as the number of elements is increased-we are ensured of monotonic convergence of the solution when compatible and complete displacement functions are used. :txamples of this convergence are given in References [1J and [11], and in Table 7-2 for the beam with loading shown of Exact solution Number of elements i~ ~----------------~----~ '\ ~ is "c ompall" Ie bdlsplatement ' fonnulation Figure :'J -13 Convergence of a finite element solution based on the compatible displacement formulation Table 7-2 Comparison of results for different numbers of elements Case Number of Nodes Number of Elements 21 39 12 2 24 I 45 85 105 32 6480 3 1.5 1.2 2 3 4 5 Exact solution [2J Aspect Ratio Vertical Displacement, v (in.) Point A -0.740 -0.980 -0.875 -1.078 -1.100 -1.152 368 £ 7 Practical Consiclerations in Modeling; Interpreting Results in Figure 7-1(a). All elements in the table are rectangular. The results in Table 7-2 indicate the influence of the numbe~ of elements (or the number of degrees of freedom as measured by the number of nodes) on the convergence toward a common solution, in this case the exact one. We again observe the influence of the aspect ratio. The higher the aspect ratio, even with a larger num~r of degrees of freedom, the worse the answer, as indicated by comparing cases 2 and 3. :l 7.4 Interpretation of Stresses~ In th~ stiffness or displacement formulation of the finite element method used throughout this text, the primary quantities determined are the interelement nodal displacements of the assemblage. The secondary quantities, such as strain and stress in an element, are then obtained through use of {t} = [B]{d} and {u} = [D][B]{d}. For elements using linear..rusplacement models, such as the bar and the constant-strain triangle, [BJ is constant, and since we assume [DJ to be constant, the stresses are constant over the element. In this easel it is common practice to assign the stress to the centroid of the element with acceptable results. However. as il1ustrated in Section 3.11 for the axial member, stresses are not predicted as accurately as the displacements (see Figures 3-32 and 3-33). For example, remember the constant-strain or constant-stress element has been used in mode1· ing the beam in Figure 7-1. Therefore, the stress in each element is assumed constant. Figure 7-14 compares the exact beam theory solution for bending stress through the beam depth at the centroidallocation of the elements next to the wall with the finite element solution of case 4 in 'Table 7-2. This finite element model consists of four elements through the beam depth. Therefore) only four stress values are o.btained y(in.) 4 t - - - - - - i i > 174.4 t2h 130.8 2 Finite element solution::; e Exact solution 39 43.6 --------.1'--'--100 ........1-'5-0-200....1.-- (T'~ (ksi) -39 -2 e-122 -3 "'------t- -4 -174.4 Figure 7-14 Com parison of the finite element solution and the exact solution of bending stress through a beam cross section 7.5 Static Condensation A.. 369 through the depth. Again, the best approximation of the stress appears to occur at the midpoint of each element, since the derivative of displacement is better predicted between the nodes than at the nodes. . For higher-order elements, such as the linear-strain triangle of Chapter 8, [BI} and hence the stresses, are functions of the coordinates. The common practice is then to evaluate directly the stresses at the centroid of the element. An alternative procedure sometimes is to use an average (possibly weighted) value of the stresses evaluated at each node of the element. This averaging method is often based on evaluating the stresses at the· Gauss poInts located within the element (described in Chapter 10) and then interpolating to the e~ment nodes using the shape functions of the specific element. Then these stresses in all elements at a common node afe averaged to represent the stress at the node. This averaging process is called smoothing. Figure 7-11 shows a maximum. principal stress "fringe carpet" (dithered) contour plot obtained by smoothing. Smoothing results in a pleasing, continuous plot which may not indicate some serious problems with the model and the results. You should always view the unsmoothed contour plots as well. Highly discontinuous contours between elements in a region of an unsmoothed plot indicate modeling probl~ms and typically require additional refinement of the element mesh in the suspect region. . 'If the discontinuities in an unsmoothed contour plot are small or are in regions of little consequence) a smoothed contour plot can normally be used with a high degree of confidence in the results. There are, however, exceptions when smoothing leads to erroneous results. For instance, if the thickness or material stiffness changes significantly between adjacent elements, the stresses will no~ally be different from one element to the next. Smoothing will likely hide the actual results. Also, for shrinkfit problems involving one cylinder being expanded enough by heating to slip over the .smaller one, the circumferential stress between the mating cylinders is normally quite different [16]. The computer program examples in Section 7.7 show additional results, such as displaced models, along with line contour stress plots and smoothed stress plots. The stresses to be plotted can be von Mises (used in the maximum distortion energy theory to predict failure of ductile materials subjected to static loading as described in Sec.tion 6.5); Tresca (used in the Tresca or maximum shear stress theory also to predict failUI;e of ductile materials subjected to static loading) [14, 16], and maximum and minimum principal stresses. "" 7.5 Static Condensation We will now consider the concept of static condensation because this concept is used in developing the stiffness matrix of a quadrilateral element in many computer programs. Consider the basic quadrilateral eiement with external nodes 1-4 shown in Figure 7-15. An imaginary node 5 is temporarily introduced at the intersection of the diagonals of the quadrilateral to create four triangles. We then superimpose the stiffness matrices of the four triangles to create the stiffness matrix of the quadrilateral 370 .. Pract~cal 7 Considerations in Modeling; Interpreting Results y," I 3 CD // " / ® 0)<, /5, .1/ CD " 4 " Figure 7-15 Quadrilateral element with an internal node J 2 I ' - - - - - - - - x. u element, where the internal imaginary node 5 degrees of freedom are said to be condensed out so as never to enter the final equations. Hence, only the degrees of freedom associated with the four actual external comer nodes enter the equations. 'We begin the static condensation procedure by partitioning the equilibrium equations as {7.5.l ) where di is the vector of internal displacements corresponding to the imaginary' internal node (node 5 in Figure 7-15), fi is the vector of loads at the internal node, and do and Fa are the actual nodal degrees of freedom and loads, respectively, at the actual nodes. Rewriting Eq. (7.5.1), we have [kll}{du } + [k!2l{di } = {Fa} + [k22J{dr} = [k21 ]{da } {F;} (7.5.2) (7.5.3) Solving for {di } in Eq. (7.5.3), we obtain {di } = -[k22 r I [k21 l{da } + [k22r l {Fi} (7.5.4) Substituting Eq. (7.5.4) into Eq_ (7.5.2), we obtain the condensed equilibriwn equation [kd{d:J} where = {Pc} (7.5.5) [krJ = [kill - [k I2 J[k:rd- 1[k2d (7.5.6) [kI211k22rl {Fi} (7.5.7) {f;} = {F lI } and [k/.J and {F;.} are called the condensed stiffness matrix and the condensed load vec 101', respectively. Equation (7.5.5) can now be solved for the actual comer node displacements in the usual manner of solving simultaneous linear equations. Both constant-strain triangular (CST) and constant-strain quadrilateral elements are used to analyze plane stress/plane strain problems. The quadrilateral element has the stiffness of four CST elements~ An advantage of the fOUI..csr quadrilateral is that the solution becomes less dependent on the skew of the subdivision mesh, as shown in u 7.5 Static Condensation ... 171 Figure 7-16. Here skew means the directional stiffness bias that can be built into a model through certain discretization patterns, since the stiffness matrix of an element is a function of its nodal coordinates) as indicted by Eq. (6.2.52). The four-CST mesh of FigiJre 7-16(c) represents a reduction in the skew effect over the meshes of Figure 7-16(a) and (b). Figure 7-16(b) is generally worse than Figure 7-16(a) because the use of long, narrow triangles results in an element stiffness matrix that is stiffer along the narrow direction of the triangle.. The resulting stiffness matrix of the quadrilateral element will be an 8 x 8 matrix consisting of the stiffnesses of four triangles, as was shown in Figure 7-1 S. The stiff~ ness matrix is first assembled according to the usual direct stiffness method. Then we apply static condensatiori as outlined in Eqs. (7.5.1)-(7.5.7) to remove the internal node 5 degrees of freedom . .The stiffness matrix of a typical triangular element (labeled dement 1 in Figure 7-15) with nodes 1,2, and 5 is given in general form by [k(1)j k(J) -12 k(l) -15 k(l) -II = k(l) [ -21 k(l} -22 1 k(l) (7.5.8) -25 k~~) k;~ k;~) where the superscript in parentheses again refers to the element number, and each submatrix [kV)] is of order 2 x 2. The stiffness matrix of the quadrilateral, assembled using Eq. (7.5.8) along with similar stiffness matrices for elements 2-4 of Figure 7-15, is given by the following (before static condensation is used): (UI. VI) V4) (us, vs) (U2' V2) (U3) V3) [k~~] [0] [ki~l [ki~)] + rki~)J [k~;)J [0] [k~~)] + [k~~l rk~!)] [k~~)] + [k~!)] (U4) (kii)] + [ki~)l fk~;)] [k;~)] + [J4~J [OJ [k~;)l [k~~)] [k]= + [k~~)I [k~)] rkl~)l [OJ [k~~)J + [ki~}] + [ki~)] [k~)l [kg)l [kWl [k~!ll + + + + ([ki~)l + [k~~)]) + [k1i1j [kg) 1 [kg)] [k~)l ([k~~)] + [k;~)J) '] (7.5.9) 372: .A.., 7 Practical Considerations in Modeling; Interpreting Results . (a) (c) _ (b) Figure 7-16 Skew effects in finite element modeling where the orders of the degrees of freedom are shown above the columns of the stiffness matrix and the partitioning scheme used in static condensation is indicated by the dotted lines. Before- static condensation is applied. the stiffness matrix is of order 10 x 10. Example 7.1 Consider the quadrilateral with internal node 5 and dimensions as shown in Figure 7-17 to illustraJ:e the application of-static condensation. 3 4 ® CD I ~ 5 0'_ ....-24_ 20'_ looe-",..- ~ '7St.~ (e) Figure 7--19 Bicycle wrench (a) Outline drawing of wrench, (b) meshed model of wrench, (c) boundary conditions and selecting surface where surface traction will be applied, (d) checked model showing the boundary conditions. and surface traction, and (e) von Mises stress plot (compliments. of Angela Moe) 7.7 Computer Program Assisted Step-by-Step Solution, Other Models and Results .. 317 --~ " ", ........... \IlIIIJI< ·1 41l'J13e.Cl1 lbVC....2) (b) Figure 7-20 (a) Conneaing rod subjected to tensile loading and (b) resulting principal stress throughout the rod strength of the ASTM A-S14 steel is 690 MPa. Therefore, the wrench is safe from yielding. Additional trials can be made if the factor of safety is satisfied and if the maximum deflection appears to be satisfactory. Figure 7-20{a) shows a finite element model of a steel connecting rod that is fixed on its left edge and loading around· tlie right inner edge of the hole with a total force of 3000 lb. For more details, including the geometry of this rod, see Figure P7-11 at the end of this chapter. Figure 7-20{b) shows the resulting maximum 378 .... 7 Practical Considerations in Modeling; Interpreting Results -:;;....;;.:.;..;.......:-.;;...;.....:..;;.;.;;...-....'li..::;:......~..;...;;.::.;..--.;..=..O';'...:...::;...:....::.:..-~~:-;..:.;..;..(...~----.. ..;-'-_\...:-.---- -:~' ~, ~ Figure 7-21 von Mises stress. plot of overload protection device principal stress plot. The largest principal stress of 12051 psi occurs at the top and bottom inside edge of the hole. Figure 7-21 shows a finite element model along with the von Mises stress plot of an overload protection device (see Problem 7-30 for details of this problem). The upper member of the device was modeled. Node S at the shear pin location was constrained from vertical motion and a node at the roller E was constrained in the horizontal direction. An equilibrium load was applied at B along line BD. The magnitude of this load was calculated as one that just makes the shear stress reach 40 MPa in the pin at S. The largest von Mises stress of 178 MPa occurs at the inner edge of the cutout section. Figure 7-22 shows the shrink plot of a finite element analysis of a tapered plate with a hole in it, subjected to tensile loading along the right edge. The left edge was fixed. For details of this problem see Problem 7-23. The shrink. plot separates the elements for a clear look at the model. The lugest principal stress of 19.'0 e6' Pa (19.0 MPa) occurs at the edge of the hole, whereas the second largest principal stress of 17.95 e6 Pa (17.95 MPa) occurs at the elbow between the smallest cross section . and where the taper begins. Figure 7-23 shows the shrink fit plot of the maximum principal stresses in an overpass subjected to vertiCal loading on the top edge. The largest principal stress of 56162 Iblft2 (390 psi) occUrs at the top inside edge. For more details of this problem see Problem 7.20. 7.7 Computer Program Assisted Step-by-Step Solution, Other Models and Results ; .. 379 '~ 1.1121DtkoOO1 1~_ -,~ --- '.'H~ 1151_ 7111_ 5:I'OGn>4 ~ ~ .3.t121mM1• .oog l[i:liii!ii!iili~l~ =!~ Figure 7-22 Shrink fit plot of principal stresses in a tapered plate with hole Finally, Figure 7-24(a) shows a finite element discretized model of a steel spur gear for stress analysis. The auto meshing feature resulted in very small elements at the base of the tooth. The applied load of 164.8 Ib and the fixed nodes around the inner hole of the gear shown. Figure 7-24(b) shows an enlarged von Mises stress are Su.s:s von Mists Ib\l(ft"2) 56161.S6 50571.37 44lSO.SS 39390.4 33799.91 28209.42 22618.93 17028.44 11437.05 5847.4136 258.9776 Figure 7-23 Shrink fit plot of principal stresses in overpass (Compliments of David Walgrave) 380 ... 7 Practical Considerations in Modeling; Interpreting Results (a) (b) Figure 7-24 (a) Finite element model of a spur gear and (b) von Mises stress plot (Compliments of Bruce Figi) plot near the root of the tooth with the applied load acting on it. Notice that the largest stress of 4315 psi occurs at the left root of the tooth. The gear model has 27761 plane stress elements. References ... 381 ;1 References [I] Desai. C. S., and Abel, J. F., 'Introdu.ction to the Finite Element Method, Van Nostrand Reinhold, New York, J972. [2] Timoshenko, S., and Goodier) 1., Theory of Elasticity, 3rd ed., McGraw-Hili, New York, 1970. [3] Glockner, P. G., "Symmetry in Structural Mechanics," Journal of rhe Structural Division, American Society of Civil Engineers, Vol. 99, No. STI, pp. 71-89, 1973. [4J Yamada, Y., "Dynamic Analysis of Civil Engineering Structures," Reeent Advances in Malrix MethQdsofStructural Analysis and Design, R. H. Gallagher, Y. Yamada, and J. T. Oden, eds., University of Alabama Press, Tuscaloosa, AL, pp. 487-512, 1970. [5] Koswara, H., A Finite Element Analysis of Underground Shelter Subjecled to Ground Shock Load, M. S. Thesis, Rose·Hulman Institute of Technology, Terre Haute, IN, 1983. 16] Duruop, P., Duncan, J. M., and Seed, H. B., "Finite Element Analyses of Slopes in Soil," Journal of the Soil Meclumics and Fozmdations Division, Proceedings of the American Society of Civil Engineers, Vol. 96, No. SM2, March 1970. [7] Cook, R. D., Malkus, D. S., Plesha, M. E., and Witt, R. J., Concepts and Applications of Finite Element Analysis, 4th ed., Wiley, New York, 2002. [8) Taylor, R. L., Beresford, P. J., and Wilson, E. L., "A Nonconforming Element for Stress Analysis," InterM tionar Journal for Numerical Methods in Engineering, Vol. to, No.6, pp. 1211-1219, 1976. !9J Melosh, R. J., "Basis for Derivation of Matrices for the Direct Stiffness Method," Journal of the American Institute of Aeronautics and Astronautics, Vol. I, No.7, pp. 1631-1637. July 1963. PO] Fraeijes de Veubeke, B., "Upper and Lower Bounds in Matrix Structural Analysis," Matrix Methods of Structural Amziysis, AGARDograpb 72, B. Fraeijes de Veubeke, ed., MacmiI1an, New York, 1964. fIll Dunder, V., and Ridlon, S.,' "Practical Applications of Finite Element Method," Journal of the Structural Division, American Socic;;ty of Civil Engineers, No. STl, pp. 9-21, 1978. [IlJ Linear Stress and Dynamics Reference Division, Docutech On-line DOCUD1entation, Algor, Inc., Pittsburgh, PA t 5238. [13J Bettess, P., "More on Infinite Elements," Internatiomzl Journal for Numerical Method.,: in Engineering, Vol. 15, pp. 1613-1626~ 1980. [14] Gere,1. M., Mechanics of Materials, 5th ed., Brooks/Cole Publishers, Pacific Grove, CA, 2001. [IS} Superdraw Reference Division, Docutech On·1ine Documentation, Algor, Inc., Pittsburgh, PA ]5238. [16] Cook, R. D., and Young, W. c., Advanced Mechanics of Materials, Macmillan, New York,1985. (17] Cook, R. D., Finite Element Modeling/or Stress Analysis> Wiley, New York, 1995. [18] Kuro\h-ski, P., "Easily Made Errors Mar FEA Results/' Machine Design, Sept. 13,2001. [19} Huebner, K. H., Dewirst, D. L., Sr:n.ith, D. E., and Byrom, T. G., The Finite Element Methodfor Engineers, Wiley, New York, 2001. [20] Demkowicz, L., Devloo, P., and Oden, J. T., "On an h·Type Mesh-Refinement Strategy Based on Minimization of Interpolation Errors," Comput: Methods Appl. Meck Eng., Vol. 53, 1985, pp. 67-89. [21] Lohner, R., Morgan, K., and Zienkiewicz, O. C., "An Adaptive Finite Element Procedure for Compressible High Speed Flows," Comput. Methods Appl. Meeh. Eng., Vol. 51, 1985, pp. 441-465. 382 .& 7 Practical Considerations in Modeling; Interpreting Results [22} [23] [24] t2S} [26] [27} [281 [29J 4. Lohner~ R., "An Adaptive Finite Element Scheme for Transient Problems in CFD," Comput. Methods Appl. Meeh. End, Vol. 61, 1987, pp. 323-338. Ramakrisbnan, R., Bey, K. S., and Thornton, E. A., "Adaptive Quadrilateral and Triangular Finite Element Scheme for Compressible Flows," AIAA J.; Vol. 28, No.1, 1990, pp.51-59. Peano, A. G., "Hierarchies of Conforming Finite Elements for Plane Elasticity and Plate Bending," Comput. Match. Appl, Vol. 2, 1976, pp. 211-224. Szabo, B. A., ':Some Recent Developments in Finite Element Analysis," Comput. Match. Appl., Vol. 5, 1979, pp. 99-115. Peano, A. G.~ Pasin~ A., Ric:cioni., R. and Sardena, L., "Adaptive Approximation in Finite Element Structural Analysis," Comput. Struct., Vol. 10, 1979, pp. 332-342. Zienkiewicz, O. c., Gago, J. P. de S. R, and Kelly, D. W., "The Hierarchical Concept in Finite Element Analysis," Comput. StrucL, Vol. 16, No. 1-4, 1983, pp. 53-65. Szaoo, B., A., "Mesh Design for the p-Version of the Finite Element Method," Comput. Meth.ods Appl Meek Eng.) Vol. 55, 1986, pp. 181-197. Toogood, Roger, ProlMECHANICA, Structural Tutorial, SDC Publications, 2001. Problems 7.1 For the finite element mesh shown in Figure P7-I, comment on the goodness of the mesh. Indicate the mistakes in the model. Explain and show how to correct them. ~IIIII III Figure P7-1 /: C ............ ../" ~ . A ./ D B Figure P7-2 7.2 Comment on the mesh sizing in Figure P7-2. Is it reasonable? lfnot, explain why not. 7.3 What happens if the material property v = 0.5 in the plane strain caSe? Is this possible? Explain. j . 7.4 Under what conditions is the structure in Figure P7-4 a plane strain problem? Under what conditions is the structure a pJane stress problem? Problems !--IOin.-t % T 1000 Ib/in. lOin. / =--r ... 383 Figure P7-4 10 in. '~_~_----J -.L I---- in.---1 20 7.5 When do problems occur using the smoothing (averaging of stress at the nodes from elements connected to the node) method for obtaining stress results? 7.6 What thickness do you think is used in computer programs for plane strain problems? 7.7 Which one of the CST models shown below is expected to give the best results for a cantilever beam subjected to an end shear load? Why? ~] I_ I I I I 14@."=4bL 6 @ 2" = 12" ~ 1111111111 nJ m 12 @ 1" = 12" - I Ifo/It.. - - - - - ' " " I (a) (b) ~-I-------li :: !(c) 6" -I- 6" -I (d) Figure P7-7 7.8 Show that Eq. (7.5.13) is obtajned by static condensation of Eq. (7.5.12). Solve the following problems using a computer program.. In some of these problems, we suggest that students be assigned separate parts (or models) to facilitate parametric studies. ' • 7.9 Determine the free-end displacements and the element stresses for the plate discretized into four triangular elements and subjected to the tensile forces ~hown in Figure P7-9. Compare your results to the soIutio~ given in Section 6.5 Why are these results dif.;. ferent? Let E = 30 x 106 psi, V = 0.30, and I = 1 in. 384 • 7 Practical Considerations in Modeling; Interpreting Results C>