Linear elasticity equation

Finite Element Formulation for Vector Field Problems Linear Elasticity finite element analysis of structures through unified formulation
GregDeamons Profile Pic
GregDeamons,New Zealand,Professional
Published Date:03-08-2017
Your Website URL(Optional)
9 Finite Element Formulation for VectorFieldProblems – Linear Elasticity The discipline underlying linear stress analysis is the theory of elasticity. Both linear and nonlinear elasticity have been studied extensively over the past three centuries, beginning with Hooke, a contem- porary of Newton. Hooke formulated what has come to be known as Hooke’s law, the stress–strain relation for linear materials. Linear elasticity is used for most industrial stress analyses, as under operating conditions most products are not expected to undergo material or geometric nonlinearities. Linear elasticity also deals with many important phenomena relevant to materials science, such as the stress and strain fields around cracks and dislocations. These are not considered in this course. We start by presenting the basic assumptions and governing equations for linear elasticity in Section 9.1, followed by the exposition of strong and weak forms in Section 9.2. Finite element formulation for linear elasticity is then given in Section 9.3. Finite element solutions for linear elasticity problems in 2D concludes this chapter. 9.1 LINEAR ELASTICITY The theory of linear elasticity hinges on the following four assumptions: 1. deformations are small; 2. the behavior of the material is linear; 3. dynamic effects are neglected; 4. no gaps or overlaps occur during the deformation of the solid. In the following, we discuss each of these assumptions. The first assumption is also made in any strength of materials course that is taught at the undergraduate level. This assumption arises because in linear stress analysis, the second-order terms in the strain– displacement equations are neglected and the body is treated as if the shape did not change under the influence of the loads. The absence of change in shape is a more useful criterion for deciding as to when linear analysis is appropriate: when the application of the forces does not significantly change the configuration of the solid or structure, then linear stress analysis is applicable. For structures that are large enough so that their behavior can readily be observed by the naked eye, this assumption implies that A First Course in Finite Elements J. Fish and T. Belytschko 2007 John Wiley & Sons, Ltd ISBNs: 0 470 85275 5 (cased) 0 470 85276 3 (Pbk)216 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY the deformations of the solid should not be visible. For example, when a car passes over a bridge, the deformations of the bridge are invisible (at least we hope so). Similarly, wind loads on a high-rise building, although often felt by the occupants, result in invisible deformations. The deformations of an engine block due to the detonations in the cylinders are also invisible. On the other hand, the deformation of a blank in a punch press is readily visible, so this problem is not amenable to linear analysis. Other examples that require nonlinear analysis are a. the deformations of a car in a crash; b. the failure of an earth embankment; c. the deformations of skin during a massage. 2 As a rough rule, the deformations should be of the order of 10 of the dimensions of a body to apply linear stress analysis. As we will see later, this implies that the terms that are quadratic in the deformations are of 2 the order of 10 of the strains, and consequently, the errors due to the assumption of linearity are of the order of 1%. Many situations are just barely linear, and the analyst must exercise significant judgment as to whether a linear analysis should be trusted. For example, the deformations of a diving board under a diver are quite visible, yet a linear analysis often suffices. Sometimes these decisions are driven by practicality. For example, you have probably seen the large motions of the wingtip of a Boeing 747 on takeoff. Would a linear analysis be adequate? It turns out that the design of the aircraft is still primarily analyzed by linear methods, because the errors due to the assumption of linearity are small and thousands of loadings need to be considered, and this becomes much more complex with nonlinear analysis. The linearity of material behavior is also a matter of judgment. Many metals exhibit a relationship betweenstressandstrainthatdeviatesfrom linearitybyonlyafewpercentuntiltheonset ofplasticyielding. Until the yield point, a linear stress–strain law very accurately reproduces the behavior of the material. Beyond the yield point, a linear analysis is useless. On the contrary, materials such as concrete and soils are often nonlinear even for small strains, but their behavior can be fit by an average linear stress–strain law. The assumption of static behavior corresponds to assuming that the accelerations sustained during the loading are small. This statement by itself gives no meaningful criterion, as one can immediately ask, ‘small compared to what?’ There are several ways to answer this question. One way is to consider the d’Alem d’Alembert forcef due to the acceleration, which is given by d’Alem jf j¼jMaj; whereM is the mass of the body anda is the acceleration; we have put absolute value signs on both sides of the equation because we are only interested in magnitudes. If the d’Alembert forces are small compared to the loads, then dynamic effects are also small. The dynamic effect can be viewed as the overshoot that you will see on a floor scale if you jump on it compared to stepping on it slowly. An easier way to judge the appropriateness of a static analysis, i.e. neglecting dynamic effects, is to comparethetimeofloadapplicationtothelowestperiodofthesolidorstructure.Thelowestperiodisthetime for a structure to complete one cycle of vibration whenvibrating freely. Ifthe time inwhich the load isapplied is large compared to the period associated with the lowest frequency, then static analysis is applicable. The fourth assumption states that as the solid deforms, it does not crack or undergo any interpenetration of material; in short, no gaps or overlaps develop in the body. Interpenetration of material is generally not possible, unless one material is liquefied or vaporized, so this part of the assumption is just common sense. The first part of the assumption states that the material does not crack or fail in some other way. Obviously, materials do fail, but linear stress analysis is then not appropriate; instead, special nonlinear finite element methods that account for cracking must be used. The last assumption can be interpreted in terms of continuity. It states that the displacement field is smooth. The order of smoothness that is required is something we have already learned and is associatedLINEAR ELASTICITY 217 with the requirements of the integrability of theweak form, but physically it can be justifiedby requiring the deformation to be such that there are no gaps or overlaps. The requirements of a linear stress analysis solution are closely related to the assumptions. The requirements are a. the body must be in equilibrium; b. it must satisfy the stress–strain law; c. the deformation must be smooth. In addition to the above, in order to write a stress–strain law, we need a measure of the strain that expresses the strain in terms of the deformation, which is called the strain–displacement equation. Equilibrium requires that the sum of the forces at any point of the solid must vanish. The other two requirements have already been discussed. 9.1.1 Kinematics The displacement vector in two dimensions is a vector with two components. We will use a Cartesian coordinate system, so the components ofthe displacement are thex-componentand they-component. It can be written in matrix and vector forms as  u x u¼ ; u¼ u iþu j; ð9:1Þ x y u y where the subscript indicates the component. Figures 9.1(a) and (b), depict the deformation of a control volumexy in thex andy directions respectively. The combined deformation is given in Figure 9.1(c). Under the assumption of small displacement gradients, we can use three independent variables to describe the deformation of a control volume. These variables correspond to the strains. The extensional strains aree ande ; sometimes the repeated subscripts are dropped and the extensional xx yy strains are written as e and e . The expressions for these strains can be derived exactly like the one- x y dimensional extensional strain. The extensional strains e and e are the changes in the lengths of the x y infinitesimal line segments in thex andy directions,x andy, respectively, divided by the original lengths of the line segments. Based on this definition, we obtain the following relations for the extensional strains: u ðxþx; yÞu ðx; yÞ u x x x e ¼ lim ¼ ; xx x0 x x ð9:2Þ u ðx; yþyÞu ðx; yÞ u y y y e ¼ lim ¼ ; yy y0 y y The shear strain,g , measures the change in angle between the unit vectors in thex andy directions in units xy of radians: u ðx; yþyÞu ðx; yÞ u ðxþx; yÞu ðx; yÞ x x y y g ¼ lim þ lim xy y0 y x0 x ð9:3Þ u u y x ¼ þ ¼ a þ a : 1 2 x y wherea are showninFigure9.1.Twoforms ofthe shearstrainappear commonlyinfiniteelement software: i the engineering shear strain g given above and the tensor shear strain e ¼ð1=2Þg . xy xy xy218 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY y y α 1 u (x + ∆ x,y + ∆ y) u (x + ∆ x,y + ∆ y) x y u (x,y + ∆ y) y α 2 α 2 ∆ y ∆ y α 1 u (x,y) y ∆ x x x ∆ x u (x,y) u (x + ∆ x,y) x x (a) (b) y u(x,y + ∆ y) u(x + ∆ x,y + ∆ y) α 2 u(x + ∆ x,y) α 1 u(x ,y) x (c) Figure 9.1 Deformation of a control volume: (a) deformation in x due to ru ; (b) deformation in y due to ru ; x y (c) deformation inx andy. Note that ifa ¼a , the shear strain vanishes. The resulting deformation is depicted in Figure 9.2. It 1 2 can be seen that the control volume undergoes axial elongations in addition to the rotation. The rotation of the control volume in two dimensions, denoted by , is computed by xy  1 u u 1 x y ¼  ¼ ða  a Þ: ð9:4Þ xy 2 1 2 y x 2 For an infinitesimal displacement field, uðx;yÞ, the rotation is very small, and therefore it does not affect xy the stress field. y u(x,y + ∆ y) u(x + ∆ x,y + ∆ y) α 2 u(x,y) α 1 x u(x + ∆ x,y) Figure 9.2 Axial strains and rotation of a control volume.LINEAR ELASTICITY 219 In finite element methods, the strains are usually arranged in a column matrixe, as shown below: T e¼½e e g  : ð9:5Þ xx yy xy Equations (9.2)–(9.3) can be written in terms of the displacements as a single matrix equation: 2 3  e xx u x 4 5 e e¼ ¼= u ¼ = ; ð9:6Þ yy S S u y g xy where= is a symmetric gradient matrix operator S 2 3 =x 0 4 5 = ¼ 0 =y : ð9:7Þ S =y =x 9.1.2 Stress and Traction Stresses in two dimensions correspond to the forces per unit area acting on the planes normal to thex ory axes (these are called tractions). The traction on the plane with the normal vectorn aligned along thex-axis is denoted by  and its vector form is  ¼  iþ j. Likewise, the traction with the outer normal unit x x xx xy vector n aligned along they-axis is denoted by  and its corresponding components are  ¼  iþ j. y y yx yy Wewill refer to  and  asstressvectors acting on the planes normal to thex andy directions, respectively. x y The stress state ina two-dimensional body isdescribed by two normalstresses and and shear stresses xx yy  and  as illustrated in Figure 9.3. From moment equilibrium in a unit square, it can be shown that xy yx  ¼  , so these stresses are identical. xy yx Figure 9.3 depicts stress components acting on two planes, the normals pointing in the positivex andy directions. Positive stress components act in the positive direction on a positive face. The first subscript on the stress corresponds to the direction of the normal to the plane; the second subscript denotes the direction of the force. The normal stresses are often written with a single subscript as and . x y Stresses can be arranged in a matrix form similarly to strains: T r ¼½   : ð9:8Þ xx yy xy Occasionally, it is convenient to arrange stress components in a 2 2 symmetric matrixs as    xx xy s¼ : ð9:9Þ   xy yy σ yy σ y y σ yx σ x −σ σ x xx σ xy −σ y x Figure 9.3 Stress components.220 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY y t n j y n dy dΓ n i x −σ x dx x −σ y Figure 9.4 Relationship between stress and traction. The stress vectors  and  can be conveniently used to obtain the tractions on any surface of the body. x y Tractions, like stresses, are forces per unit area, but they are associated with a specific surface, whereas the stresses provideinformation abouttractionsonanysurfaceat apoint. Therelationship betweenstresses and tractions is written in terms of the unit normal to the surface n, as illustrated in Figure 9.4. Consider the triangular body shown in Figure 9.4. The thickness of the triangle is taken to be unity. At the surfacewiththe unit normal vector n,the tractionvectorist.On the planes normaltothe coordinate axes,the traction vectors are  and  . The components of the unit normal vector are n given by x y n¼ n iþn j: x y The force equilibrium of the triangular body shown in Figure 9.4 requires that t d  dy  dx¼ 0: x y Dividing the above equation by d and noting that dy¼ n d and dx¼ n d, we obtain x y t  n   n ¼ 0: x x y y Multiplying the above by unit vectorsi andj yields, respectively, t ¼  n þ n ¼   n; x xx x xy y x ð9:10Þ t ¼  n þ n ¼   n; y xy x yy y y where we have used the relationst ¼ ti,t ¼ tj, ¼  i, ¼  i, ¼  j and ¼  j. x y xx x xy y xy x yy y Equation (9.10) can be written in the matrix form as t¼ sn: ð9:11Þ 9.1.3 Equilibrium Consider an arbitrarily shaped body shown in Figure 9.5 of unit thickness; the body force and the surface traction are assumed to be acting in thexy-plane.LINEAR ELASTICITY 221 ∆ y (, xy + ) σ t y 2 n b Ox (,y) ∆ y x ∆ x ∆ y (, xy + ) −σ(, xy − ) σ x x 2 ∆ x Γ 2 ∆ x Ω y (, xy ) −σ − y 2 (a) (b) Figure 9.5 Problem definition: (a) domain of the unit-thickness plate and (b) traction vectors acting on the infinite- simal element. The forces acting on the body are the traction vectort along the boundary and the body forceb per unit volume. The body force and the tractionvectorsarewritten asb¼ b iþb j and t¼ t iþt j, respectively. x y x y Examples of the body forces are gravity and magnetic forces. Thermal stresses also manifest themselves as body forces. Next, consider the equilibrium of the infinitesimal domain of unit thickness depicted in Figure 9.5(b). For a static problem (no dynamic effects), the equilibrium equation on the infinitesimal domain is given by   x x   x ;y yþ  xþ ;y y x x 2 2   y y   x;y xþ  x;yþ xþbðx;yÞxy¼ 0: y y 2 2 Dividing the above byxy, taking the limit asx 0,y 0 and recalling the definition of partial derivatives,  x x  xþ ;y   x ;y x x  2 2 x lim ¼ ; x0 x x  y y  x;yþ   x;y y y 2 2  y lim ¼ : y0 y y Combining the above two equations yields the equilibrium equation:   x y þ þb¼ 0: ð9:12Þ x y Multiplying (9.12) by unit vectorsi andj gives two equilibrium equations:   xx xy þ þb ¼ 0; x x y   yx yy þ þb ¼ 0; ð9:13Þ y x y or in the vector form: r  þb ¼ 0; r  þb ¼ 0: ð9:14Þ x x y y222 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY The equilibrium equations will also be considered in the matrix form. If you consider the transpose of the symmetric gradient operator given in (9.7) and the column matrix form of the stress: 2 3 2 3 0  xx 6 7 T x y 6 7 4 5 = ¼ ; r¼  ; yy S 4 5 0  xy y x then the matrix form of equilibrium equations (9.13) can be written as T = rþ b¼ 0: ð9:15Þ S The fact that the equilibrium equation (9.15) is the transpose of the strain–displacement equation (9.6) is an interesting feature that characterizes what are called self-adjoint (or symmetric) systems of partial differential equations. The heat conduction (or diffusion) equations are similarly self-adjoint. The self- adjointness of these partial differential equations is the underlying reason for the symmetry of the discrete equations, i.e. the stiffness matrix and the conductance matrix. 9.1.4 Constitutive Equation Now let us consider the relation between stresses and strains, which is called the constitutive equation. Examples of constitutive equations are elasticity, plasticity, viscoelasticity, viscoplasticity and creep. Here, we focus on the simplest constitutive theory, linear elasticity. Recall that in one dimension, a linear elastic material is governed by Hooke’s law¼ Ee, where the material constantE is Young’s modulus. In two dimensions, the most general linear relation between the stress and strain matrices can be written as r¼ De; ð9:16Þ where D is a 3 3 matrix. This expression is called the generalized Hooke’s law. It is always a symmetric, positive-definite matrix; these two properties are due to energy considerations, which will not be discussed here but can be found in any text on continuum mechanics or elasticity. In two-dimensional problems, the matrix D depends on whether one assumes a plane stress or a plane strain condition. These assumptions determine how the model is simplified from a three-dimensional physical body to a two-dimensional model. A plane strain model assumes that the body is thick relative to thexy-plane in which the model is constructed. Consequently, the strain normal to the plane,e , is zero and z the shear strains that involve angles normal to the plane,g andg , are assumed to vanish. A plane stress xz yz model is appropriate when the object is thin relative to the dimensions in the xy-plane. In that case, we assume that no loads are applied on thez-faces of the body and that the stress normal to thexy-plane, ,is zz assumed to vanish. The physical arguments for these assumptions are as follows. If a body is thin, as the stress must vanish on the outside surfaces, there is no mechanism for developing a significant nonzero zz stress  . On the other hand, when a body is thick, significant stresses can develop on the z-faces, in zz particular the normal stress can be quite large. zz The D matrix depends on the symmetry properties of the material. An isotropic material is a material whose stress–strain law is independent of the coordinate system, which means that regardless of the orientation of the coordinate system, the elasticity matrix is the same. Many materials, such as most steels, aluminums, soil and concrete, are modeled as isotropic, even though manufacturing processes, such as sheet metal forming, may induce some anisotropy.STRONG AND WEAK FORMS 223 For an isotropic material, the D matrix is given by Planestress: 2 3 1  0 E 4 5 D¼ 10 : 2 1 00 ð1Þ=2 Planestrain: 2 3 1 0 E 4 5 D¼  1 0 : ð1þÞð1 2Þ 00 ð1 2Þ=2 As can be seen from the above, for an isotropic material, the Hookean matrix D has two independent material constants: Young’s modulusE and Poisson’s ratio. Note that for plane strain, as 0:5, the Hookean matrix becomes infinite. A Poisson’s ratio of 0.5 corresponds to an incompressible material. This behavior of the Hookean matrix as the material tends toward incompressibility and other features of finite elements make the analysis for incompressible and nearly incompressible materials more difficult than for compressible materials. Therefore, special elements must be used for incompressible materials. These difficulties do not occur for plane stress problems, but they do occur in three dimensions. The Hookean matrix for an isotropic material can also be written in terms of alternative material constants, such as the bulk modulusK ¼ E=3ð1Þ and the shear modulusG¼ E=2ð1þÞ. In some circumstances, a two-dimensional model is appropriate but the standard plane stress or plane strain assumptions are not appropriate because although the z components of the stress or strain are constant, theyare nonzero. This is called a state of generalized plane stress orgeneralized plane strain when  ore are constant, respectively. zz zz 9.2 STRONG AND WEAK FORMS Let us summarize the relations established so far for 2D linear elasticity. Equilibriumequation: T = rþ b¼ 0; or r  þb ¼ 0 and r  þb ¼ 0: ð9:17Þ x x y y S Kinematicsequation (strain–displacementrelation): e¼= u: S Constitutiveequation (stress–strainrelation): r¼ De: As in one dimension, we consider two types of boundary conditions: The portion of the boundary where the traction is prescribed is denoted by , and the portion of the boundary where the displacement is prescribed t is denoted by . The traction boundary condition is written as u    sn¼ t on  ; or   n¼ t and   n¼ t on  : ð9:18Þ t x x y y t224 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY The displacement boundary condition is written as u¼ u on  ; or u¼ u on  : ð9:19Þ u u The displacement boundary condition is an essential boundary condition, i.e. it must be satisfied by the displacement field. The traction boundary condition is a natural boundary condition. As before, the displacement and traction both cannot be prescribed on any portion of the boundary, so  \ ¼ 0: u t However, on any portion of the boundary, either the displacement or the traction must be prescribed, so   ¼ : u t We summarize the strong form for the linear elasticity problem in 2D in Box 9.1 in the mixed vector–matrix notation, relevant for derivation of the weak form. Box 9.1. Strong form for linear elasticity ðaÞ r  þb ¼ 0 and r  þb ¼0on ; x x y y ðbÞ r¼ D= u; S ð9:20Þ   ðcÞ   n¼ t and   n¼ t on  ; x x y y t ðdÞ u¼ u on  : u To obtain the weak form, we first define the admissible weight functions and trial solutions as in Section 3.5.2. We then premultiply the equilibrium equations in thex andy directions (9.20a) and the two natural boundary conditions (9.20c) by the corresponding weight functions and integrate over the corresponding domains, which gives Z Z ðaÞ w r  dþ w b d¼ 0 8w 2 U ; x x x x x 0   Z Z ðbÞ w r  dþ w b d¼ 0 8w 2 U ; y y y y y 0   Z ð9:21Þ  ðcÞ w ðt    nÞ d¼ 0 8w 2 U ; x x x x 0  t Z  ðdÞ w ðt    nÞ d¼ 0 8w 2 U ; y y y y 0  t where  w x w¼ ; w¼ w iþw j: x y w yFINITE ELEMENT DISCRETIZATION 225 Green’s theorem is applied (see Chapter 6) to the first term in equations (9.21a) and (9.21b), which yields Z I Z w r  d¼ w   n d rw   d; x x x x x x    Z I Z ð9:22Þ w r  d¼ w   n d rw   d: y y y y y y    Adding the two equations in (9.22) and recalling that the weight functionsw andw vanish on yields x y u Z I Z ð rw   þ rw   Þ d¼ ðw   nþw   nÞ dþ ðw b þw b Þ d: ð9:23Þ x x y y x x y y x x y y    t Substituting (9.21c) and (9.21d) into (9.23) and writing the RHS in (9.23) in the vector form gives Z I Z ð rw   þ rw   Þ d¼ wt dþ wb d: ð9:24Þ x x y y    t Expanding the integrand on the LHS of (9.24) yields w w w w x x y y rw   þ rw   ¼  þ  þ  þ  x x y y xx xy xy yy x y x y 2 3    xx ð9:25Þ w w w w x y x y 6 7 T ¼ þ 4 5¼ð= wÞ r: yy S x y y x  xy Inserting (9.25) into (9.24) and writing the RHS of (9.24) in the matrix form gives Z Z Z T T T  ð= wÞ r d¼ w t dþ w b d 8w2 U : S 0  t  After the substitution of (9.20b) forr the weak form in two dimensions can be written as follows: Find u2 U such that Z Z Z T T T  ð= wÞ D= u d¼ w t dþ w b d 8w2 U ; S S 0    ð9:26Þ t 1 1 where U¼fuju2 H ; u¼ u on  g; U ¼fwjw2 H ;w¼0on  g: u 0 u 9.3 FINITE ELEMENT DISCRETIZATION Consider a problem domain with boundary discretized with two-dimensional elements (triangles or quadrilaterals) as shown in Figure 9.6.; the total number of elements is denoted byn . el T Thex andy components of the displacement field u¼½u u  are generally approximated by the same x y shape functions, although in principle different shape functions could be used for each of the components.226 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY y τ n = t on Γ t Ω on Γ x u = u u Figure 9.6 Finite element mesh in two dimensions. There are two degrees-of-freedom per node corresponding to the two components of the global displacements, so the nodal displacement matrix is: h T d¼ u u u u ... u u  ynnp x1 y1 x2 y2 xn np wheren is the number of nodes in the finite element mesh. The displacement field in the finite element is np written in terms of the shape functions, which from Chapter 7 we know depend on the type of element and the number of nodes. The finite element approximation of the trial solution and weight function on each element can be expressed by: e e e e uðx;yÞ u ðx;yÞ¼ NðÞ x;y d ðx;yÞ2  ð9:27Þ T eT eT e T e wðÞ x;ywðÞ x;y¼ w NðÞ x; y ðx; yÞ2 e where element shape function matrix N in Eq. (9.27) is given as  e e e N 0 N 0 ... N 0 e 1 2 n en N ¼ e e e 0 N 0 N ... 0 N 1 2 n en  T e e e e e e e u u u u ... u u and d ¼ are the element nodal displacements and x1 y1 x2 y2 xn yn en en  T e e e e e e e w w w w ... w w w ¼  are the element nodal values of weight functions. x1 y1 x2 y2 xnen ynen 0 Recall from Chapter 6 that the finite element approximation is C continuous, i.e. it is smooth over element domains but have kinks at the element boundaries. Therefore, the integral over in the weak form e (9.26) is computed as a sum of integrals over element domains () Z Z Z nel X eT e e eT eT  = w D = u d w td w bd ¼ 0 ð9:28Þ S S e e e    e¼1 t Next we express the strains in terms of the element shape functions and the nodal displacements. Recall the strain-displacement equations (9.6) expressed in terms of the symmetric gradient operator. Applying the e symmetric gradient operator to N gives 2 3 e xx e e e e e e 4 5 e e¼  e ¼= u ¼= N d ¼ B d ; ð9:29Þ yy S S g xyFINITE ELEMENT DISCRETIZATION 227 e where the strain–displacement matrix B is defined as 2 3 e e e N N N n 1 2 en 0 0  0 6 7 x x x 6 7 6 7 e e e 6 N 7 N N n e e 1 2 en 6 7 B = N ¼ 0 0  0 : S 6 7 y y y 6 7 6 7 e e e e e e 4 5 N N N N N N n n 1 1 2 2 en en  y x y x y x The derivatives of weight functions are: T T e e e eT eT ð9:30Þ ð= w Þ ¼ðB w Þ ¼ w B : S e e eT eT T Substituting (9.30), (9.29) and (9.27) into (9.28) and recalling that d ¼ L d, w ¼ w L yields 8 9 2 3 Z Z Z n = el X 6 7 T eT eT e e e eT eT  w L B D B d L d N td N b d ¼ 0 8w : ð9:31Þ 4 5 F : ; e¼1 e e e    t In the above, we have replaced the arbitrary weight functionswðx;yÞ by arbitrary parameters w . w is the F F portion of wcorresponding to nodesthat are noton an essential boundary.Following the derivation outlined in Chapters 5 and 8, the element matrices are given as follows: Elementstiffnessmatrix: Z e eT e e K ¼ B D B d: ð9:32Þ e  Elementexternalforcematrix: Z Z e eT eT  f ¼ N b dþ N t d; ð9:33Þ e e   t fflfflfflfflfflfflfflfflzfflfflfflfflfflfflfflffl fflfflfflfflfflfflfflzfflfflfflfflfflfflffl e e f f   e e where f and f in (9.33) are the body and boundary force matrices.   The weak form can then be written as 20 1 0 13 6B C B C7 n n el el X X 6B C B C7 T eT e e eT e 6B C B C7 w L K L d L f ¼ 0 8w : ð9:34Þ F 6B C B C7 4 A A5 e¼1 e¼1 fflfflfflfflfflfflfflfflfflzfflfflfflfflfflfflfflfflffl fflfflfflfflfflzfflfflfflfflffl K f228 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY Using (9.32) and (9.33) and the assembly operations (5.13) and (5.14), the system (9.34) reduces to T w ðKd fÞ¼ 0 8w : ð9:35Þ F Equation (9.35) can be written as T w r¼ 0 8w : ð9:36Þ F Partitioning Equation (9.36) into E- and F-nodes gives T T w r þ w r ¼ 0 8w : F E F F E Asw ¼ 0andw isarbitrary,it followsthatr ¼ 0.Consequently,theaboveequationcanbeconveniently F E F rewritten as   K K f þ r E EF d E E E ¼ ; ð9:37Þ T K K d f F EF F F where K , K and K are partitioned to be congruent with the partition of d and f. Equation (9.37) is E F EF solved using the two-step partition approach discussed in Chapter 5. 9.4 THREE-NODE TRIANGULAR ELEMENT The triangular three-node element is illustrated in Figure 9.7. It is a linear displacement element. The strains are constant in the element. The nodes must be numbered counterclockwise as shown in the figure. e Each node has two degrees of freedom, so the column matrix d consists of six terms: e e e e e e e T d ¼½u ;u ;u ;u ;u ;u  : ð9:38Þ x1 y1 x2 y2 x3 y3 The displacement field in the element can then be expressed in the form of  e e e e u N 0 N 0 N 0 x 1 2 3 e ¼ d : e e e u 0 N 0 N 0 N y 1 2 3 e u y1 1 e u x 1 e u y 2 e 2 u e y 3 u x 2 e u x 3 3 Figure 9.7 A single triangular finite element.THREE-NODE TRIANGULAR ELEMENT 229 Applying the symmetric gradient operator (9.6) gives 2 3 2 3 e e e e e N 0 N 0 N 0 xx 1;x 2;x 3;x e e e e 4 5 4 5 0 N 0 N 0 N e ¼ d ; ð9:39Þ yy 1;y 2;y 3;y e e e e e e g N N N N N N xy 1;y 1;x 2;y 2;x 3;y 3;x e e N N e I e I whereN ¼ andN ¼ . Using the relations given in Chapter 7, it follows that I;x x I;y y 2 3 2 3 e e e e e y 0 y 0 y 0 xx 23 31 12 1 e e e e e 4 5 4 5 e 0 x 0 x 0 x e ¼ yy ¼ d ; ð9:40Þ 32 13 21 e 2A e e e e e e g x y x y x y xy 32 23 13 31 21 12 e e e e e wherex ¼ x x , which defines the B matrix for the element. It can be seen that as expected, the B IJ I J matrix is not a function ofx ory, i.e. the strain is constant in the element. The stiffness matrix is given by (9.32): Z e eT e e K ¼ B D B d: e  In most cases, for a low-order element such as this, the material properties are assumed constant in the element. Consequently, the integrand is a constant, and for an element of unit thickness, we have e e eT e e K ¼ A B D B : The stiffness matrix is 6 6, and it is quite large for manual computations, so it is usually evaluated by computer. 9.4.1 Element Body Force Matrix The element body force matrix is given by (9.21): Z e eT f ¼ N b d: ð9:41Þ   There are two ways of evaluating this matrix: (i) by direct numerical integration, and (ii) by interpolating b, usually with a linear function, and integrating the result in the closed form. Note that in direct integration, interpolation is still often required as the body forces may only be given at discrete points and interpolation is required to evaluate the integral. Evaluation of the matrix in the closed form is extremely difficult unless triangular coordinates are used, so we will use them here. We interpolate the body force in the element by the linear shape functions in the triangular coordinates as   3 X b b x 3T xI b¼ ¼ N ; ð9:42Þ I b b y yI I¼1230 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY whereb andb are thex andy components of the body force at nodeI. Substituting (9.42) into (9.41), we xI yI obtain 2 3 2 3 3T N 0 2b þb þb x1 x2 x3 1 3T 6 7 0 N 6 2b þb þb 7 1 y1 y2 y3 6 7  Z 3 e6 7 3T X 6 7 b A N 0 6 b þ 2b þb 7 e 3T xI x1 x2 x3 2 6 7 f ¼ N d¼ 6 7: ð9:43Þ  3T I 6 7 b þ 2b þb 0 N b y1 y2 y3 yI 126 7 6 2 7 I¼1 4 5 e  3T 4 5 b þb þ 2b N 0 x1 x2 x3 3 3T b þb þ 2b 0 N y1 y2 y3 3 The last step was performed by using the integration formulas as given in Section 7.8.2. 9.4.2 Boundary Force Matrix The boundary force matrix is given by Z e eT  f ¼ N t d: ð9:44Þ  e  t Asfor the bodyforces, they can be evaluated by direct integrationor by interpolation.Weillustrate the latter approach for a linearly interpolated traction. To simplify the explanation, consider the triangular element shown in Figure 9.8; in the figure, the traction is applied to the edge joining nodes 1 and 2, but the results are easily applied to any node numbers. e We know from the Kronecker delta property of shape functions thatN vanishes at nodes 1 and 2, and as the 3 2L 2L shape function is linear along the edge, it vanishes along the entire edge. Furthermore,N andN are 1 2 linear along the edge and can be written in terms of the edge parameter as 2L 2L N ¼ 1; N ¼ : 1 2 The integral (9.44) then becomes 2 3 1 0 6 7 01 1 6 7 Z 6 7  0 t ð1Þþt  e x1 x2 6 7 f ¼ l d;  6 7 0  t ð1Þþt  y1 y2 6 7 0 4 5 00 00 e e y u , f y2 y2 2 e e u , f x2 x2 e e u , f x1 y1 element e e e u , f 1 x1 x1 e e u , f 3 y3 y3 e e u , f x3 x3 x Figure 9.8 Triangular three-node element showing nodal displacements and nodal forces (they are shown as collinear that usually are not).GENERALIZATION OF BOUNDARY CONDITIONS 231 where we have used d¼ l d and changed the limits of integration to 0 to 1 (l is the length of the edge). Note that we have used a linear interpolation of the two traction components. The above is easily integrated in the closed form, giving 2 3 2t þt x1 x2 6 7 2t þt y1 y2 6 7 6 7 l t þ 2t e x1 x2 6 7 f ¼ :  6 7 t þ 2t 6 y1 y2 6 7 4 5 0 0 Thus, there are no nodal forces on node 3 due to the tractions on the edge connecting nodes 1 and 2. The nodal force at node 1 (or 2) is more heavily weighted by the traction at node 1 (or 2). For a constant traction,   t ¼ t ¼ t andt ¼ t ¼ t , we obtain x1 x2 x y1 y2 y 2 3  t x 6 7  t y 6 7 6 7 l  t x e 6 7 f ¼ ;  6 7  t 2 y 6 7 4 5 0 0 which shows that the total forces (the thickness is unity) are split equally among the two nodes. 9.5 GENERALIZATION OF BOUNDARY CONDITIONS Although we have subdivided the boundary into prescribed displacement and prescribed traction bound- aries, in fact, one has substantially more versatility in stress analysis: on any portion of the outside surface, anycomponent of the traction or the displacement can be prescribed. To specify this mathematically, we denote the portion of the surface on which theith component of the traction is prescribed by (thei¼ 1 t i component is the x-component, the i¼ 2 component is they-component). Similarly, the portion of the boundary on which theith component of the displacement is prescribed is denoted by . The boundary u i conditions are then written as   n¼ t on  ; x x tx    n¼ t on  ; y y ty u ¼ u on  ; x x ux u ¼ u on  : y y uy This weak form can be derived by an appropriate choice ofw andw on the boundary. Note that the same x y component of traction and displacement cannot be prescribed on any part of the boundary, so  \ ¼ 0;  \ ¼ 0: ux tx uy ty Furthermore, for each component, either the traction or the displacement can be prescribed, so   ¼ ;   ¼ : ux tx uy ty232 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY Observe that these boundary conditions conform to the rule that any two variables that are conjugate in work cannot be prescribed. Thus,u andt are conjugate in work in the sense that an increment of work is x x given by dW ¼ t du , whereast andu are not conjugate in work, so they can be prescribed on any portion x x x y of the boundary. Example 9.1: Illustration of boundary conditions In the following, we describe how to specify boundary conditions for various problems. We start with some simple idealized problems and then proceed to situations that are more realistic. For the latter, choosing appropriate boundary conditions is often an art. Consider the plate with a hole shown in Figure 9.9 with loads applied at the top and bottom. Sides AD and BC are traction free, and nothing needs to be done in a finite element model to enforce a homogeneous (zero) natural boundary condition. Sides CD and AB are also natural boundaries, but the tractions must be incorporated in the equations through the boundary force matrix f . However, these boundary conditions do not suffice to render the system solvable, as  these boundary conditions admit rigid body motion, so there are an infinite number of solutions and K is singular. To eliminate rigid body motion, at least three nodal displacement components must be specified so that translation and rotation of the body is prevented (corresponding to translations in the x and y directions and rotation about the z-axis). One way to make K regular (nonsingular) is to let u ¼ u ¼ u ¼ 0: xA yA yB Note that if you replaceu ¼ 0byu ¼ 0, K is still singular as rotation is not prevented. The above yB xB conditions prevent both rigid body translation and rotation. Another way to model this problem is to use symmetry, resulting in the model shown in Figure 9.9(b). The lines of symmetry are FG and HK. Along a line of symmetry, the displacement component normal to the line (or plane) of symmetry must vanish. Otherwise, as the displacement fields in symmetric subdomains, i.e. and in Figure 9.9(c), are mirror images, so a A B y y F D C Line of symmetry G x H K (b) Line of x symmetry F B Ω Overlap A Gap B Ω A G (c) (a) Figure 9.9 Plate with a hole: (a) a model of complete problem; (b) a model of symmetric portion; (c) an illustration of why displacements normal to a line of symmetry must vanish.GENERALIZATION OF BOUNDARY CONDITIONS 233 y A H E D G F x B C (b) (a) Figure 9.10 A bracket and its model. nonzero normal displacement along the line (or plane) of symmetry results in either gaps or overlaps, which violates compatibility. The other symmetry condition is that the shear on the line of symmetry must vanish. To summarize, for Figure 9.9(b), u ¼ 0 and t ¼  ¼ 0on FG; x y xy u ¼ 0 and t ¼  ¼ 0on HK: y x xy As the above traction (natural) boundary conditions are homogeneous, they are naturally satisfied if we do not constrain the corresponding displacement. Figure 9.10 shows a bracket and a simplified model, which is aimed at finding the maximum stress in the bracket. In many cases, it would be desirable to model the bolt and vertical rod, but this would entail substantially more computational effort and the use of contact interfaces, which arenonlinear. Therefore we model them with prescribed displacements and applied loads. The boundary conditions are as follows: 1. along AB,u ¼ 0 and at one nodeu ¼ 0; x y 2. the remaining surfaces are all traction free, i.e.t ¼ t ¼ 0, except on the segment FG. x y Note that the frictional force along AB is not modeled; friction isnonlinear and the effect of the frictional forces would be small. Fillets are also not modeled. Example 9.2: Quadrilateral element Consider a linear elasticity problem on the trapezoidal panel domain as shown in Figure 9.11. The  vertical left edge is fixed. The bottom and the right vertical edges are traction free, i.e. t¼ 0. 1  Traction t ¼20 N m is applied on the top horizontal edge. Material properties are Young’s y 7 modulus E¼ 3 10 Pa and Poisson’s ratio  ¼ 0:3. Plane stress conditions are considered. The problem is discretized using one quadrilateral element. The finite element mesh and nodal coordi- nates in meters are shown in Figure 9.12. The constitutive matrix D is 2 3 2 3 1  0 10:30 E 6 7 7 10 4 5 D¼ 4 5¼ 3:3 10 0:31 0 : 2 1 1 00 00 0:35 2234 FINITE ELEMENT FORMULATION FOR VECTOR FIELD PROBLEMS – LINEAR ELASTICITY t = −20 y 2 m 0.5 m 1 m t on Γ t Figure 9.11 Problem definition of Example 9.2. Figure 9.12 Finite element mesh for Example 9.2. The coordinate matrix is 2 3 2 3 e e x y 01 1 1 e e 6 7 6 7 x y 00 e e 2 2 6 7 6 7 ½ x y¼ ¼ : e e 4 5 4 5 x y 20:5 3 3 e e x y 21 4 4 The Lagrangian shape functions in the parent element are   1 4Q 2 4 N ð;Þ¼ ¼ ð1Þð1Þ; 1     4 1 2 1 4   1 1 4 4Q N ð;Þ¼ ¼ ð1þÞð1Þ; 2     4 2 1 1 4   1 4Q 1 1 N ð;Þ¼ ¼ ð1þÞð1þÞ; 3     4 2 1 4 1   1 4Q 2 1 N ð;Þ¼ ¼ ð1Þð1þÞ; 4     4 1 2 4 1

Advise: Why You Wasting Money in Costly SEO Tools, Use World's Best Free SEO Tool Ubersuggest.