# Hutton - Fundamentals of Finite Element Analysis

Посмотреть архив целикомПросмотреть файл в отдельном окне: bb806d45c329fd25e619128f03fb6a12.pdf

Hutton: Fundamentals of

Finite Element Analysis

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

PREFACE

F

undamentals of Finite Element Analysis is intended to be the text for a

senior-level ?nite element course in engineering programs. The most

appropriate major programs are civil engineering, engineering mechanics, and mechanical engineering. The ?nite element method is such a widely used

analysis-and-design technique that it is essential that undergraduate engineering

students have a basic knowledge of the theory and applications of the technique.

Toward that objective, I developed and taught an undergraduate “special topics”

course on the ?nite element method at Washington State University in the summer of 1992. The course was composed of approximately two-thirds theory and

one-third use of commercial software in solving ?nite element problems. Since

that time, the course has become a regularly offered technical elective in the

mechanical engineering program and is generally in high demand. During

the developmental process for the course, I was never satis?ed with any text that

was used, and we tried many. I found the available texts to be at one extreme or

the other; namely, essentially no theory and all software application, or all theory

and no software application. The former approach, in my opinion, represents

training in using computer programs, while the latter represents graduate-level

study. I have written this text to seek a middle ground.

Pedagogically, I believe that training undergraduate engineering students to

use a particular software package without providing knowledge of the underlying

theory is a disservice to the student and can be dangerous for their future employers. While I am acutely aware that most engineering programs have a speci?c

?nite element software package available for student use, I do not believe that the

text the students use should be tied only to that software. Therefore, I have written this text to be software-independent. I emphasize the basic theory of the ?nite

element method, in a context that can be understood by undergraduate engineering students, and leave the software-speci?c portions to the instructor.

As the text is intended for an undergraduate course, the prerequisites required

are statics, dynamics, mechanics of materials, and calculus through ordinary differential equations. Of necessity, partial differential equations are introduced

but in a manner that should be understood based on the stated prerequisites.

Applications of the ?nite element method to heat transfer and ?uid mechanics are

included, but the necessary derivations are such that previous coursework in

those topics is not required. Many students will have taken heat transfer and ?uid

mechanics courses, and the instructor can expand the topics based on the students’ background.

Chapter 1 is a general introduction to the ?nite element method and includes a description of the basic concept of dividing a domain into ?nite-size

subdomains. The ?nite difference method is introduced for comparison to the

xi

Hutton: Fundamentals of

Finite Element Analysis

xii

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

?nite element method. A general procedure in the sequence of model de?nition,

solution, and interpretation of results is discussed and related to the generally

accepted terms of preprocessing, solution, and postprocessing. A brief history of

the ?nite element method is included, as are a few examples illustrating application of the method.

Chapter 2 introduces the concept of a ?nite element stiffness matrix and

associated displacement equation, in terms of interpolation functions, using the

linear spring as a ?nite element. The linear spring is known to most undergraduate students so the mechanics should not be new. However, representation of

the spring as a ?nite element is new but provides a simple, concise example of

the ?nite element method. The premise of spring element formulation is extended to the bar element, and energy methods are introduced. The ?rst theorem

of Castigliano is applied, as is the principle of minimum potential energy.

Castigliano’s theorem is a simple method to introduce the undergraduate student

to minimum principles without use of variational calculus.

Chapter 3 uses the bar element of Chapter 2 to illustrate assembly of global

equilibrium equations for a structure composed of many ?nite elements. Transformation from element coordinates to global coordinates is developed and

illustrated with both two- and three-dimensional examples. The direct stiffness

method is utilized and two methods for global matrix assembly are presented.

Application of boundary conditions and solution of the resultant constraint equations is discussed. Use of the basic displacement solution to obtain element strain

and stress is shown as a postprocessing operation.

Chapter 4 introduces the beam/?exure element as a bridge to continuity

requirements for higher-order elements. Slope continuity is introduced and this

requires an adjustment to the assumed interpolation functions to insure continuity.

Nodal load vectors are discussed in the context of discrete and distributed loads,

using the method of work equivalence.

Chapters 2, 3, and 4 introduce the basic procedures of ?nite-element modeling in the context of simple structural elements that should be well-known to the

student from the prerequisite mechanics of materials course. Thus the emphasis

in the early part of the course in which the text is used can be on the ?nite element method without introduction of new physical concepts. The bar and beam

elements can be used to give the student practical truss and frame problems for

solution using available ?nite element software. If the instructor is so inclined,

the bar and beam elements (in the two-dimensional context) also provide a relatively simple framework for student development of ?nite element software

using basic programming languages.

Chapter 5 is the springboard to more advanced concepts of ?nite element

analysis. The method of weighted residuals is introduced as the fundamental

technique used in the remainder of the text. The Galerkin method is utilized

exclusively since I have found this method is both understandable for undergraduate students and is amenable to a wide range of engineering problems. The

material in this chapter repeats the bar and beam developments and extends the

?nite element concept to one-dimensional heat transfer. Application to the bar

Hutton: Fundamentals of

Finite Element Analysis

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

and beam elements illustrates that the method is in agreement with the basic mechanics approach of Chapters 2–4. Introduction of heat transfer exposes the student to additional applications of the ?nite element method that are, most likely,

new to the student.

Chapter 6 is a stand-alone description of the requirements of interpolation

functions used in developing ?nite element models for any physical problem.

Continuity and completeness requirements are delineated. Natural (serendipity)

coordinates, triangular coordinates, and volume coordinates are de?ned and used

to develop interpolation functions for several element types in two- and threedimensions. The concept of isoparametric mapping is introduced in the context of

the plane quadrilateral element. As a precursor to following chapters, numerical

integration using Gaussian quadrature is covered and several examples included.

The use of two-dimensional elements to model three-dimensional axisymmetric

problems is included.

Chapter 7 uses Galerkin’s ?nite element method to develop the ?nite element equations for several commonly encountered situations in heat transfer.

One-, two- and three-dimensional formulations are discussed for conduction and

convection. Radiation is not included, as that phenomenon introduces a nonlinearity that undergraduate students are not prepared to deal with at the intended

level of the text. Heat transfer with mass transport is included. The ?nite difference method in conjunction with the ?nite element method is utilized to present

methods of solving time-dependent heat transfer problems.

Chapter 8 introduces ?nite element applications to ?uid mechanics. The

general equations governing ?uid ?ow are so complex and nonlinear that the

topic is introduced via ideal ?ow. The stream function and velocity potential

function are illustrated and the applicable restrictions noted. Example problems

are included that note the analogy with heat transfer and use heat transfer ?nite

element solutions to solve ideal ?ow problems. A brief discussion of viscous

?ow shows the nonlinearities that arise when nonideal ?ows are considered.

Chapter 9 applies the ?nite element method to problems in solid mechanics

with the proviso that the material response is linearly elastic and small de?ection.

Both plane stress and plane strain are de?ned and the ?nite element formulations

developed for each case. General three-dimensional states of stress and axisymmetric stress are included. A model for torsion of noncircular sections is developed using the Prandtl stress function. The purpose of the torsion section is to

make the student aware that all torsionally loaded objects are not circular and the

analysis methods must be adjusted to suit geometry.

Chapter 10 introduces the concept of dynamic motion of structures. It is not

presumed that the student has taken a course in mechanical vibrations; as a result, this chapter includes a primer on basic vibration theory. Most of this material is drawn from my previously published text Applied Mechanical Vibrations.

The concept of the mass or inertia matrix is developed by examples of simple

spring-mass systems and then extended to continuous bodies. Both lumped and

consistent mass matrices are de?ned and used in examples. Modal analysis is the

basic method espoused for dynamic response; hence, a considerable amount of

xiii

Hutton: Fundamentals of

Finite Element Analysis

xiv

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

text material is devoted to determination of natural modes, orthogonality, and

modal superposition. Combination of ?nite difference and ?nite element methods for solving transient dynamic structural problems is included.

The appendices are included in order to provide the student with material

that might be new or may be “rusty” in the student’s mind.

Appendix A is a review of matrix algebra and should be known to the student from a course in linear algebra.

Appendix B states the general three-dimensional constitutive relations for

a homogeneous, isotropic, elastic material. I have found over the years that undergraduate engineering students do not have a ?rm grasp of these relations. In

general, the student has been exposed to so many special cases that the threedimensional equations are not truly understood.

Appendix C covers three methods for solving linear algebraic equations.

Some students may use this material as an outline for programming solution

methods. I include the appendix only so the reader is aware of the algorithms underlying the software he/she will use in solving ?nite element problems.

Appendix D describes the basic computational capabilities of the FEPC

software. The FEPC (FEP?nite element program for the PCpersonal computer)

was developed by the late Dr. Charles Knight of Virginia Polytechnic Institute

and State University and is used in conjunction with this text with permission of

his estate. Dr. Knight’s programs allow analysis of two-dimensional programs

using bar, beam, and plane stress elements. The appendix describes in general

terms the capabilities and limitations of the software. The FEPC program is

available to the student at www.mhhe.com/hutton.

Appendix E includes problems for several chapters of the text that should be

solved via commercial ?nite element software. Whether the instructor has available ANSYS, ALGOR, COSMOS, etc., these problems are oriented to systems

having many degrees of freedom and not amenable to hand calculation. Additional problems of this sort will be added to the website on a continuing basis.

The textbook features a Web site (www.mhhe.com/hutton) with ?nite element analysis links and the FEPC program. At this site, instructors will have

access to PowerPoint images and an Instructors’ Solutions Manual. Instructors

can access these tools by contacting their local McGraw-Hill sales representative

for password information.

I thank Raghu Agarwal, Rong Y. Chen, Nels Madsen, Robert L. Rankin,

Joseph J. Rencis, Stephen R. Swanson, and Lonny L. Thompson, who reviewed

some or all of the manuscript and provided constructive suggestions and criticisms that have helped improve the book.

I am grateful to all the staff at McGraw-Hill who have labored to make this

project a reality. I especially acknowledge the patient encouragement and professionalism of Jonathan Plant, Senior Editor, Lisa Kalner Williams, Developmental Editor, and Kay Brimeyer, Senior Project Manager.

David V. Hutton

Pullman, WA

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

C H A P T E R

1

Basic Concepts of the

Finite Element Method

1.1 INTRODUCTION

The ?nite element method (FEM), sometimes referred to as ?nite element

analysis (FEA), is a computational technique used to obtain approximate solutions of boundary value problems in engineering. Simply stated, a boundary

value problem is a mathematical problem in which one or more dependent variables must satisfy a differential equation everywhere within a known domain of

independent variables and satisfy speci?c conditions on the boundary of the

domain. Boundary value problems are also sometimes called ?eld problems. The

?eld is the domain of interest and most often represents a physical structure.

The ?eld variables are the dependent variables of interest governed by the differential equation. The boundary conditions are the speci?ed values of the ?eld

variables (or related variables such as derivatives) on the boundaries of the ?eld.

Depending on the type of physical problem being analyzed, the ?eld variables

may include physical displacement, temperature, heat ?ux, and ?uid velocity to

name only a few.

1.2 HOW DOES THE FINITE ELEMENT

METHOD WORK?

The general techniques and terminology of ?nite element analysis will be introduced with reference to Figure 1.1. The ?gure depicts a volume of some material

or materials having known physical properties. The volume represents the

domain of a boundary value problem to be solved. For simplicity, at this point,

we assume a two-dimensional case with a single ?eld variable (x, y) to be

determined at every point P(x, y) such that a known governing equation (or equations) is satis?ed exactly at every such point. Note that this implies an exact

1

Hutton: Fundamentals of

Finite Element Analysis

2

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

3

P(x, y)

2

1

(a)

(b)

(c)

Figure 1.1

(a) A general two-dimensional domain of ?eld variable (x, y).

(b) A three-node ?nite element de?ned in the domain. (c) Additional

elements showing a partial ?nite element mesh of the domain.

mathematical solution is obtained; that is, the solution is a closed-form algebraic

expression of the independent variables. In practical problems, the domain may

be geometrically complex as is, often, the governing equation and the likelihood

of obtaining an exact closed-form solution is very low. Therefore, approximate

solutions based on numerical techniques and digital computation are most

often obtained in engineering analyses of complex problems. Finite element

analysis is a powerful technique for obtaining such approximate solutions with

good accuracy.

A small triangular element that encloses a ?nite-sized subdomain of the area

of interest is shown in Figure 1.1b. That this element is not a differential element

of size dx ? dy makes this a ?nite element. As we treat this example as a twodimensional problem, it is assumed that the thickness in the z direction is constant and z dependency is not indicated in the differential equation. The vertices

of the triangular element are numbered to indicate that these points are nodes. A

node is a speci?c point in the ?nite element at which the value of the ?eld variable is to be explicitly calculated. Exterior nodes are located on the boundaries

of the ?nite element and may be used to connect an element to adjacent ?nite

elements. Nodes that do not lie on element boundaries are interior nodes and

cannot be connected to any other element. The triangular element of Figure 1.1b

has only exterior nodes.

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

If the values of the ?eld variable are computed only at nodes, how are values

obtained at other points within a ?nite element? The answer contains the crux of

the ?nite element method: The values of the ?eld variable computed at the nodes

are used to approximate the values at nonnodal points (that is, in the element

interior) by interpolation of the nodal values. For the three-node triangle example, the nodes are all exterior and, at any other point within the element, the ?eld

variable is described by the approximate relation

(x , y) = N 1 (x , y)1 + N 2 (x , y)2 + N 3 (x , y)3

(1.1)

where 1 , 2 , and 3 are the values of the ?eld variable at the nodes, and N1, N2,

and N3 are the interpolation functions, also known as shape functions or blending functions. In the ?nite element approach, the nodal values of the ?eld variable are treated as unknown constants that are to be determined. The interpolation functions are most often polynomial forms of the independent variables,

derived to satisfy certain required conditions at the nodes. These conditions are

discussed in detail in subsequent chapters. The major point to be made here is

that the interpolation functions are predetermined, known functions of the independent variables; and these functions describe the variation of the ?eld variable

within the ?nite element.

The triangular element described by Equation 1.1 is said to have 3 degrees

of freedom, as three nodal values of the ?eld variable are required to describe

the ?eld variable everywhere in the element. This would be the case if the ?eld

variable represents a scalar ?eld, such as temperature in a heat transfer problem

(Chapter 7). If the domain of Figure 1.1 represents a thin, solid body subjected to

plane stress (Chapter 9), the ?eld variable becomes the displacement vector and

the values of two components must be computed at each node. In the latter case,

the three-node triangular element has 6 degrees of freedom. In general, the number of degrees of freedom associated with a ?nite element is equal to the product

of the number of nodes and the number of values of the ?eld variable (and possibly its derivatives) that must be computed at each node.

How does this element-based approach work over the entire domain of interest? As depicted in Figure 1.1c, every element is connected at its exterior

nodes to other elements. The ?nite element equations are formulated such that, at

the nodal connections, the value of the ?eld variable at any connection is the

same for each element connected to the node. Thus, continuity of the ?eld variable at the nodes is ensured. In fact, ?nite element formulations are such that

continuity of the ?eld variable across interelement boundaries is also ensured.

This feature avoids the physically unacceptable possibility of gaps or voids occurring in the domain. In structural problems, such gaps would represent physical separation of the material. In heat transfer, a “gap” would manifest itself in

the form of different temperatures at the same physical point.

Although continuity of the ?eld variable from element to element is inherent

to the ?nite element formulation, interelement continuity of gradients (i.e., derivatives) of the ?eld variable does not generally exist. This is a critical observation. In most cases, such derivatives are of more interest than are ?eld variable

values. For example, in structural problems, the ?eld variable is displacement but

3

Hutton: Fundamentals of

Finite Element Analysis

4

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

the true interest is more often in strain and stress. As strain is de?ned in terms of

?rst derivatives of displacement components, strain is not continuous across element boundaries. However, the magnitudes of discontinuities of derivatives can

be used to assess solution accuracy and convergence as the number of elements

is increased, as is illustrated by the following example.

1.2.1 Comparison of Finite Element and Exact Solutions

The process of representing a physical domain with ?nite elements is referred to

as meshing, and the resulting set of elements is known as the ?nite element mesh.

As most of the commonly used element geometries have straight sides, it is generally impossible to include the entire physical domain in the element mesh if the

domain includes curved boundaries. Such a situation is shown in Figure 1.2a,

where a curved-boundary domain is meshed (quite coarsely) using square elements. A re?ned mesh for the same domain is shown in Figure 1.2b, using

smaller, more numerous elements of the same type. Note that the re?ned mesh

includes signi?cantly more of the physical domain in the ?nite element representation and the curved boundaries are more closely approximated. (Triangular

elements could approximate the boundaries even better.)

If the interpolation functions satisfy certain mathematical requirements

(Chapter 6), a ?nite element solution for a particular problem converges to the

exact solution of the problem. That is, as the number of elements is increased and

the physical dimensions of the elements are decreased, the ?nite element solution

changes incrementally. The incremental changes decrease with the mesh re?nement process and approach the exact solution asymptotically. To illustrate

convergence, we consider a relatively simple problem that has a known solution.

Figure 1.3a depicts a tapered, solid cylinder ?xed at one end and subjected to

a tensile load at the other end. Assuming the displacement at the point of load

application to be of interest, a ?rst approximation is obtained by considering

the cylinder to be uniform, having a cross-sectional area equal to the average area

(a)

(b)

Figure 1.2

(a) Arbitrary curved-boundary domain modeled using square elements. Stippled

areas are not included in the model. A total of 41 elements is shown. (b) Re?ned

?nite element mesh showing reduction of the area not included in the model. A

total of 192 elements is shown.

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

of the cylinder (Figure 1.3b). The uniform bar is a link or bar ?nite element

(Chapter 2), so our ?rst approximation is a one-element, ?nite element model.

The solution is obtained using the strength of materials theory. Next, we model

the tapered cylinder as two uniform bars in series, as in Figure 1.3c. In the twoelement model, each element is of length equal to half the total length of the

cylinder and has a cross-sectional area equal to the average area of the corresponding half-length of the cylinder. The mesh re?nement is continued using a

four-element model, as in Figure 1.3d, and so on. For this simple problem, the

displacement of the end of the cylinder for each of the ?nite element models is as

shown in Figure 1.4a, where the dashed line represents the known solution. Convergence of the ?nite element solutions to the exact solution is clearly indicated.

ro

x

r

L

A

Ao AL

2

rL

F

(a)

(b)

Element 1

Element 2

(c)

(d)

Figure 1.3

(a) Tapered circular cylinder subjected to tensile loading:

r(x) r0 (x/L)(r0 rL). (b) Tapered cylinder as a single axial

(bar) element using an average area. Actual tapered cylinder

is shown as dashed lines. (c) Tapered cylinder modeled as

two, equal-length, ?nite elements. The area of each element

is average over the respective tapered cylinder length.

(d) Tapered circular cylinder modeled as four, equal-length

?nite elements. The areas are average over the respective

length of cylinder (element length L?4).

5

Hutton: Fundamentals of

Finite Element Analysis

6

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

Exact

Four elements

?

( Lx )

?(x L)

Exact

1

2

3

Number of elements

(a)

4

0.25

0.5

0.75

1.0

x

L

(b)

Figure 1.4

(a) Displacement at x L for tapered cylinder in tension of Figure 1.3. (b) Comparison of the exact solution

and the four-element solution for a tapered cylinder in tension.

On the other hand, if we plot displacement as a function of position along the

length of the cylinder, we can observe convergence as well as the approximate

nature of the ?nite element solutions. Figure 1.4b depicts the exact strength of

materials solution and the displacement solution for the four-element models.

We note that the displacement variation in each element is a linear approximation

to the true nonlinear solution. The linear variation is directly attributable to the

fact that the interpolation functions for a bar element are linear. Second, we note

that, as the mesh is re?ned, the displacement solution converges to the nonlinear

solution at every point in the solution domain.

The previous paragraph discussed convergence of the displacement of the

tapered cylinder. As will be seen in Chapter 2, displacement is the primary ?eld

variable in structural problems. In most structural problems, however, we are

interested primarily in stresses induced by speci?ed loadings. The stresses must

be computed via the appropriate stress-strain relations, and the strain components are derived from the displacement ?eld solution. Hence, strains and

stresses are referred to as derived variables. For example, if we plot the element

stresses for the tapered cylinder example just cited for the exact solution as well

as the ?nite element solutions for two- and four-element models as depicted in

Figure 1.5, we observe that the stresses are constant in each element and represent a discontinuous solution of the problem in terms of stresses and strains. We

also note that, as the number of elements increases, the jump discontinuities in

stress decrease in magnitude. This phenomenon is characteristic of the ?nite element method. The formulation of the ?nite element method for a given problem

is such that the primary ?eld variable is continuous from element to element but

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

4.5

4.0

Exact

Two elements

Four elements

3.5

?

?0

3.0

2.5

2.0

1.5

1.0

0.25

0.5

x

L

0.75

1.0

Figure 1.5

Comparison of the computed axial stress value in a

tapered cylinder: 0 F?A0.

the derived variables are not necessarily continuous. In the limiting process of

mesh re?nement, the derived variables become closer and closer to continuity.

Our example shows how the ?nite element solution converges to a known

exact solution (the exactness of the solution in this case is that of strength of

materials theory). If we know the exact solution, we would not be applying the

?nite element method! So how do we assess the accuracy of a ?nite element solution for a problem with an unknown solution? The answer to this question is not

simple. If we did not have the dashed line in Figure 1.3 representing the exact

solution, we could still discern convergence to a solution. Convergence of a

numerical method (such as the ?nite element method) is by no means assurance

that the convergence is to the correct solution. A person using the ?nite element

analysis technique must examine the solution analytically in terms of (1) numerical convergence, (2) reasonableness (does the result make sense?), (3) whether the

physical laws of the problem are satis?ed (is the structure in equilibrium? Does the

heat output balance with the heat input?), and (4) whether the discontinuities in

value of derived variables across element boundaries are reasonable. Many

such questions must be posed and examined prior to accepting the results of a ?nite

element analysis as representative of a correct solution useful for design purposes.

1.2.2 Comparison of Finite Element and Finite

Difference Methods

The ?nite difference method is another numerical technique frequently used to

obtain approximate solutions of problems governed by differential equations.

Details of the technique are discussed in Chapter 7 in the context of transient heat

7

Hutton: Fundamentals of

Finite Element Analysis

8

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

transfer. The method is also illustrated in Chapter 10 for transient dynamic analysis of structures. Here, we present the basic concepts of the ?nite difference

method for purposes of comparison.

The ?nite difference method is based on the de?nition of the derivative of a

function f (x ) that is

d f (x )

f (x + x ) ? f (x )

= lim

x?0

dx

x

(1.2)

where x is the independent variable. In the ?nite difference method, as implied

by its name, derivatives are calculated via Equation 1.2 using small, but ?nite,

values of x to obtain

d f (x )

f (x + x ) ? f (x )

?

dx

x

(1.3)

A differential equation such as

df

+x =0

dx

0?x ?1

(1.4)

f (x + x ) ? f (x )

+x =0

x

(1.5)

is expressed as

in the ?nite difference method. Equation 1.5 can be rewritten as

f (x + x ) = f (x ) ? x (x )

(1.6)

where we note that the equality must be taken as “approximately equals.” From

differential equation theory, we know that the solution of a ?rst-order differential

equation contains one constant of integration. The constant of integration must

be determined such that one given condition (a boundary condition or initial condition) is satis?ed. In the current example, we assume that the speci?ed condition

is x (0) = A = constant. If we choose an integration step x to be a small, constant value (the integration step is not required to be constant), then we can write

x i+1 = x i + x

i = 0, N

(1.7)

where N is the total number of steps required to cover the domain. Equation 1.6

is then

f i+1 = f i ? x i (x )

f0 = A

i = 0, N

(1.8)

Equation 1.8 is known as a recurrence relation and provides an approximation to

the value of the unknown function f (x) at a number of discrete points in the domain of the problem.

To illustrate, Figure 1.6a shows the exact solution f (x ) = 1 ? x 2 /2 and a

?nite difference solution obtained with x = 0.1. The ?nite difference solution is

shown at the discrete points of function evaluation only. The manner of variation

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

1

0.8

f (x)

0.6

0.4

0.2

0

0

0.2

0.6

0.4

0.8

1

x

Figure 1.6

Comparison of the exact and ?nite difference

solutions of Equation 1.4 with f0 A 1.

of the function between the calculated points is not known in the ?nite difference

method. One can, of course, linearly interpolate the values to produce an approximation to the curve of the exact solution but the manner of interpolation is

not an a priori determination in the ?nite difference method.

To contrast the ?nite difference method with the ?nite element method,

we note that, in the ?nite element method, the variation of the ?eld variable in

the physical domain is an integral part of the procedure. That is, based on the

selected interpolation functions, the variation of the ?eld variable throughout a

?nite element is speci?ed as an integral part of the problem formulation. In the

?nite difference method, this is not the case: The ?eld variable is computed at

speci?ed points only. The major rami?cation of this contrast is that derivatives

(to a certain level) can be computed in the ?nite element approach, whereas the

?nite difference method provides data only on the variable itself. In a structural

problem, for example, both methods provide displacement solutions, but the

?nite element solution can be used to directly compute strain components (?rst

derivatives). To obtain strain data in the ?nite difference method requires additional considerations not inherent to the mathematical model.

There are also certain similarities between the two methods. The integration

points in the ?nite difference method are analogous to the nodes in a ?nite

element model. The variable of interest is explicitly evaluated at such points.

Also, as the integration step (step size) in the ?nite difference method is reduced,

the solution is expected to converge to the exact solution. This is similar to the

expected convergence of a ?nite element solution as the mesh of elements is

re?ned. In both cases, the re?nement represents reduction of the mathematical

model from ?nite to in?nitesimal. And in both cases, differential equations are

reduced to algebraic equations.

9

Hutton: Fundamentals of

Finite Element Analysis

10

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

Probably the most descriptive way to contrast the two methods is to note that

the ?nite difference method models the differential equation(s) of the problem

and uses numerical integration to obtain the solution at discrete points. The ?nite

element method models the entire domain of the problem and uses known physical principles to develop algebraic equations describing the approximate solutions. Thus, the ?nite difference method models differential equations while the

?nite element method can be said to more closely model the physical problem at

hand. As will be observed in the remainder of this text, there are cases in which

a combination of ?nite element and ?nite difference methods is very useful and

ef?cient in obtaining solutions to engineering problems, particularly where dynamic (transient) effects are important.

1.3 A GENERAL PROCEDURE FOR FINITE

ELEMENT ANALYSIS

Certain steps in formulating a ?nite element analysis of a physical problem are

common to all such analyses, whether structural, heat transfer, ?uid ?ow, or

some other problem. These steps are embodied in commercial ?nite element

software packages (some are mentioned in the following paragraphs) and are

implicitly incorporated in this text, although we do not necessarily refer to the

steps explicitly in the following chapters. The steps are described as follows.

1.3.1 Preprocessing

The preprocessing step is, quite generally, described as de?ning the model and

includes

De?ne the geometric domain of the problem.

De?ne the element type(s) to be used (Chapter 6).

De?ne the material properties of the elements.

De?ne the geometric properties of the elements (length, area, and the like).

De?ne the element connectivities (mesh the model).

De?ne the physical constraints (boundary conditions).

De?ne the loadings.

The preprocessing (model de?nition) step is critical. In no case is there a better

example of the computer-related axiom “garbage in, garbage out.” A perfectly

computed ?nite element solution is of absolutely no value if it corresponds to the

wrong problem.

1.3.2 Solution

During the solution phase, ?nite element software assembles the governing algebraic equations in matrix form and computes the unknown values of the primary

?eld variable(s). The computed values are then used by back substitution to

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

1.4 Brief History of the Finite Element Method

compute additional, derived variables, such as reaction forces, element stresses,

and heat ?ow.

As it is not uncommon for a ?nite element model to be represented by tens

of thousands of equations, special solution techniques are used to reduce data

storage requirements and computation time. For static, linear problems, a wave

front solver, based on Gauss elimination (Appendix C), is commonly used. While

a complete discussion of the various algorithms is beyond the scope of this text,

the interested reader will ?nd a thorough discussion in the Bathe book [1].

1.3.3 Postprocessing

Analysis and evaluation of the solution results is referred to as postprocessing.

Postprocessor software contains sophisticated routines used for sorting, printing,

and plotting selected results from a ?nite element solution. Examples of operations that can be accomplished include

Sort element stresses in order of magnitude.

Check equilibrium.

Calculate factors of safety.

Plot deformed structural shape.

Animate dynamic model behavior.

Produce color-coded temperature plots.

While solution data can be manipulated many ways in postprocessing, the most

important objective is to apply sound engineering judgment in determining

whether the solution results are physically reasonable.

1.4 BRIEF HISTORY OF THE FINITE

ELEMENT METHOD

The mathematical roots of the ?nite element method dates back at least a half

century. Approximate methods for solving differential equations using trial solutions are even older in origin. Lord Rayleigh [2] and Ritz [3] used trial functions

(in our context, interpolation functions) to approximate solutions of differential

equations. Galerkin [4] used the same concept for solution. The drawback in the

earlier approaches, compared to the modern ?nite element method, is that the

trial functions must apply over the entire domain of the problem of concern.

While the Galerkin method provides a very strong basis for the ?nite element

method (Chapter 5), not until the 1940s, when Courant [5] introduced the concept of piecewise-continuous functions in a subdomain, did the ?nite element

method have its real start.

In the late 1940s, aircraft engineers were dealing with the invention of the jet

engine and the needs for more sophisticated analysis of airframe structures to

withstand larger loads associated with higher speeds. These engineers, without

the bene?t of modern computers, developed matrix methods of force analysis,

© The McGraw?Hill

Companies, 2004

11

Hutton: Fundamentals of

Finite Element Analysis

12

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

collectively known as the ?exibility method, in which the unknowns are the

forces and the knowns are displacements. The ?nite element method, in its most

often-used form, corresponds to the displacement method, in which the unknowns are system displacements in response to applied force systems. In this

text, we adhere exclusively to the displacement method. As will be seen as we

proceed, the term displacement is quite general in the ?nite element method and

can represent physical displacement, temperature, or ?uid velocity, for example.

The term ?nite element was ?rst used by Clough [6] in 1960 in the context of

plane stress analysis and has been in common usage since that time.

During the decades of the 1960s and 1970s, the ?nite element method was

extended to applications in plate bending, shell bending, pressure vessels, and

general three-dimensional problems in elastic structural analysis [7–11] as well

as to ?uid ?ow and heat transfer [12, 13]. Further extension of the method to

large de?ections and dynamic analysis also occurred during this time period

[14 , 15]. An excellent history of the ?nite element method and detailed bibliography is given by Noor [16].

The ?nite element method is computationally intensive, owing to the required

operations on very large matrices. In the early years, applications were performed

using mainframe computers, which, at the time, were considered to be very powerful, high-speed tools for use in engineering analysis. During the 1960s, the ?nite

element software code NASTRAN [17] was developed in conjunction with the

space exploration program of the United States. NASTRAN was the ?rst major

?nite element software code. It was, and still is, capable of hundreds of thousands

of degrees of freedom (nodal ?eld variable computations). In the years since the

development of NASTRAN, many commercial software packages have been introduced for ?nite element analysis. Among these are ANSYS [18], ALGOR [19],

and COSMOS/M [20]. In today’s computational environment, most of these

packages can be used on desktop computers and engineering workstations to

obtain solutions to large problems in static and dynamic structural analysis, heat

transfer, ?uid ?ow, electromagnetics, and seismic response. In this text, we do not

utilize or champion a particular code. Rather, we develop the fundamentals for

understanding of ?nite element analysis to enable the reader to use such software

packages with an educated understanding.

1.5 EXAMPLES OF FINITE ELEMENT

ANALYSIS

We now present, brie?y, a few examples of the types of problems that can be

analyzed via the ?nite element method. Figure 1.7 depicts a rectangular region

with a central hole. The area has been “meshed” with a ?nite element grid of twodimensional elements assumed to have a constant thickness in the z direction.

Note that the mesh of elements is irregular: The element shapes (triangles and

quadrilaterals) and sizes vary. In particular, note that around the geometric discontinuity of the hole, the elements are of smaller size. This represents not only

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.5 Examples of Finite Element Analysis

Figure 1.7

A mesh of ?nite elements over a rectangular region having a

central hole.

an improvement in geometric accuracy in the vicinity of the discontinuity but

also solution accuracy, as is discussed in subsequent chapters.

The geometry depicted in Figure 1.7 could represent the ?nite element

model of several physical problems. For plane stress analysis, the geometry

would represent a thin plate with a central hole subjected to edge loading in the

plane depicted. In this case, the ?nite element solution would be used to examine stress concentration effects in the vicinity of the hole. The element mesh

shown could also represent the case of ?uid ?ow around a circular cylinder. In

yet another application, the model shown could depict a heat transfer ?n attached to a pipe (the hole) from which heat is transferred to the ?n for dissipation to the surroundings. In each case, the formulation of the equations governing physical behavior of the elements in response to external in?uences is quite

different.

Figure 1.8a shows a truss module that was at one time considered a

building-block element for space station construction [21]. Designed to fold in

accordion fashion into a small volume for transport into orbit, the module, when

deployed, extends to overall dimensions 1.4 m ? 1.4 m ? 2.8 m. By attaching

such modules end-to-end, a truss of essentially any length could be obtained.

The structure was analyzed via the ?nite element method to determine the

vibration characteristics as the number of modules, thus overall length, was

varied. As the connections between the various structural members are pin or

ball-and-socket joints, a simple axial tension-compression element (Chapter 2)

was used in the model. The ?nite element model of one module was composed

of 33 elements. A sample vibration shape of a ?ve-module truss is shown in

Figure 1.8b.

The truss example just described involves a rather large structure modeled

by a small number of relatively large ?nite elements. In contrast, Figure 1.9

shows the ?nite element model of a very thin tube designed for use in heat

13

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

XG

YG

ZG

(a)

(b)

Figure 1.8

(a) Deployable truss module showing details of folding joints.

(b) A sample vibration-mode shape of a ?ve-module truss as obtained

via ?nite element analysis. (Courtesy: AIAA)

14

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.5 Examples of Finite Element Analysis

0.00197

Z

X

0.488

0.25

Figure 1.9

Finite element model of a thin-walled

heat exchanger tube.

transfer in a spacecraft application. The tube has inside diameter of 0.976 in. and

wall thickness 0.00197 in. and overall length 36 in. Materials considered for

construction of the tube were copper and titanium alloys. Owing to the wall

thickness, prototype tubes were found to be very fragile and dif?cult to handle

without damage. The objectives of the ?nite element analysis were to examine

the bending, torsional, and buckling loads allowable. The ?gure shows the ?nite

element mesh used to model a section of the tube only 0.25 in. in length. This

model contains 1920 three-dimensional solid elements, each having eight nodes

with 3 degrees of freedom at each node. Such a large number of elements was

required for a small structure in consideration of computational accuracy. The

concern here was the so-called aspect ratio of the elements, as is de?ned and

discussed in subsequent chapters.

As a ?nal example, Figure 1.10a represents the ?nite element model of the

main load-carrying component of a prosthetic device. The device is intended to

be a hand attachment to an arti?cial arm. In use, the hand would allow a lower

arm amputee to engage in weight lifting as part of a physical ?tness program.

The ?nite element model was used to determine the stress distribution in the

component in terms of the range of weight loading anticipated, so as to properly

size the component and select the material. Figure 1.10b shows a prototype of the

completed hand design.

15

Hutton: Fundamentals of

Finite Element Analysis

16

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

(a)

(b)

Figure 1.10

(a) A ?nite element model of a prosthetic hand for weightlifting. (b) Completed

prototype of a prosthetic hand, attached to a bar.

(Courtesy of Payam Sadat. All rights reserved.)

1.6 OBJECTIVES OF THE TEXT

I wrote Fundamentals of Finite Element Analysis for use in senior-level ?nite

element courses in engineering programs. The majority of available textbooks

on the ?nite element method are written for graduate-level courses. These

texts are heavy on the theory of ?nite element analysis and rely on mathematical

techniques (notably, variational calculus) that are not usually in the repertoire of

undergraduate engineering students. Knowledge of advanced mathematical techniques is not required for successful use of this text. The prerequisite study is

based on the undergraduate coursework common to most engineering programs:

linear algebra, calculus through differential equations, and the usual series of

statics, dynamics, and mechanics of materials. Although not required, prior study

of ?uid mechanics and heat transfer is helpful. Given this assumed background,

the ?nite element method is developed on the basis of physical laws (equilibrium, conservation of mass, and the like), the principle of minimum potential energy (Chapter 2), and Galerkin’s ?nite element method (introduced and developed in Chapter 5).

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

References

As the reader progresses through the text, he or she will discern that we

cover a signi?cant amount of ?nite element theory in addition to application

examples. Given the availability of many powerful and sophisticated ?nite

element software packages, why study the theory? The ?nite element method is

a tool, and like any other tool, using it without proper instruction can be quite

dangerous. My premise is that the proper instruction in this context includes

understanding the basic theory underlying formulation of ?nite element models

of physical problems. As stated previously, critical analysis of the results of a

?nite element model computation is essential, since those results may eventually

become the basis for design. Knowledge of the theory is necessary for both

proper modeling and evaluation of computational results.

REFERENCES

1.

Bathe, K-J. Finite Element Procedures. Englewood Cliffs, NJ: Prentice-Hall,

1996.

2. Lord Rayleigh. “On the Theory of Resonance.” Transactions of the Royal Society

(London) A161 (1870).

3. Ritz, W. “Uber eine neue Methode zur Losung gewissen Variations-Probleme der

mathematischen Physik.” J. Reine Angew. Math. 135 (1909).

4. Galerkin, B. G. “Series Solution of Some Problems of Elastic Equilibrium of Rods

and Plates” [in Russian]. Vestn. Inzh. Tekh. 19 (1915).

5. Courant, R. “Variational Methods for the Solution of Problems of Equilibrium and

Vibrations.” Bulletin of the American Mathematical Society 49 (1943).

6. Clough, R. W. “The Finite Element Method in Plane Stress Analysis.”

Proceedings, American Society of Civil Engineers, Second Conference on

Electronic Computation, Pittsburgh, 1960.

7. Melosh, R. J. “A Stiffness Method for the Analysis of Thin Plates in Bending.”

Journal of Aerospace Sciences 28, no. 1 (1961).

8. Grafton, P. E., and D. R. Strome. “Analysis of Axisymmetric Shells by the Direct

Stiffness Method.” Journal of the American Institute of Aeronautics and

Astronautics 1, no. 10 (1963).

9. Gallagher, R. H. “Analysis of Plate and Shell Structures.” Proceedings,

Symposium on the Application of Finite Element Methods in Civil Engineering,

Vanderbilt University, Nashville, 1969.

10. Wilson, E. L. “Structural Analysis of Axisymmetric Solids.” Journal of the

American Institute of Aeronautics and Astronautics 3, (1965).

11. Melosh, R. J. “Structural Analysis of Solids.” Journal of the Structural Division,

Proceedings of the American Society of Civil Engineers, August 1963.

12. Martin, H. C. “Finite Element Analysis of Fluid Flows.” Proceedings of the

Second Conference on Matrix Methods in Structural Mechanics, Wright-Patterson

Air Force Base, Kilborn, Ohio, October 1968.

13. Wilson, E. L., and R. E. Nickell. “Application of the Finite Element Method to

Heat Conduction Analysis.” Nuclear Engineering Design 4 (1966).

17

Hutton: Fundamentals of

Finite Element Analysis

18

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

14. Turner, M. J., E. H. Dill, H. C. Martin, and R. J. Melosh. “Large De?ections of

Structures Subjected to Heating and External Loads.” Journal of Aeronautical

Sciences 27 (1960).

15. Archer, J. S. “Consistent Mass Matrix Formulations for Structural Analysis Using

Finite Element Techniques.” Journal of the American Institute of Aeronautics and

Astronautics 3, no. 10 (1965).

16. Noor, A. K. “Bibliography of Books and Monographs on Finite Element

Technology.” Applied Mechanics Reviews 44, no. 6 (1991).

17. MSC/NASTRAN. Lowell, MA: MacNeal-Schwindler Corp.

18. ANSYS. Houston, PA: Swanson Analysis Systems Inc.

19. ALGOR. Pittsburgh: Algor Interactive Systems.

20. COSMOS/M. Los Angeles: Structural Research and Analysis Corp.

21. Hutton, D. V. “Modal Analysis of a Deployable Truss Using the Finite Element

Method.” Journal of Spacecraft and Rockets 21, no. 5 (1984).

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

C H A P T E R

2

Stiffness Matrices, Spring

and Bar Elements

2.1 INTRODUCTION

The primary characteristics of a ?nite element are embodied in the element

stiffness matrix. For a structural ?nite element, the stiffness matrix contains the

geometric and material behavior information that indicates the resistance of

the element to deformation when subjected to loading. Such deformation may

include axial, bending, shear, and torsional effects. For ?nite elements used in

nonstructural analyses, such as ?uid ?ow and heat transfer, the term stiffness

matrix is also used, since the matrix represents the resistance of the element to

change when subjected to external in?uences.

This chapter develops the ?nite element characteristics of two relatively

simple, one-dimensional structural elements, a linearly elastic spring and an elastic tension-compression member. These are selected as introductory elements because the behavior of each is relatively well-known from the commonly studied

engineering subjects of statics and strength of materials. Thus, the “bridge” to the

?nite element method is not obscured by theories new to the engineering student.

Rather, we build on known engineering principles to introduce ?nite element

concepts. The linear spring and the tension-compression member (hereafter referred to as a bar element and also known in the ?nite element literature as a spar,

link, or truss element) are also used to introduce the concept of interpolation

functions. As mentioned brie?y in Chapter 1, the basic premise of the ?nite element method is to describe the continuous variation of the ?eld variable (in this

chapter, physical displacement) in terms of discrete values at the ?nite element

nodes. In the interior of a ?nite element, as well as along the boundaries (applicable to two- and three-dimensional problems), the ?eld variable is described via

interpolation functions (Chapter 6) that must satisfy prescribed conditions.

Finite element analysis is based, dependent on the type of problem, on several mathematic/physical principles. In the present introduction to the method,

19

20

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

we present several such principles applicable to ?nite element analysis. First, and

foremost, for spring and bar systems, we utilize the principle of static equilibrium but—and this is essential—we include deformation in the development;

that is, we are not dealing with rigid body mechanics. For extension of the ?nite

element method to more complicated elastic structural systems, we also state and

apply the ?rst theorem of Castigliano [1] and the more widely used principle of

minimum potential energy [2]. Castigliano’s ?rst theorem, in the form presented,

may be new to the reader. The ?rst theorem is the counterpart of Castigliano’s

second theorem, which is more often encountered in the study of elementary

strength of materials [3]. Both theorems relate displacements and applied forces

to the equilibrium conditions of a mechanical system in terms of mechanical

energy. The use here of Castigliano’s ?rst theorem is for the distinct purpose of

introducing the concept of minimum potential energy without resort to the higher

mathematic principles of the calculus of variations, which is beyond the mathematical level intended for this text.

2.2 LINEAR SPRING AS A FINITE ELEMENT

A linear elastic spring is a mechanical device capable of supporting axial loading

only and constructed such that, over a reasonable operating range (meaning extension or compression beyond undeformed length), the elongation or contraction of the spring is directly proportional to the applied axial load. The constant

of proportionality between deformation and load is referred to as the spring constant, spring rate, or spring stiffness [4], generally denoted as k, and has units

of force per unit length. Formulation of the linear spring as a ?nite element is

accomplished with reference to Figure 2.1a. As an elastic spring supports axial

loading only, we select an element coordinate system (also known as a local coordinate system) as an x axis oriented along the length of the spring, as shown.

The element coordinate system is embedded in the element and chosen, by geometric convenience, for simplicity in describing element behavior. The element

u2

u1

f1

1

Force, f

Hutton: Fundamentals of

Finite Element Analysis

2

k

(a)

f2

k

1

x

Deflection, ? u2 u1

(b)

Figure 2.1

(a) Linear spring element with nodes, nodal displacements, and nodal forces.

(b) Load-de?ection curve.

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

or local coordinate system is contrasted with the global coordinate system. The

global coordinate system is that system in which the behavior of a complete

structure is to be described. By complete structure is meant the assembly of

many ?nite elements (at this point, several springs) for which we desire to compute response to loading conditions. In this chapter, we deal with cases in which

the local and global coordinate systems are essentially the same except for translation of origin. In two- and three-dimensional cases, however, the distinctions

are quite different and require mathematical recti?cation of element coordinate

systems to a common basis. The common basis is the global coordinate system.

Returning attention to Figure 2.1a, the ends of the spring are the nodes and

the nodal displacements are denoted by u1 and u2 and are shown in the positive

sense. If these nodal displacements are known, the total elongation or contraction

of the spring is known as is the net force in the spring. At this point in our development, we require that forces be applied to the element only at the nodes (distributed forces are accommodated for other element types later), and these are

denoted as f1 and f2 and are also shown in the positive sense.

Assuming that both the nodal displacements are zero when the spring is undeformed, the net spring deformation is given by

= u2 ? u1

(2.1)

and the resultant axial force in the spring is

f = k = k(u 2 ? u 1 )

(2.2)

as is depicted in Figure 2.1b.

For equilibrium, f 1 + f 2 = 0 or f 1 = ? f 2 , and we can rewrite Equation 2.2

in terms of the applied nodal forces as

f 1 = ?k(u 2 ? u 1 )

(2.3a)

f 2 = k(u 2 ? u 1 )

(2.3b)

which can be expressed in matrix form (see Appendix A for a review of matrix

algebra) as

k

?k

u1

f1

=

(2.4)

?k

k

u2

f2

or

[k e ]{u} = { f }

where

[k e ] =

k

?k

?k

k

(2.5)

(2.6)

is de?ned as the element stiffness matrix in the element coordinate system (or

local system), {u} is the column matrix (vector) of nodal displacements, and {f}

is the column matrix (vector) of element nodal forces. (In subsequent chapters,

21

Hutton: Fundamentals of

Finite Element Analysis

22

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

the matrix notation is used extensively. A general matrix is designated by

brackets [ ] and a column matrix (vector) by braces { }.)

Equation 2.6 shows that the element stiffness matrix for the linear spring

element is a 2 ? 2 matrix. This corresponds to the fact that the element exhibits

two nodal displacements (or degrees of freedom) and that the two displacements

are not independent (that is, the body is continuous and elastic). Furthermore, the

matrix is symmetric. A symmetric matrix has off-diagonal terms such that k i j =

kji. Symmetry of the stiffness matrix is indicative of the fact that the body is linearly elastic and each displacement is related to the other by the same physical

phenomenon. For example, if a force F (positive, tensile) is applied at node 2

with node 1 held ?xed, the relative displacement of the two nodes is the same as

if the force is applied symmetrically (negative, tensile) at node 1 with node 2

?xed. (Counterexamples to symmetry are seen in heat transfer and ?uid ?ow

analyses in Chapters 7 and 8.) As will be seen as more complicated structural

elements are developed, this is a general result: An element exhibiting N degrees

of freedom has a corresponding N ? N, symmetric stiffness matrix.

Next consider solution of the system of equations represented by Equation 2.4. In general, the nodal forces are prescribed and the objective is to solve

for the unknown nodal displacements. Formally, the solution is represented by

u1

f1

?1

= [k e ]

(2.7)

u2

f2

where [k e ]?1 is the inverse of the element stiffness matrix. However, this inverse

matrix does not exist, since the determinant of the element stiffness matrix is

identically zero. Therefore, the element stiffness matrix is singular, and this also

proves to be a general result in most cases. The physical signi?cance of the

singular nature of the element stiffness matrix is found by reexamination of

Figure 2.1a, which shows that no displacement constraint whatever has been imposed on motion of the spring element; that is, the spring is not connected to any

physical object that would prevent or limit motion of either node. With no constraint, it is not possible to solve for the nodal displacements individually.

Instead, only the difference in nodal displacements can be determined, as this

difference represents the elongation or contraction of the spring element owing

to elastic effects. As discussed in more detail in the general formulation of interpolation functions (Chapter 6) and structural dynamics (Chapter 10), a properly

formulated ?nite element must allow for constant value of the ?eld variable. In

the example at hand, this means rigid body motion. We can see the rigid body

motion capability in terms of a single spring (element) and in the context of several connected elements. For a single, unconstrained element, if arbitrary forces

are applied at each node, the spring not only deforms axially but also undergoes

acceleration according to Newton’s second law. Hence, there exists not only

deformation but overall motion. If, in a connected system of spring elements, the

overall system response is such that nodes 1 and 2 of a particular element displace the same amount, there is no elastic deformation of the spring and therefore

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

no elastic force in the spring. This physical situation must be included in the

element formulation. The capability is indicated mathematically by singularity

of the element stiffness matrix. As the stiffness matrix is formulated on the basis

of deformation of the element, we cannot expect to compute nodal displacements

if there is no deformation of the element.

Equation 2.7 indicates the mathematical operation of inverting the stiffness

matrix to obtain solutions. In the context of an individual element, the singular

nature of an element stiffness matrix precludes this operation, as the inverse of a

singular matrix does not exist. As is illustrated profusely in the remainder of the

text, the general solution of a ?nite element problem, in a global, as opposed to

element, context, involves the solution of equations of the form of Equation 2.5. For

realistic ?nite element models, which are of huge dimension in terms of the matrix

order (N ? N) involved, computing the inverse of the stiffness matrix is a very inef?cient, time-consuming operation, which should not be undertaken except for the

very simplest of systems. Other, more-ef?cient solution techniques are available,

and these are discussed subsequently. (Many of the end-of-chapter problems

included in this text are of small order and can be ef?ciently solved via matrix inversion using “spreadsheet” software functions or software such as MATLAB.)

2.2.1 System Assembly in Global Coordinates

Derivation of the element stiffness matrix for a spring element was based on

equilibrium conditions. The same procedure can be applied to a connected system of spring elements by writing the equilibrium equation for each node. However, rather than drawing free-body diagrams of each node and formally writing

the equilibrium equations, the nodal equilibrium equations can be obtained more

ef?ciently by considering the effect of each element separately and adding the

element force contribution to each nodal equation. The process is described as

“assembly,” as we take individual stiffness components and “put them together”

to obtain the system equations. To illustrate, via a simple example, the assembly

of element characteristics into global (or system) equations, we next consider the

system of two linear spring elements connected as shown in Figure 2.2.

For generality, it is assumed that the springs have different spring constants

k1 and k2. The nodes are numbered 1, 2, and 3 as shown, with the springs sharing

node 2 as the physical connection. Note that these are global node numbers. The

global nodal displacements are identi?ed as U1, U2, and U3, where the upper case

is used to indicate that the quantities represented are global or system displacements as opposed to individual element displacements. Similarly, applied nodal

U2

U1

1

F1

1

k1

U3

2

2

F2

k2

3

F3

Figure 2.2 System of two springs with node numbers,

element numbers, nodal displacements, and nodal forces.

23

Hutton: Fundamentals of

Finite Element Analysis

24

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

u1(1)

u(1)

2

u1(2)

1

f 1(1)

1

u(2)

2

2

2

k1

f (1)

2

f (2)

2

2

3

k2

(a)

f (2)

3

(b)

f 1(1)

F1

(2)

f (1)

2 f2

1

2

f (2)

3

F2

(d)

(c)

F3

(e)

Figure 2.3 Free-body diagrams of elements and nodes for the

two-element system of Figure 2.2.

forces are F1, F2, and F3. Assuming the system of two spring elements to be

in equilibrium, we examine free-body diagrams of the springs individually (Figure 2.3a and 2.3b) and express the equilibrium conditions for each spring, using

Equation 2.4, as

(1)

(1)

u1

f 1

k1

?k 1

=

(2.8a)

(1)

(1)

?k 1

k1

u2

f 2

(2)

(2)

u1

f 2

k2

?k 2

=

(2.8b)

(2)

(2)

?k 2

k2

u

f

2

3

where the superscript is element number.

To begin “assembling” the equilibrium equations describing the behavior

of the system of two springs, the displacement compatibility conditions, which

relate element displacements to system displacements, are written as

u

(1)

1

= U1

u

(1)

2

= U2

u

(2)

1

= U2

u

(2)

2

= U3

(2.9)

The compatibility conditions state the physical fact that the springs are connected at node 2, remain connected at node 2 after deformation, and hence, must

have the same nodal displacement at node 2. Thus, element-to-element displacement continuity is enforced at nodal connections. Substituting Equations 2.9 into

Equations 2.8, we obtain

(1)

f 1

k1

?k 1

U1

=

(2.10a)

(1)

?k 1

k1

U2

f 2

and

(2)

f 2

k2

?k 2

U2

=

(2.10b)

(2)

?k 2

k2

U3

f 3

Here, we use the notation f

node i.

( j)

i

to represent the force exerted on element j at

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

Equation 2.10 is the equilibrium equations for each spring element expressed

in terms of the speci?ed global displacements. In this form, the equations clearly

show that the elements are physically connected at node 2 and have the same displacement U2 at that node. These equations are not yet amenable to direct combination, as the displacement vectors are not the same. We expand both matrix

equations to 3 ? 3 as follows (while formally expressing the facts that element 1

is not connected to node 3 and element 2 is not connected to node 1):

?

?

? f (1) ?

U1

k1 ?k1 0

? 1 ?

U2 = f (1)

?k1 k1 0

(2.11)

?

?

? 2 ?

0

0

0

0

0

?

?

? 0 ?

0

0

0

0

?

?

(2)

0 k2 ?k2

U2 = f 2

(2.12)

?

?

? (2) ?

0 ?k2 k2

U3

f

3

The addition of Equations 2.11 and 2.12 yields

k1

?k1

0

?k1

k1 + k2

?k2

0

?k2

k2

U1

U2

U3

=

?

?

?

?

?

f (1)

1

(2)

f (1)

2 + f 2

f (2)

3

?

?

?

?

?

(2.13)

Next, we refer to the free-body diagrams of each of the three nodes depicted in

Figure 2.3c, 2.3d, and 2.3e. The equilibrium conditions for nodes 1, 2, and 3

show that

f

(1)

1

= F1

f

(1)

2

+ f

(2)

2

= F2

f

(2)

3

= F3

respectively. Substituting into Equation 2.13, we obtain the ?nal result:

?k 1

0

k1

U1

F1

?k 1 k 1 + k 2 ?k 2

U 2 = F2

0

?k 2

k2

U3

F3

(2.14)

(2.15)

which is of the form [K ]{U} = {F}, similar to Equation 2.5. However, Equation 2.15 represents the equations governing the system composed of two connected spring elements. By direct consideration of the equilibrium conditions,

we obtain the system stiffness matrix [K ] (note use of upper case) as

k1

?k 1

0

[K ] = ?k 1 k 1 + k 2 ?k 2

(2.16)

0

?k 2

k2

Note that the system stiffness matrix is (1) symmetric, as is the case with all linear systems referred to orthogonal coordinate systems; (2) singular, since no

constraints are applied to prevent rigid body motion of the system; and (3) the

system matrix is simply a superposition of the individual element stiffness

matrices with proper assignment of element nodal displacements and associated

stiffness coef?cients to system nodal displacements. The superposition procedure is formalized in the context of frame structures in the following paragraphs.

25

Hutton: Fundamentals of

Finite Element Analysis

26

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

EXAMPLE 2.1

Consider the two element system depicted in Figure 2.2 given that

Node 1 is attached to a ?xed support, yielding the displacement constraint U1 = 0.

k1 = 50 lb./in., k2 = 75 lb./in., F2 = F3 = 75 lb.

for these conditions determine nodal displacements U2 and U3.

? Solution

Substituting the speci?ed values into Equation 2.15 yields

?? ? ? ?

50 ?50

0

? 0 ? ? F1 ?

? ?50 125 ?75 ? U2 = 75

? ? ? ?

75

U3

0

?75 75

?

and we note that, owing to the constraint of zero displacement at node 1, nodal force F1

becomes an unknown reaction force. Formally, the ?rst algebraic equation represented in

this matrix equation becomes

?50U 2 = F1

and this is known as a constraint equation, as it represents the equilibrium condition

of a node at which the displacement is constrained. The second and third equations

become

125

?75

?75

75

U2

U3

=

75

75

which can be solved to obtain U2 = 3 in. and U3 = 4 in. Note that the matrix equations

governing the unknown displacements are obtained by simply striking out the ?rst row

and column of the 3 ? 3 matrix system, since the constrained displacement is zero.

Hence, the constraint does not affect the values of the active displacements (we use the

term active to refer to displacements that are unknown and must be computed). Substitution of the calculated values of U2 and U3 into the constraint equation yields the value

F1 = ?150 lb., which value is clearly in equilibrium with the applied nodal forces of

75 lb. each. We also illustrate element equilibrium by writing the equations for each

element as

50

?50

75

?75

(1)

f 1

?150

0

=

lb.

=

(1)

150

3

f 2

(2)

f 2

?75

?75

3

=

lb.

=

(2)

75

75

4

f 3

?50

50

for element 1

for element 2

Example 2.1 illustrates the general procedure for solution of ?nite element models: Formulate the system equilibrium equations, apply the speci?ed constraint

conditions, solve the reduced set of equations for the “active” displacements, and

substitute the computed displacements into the constraint equations to obtain the

unknown reactions. While not directly applicable for the spring element, for

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

27

more general ?nite element formulations, the computed displacements are also

substituted into the strain relations to obtain element strains, and the strains are,

in turn, substituted into the applicable stress-strain equations to obtain element

stress values.

EXAMPLE 2.2

Figure 2.4a depicts a system of three linearly elastic springs supporting three equal

weights W suspended in a vertical plane. Treating the springs as ?nite elements, determine the vertical displacement of each weight.

? Solution

To treat this as a ?nite element problem, we assign node and element numbers as shown

in Figure 2.4b and ignore, for the moment, that displacement U1 is known to be zero by

the ?xed support constraint. Per Equation 2.6, the stiffness matrix of each element is

(preprocessing)

k

(1)

=

k (2) =

k (3) =

3k

3k

?3k

?3k

3k

2k

?2k

?2k

2k

?k

k

k

?k

1

U1

3k

1

2k

2

k

3

W

2

U2

2k

W

3

U3

k

W

4

U4

(a)

(b)

Figure 2.4 Example 2.2: elastic

spring supporting weights.

Hutton: Fundamentals of

Finite Element Analysis

28

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

The element-to-global displacement relations are

u

(1)

1

= U1

u

(1)

2

=u

(2)

1

= U2

u

(2)

2

=u

(3)

1

= U3

u

(3)

2

= U4

Proceeding as in the previous example, we then write the individual element equations as

?? ? ?

?

U

3k ?3k 0 0 ?

? 1?

? ?

?

? ?3k 3k 0 0 ?? U2 ? ?

?

?

=

? 0

U ? ?

0

0 0 ??

?

? 3?

? ?

?

?

U4

0

0

0 0

?? ? ?

?

?

U

0

0

0

0 ?

? 1?

? ?

?

? 0 2k ?2k 0 ?? U2 ? ?

?

?

? 0 ?2k 2k 0 ?? U3 ? = ?

?

? ?

? ?

?

?

U4

0

0

0

0

?

?? ? ?

?

U

0 0 0

0 ?

? 1?

? ?

?

?0 0 0

?? U2 ? ?

0

?

?

=

? 0 0 k ?k ?? U3 ? ?

?

? ?

? ?

?

?

U4

0 0 ?k k

?

?

f (1)

1 ?

?

?

?

f (1)

2

?

0 ?

?

?

0

?

0 ?

(2) ?

?

f 1 ?

?

f (2)

2 ?

?

?

0

?

0 ?

?

?

0 ?

f (3)

1 ?

?

?

(3) ?

f 2

(1)

(2)

(3)

Adding Equations 1–3, we obtain

?? ? ? ?

U

F

3 ?3 0

0 ?

?

? ?

? 1?

? 1?

? ?3 5 ?2 0 ?? U2 ? ? W ?

?

k?

=

? 0 ?2 3 ?1 ?? U3 ? ? W ?

?

?

? ?

? ?

? ?

U4

W

0

0 ?1 1

?

(4)

where we utilize the fact that the sum of the element forces at each node must equal the

applied force at that node and, at node 1, the force is an unknown reaction.

Applying the displacement constraint U1 = 0 (this is also preprocessing), we obtain

?3kU 2 = F1

(5)

as the constraint equation and the matrix equation

?? ? ? ?

5 ?2 0 ? U2 ? ? W ?

k ? ?2 3 ?1 ? U3 = W

? ? ? ?

W

U4

0 ?1 1

?

(6)

for the active displacements. Again note that Equation 6 is obtained by eliminating the

constraint equation from 4 corresponding to the prescribed zero displacement.

Simultaneous solution (the solution step) of the algebraic equations represented by

Equation 6 yields the displacements as

U2 =

W

k

U3 =

2W

k

and Equation 5 gives the reaction force as

F1 = ?3W

(This is postprocessing.)

U4 =

3W

k

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

29

Note that the solution is exactly that which would be obtained by the usual statics

equations. Also note the general procedure as follows:

Formulate the individual element stiffness matrices.

Write the element to global displacement relations.

Assemble the global equilibrium equation in matrix form.

Reduce the matrix equations according to speci?ed constraints.

Solve the system of equations for the unknown nodal displacements (primary

variables).

Solve for the reaction forces (secondary variable) by back-substitution.

EXAMPLE 2.3

Figure 2.5 depicts a system of three linear spring elements connected as shown. The node

and element numbers are as indicated. Node 1 is ?xed to prevent motion, and node 3 is

given a speci?ed displacement as shown. Forces F2 = ? F and F4 = 2F are applied at

nodes 2 and 4. Determine the displacement of each node and the force required at node 3

for the speci?ed conditions.

? Solution

This example includes a nonhomogeneous boundary condition. In previous examples, the

boundary conditions were represented by zero displacements. In this example, we have

both a zero (homogeneous) and a speci?ed nonzero (nonhomogeneous) displacement

condition. The algebraic treatment must be different as follows. The system equilibrium

equations are expressed in matrix form (Problem 2.6) as

?

?k

4k

?3k

0

k

? ?k

?

? 0

0

0

?3k

5k

?2k

?

?? ? ? ? ?

0

U1 ?

F1 ?

F1 ?

?

?

?

?

?

?

?

?

?

? ? ? ? ?

?

0 ?

? U2 = F2 = ?F

?

?2k ?

U ? ?

F ? ?

F ?

?

? 3?

? ?

? ?

?

? 3?

? 3 ?

2k

U4

F4

2F

Substituting the speci?ed conditions U 1 = 0 and U 3 = results in the system of

equations

?

k

? ?k

?

? 0

0

?k

4k

?3k

0

k

?

?? ? ?

0

0 ?

F1 ?

?

?

?

?

?

?

? ? ?

?

0 ?

? U2 = ?F

?

? F3 ?

?2k ?

? ?

?

?

?

? ?

?

?

2k

2F

U4

F2 F

1

1

0

?3k

5k

?2k

3

3

2

3k

2

4

F4 2F

2k

?

Figure 2.5 Example 2.3: Three-element system with speci?ed

nonzero displacement at node 3.

Hutton: Fundamentals of

Finite Element Analysis

30

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Since U 1 = 0 , we remove the ?rst row and column to obtain

?

4k

? ?3k

0

?3k

5k

?2k

?

?? ? ?

0

? U2 ? ? ?F ?

=

F

?2k ?

? ? ? 3 ?

2F

U4

2k

as the system of equations governing displacements U2 and U4 and the unknown nodal

force F3. This last set of equations clearly shows that we cannot simply strike out the row

and column corresponding to the nonzero speci?ed displacement because it appears in

the equations governing the active displacements. To illustrate a general procedure, we

rewrite the last matrix equation as

?

5k

? ?3k

?2k

?3k

4k

0

?

?? ? ?

?2k ? ? ? F3 ?

0 ? U2 = ?F

?

? ? ?

U4

2F

2k

Next, we formally partition the stiffness matrix and write

?

5k

? ?3k

?2k

?3k

4k

0

?? ?

?2k ? ?

{}

{F }

[K ] [K U ]

=

0 ? U2 =

? ?

{FU }

[K U ] [K UU ] {U }

U4

2k

with

[K ] = [5k]

[K U ] = [?3k

?2k]

?3k

[K U ] = [K U ] T =

?2k

4k 0

[K UU ] =

0 2k

{} = {}

U2

{U } =

U4

{F } = {F3 }

?F

{FU } =

2F

From the second “row” of the partitioned matrix equations, we have

[K U ]{} + [K UU ]{U } = {FU }

and this can be solved for the unknown displacements to obtain

{U } = [K UU ]?1 ({F } ? [K U ]{})

provided that [K UU ]?1 exists. Since the constraints have been applied correctly, this

inverse does exist and is given by

?

[K UU ]?1

1

? 4k

=?

?

0

?

0 ?

?

1 ?

2k

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

Substituting, we obtain the unknown displacements as

?

?

?

3 ?

F

?

?

?

0 ?

?

?

+

?

?

4k

4 ?

? ?F + 3k

=

?

1 ? 2F + 2k

?

?

F

?

?

?

+ ?

?

?

2k

k

?

1

?

U2

? 4k

=?

{U } =

U4

?

0

The required force at node 3 is obtained by substitution of the displacement into the upper

partition to obtain

5

3

F3 = ? F + k

4

4

Finally, the reaction force at node 1 is

F

3

? k

4

4

F1 = ?kU 2 =

As a check on the results, we substitute the computed and prescribed displacements into

the individual element equations to insure that equilibrium is satis?ed.

Element 1

?k

k

k

?k

0

U2

=

?kU2

kU2

=

?

?

? f (1) ?

1

? f (1) ?

2

which shows that the nodal forces on element 1 are equal and opposite as required for

equilibrium.

Element 2

3k

?3k

?3k

3k

U2

U3

?

?

3

F

?3k ? ? + ?

4k

4

?

3k ?

?

?

3

3F

?

(2)

?

?

? k ?

?

??

f 2

4k

4

=

=

?

?

?

f (2)

?

? 3F + 3 k ?

3

4k

4

3k

=

?3k

which also veri?es equilibrium.

Element 3

2k

?2k

?2k

2k

U3

U4

=

2k

?2k

?2k

2k

Therefore element 3 is in equilibrium as well.

F

+

k

=

?2F

2F

=

f

f

(3)

3

(3)

4

2.3 ELASTIC BAR, SPAR/LINK/TRUSS ELEMENT

While the linear elastic spring serves to introduce the concept of the stiffness matrix, the usefulness of such an element in ?nite element analysis is rather limited.

Certainly, springs are used in machinery in many cases and the availability of a

?nite element representation of a linear spring is quite useful in such cases. The

31

Hutton: Fundamentals of

Finite Element Analysis

32

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

spring element is also often used to represent the elastic nature of supports for

more complicated systems. A more generally applicable, yet similar, element is

an elastic bar subjected to axial forces only. This element, which we simply call

a bar element, is particularly useful in the analysis of both two- and threedimensional frame or truss structures. Formulation of the ?nite element characteristics of an elastic bar element is based on the following assumptions:

1.

2.

3.

4.

The bar is geometrically straight.

The material obeys Hooke’s law.

Forces are applied only at the ends of the bar.

The bar supports axial loading only; bending, torsion, and shear are not

transmitted to the element via the nature of its connections to other

elements.

The last assumption, while quite restrictive, is not impractical; this condition is

satis?ed if the bar is connected to other structural members via pins (2-D) or balland-socket joints (3-D). Assumptions 1 and 4, in combination, show that this is

inherently a one-dimensional element, meaning that the elastic displacement of

any point along the bar can be expressed in terms of a single independent variable. As will be seen, however, the bar element can be used in modeling both

two- and three-dimensional structures. The reader will recognize this element

as the familiar two-force member of elementary statics, meaning, for equilibrium, the forces exerted on the ends of the element must be colinear, equal in

magnitude, and opposite in sense.

Figure 2.6 depicts an elastic bar of length L to which is af?xed a uniaxial

coordinate system x with its origin arbitrarily placed at the left end. This is the

element coordinate system or reference frame. Denoting axial displacement at

any position along the length of the bar as u(x), we de?ne nodes 1 and 2 at each

end as shown and introduce the nodal displacements u 1 = u(x = 0) and

u 2 = u(x = L ) . Thus, we have the continuous ?eld variable u(x), which is to be

expressed (approximately) in terms of two nodal variables u 1 and u 2 . To accomplish this discretization, we assume the existence of interpolation functions

N 1 (x ) and N 2 (x ) (also known as shape or blending functions) such that

u(x ) = N 1 (x )u 1 + N 2 (x )u 2

(2.17)

u1

u2

1

2

x

u(x)

L

Figure 2.6 A bar (or truss) element with element

coordinate system and nodal displacement

notation.

x

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

(It must be emphasized that, although an equality is indicated by Equation 2.17,

the relation, for ?nite elements in general, is an approximation. For the bar element, the relation, in fact, is exact.) To determine the interpolation functions, we

require that the boundary values of u(x ) (the nodal displacements) be identically

satis?ed by the discretization such that

u(x = 0) = u 1

u(x = L ) = u 2

(2.18)

Equations 2.17 and 2.18 lead to the following boundary (nodal) conditions:

N 1 (0) = 1

N 2 (0) = 0

(2.19)

N1( L ) = 0

N2( L ) = 1

(2.20)

which must be satis?ed by the interpolation functions. It is required that the displacement expression, Equation 2.17, satisfy the end (nodal) conditions identically, since the nodes will be the connection points between elements and the

displacement continuity conditions are enforced at those connections. As we

have two conditions that must be satis?ed by each of two one-dimensional functions, the simplest forms for the interpolation functions are polynomial forms:

N 1 (x ) = a0 + a1 x

(2.21)

N 2 (x ) = b0 + b1 x

(2.22)

where the polynomial coef?cients are to be determined via satisfaction of the

boundary (nodal) conditions. We note here that any number of mathematical

forms of the interpolation functions could be assumed while satisfying the

required conditions. The reasons for the linear form is explained in detail in

Chapter 6.

Application of conditions represented by Equation 2.19 yields a0 = 1,

b0 = 0 while Equation 2.20 results in a1 = ?(1/L ) and b1 = x /L . Therefore,

the interpolation functions are

N 1 (x ) = 1 ? x /L

(2.23)

N 2 (x ) = x /L

(2.24)

and the continuous displacement function is represented by the discretization

u(x ) = (1 ? x /L )u 1 + (x /L )u 2

(2.25)

As will be found most convenient subsequently, Equation 2.25 can be expressed

in matrix form as

u

u(x ) = [N 1 (x ) N 2 (x )] 1 = [N ] {u}

(2.26)

u2

where [N ] is the row matrix of interpolation functions and {u} is the column

matrix (vector) of nodal displacements.

Having expressed the displacement ?eld in terms of the nodal variables, the

task remains to determine the relation between the nodal displacements and

applied forces to obtain the stiffness matrix for the bar element. Recall from

33

Hutton: Fundamentals of

Finite Element Analysis

34

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

elementary strength of materials that the de?ection of an elastic bar of length L

and uniform cross-sectional area A when subjected to axial load P is given by

PL

=

(2.27)

AE

where E is the modulus of elasticity of the material. Using Equation 2.27, we

obtain the equivalent spring constant of an elastic bar as

P

AE

k=

=

(2.28)

L

and could, by analogy with the linear elastic spring, immediately write the stiffness matrix as Equation 2.6. While the result is exactly correct, we take a more

general approach to illustrate the procedures to be used with more complicated

element formulations.

Ultimately, we wish to compute the nodal displacements given some loading

condition on the element. To obtain the necessary equilibrium equations relating

the displacements to applied forces, we proceed from displacement to strain,

strain to stress, and stress to loading, as follows. In uniaxial loading, as in the bar

element, we need consider only the normal strain component, de?ned as

du

?x =

(2.29)

dx

which, when applied to Equation 2.25, gives

u2 ? u1

?x =

(2.30)

L

which shows that the spar element is a constant strain element. This is in accord

with strength of materials theory: The element has constant cross-sectional area

and is subjected to constant forces at the end points, so the strain does not vary

along the length. The axial stress, by Hooke’s law, is then

u2 ? u1

x = E ? x = E

(2.31)

L

and the associated axial force is

AE

P = x A =

(u 2 ? u 1 )

(2.32)

L

Taking care to observe the correct algebraic sign convention, Equation 2.32 is

now used to relate the applied nodal forces f 1 and f 2 to the nodal displacements

u 1 and u 2 . Observing that, if Equation 2.32 has a positive sign, the element is in

tension and nodal force f 2 must be in the positive coordinate direction while

nodal force f 1 must be equal and opposite for equilibrium; therefore,

AE

f1 = ?

(u 2 ? u 1 )

(2.33)

L

AE

f2 =

(u 2 ? u 1 )

(2.34)

L

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

Equations 2.33 and 2.34 are expressed in matrix form as

AE 1 ?1

u1

f1

=

?1

1

u

f2

L

2

35

(2.35)

Comparison of Equation 2.35 to Equation 2.4 shows that the stiffness matrix for

the bar element is given by

AE 1 ?1

[k e ] =

(2.36)

L ?1 1

As is the case with the linear spring, we observe that the element stiffness matrix

for the bar element is symmetric, singular, and of order 2 ? 2 in correspondence

with two nodal displacements or degrees of freedom. It must be emphasized that

the stiffness matrix given by Equation 2.36 is expressed in the element coordinate system, which in this case is one-dimensional. Application of this element

formulation to analysis of two- and three-dimensional structures is considered in

the next chapter.

EXAMPLE 2.4

Figure 2.7a depicts a tapered elastic bar subjected to an applied tensile load P at one end

and attached to a ?xed support at the other end. The cross-sectional area varies linearly

from A 0 at the ?xed support at x = 0 to A 0 /2 at x = L . Calculate the displacement of the

end of the bar (a) by modeling the bar as a single element having cross-sectional area

equal to the area of the actual bar at its midpoint along the length, (b) using two bar

elements of equal length and similarly evaluating the area at the midpoint of each, and

(c) using integration to obtain the exact solution.

? Solution

(a) For a single element, the cross-sectional area is 3A 0 /4 and the element “spring

constant” is

k=

AE

3A 0 E

=

L

4L

and the element equations are

3A 0 E

4L

1

?1

?1

?1

U1

U2

=

F1

P

The element and nodal displacements are as shown in Figure 2.7b. Applying the

constraint condition U 1 = 0 , we ?nd

U2 =

4PL

PL

= 1.333

3A 0 E

A0 E

as the displacement at x = L .

(b) Two elements of equal length L /2 with associated nodal displacements are

depicted in Figure 2.7c. For element 1, A 1 = 7A 0 /8 so

k1 =

7A 0 E

A1E

7A 0 E

=

=

L1

8( L /2)

4L

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

1

x

7

A

8 o

L

(

A(x) Ao 1

x

2L

A

)

3

A

4 o

5

A

8 o

2

P

(b)

(a)

(c)

1.4

Exact

One element

Two elements

1.2

P

Ao

x

?(x)

0.8

?Y

1.0

0.6

0.4

0.2

0

0

0.1

0.2

0.3

0.4

P

0.5

x

L

0.6

0.7

0.8

0.9

1.0

(e)

(d)

2

Exact

One element

Two elements

1.8

(

?

P

?o ?o Ao

)

1.6

1.4

1.2

1

0.8

0.6

0.4

0.2

0

0

0.1

0.2

0.3

0.4

0.5

x

L

0.6

0.7

0.8

0.9

1.0

(f)

Figure 2.7

(a) Tapered axial bar, (b) one-element model, (c) two-element model, (d) free-body diagram

for an exact solution, (e) displacement solutions, (f) stress solutions.

36

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

while for element 2, we have

5A 0

8

A1 =

A2 E

5A 0 E

5A 0 E

=

=

L2

8( L /2)

4L

and k 2 =

Since no load is applied at the center of the bar, the equilibrium equations for the

system of two elements is

?

?? ? ? ?

0 ? U1 ? ? F1 ?

?k2 ? U2 =

0

? ? ? ?

P

k2

U3

?k1

k1 + k2

?k2

k1

? ?k1

0

Applying the constraint condition U 1 = 0 results in

k1 + k2

?k 2

?k 2

k2

U2

U3

=

0

P

Adding the two equations gives

P

4PL

=

k1

7A 0 E

U2 =

and substituting this result into the ?rst equation results in

k1 + k2

48PL

PL

=

= 1.371

k2

35A 0 E

A0 E

U3 =

(c) To obtain the exact solution, we refer to Figure 2.7d, which is a free-body diagram of

a section of the bar between an arbitrary position x and the end x = L . For equilibrium,

x A = P

x

A = A(x ) = A 0 1 ?

2L

and since

the axial stress variation along the length of the bar is described by

x =

A0

P

x

1?

2L

Therefore, the axial strain is

?x =

x

=

E

E A0

P

x

1?

2L

Since the bar is ?xed at x = 0 , the displacement at x = L is given by

L

=

0

=

P

? x dx =

EA 0

L

0

dx

x

1?

2L

L

2 PL

2 PL

2 PL

PL

[?ln(2L ? x )]0 =

[ln(2L ) ? ln L ] =

ln 2 = 1.386

E A0

E A0

E A0

A0 E

37

Hutton: Fundamentals of

Finite Element Analysis

38

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Comparison of the results of parts b and c reveals that the two element solution

exhibits an error of only about 1 percent in comparison to the exact solution from

strength of materials theory. Figure 2.7e shows the displacement variation along the

length for the three solutions. It is extremely important to note, however, that the

computed axial stress for the ?nite element solutions varies signi?cantly from that of

the exact solution. The axial stress for the two-element solution is shown in Figure 2.7f, along with the calculated stress from the exact solution. Note particularly

the discontinuity of calculated stress values for the two elements at the connecting

node. This is typical of the derived, or secondary, variables, such as stress and strain,

as computed in the ?nite element method. As more and more smaller elements are

used in the model, the values of such discontinuities decrease, indicating solution

convergence. In structural analyses, the ?nite element user is most often more interested in stresses than displacements, hence it is essential that convergence of the

secondary variables be monitored.

2.4 STRAIN ENERGY, CASTIGLIANO’S

FIRST THEOREM

When external forces are applied to a body, the mechanical work done by those

forces is converted, in general, into a combination of kinetic and potential energies. In the case of an elastic body constrained to prevent motion, all the work

is stored in the body as elastic potential energy, which is also commonly

referred to as strain energy. Here, strain energy is denoted U e and mechanical

work W. From elementary statics, the mechanical work performed by a force F

as its point of application moves along a path from position 1 to position 2 is

de?ned as

2

r

W = F · d

(2.37)

1

where

d

r = dx i + dy j + dz k

(2.38)

is a differential vector along the path of motion. In Cartesian coordinates, work

is given by

x2

y2

z2

W = Fx dx + Fy dy + Fz dz

x1

y1

(2.39)

z1

where Fx , Fy , and Fz are the Cartesian components of the force vector.

For linearly elastic deformations, de?ection is directly proportional to applied force as, for example, depicted in Figure 2.8 for a linear spring. The slope

of the force-de?ection line is the spring constant such that F = k. Therefore,

the work required to deform such a spring by an arbitrary amount 0 from its

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

Force, F

2.4 Strain Energy, Castigliano’s First Theorem

k

1

Deflection, ?

Figure 2.8 Force-de?ection

relation for a linear elastic

spring.

free length is

0

0

1

W = F d = k d = k02 = U e

2

0

(2.40)

0

and we observe that the work and resulting elastic potential energy are quadratic

functions of displacement and have the units of force-length. This is a general

result for linearly elastic systems, as will be seen in many examples throughout

this text.

Utilizing Equation 2.28, the strain energy for an axially loaded elastic bar

?xed at one end can immediately be written as

Ue =

1 2

1 AE 2

k =

2

2 L

(2.41)

However, for a more general purpose, this result is converted to a different form

(applicable to a bar element only) as follows:

2

1 2

1 AE PL

1 P

P

1

U e = k =

=

AL = ?V

(2.42)

2

2 L

AE

2 A

AE

2

where V is the total volume of deformed material and the quantity 12 ? is strain

energy per unit volume, also known as strain energy density. In Equation 2.42,

stress and strain values are those corresponding to the ?nal value of applied

force. The factor 12 arises from the linear relation between stress and strain as the

load is applied from zero to the ?nal value P. In general, for uniaxial loading, the

strain energy per unit volume u e is de?ned by

?

u e = d?

(2.43)

0

which is extended to more general states of stress in subsequent chapters. We note

that Equation 2.43 represents the area under the elastic stress-strain diagram.

39

Hutton: Fundamentals of

Finite Element Analysis

40

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Presently, we will use the work-strain energy relation to obtain the governing equations for the bar element using the following theorem.

Castigliano’s First Theorem [1]

For an elastic system in equilibrium, the partial derivative of total strain energy

with respect to de?ection at a point is equal to the applied force in the direction

of the de?ection at that point.

Consider an elastic body subjected to N forces Fj for which the total strain

energy is expressed as

j

N

Ue = W =

Fj dj

(2.44)

j=1

0

where j is the de?ection at the point of application of force Fj in the direction of

the line of action of the force. If all points of load application are ?xed except

one, say, i, and that point is made to de?ect an in?nitesimal amount i by an

incremental in?nitesimal force Fi , the change in strain energy is

i

U e = W = Fi i + Fi di

(2.45)

0

where it is assumed that the original force Fi is constant during the in?nitesimal

change. The integral term in Equation 2.45 involves the product of in?nitesimal

quantities and can be neglected to obtain

U e

= Fi

i

(2.46)

which in the limit as i approaches zero becomes

?U

= Fi

? i

(2.47)

The ?rst theorem of Castigliano is a powerful tool for ?nite element formulation, as is now illustrated for the bar element. Combining Equations 2.30, 2.31,

and 2.43, total strain energy for the bar element is given by

2

1

1

u2 ? u1

U e = x ? x V = E

AL

(2.48)

2

2

L

Applying Castigliano’s theorem with respect to each displacement yields

?U e

AE

=

(u 1 ? u 2 ) = f 1

? u1

L

(2.49)

?U e

AE

=

(u 2 ? u 1 ) = f 2

? u2

L

(2.50)

which are observed to be identical to Equations 2.33 and 2.34.

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.4 Strain Energy, Castigliano’s First Theorem

41

The ?rst theorem of Castigliano is also applicable to rotational displacements. In the case of rotation, the partial derivative of strain energy with respect

to a rotational displacement is equal to the moment/torque applied at the point of

concern in the sense of the rotation. The following example illustrates the application in terms of a simple torsional member.

EXAMPLE 2.5

A solid circular shaft of radius R and length L is subjected to constant torque T. The shaft

is ?xed at one end, as shown in Figure 2.9. Formulate the elastic strain energy in terms of

the angle of twist at x = L and show that Castigliano’s ?rst theorem gives the correct

expression for the applied torque.

? Solution

From strength of materials theory, the shear stress at any cross section along the length of

the member is given by

=

Tr

J

where r is radial distance from the axis of the member and J is polar moment of inertia of

the cross section. For elastic behavior, we have

Tr

=

G

JG

=

where G is the shear modulus of the material, and the strain energy is then

1

Ue =

2

1

dV =

2

V

=

?

?

L

Tr

Tr

?

dA?dx

J

JG

A

0

T2

2J 2 G

L

r 2 dA dx =

T2L

2JG

0 A

where we have used the de?nition of the polar moment of inertia

J =

r2 dA

A

R

L

T

Figure 2.9 Example 2.5:

Circular cylinder subjected to

torsion.

Hutton: Fundamentals of

Finite Element Analysis

42

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Again invoking the strength of materials results, the angle of twist at the end of the

Finite Element Analysis

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

PREFACE

F

undamentals of Finite Element Analysis is intended to be the text for a

senior-level ?nite element course in engineering programs. The most

appropriate major programs are civil engineering, engineering mechanics, and mechanical engineering. The ?nite element method is such a widely used

analysis-and-design technique that it is essential that undergraduate engineering

students have a basic knowledge of the theory and applications of the technique.

Toward that objective, I developed and taught an undergraduate “special topics”

course on the ?nite element method at Washington State University in the summer of 1992. The course was composed of approximately two-thirds theory and

one-third use of commercial software in solving ?nite element problems. Since

that time, the course has become a regularly offered technical elective in the

mechanical engineering program and is generally in high demand. During

the developmental process for the course, I was never satis?ed with any text that

was used, and we tried many. I found the available texts to be at one extreme or

the other; namely, essentially no theory and all software application, or all theory

and no software application. The former approach, in my opinion, represents

training in using computer programs, while the latter represents graduate-level

study. I have written this text to seek a middle ground.

Pedagogically, I believe that training undergraduate engineering students to

use a particular software package without providing knowledge of the underlying

theory is a disservice to the student and can be dangerous for their future employers. While I am acutely aware that most engineering programs have a speci?c

?nite element software package available for student use, I do not believe that the

text the students use should be tied only to that software. Therefore, I have written this text to be software-independent. I emphasize the basic theory of the ?nite

element method, in a context that can be understood by undergraduate engineering students, and leave the software-speci?c portions to the instructor.

As the text is intended for an undergraduate course, the prerequisites required

are statics, dynamics, mechanics of materials, and calculus through ordinary differential equations. Of necessity, partial differential equations are introduced

but in a manner that should be understood based on the stated prerequisites.

Applications of the ?nite element method to heat transfer and ?uid mechanics are

included, but the necessary derivations are such that previous coursework in

those topics is not required. Many students will have taken heat transfer and ?uid

mechanics courses, and the instructor can expand the topics based on the students’ background.

Chapter 1 is a general introduction to the ?nite element method and includes a description of the basic concept of dividing a domain into ?nite-size

subdomains. The ?nite difference method is introduced for comparison to the

xi

Hutton: Fundamentals of

Finite Element Analysis

xii

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

?nite element method. A general procedure in the sequence of model de?nition,

solution, and interpretation of results is discussed and related to the generally

accepted terms of preprocessing, solution, and postprocessing. A brief history of

the ?nite element method is included, as are a few examples illustrating application of the method.

Chapter 2 introduces the concept of a ?nite element stiffness matrix and

associated displacement equation, in terms of interpolation functions, using the

linear spring as a ?nite element. The linear spring is known to most undergraduate students so the mechanics should not be new. However, representation of

the spring as a ?nite element is new but provides a simple, concise example of

the ?nite element method. The premise of spring element formulation is extended to the bar element, and energy methods are introduced. The ?rst theorem

of Castigliano is applied, as is the principle of minimum potential energy.

Castigliano’s theorem is a simple method to introduce the undergraduate student

to minimum principles without use of variational calculus.

Chapter 3 uses the bar element of Chapter 2 to illustrate assembly of global

equilibrium equations for a structure composed of many ?nite elements. Transformation from element coordinates to global coordinates is developed and

illustrated with both two- and three-dimensional examples. The direct stiffness

method is utilized and two methods for global matrix assembly are presented.

Application of boundary conditions and solution of the resultant constraint equations is discussed. Use of the basic displacement solution to obtain element strain

and stress is shown as a postprocessing operation.

Chapter 4 introduces the beam/?exure element as a bridge to continuity

requirements for higher-order elements. Slope continuity is introduced and this

requires an adjustment to the assumed interpolation functions to insure continuity.

Nodal load vectors are discussed in the context of discrete and distributed loads,

using the method of work equivalence.

Chapters 2, 3, and 4 introduce the basic procedures of ?nite-element modeling in the context of simple structural elements that should be well-known to the

student from the prerequisite mechanics of materials course. Thus the emphasis

in the early part of the course in which the text is used can be on the ?nite element method without introduction of new physical concepts. The bar and beam

elements can be used to give the student practical truss and frame problems for

solution using available ?nite element software. If the instructor is so inclined,

the bar and beam elements (in the two-dimensional context) also provide a relatively simple framework for student development of ?nite element software

using basic programming languages.

Chapter 5 is the springboard to more advanced concepts of ?nite element

analysis. The method of weighted residuals is introduced as the fundamental

technique used in the remainder of the text. The Galerkin method is utilized

exclusively since I have found this method is both understandable for undergraduate students and is amenable to a wide range of engineering problems. The

material in this chapter repeats the bar and beam developments and extends the

?nite element concept to one-dimensional heat transfer. Application to the bar

Hutton: Fundamentals of

Finite Element Analysis

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

and beam elements illustrates that the method is in agreement with the basic mechanics approach of Chapters 2–4. Introduction of heat transfer exposes the student to additional applications of the ?nite element method that are, most likely,

new to the student.

Chapter 6 is a stand-alone description of the requirements of interpolation

functions used in developing ?nite element models for any physical problem.

Continuity and completeness requirements are delineated. Natural (serendipity)

coordinates, triangular coordinates, and volume coordinates are de?ned and used

to develop interpolation functions for several element types in two- and threedimensions. The concept of isoparametric mapping is introduced in the context of

the plane quadrilateral element. As a precursor to following chapters, numerical

integration using Gaussian quadrature is covered and several examples included.

The use of two-dimensional elements to model three-dimensional axisymmetric

problems is included.

Chapter 7 uses Galerkin’s ?nite element method to develop the ?nite element equations for several commonly encountered situations in heat transfer.

One-, two- and three-dimensional formulations are discussed for conduction and

convection. Radiation is not included, as that phenomenon introduces a nonlinearity that undergraduate students are not prepared to deal with at the intended

level of the text. Heat transfer with mass transport is included. The ?nite difference method in conjunction with the ?nite element method is utilized to present

methods of solving time-dependent heat transfer problems.

Chapter 8 introduces ?nite element applications to ?uid mechanics. The

general equations governing ?uid ?ow are so complex and nonlinear that the

topic is introduced via ideal ?ow. The stream function and velocity potential

function are illustrated and the applicable restrictions noted. Example problems

are included that note the analogy with heat transfer and use heat transfer ?nite

element solutions to solve ideal ?ow problems. A brief discussion of viscous

?ow shows the nonlinearities that arise when nonideal ?ows are considered.

Chapter 9 applies the ?nite element method to problems in solid mechanics

with the proviso that the material response is linearly elastic and small de?ection.

Both plane stress and plane strain are de?ned and the ?nite element formulations

developed for each case. General three-dimensional states of stress and axisymmetric stress are included. A model for torsion of noncircular sections is developed using the Prandtl stress function. The purpose of the torsion section is to

make the student aware that all torsionally loaded objects are not circular and the

analysis methods must be adjusted to suit geometry.

Chapter 10 introduces the concept of dynamic motion of structures. It is not

presumed that the student has taken a course in mechanical vibrations; as a result, this chapter includes a primer on basic vibration theory. Most of this material is drawn from my previously published text Applied Mechanical Vibrations.

The concept of the mass or inertia matrix is developed by examples of simple

spring-mass systems and then extended to continuous bodies. Both lumped and

consistent mass matrices are de?ned and used in examples. Modal analysis is the

basic method espoused for dynamic response; hence, a considerable amount of

xiii

Hutton: Fundamentals of

Finite Element Analysis

xiv

Front Matter

Preface

© The McGraw?Hill

Companies, 2004

Preface

text material is devoted to determination of natural modes, orthogonality, and

modal superposition. Combination of ?nite difference and ?nite element methods for solving transient dynamic structural problems is included.

The appendices are included in order to provide the student with material

that might be new or may be “rusty” in the student’s mind.

Appendix A is a review of matrix algebra and should be known to the student from a course in linear algebra.

Appendix B states the general three-dimensional constitutive relations for

a homogeneous, isotropic, elastic material. I have found over the years that undergraduate engineering students do not have a ?rm grasp of these relations. In

general, the student has been exposed to so many special cases that the threedimensional equations are not truly understood.

Appendix C covers three methods for solving linear algebraic equations.

Some students may use this material as an outline for programming solution

methods. I include the appendix only so the reader is aware of the algorithms underlying the software he/she will use in solving ?nite element problems.

Appendix D describes the basic computational capabilities of the FEPC

software. The FEPC (FEP?nite element program for the PCpersonal computer)

was developed by the late Dr. Charles Knight of Virginia Polytechnic Institute

and State University and is used in conjunction with this text with permission of

his estate. Dr. Knight’s programs allow analysis of two-dimensional programs

using bar, beam, and plane stress elements. The appendix describes in general

terms the capabilities and limitations of the software. The FEPC program is

available to the student at www.mhhe.com/hutton.

Appendix E includes problems for several chapters of the text that should be

solved via commercial ?nite element software. Whether the instructor has available ANSYS, ALGOR, COSMOS, etc., these problems are oriented to systems

having many degrees of freedom and not amenable to hand calculation. Additional problems of this sort will be added to the website on a continuing basis.

The textbook features a Web site (www.mhhe.com/hutton) with ?nite element analysis links and the FEPC program. At this site, instructors will have

access to PowerPoint images and an Instructors’ Solutions Manual. Instructors

can access these tools by contacting their local McGraw-Hill sales representative

for password information.

I thank Raghu Agarwal, Rong Y. Chen, Nels Madsen, Robert L. Rankin,

Joseph J. Rencis, Stephen R. Swanson, and Lonny L. Thompson, who reviewed

some or all of the manuscript and provided constructive suggestions and criticisms that have helped improve the book.

I am grateful to all the staff at McGraw-Hill who have labored to make this

project a reality. I especially acknowledge the patient encouragement and professionalism of Jonathan Plant, Senior Editor, Lisa Kalner Williams, Developmental Editor, and Kay Brimeyer, Senior Project Manager.

David V. Hutton

Pullman, WA

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

C H A P T E R

1

Basic Concepts of the

Finite Element Method

1.1 INTRODUCTION

The ?nite element method (FEM), sometimes referred to as ?nite element

analysis (FEA), is a computational technique used to obtain approximate solutions of boundary value problems in engineering. Simply stated, a boundary

value problem is a mathematical problem in which one or more dependent variables must satisfy a differential equation everywhere within a known domain of

independent variables and satisfy speci?c conditions on the boundary of the

domain. Boundary value problems are also sometimes called ?eld problems. The

?eld is the domain of interest and most often represents a physical structure.

The ?eld variables are the dependent variables of interest governed by the differential equation. The boundary conditions are the speci?ed values of the ?eld

variables (or related variables such as derivatives) on the boundaries of the ?eld.

Depending on the type of physical problem being analyzed, the ?eld variables

may include physical displacement, temperature, heat ?ux, and ?uid velocity to

name only a few.

1.2 HOW DOES THE FINITE ELEMENT

METHOD WORK?

The general techniques and terminology of ?nite element analysis will be introduced with reference to Figure 1.1. The ?gure depicts a volume of some material

or materials having known physical properties. The volume represents the

domain of a boundary value problem to be solved. For simplicity, at this point,

we assume a two-dimensional case with a single ?eld variable (x, y) to be

determined at every point P(x, y) such that a known governing equation (or equations) is satis?ed exactly at every such point. Note that this implies an exact

1

Hutton: Fundamentals of

Finite Element Analysis

2

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

3

P(x, y)

2

1

(a)

(b)

(c)

Figure 1.1

(a) A general two-dimensional domain of ?eld variable (x, y).

(b) A three-node ?nite element de?ned in the domain. (c) Additional

elements showing a partial ?nite element mesh of the domain.

mathematical solution is obtained; that is, the solution is a closed-form algebraic

expression of the independent variables. In practical problems, the domain may

be geometrically complex as is, often, the governing equation and the likelihood

of obtaining an exact closed-form solution is very low. Therefore, approximate

solutions based on numerical techniques and digital computation are most

often obtained in engineering analyses of complex problems. Finite element

analysis is a powerful technique for obtaining such approximate solutions with

good accuracy.

A small triangular element that encloses a ?nite-sized subdomain of the area

of interest is shown in Figure 1.1b. That this element is not a differential element

of size dx ? dy makes this a ?nite element. As we treat this example as a twodimensional problem, it is assumed that the thickness in the z direction is constant and z dependency is not indicated in the differential equation. The vertices

of the triangular element are numbered to indicate that these points are nodes. A

node is a speci?c point in the ?nite element at which the value of the ?eld variable is to be explicitly calculated. Exterior nodes are located on the boundaries

of the ?nite element and may be used to connect an element to adjacent ?nite

elements. Nodes that do not lie on element boundaries are interior nodes and

cannot be connected to any other element. The triangular element of Figure 1.1b

has only exterior nodes.

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

If the values of the ?eld variable are computed only at nodes, how are values

obtained at other points within a ?nite element? The answer contains the crux of

the ?nite element method: The values of the ?eld variable computed at the nodes

are used to approximate the values at nonnodal points (that is, in the element

interior) by interpolation of the nodal values. For the three-node triangle example, the nodes are all exterior and, at any other point within the element, the ?eld

variable is described by the approximate relation

(x , y) = N 1 (x , y)1 + N 2 (x , y)2 + N 3 (x , y)3

(1.1)

where 1 , 2 , and 3 are the values of the ?eld variable at the nodes, and N1, N2,

and N3 are the interpolation functions, also known as shape functions or blending functions. In the ?nite element approach, the nodal values of the ?eld variable are treated as unknown constants that are to be determined. The interpolation functions are most often polynomial forms of the independent variables,

derived to satisfy certain required conditions at the nodes. These conditions are

discussed in detail in subsequent chapters. The major point to be made here is

that the interpolation functions are predetermined, known functions of the independent variables; and these functions describe the variation of the ?eld variable

within the ?nite element.

The triangular element described by Equation 1.1 is said to have 3 degrees

of freedom, as three nodal values of the ?eld variable are required to describe

the ?eld variable everywhere in the element. This would be the case if the ?eld

variable represents a scalar ?eld, such as temperature in a heat transfer problem

(Chapter 7). If the domain of Figure 1.1 represents a thin, solid body subjected to

plane stress (Chapter 9), the ?eld variable becomes the displacement vector and

the values of two components must be computed at each node. In the latter case,

the three-node triangular element has 6 degrees of freedom. In general, the number of degrees of freedom associated with a ?nite element is equal to the product

of the number of nodes and the number of values of the ?eld variable (and possibly its derivatives) that must be computed at each node.

How does this element-based approach work over the entire domain of interest? As depicted in Figure 1.1c, every element is connected at its exterior

nodes to other elements. The ?nite element equations are formulated such that, at

the nodal connections, the value of the ?eld variable at any connection is the

same for each element connected to the node. Thus, continuity of the ?eld variable at the nodes is ensured. In fact, ?nite element formulations are such that

continuity of the ?eld variable across interelement boundaries is also ensured.

This feature avoids the physically unacceptable possibility of gaps or voids occurring in the domain. In structural problems, such gaps would represent physical separation of the material. In heat transfer, a “gap” would manifest itself in

the form of different temperatures at the same physical point.

Although continuity of the ?eld variable from element to element is inherent

to the ?nite element formulation, interelement continuity of gradients (i.e., derivatives) of the ?eld variable does not generally exist. This is a critical observation. In most cases, such derivatives are of more interest than are ?eld variable

values. For example, in structural problems, the ?eld variable is displacement but

3

Hutton: Fundamentals of

Finite Element Analysis

4

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

the true interest is more often in strain and stress. As strain is de?ned in terms of

?rst derivatives of displacement components, strain is not continuous across element boundaries. However, the magnitudes of discontinuities of derivatives can

be used to assess solution accuracy and convergence as the number of elements

is increased, as is illustrated by the following example.

1.2.1 Comparison of Finite Element and Exact Solutions

The process of representing a physical domain with ?nite elements is referred to

as meshing, and the resulting set of elements is known as the ?nite element mesh.

As most of the commonly used element geometries have straight sides, it is generally impossible to include the entire physical domain in the element mesh if the

domain includes curved boundaries. Such a situation is shown in Figure 1.2a,

where a curved-boundary domain is meshed (quite coarsely) using square elements. A re?ned mesh for the same domain is shown in Figure 1.2b, using

smaller, more numerous elements of the same type. Note that the re?ned mesh

includes signi?cantly more of the physical domain in the ?nite element representation and the curved boundaries are more closely approximated. (Triangular

elements could approximate the boundaries even better.)

If the interpolation functions satisfy certain mathematical requirements

(Chapter 6), a ?nite element solution for a particular problem converges to the

exact solution of the problem. That is, as the number of elements is increased and

the physical dimensions of the elements are decreased, the ?nite element solution

changes incrementally. The incremental changes decrease with the mesh re?nement process and approach the exact solution asymptotically. To illustrate

convergence, we consider a relatively simple problem that has a known solution.

Figure 1.3a depicts a tapered, solid cylinder ?xed at one end and subjected to

a tensile load at the other end. Assuming the displacement at the point of load

application to be of interest, a ?rst approximation is obtained by considering

the cylinder to be uniform, having a cross-sectional area equal to the average area

(a)

(b)

Figure 1.2

(a) Arbitrary curved-boundary domain modeled using square elements. Stippled

areas are not included in the model. A total of 41 elements is shown. (b) Re?ned

?nite element mesh showing reduction of the area not included in the model. A

total of 192 elements is shown.

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

of the cylinder (Figure 1.3b). The uniform bar is a link or bar ?nite element

(Chapter 2), so our ?rst approximation is a one-element, ?nite element model.

The solution is obtained using the strength of materials theory. Next, we model

the tapered cylinder as two uniform bars in series, as in Figure 1.3c. In the twoelement model, each element is of length equal to half the total length of the

cylinder and has a cross-sectional area equal to the average area of the corresponding half-length of the cylinder. The mesh re?nement is continued using a

four-element model, as in Figure 1.3d, and so on. For this simple problem, the

displacement of the end of the cylinder for each of the ?nite element models is as

shown in Figure 1.4a, where the dashed line represents the known solution. Convergence of the ?nite element solutions to the exact solution is clearly indicated.

ro

x

r

L

A

Ao AL

2

rL

F

(a)

(b)

Element 1

Element 2

(c)

(d)

Figure 1.3

(a) Tapered circular cylinder subjected to tensile loading:

r(x) r0 (x/L)(r0 rL). (b) Tapered cylinder as a single axial

(bar) element using an average area. Actual tapered cylinder

is shown as dashed lines. (c) Tapered cylinder modeled as

two, equal-length, ?nite elements. The area of each element

is average over the respective tapered cylinder length.

(d) Tapered circular cylinder modeled as four, equal-length

?nite elements. The areas are average over the respective

length of cylinder (element length L?4).

5

Hutton: Fundamentals of

Finite Element Analysis

6

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

Exact

Four elements

?

( Lx )

?(x L)

Exact

1

2

3

Number of elements

(a)

4

0.25

0.5

0.75

1.0

x

L

(b)

Figure 1.4

(a) Displacement at x L for tapered cylinder in tension of Figure 1.3. (b) Comparison of the exact solution

and the four-element solution for a tapered cylinder in tension.

On the other hand, if we plot displacement as a function of position along the

length of the cylinder, we can observe convergence as well as the approximate

nature of the ?nite element solutions. Figure 1.4b depicts the exact strength of

materials solution and the displacement solution for the four-element models.

We note that the displacement variation in each element is a linear approximation

to the true nonlinear solution. The linear variation is directly attributable to the

fact that the interpolation functions for a bar element are linear. Second, we note

that, as the mesh is re?ned, the displacement solution converges to the nonlinear

solution at every point in the solution domain.

The previous paragraph discussed convergence of the displacement of the

tapered cylinder. As will be seen in Chapter 2, displacement is the primary ?eld

variable in structural problems. In most structural problems, however, we are

interested primarily in stresses induced by speci?ed loadings. The stresses must

be computed via the appropriate stress-strain relations, and the strain components are derived from the displacement ?eld solution. Hence, strains and

stresses are referred to as derived variables. For example, if we plot the element

stresses for the tapered cylinder example just cited for the exact solution as well

as the ?nite element solutions for two- and four-element models as depicted in

Figure 1.5, we observe that the stresses are constant in each element and represent a discontinuous solution of the problem in terms of stresses and strains. We

also note that, as the number of elements increases, the jump discontinuities in

stress decrease in magnitude. This phenomenon is characteristic of the ?nite element method. The formulation of the ?nite element method for a given problem

is such that the primary ?eld variable is continuous from element to element but

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

4.5

4.0

Exact

Two elements

Four elements

3.5

?

?0

3.0

2.5

2.0

1.5

1.0

0.25

0.5

x

L

0.75

1.0

Figure 1.5

Comparison of the computed axial stress value in a

tapered cylinder: 0 F?A0.

the derived variables are not necessarily continuous. In the limiting process of

mesh re?nement, the derived variables become closer and closer to continuity.

Our example shows how the ?nite element solution converges to a known

exact solution (the exactness of the solution in this case is that of strength of

materials theory). If we know the exact solution, we would not be applying the

?nite element method! So how do we assess the accuracy of a ?nite element solution for a problem with an unknown solution? The answer to this question is not

simple. If we did not have the dashed line in Figure 1.3 representing the exact

solution, we could still discern convergence to a solution. Convergence of a

numerical method (such as the ?nite element method) is by no means assurance

that the convergence is to the correct solution. A person using the ?nite element

analysis technique must examine the solution analytically in terms of (1) numerical convergence, (2) reasonableness (does the result make sense?), (3) whether the

physical laws of the problem are satis?ed (is the structure in equilibrium? Does the

heat output balance with the heat input?), and (4) whether the discontinuities in

value of derived variables across element boundaries are reasonable. Many

such questions must be posed and examined prior to accepting the results of a ?nite

element analysis as representative of a correct solution useful for design purposes.

1.2.2 Comparison of Finite Element and Finite

Difference Methods

The ?nite difference method is another numerical technique frequently used to

obtain approximate solutions of problems governed by differential equations.

Details of the technique are discussed in Chapter 7 in the context of transient heat

7

Hutton: Fundamentals of

Finite Element Analysis

8

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

transfer. The method is also illustrated in Chapter 10 for transient dynamic analysis of structures. Here, we present the basic concepts of the ?nite difference

method for purposes of comparison.

The ?nite difference method is based on the de?nition of the derivative of a

function f (x ) that is

d f (x )

f (x + x ) ? f (x )

= lim

x?0

dx

x

(1.2)

where x is the independent variable. In the ?nite difference method, as implied

by its name, derivatives are calculated via Equation 1.2 using small, but ?nite,

values of x to obtain

d f (x )

f (x + x ) ? f (x )

?

dx

x

(1.3)

A differential equation such as

df

+x =0

dx

0?x ?1

(1.4)

f (x + x ) ? f (x )

+x =0

x

(1.5)

is expressed as

in the ?nite difference method. Equation 1.5 can be rewritten as

f (x + x ) = f (x ) ? x (x )

(1.6)

where we note that the equality must be taken as “approximately equals.” From

differential equation theory, we know that the solution of a ?rst-order differential

equation contains one constant of integration. The constant of integration must

be determined such that one given condition (a boundary condition or initial condition) is satis?ed. In the current example, we assume that the speci?ed condition

is x (0) = A = constant. If we choose an integration step x to be a small, constant value (the integration step is not required to be constant), then we can write

x i+1 = x i + x

i = 0, N

(1.7)

where N is the total number of steps required to cover the domain. Equation 1.6

is then

f i+1 = f i ? x i (x )

f0 = A

i = 0, N

(1.8)

Equation 1.8 is known as a recurrence relation and provides an approximation to

the value of the unknown function f (x) at a number of discrete points in the domain of the problem.

To illustrate, Figure 1.6a shows the exact solution f (x ) = 1 ? x 2 /2 and a

?nite difference solution obtained with x = 0.1. The ?nite difference solution is

shown at the discrete points of function evaluation only. The manner of variation

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.2 How Does the Finite Element Method Work?

1

0.8

f (x)

0.6

0.4

0.2

0

0

0.2

0.6

0.4

0.8

1

x

Figure 1.6

Comparison of the exact and ?nite difference

solutions of Equation 1.4 with f0 A 1.

of the function between the calculated points is not known in the ?nite difference

method. One can, of course, linearly interpolate the values to produce an approximation to the curve of the exact solution but the manner of interpolation is

not an a priori determination in the ?nite difference method.

To contrast the ?nite difference method with the ?nite element method,

we note that, in the ?nite element method, the variation of the ?eld variable in

the physical domain is an integral part of the procedure. That is, based on the

selected interpolation functions, the variation of the ?eld variable throughout a

?nite element is speci?ed as an integral part of the problem formulation. In the

?nite difference method, this is not the case: The ?eld variable is computed at

speci?ed points only. The major rami?cation of this contrast is that derivatives

(to a certain level) can be computed in the ?nite element approach, whereas the

?nite difference method provides data only on the variable itself. In a structural

problem, for example, both methods provide displacement solutions, but the

?nite element solution can be used to directly compute strain components (?rst

derivatives). To obtain strain data in the ?nite difference method requires additional considerations not inherent to the mathematical model.

There are also certain similarities between the two methods. The integration

points in the ?nite difference method are analogous to the nodes in a ?nite

element model. The variable of interest is explicitly evaluated at such points.

Also, as the integration step (step size) in the ?nite difference method is reduced,

the solution is expected to converge to the exact solution. This is similar to the

expected convergence of a ?nite element solution as the mesh of elements is

re?ned. In both cases, the re?nement represents reduction of the mathematical

model from ?nite to in?nitesimal. And in both cases, differential equations are

reduced to algebraic equations.

9

Hutton: Fundamentals of

Finite Element Analysis

10

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

Probably the most descriptive way to contrast the two methods is to note that

the ?nite difference method models the differential equation(s) of the problem

and uses numerical integration to obtain the solution at discrete points. The ?nite

element method models the entire domain of the problem and uses known physical principles to develop algebraic equations describing the approximate solutions. Thus, the ?nite difference method models differential equations while the

?nite element method can be said to more closely model the physical problem at

hand. As will be observed in the remainder of this text, there are cases in which

a combination of ?nite element and ?nite difference methods is very useful and

ef?cient in obtaining solutions to engineering problems, particularly where dynamic (transient) effects are important.

1.3 A GENERAL PROCEDURE FOR FINITE

ELEMENT ANALYSIS

Certain steps in formulating a ?nite element analysis of a physical problem are

common to all such analyses, whether structural, heat transfer, ?uid ?ow, or

some other problem. These steps are embodied in commercial ?nite element

software packages (some are mentioned in the following paragraphs) and are

implicitly incorporated in this text, although we do not necessarily refer to the

steps explicitly in the following chapters. The steps are described as follows.

1.3.1 Preprocessing

The preprocessing step is, quite generally, described as de?ning the model and

includes

De?ne the geometric domain of the problem.

De?ne the element type(s) to be used (Chapter 6).

De?ne the material properties of the elements.

De?ne the geometric properties of the elements (length, area, and the like).

De?ne the element connectivities (mesh the model).

De?ne the physical constraints (boundary conditions).

De?ne the loadings.

The preprocessing (model de?nition) step is critical. In no case is there a better

example of the computer-related axiom “garbage in, garbage out.” A perfectly

computed ?nite element solution is of absolutely no value if it corresponds to the

wrong problem.

1.3.2 Solution

During the solution phase, ?nite element software assembles the governing algebraic equations in matrix form and computes the unknown values of the primary

?eld variable(s). The computed values are then used by back substitution to

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

1.4 Brief History of the Finite Element Method

compute additional, derived variables, such as reaction forces, element stresses,

and heat ?ow.

As it is not uncommon for a ?nite element model to be represented by tens

of thousands of equations, special solution techniques are used to reduce data

storage requirements and computation time. For static, linear problems, a wave

front solver, based on Gauss elimination (Appendix C), is commonly used. While

a complete discussion of the various algorithms is beyond the scope of this text,

the interested reader will ?nd a thorough discussion in the Bathe book [1].

1.3.3 Postprocessing

Analysis and evaluation of the solution results is referred to as postprocessing.

Postprocessor software contains sophisticated routines used for sorting, printing,

and plotting selected results from a ?nite element solution. Examples of operations that can be accomplished include

Sort element stresses in order of magnitude.

Check equilibrium.

Calculate factors of safety.

Plot deformed structural shape.

Animate dynamic model behavior.

Produce color-coded temperature plots.

While solution data can be manipulated many ways in postprocessing, the most

important objective is to apply sound engineering judgment in determining

whether the solution results are physically reasonable.

1.4 BRIEF HISTORY OF THE FINITE

ELEMENT METHOD

The mathematical roots of the ?nite element method dates back at least a half

century. Approximate methods for solving differential equations using trial solutions are even older in origin. Lord Rayleigh [2] and Ritz [3] used trial functions

(in our context, interpolation functions) to approximate solutions of differential

equations. Galerkin [4] used the same concept for solution. The drawback in the

earlier approaches, compared to the modern ?nite element method, is that the

trial functions must apply over the entire domain of the problem of concern.

While the Galerkin method provides a very strong basis for the ?nite element

method (Chapter 5), not until the 1940s, when Courant [5] introduced the concept of piecewise-continuous functions in a subdomain, did the ?nite element

method have its real start.

In the late 1940s, aircraft engineers were dealing with the invention of the jet

engine and the needs for more sophisticated analysis of airframe structures to

withstand larger loads associated with higher speeds. These engineers, without

the bene?t of modern computers, developed matrix methods of force analysis,

© The McGraw?Hill

Companies, 2004

11

Hutton: Fundamentals of

Finite Element Analysis

12

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

collectively known as the ?exibility method, in which the unknowns are the

forces and the knowns are displacements. The ?nite element method, in its most

often-used form, corresponds to the displacement method, in which the unknowns are system displacements in response to applied force systems. In this

text, we adhere exclusively to the displacement method. As will be seen as we

proceed, the term displacement is quite general in the ?nite element method and

can represent physical displacement, temperature, or ?uid velocity, for example.

The term ?nite element was ?rst used by Clough [6] in 1960 in the context of

plane stress analysis and has been in common usage since that time.

During the decades of the 1960s and 1970s, the ?nite element method was

extended to applications in plate bending, shell bending, pressure vessels, and

general three-dimensional problems in elastic structural analysis [7–11] as well

as to ?uid ?ow and heat transfer [12, 13]. Further extension of the method to

large de?ections and dynamic analysis also occurred during this time period

[14 , 15]. An excellent history of the ?nite element method and detailed bibliography is given by Noor [16].

The ?nite element method is computationally intensive, owing to the required

operations on very large matrices. In the early years, applications were performed

using mainframe computers, which, at the time, were considered to be very powerful, high-speed tools for use in engineering analysis. During the 1960s, the ?nite

element software code NASTRAN [17] was developed in conjunction with the

space exploration program of the United States. NASTRAN was the ?rst major

?nite element software code. It was, and still is, capable of hundreds of thousands

of degrees of freedom (nodal ?eld variable computations). In the years since the

development of NASTRAN, many commercial software packages have been introduced for ?nite element analysis. Among these are ANSYS [18], ALGOR [19],

and COSMOS/M [20]. In today’s computational environment, most of these

packages can be used on desktop computers and engineering workstations to

obtain solutions to large problems in static and dynamic structural analysis, heat

transfer, ?uid ?ow, electromagnetics, and seismic response. In this text, we do not

utilize or champion a particular code. Rather, we develop the fundamentals for

understanding of ?nite element analysis to enable the reader to use such software

packages with an educated understanding.

1.5 EXAMPLES OF FINITE ELEMENT

ANALYSIS

We now present, brie?y, a few examples of the types of problems that can be

analyzed via the ?nite element method. Figure 1.7 depicts a rectangular region

with a central hole. The area has been “meshed” with a ?nite element grid of twodimensional elements assumed to have a constant thickness in the z direction.

Note that the mesh of elements is irregular: The element shapes (triangles and

quadrilaterals) and sizes vary. In particular, note that around the geometric discontinuity of the hole, the elements are of smaller size. This represents not only

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.5 Examples of Finite Element Analysis

Figure 1.7

A mesh of ?nite elements over a rectangular region having a

central hole.

an improvement in geometric accuracy in the vicinity of the discontinuity but

also solution accuracy, as is discussed in subsequent chapters.

The geometry depicted in Figure 1.7 could represent the ?nite element

model of several physical problems. For plane stress analysis, the geometry

would represent a thin plate with a central hole subjected to edge loading in the

plane depicted. In this case, the ?nite element solution would be used to examine stress concentration effects in the vicinity of the hole. The element mesh

shown could also represent the case of ?uid ?ow around a circular cylinder. In

yet another application, the model shown could depict a heat transfer ?n attached to a pipe (the hole) from which heat is transferred to the ?n for dissipation to the surroundings. In each case, the formulation of the equations governing physical behavior of the elements in response to external in?uences is quite

different.

Figure 1.8a shows a truss module that was at one time considered a

building-block element for space station construction [21]. Designed to fold in

accordion fashion into a small volume for transport into orbit, the module, when

deployed, extends to overall dimensions 1.4 m ? 1.4 m ? 2.8 m. By attaching

such modules end-to-end, a truss of essentially any length could be obtained.

The structure was analyzed via the ?nite element method to determine the

vibration characteristics as the number of modules, thus overall length, was

varied. As the connections between the various structural members are pin or

ball-and-socket joints, a simple axial tension-compression element (Chapter 2)

was used in the model. The ?nite element model of one module was composed

of 33 elements. A sample vibration shape of a ?ve-module truss is shown in

Figure 1.8b.

The truss example just described involves a rather large structure modeled

by a small number of relatively large ?nite elements. In contrast, Figure 1.9

shows the ?nite element model of a very thin tube designed for use in heat

13

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

XG

YG

ZG

(a)

(b)

Figure 1.8

(a) Deployable truss module showing details of folding joints.

(b) A sample vibration-mode shape of a ?ve-module truss as obtained

via ?nite element analysis. (Courtesy: AIAA)

14

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

1.5 Examples of Finite Element Analysis

0.00197

Z

X

0.488

0.25

Figure 1.9

Finite element model of a thin-walled

heat exchanger tube.

transfer in a spacecraft application. The tube has inside diameter of 0.976 in. and

wall thickness 0.00197 in. and overall length 36 in. Materials considered for

construction of the tube were copper and titanium alloys. Owing to the wall

thickness, prototype tubes were found to be very fragile and dif?cult to handle

without damage. The objectives of the ?nite element analysis were to examine

the bending, torsional, and buckling loads allowable. The ?gure shows the ?nite

element mesh used to model a section of the tube only 0.25 in. in length. This

model contains 1920 three-dimensional solid elements, each having eight nodes

with 3 degrees of freedom at each node. Such a large number of elements was

required for a small structure in consideration of computational accuracy. The

concern here was the so-called aspect ratio of the elements, as is de?ned and

discussed in subsequent chapters.

As a ?nal example, Figure 1.10a represents the ?nite element model of the

main load-carrying component of a prosthetic device. The device is intended to

be a hand attachment to an arti?cial arm. In use, the hand would allow a lower

arm amputee to engage in weight lifting as part of a physical ?tness program.

The ?nite element model was used to determine the stress distribution in the

component in terms of the range of weight loading anticipated, so as to properly

size the component and select the material. Figure 1.10b shows a prototype of the

completed hand design.

15

Hutton: Fundamentals of

Finite Element Analysis

16

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

(a)

(b)

Figure 1.10

(a) A ?nite element model of a prosthetic hand for weightlifting. (b) Completed

prototype of a prosthetic hand, attached to a bar.

(Courtesy of Payam Sadat. All rights reserved.)

1.6 OBJECTIVES OF THE TEXT

I wrote Fundamentals of Finite Element Analysis for use in senior-level ?nite

element courses in engineering programs. The majority of available textbooks

on the ?nite element method are written for graduate-level courses. These

texts are heavy on the theory of ?nite element analysis and rely on mathematical

techniques (notably, variational calculus) that are not usually in the repertoire of

undergraduate engineering students. Knowledge of advanced mathematical techniques is not required for successful use of this text. The prerequisite study is

based on the undergraduate coursework common to most engineering programs:

linear algebra, calculus through differential equations, and the usual series of

statics, dynamics, and mechanics of materials. Although not required, prior study

of ?uid mechanics and heat transfer is helpful. Given this assumed background,

the ?nite element method is developed on the basis of physical laws (equilibrium, conservation of mass, and the like), the principle of minimum potential energy (Chapter 2), and Galerkin’s ?nite element method (introduced and developed in Chapter 5).

Hutton: Fundamentals of

Finite Element Analysis

1. Basic Concepts of the

Finite Element Method

Text

© The McGraw?Hill

Companies, 2004

References

As the reader progresses through the text, he or she will discern that we

cover a signi?cant amount of ?nite element theory in addition to application

examples. Given the availability of many powerful and sophisticated ?nite

element software packages, why study the theory? The ?nite element method is

a tool, and like any other tool, using it without proper instruction can be quite

dangerous. My premise is that the proper instruction in this context includes

understanding the basic theory underlying formulation of ?nite element models

of physical problems. As stated previously, critical analysis of the results of a

?nite element model computation is essential, since those results may eventually

become the basis for design. Knowledge of the theory is necessary for both

proper modeling and evaluation of computational results.

REFERENCES

1.

Bathe, K-J. Finite Element Procedures. Englewood Cliffs, NJ: Prentice-Hall,

1996.

2. Lord Rayleigh. “On the Theory of Resonance.” Transactions of the Royal Society

(London) A161 (1870).

3. Ritz, W. “Uber eine neue Methode zur Losung gewissen Variations-Probleme der

mathematischen Physik.” J. Reine Angew. Math. 135 (1909).

4. Galerkin, B. G. “Series Solution of Some Problems of Elastic Equilibrium of Rods

and Plates” [in Russian]. Vestn. Inzh. Tekh. 19 (1915).

5. Courant, R. “Variational Methods for the Solution of Problems of Equilibrium and

Vibrations.” Bulletin of the American Mathematical Society 49 (1943).

6. Clough, R. W. “The Finite Element Method in Plane Stress Analysis.”

Proceedings, American Society of Civil Engineers, Second Conference on

Electronic Computation, Pittsburgh, 1960.

7. Melosh, R. J. “A Stiffness Method for the Analysis of Thin Plates in Bending.”

Journal of Aerospace Sciences 28, no. 1 (1961).

8. Grafton, P. E., and D. R. Strome. “Analysis of Axisymmetric Shells by the Direct

Stiffness Method.” Journal of the American Institute of Aeronautics and

Astronautics 1, no. 10 (1963).

9. Gallagher, R. H. “Analysis of Plate and Shell Structures.” Proceedings,

Symposium on the Application of Finite Element Methods in Civil Engineering,

Vanderbilt University, Nashville, 1969.

10. Wilson, E. L. “Structural Analysis of Axisymmetric Solids.” Journal of the

American Institute of Aeronautics and Astronautics 3, (1965).

11. Melosh, R. J. “Structural Analysis of Solids.” Journal of the Structural Division,

Proceedings of the American Society of Civil Engineers, August 1963.

12. Martin, H. C. “Finite Element Analysis of Fluid Flows.” Proceedings of the

Second Conference on Matrix Methods in Structural Mechanics, Wright-Patterson

Air Force Base, Kilborn, Ohio, October 1968.

13. Wilson, E. L., and R. E. Nickell. “Application of the Finite Element Method to

Heat Conduction Analysis.” Nuclear Engineering Design 4 (1966).

17

Hutton: Fundamentals of

Finite Element Analysis

18

1. Basic Concepts of the

Finite Element Method

CHAPTER 1

Text

© The McGraw?Hill

Companies, 2004

Basic Concepts of the Finite Element Method

14. Turner, M. J., E. H. Dill, H. C. Martin, and R. J. Melosh. “Large De?ections of

Structures Subjected to Heating and External Loads.” Journal of Aeronautical

Sciences 27 (1960).

15. Archer, J. S. “Consistent Mass Matrix Formulations for Structural Analysis Using

Finite Element Techniques.” Journal of the American Institute of Aeronautics and

Astronautics 3, no. 10 (1965).

16. Noor, A. K. “Bibliography of Books and Monographs on Finite Element

Technology.” Applied Mechanics Reviews 44, no. 6 (1991).

17. MSC/NASTRAN. Lowell, MA: MacNeal-Schwindler Corp.

18. ANSYS. Houston, PA: Swanson Analysis Systems Inc.

19. ALGOR. Pittsburgh: Algor Interactive Systems.

20. COSMOS/M. Los Angeles: Structural Research and Analysis Corp.

21. Hutton, D. V. “Modal Analysis of a Deployable Truss Using the Finite Element

Method.” Journal of Spacecraft and Rockets 21, no. 5 (1984).

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

C H A P T E R

2

Stiffness Matrices, Spring

and Bar Elements

2.1 INTRODUCTION

The primary characteristics of a ?nite element are embodied in the element

stiffness matrix. For a structural ?nite element, the stiffness matrix contains the

geometric and material behavior information that indicates the resistance of

the element to deformation when subjected to loading. Such deformation may

include axial, bending, shear, and torsional effects. For ?nite elements used in

nonstructural analyses, such as ?uid ?ow and heat transfer, the term stiffness

matrix is also used, since the matrix represents the resistance of the element to

change when subjected to external in?uences.

This chapter develops the ?nite element characteristics of two relatively

simple, one-dimensional structural elements, a linearly elastic spring and an elastic tension-compression member. These are selected as introductory elements because the behavior of each is relatively well-known from the commonly studied

engineering subjects of statics and strength of materials. Thus, the “bridge” to the

?nite element method is not obscured by theories new to the engineering student.

Rather, we build on known engineering principles to introduce ?nite element

concepts. The linear spring and the tension-compression member (hereafter referred to as a bar element and also known in the ?nite element literature as a spar,

link, or truss element) are also used to introduce the concept of interpolation

functions. As mentioned brie?y in Chapter 1, the basic premise of the ?nite element method is to describe the continuous variation of the ?eld variable (in this

chapter, physical displacement) in terms of discrete values at the ?nite element

nodes. In the interior of a ?nite element, as well as along the boundaries (applicable to two- and three-dimensional problems), the ?eld variable is described via

interpolation functions (Chapter 6) that must satisfy prescribed conditions.

Finite element analysis is based, dependent on the type of problem, on several mathematic/physical principles. In the present introduction to the method,

19

20

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

we present several such principles applicable to ?nite element analysis. First, and

foremost, for spring and bar systems, we utilize the principle of static equilibrium but—and this is essential—we include deformation in the development;

that is, we are not dealing with rigid body mechanics. For extension of the ?nite

element method to more complicated elastic structural systems, we also state and

apply the ?rst theorem of Castigliano [1] and the more widely used principle of

minimum potential energy [2]. Castigliano’s ?rst theorem, in the form presented,

may be new to the reader. The ?rst theorem is the counterpart of Castigliano’s

second theorem, which is more often encountered in the study of elementary

strength of materials [3]. Both theorems relate displacements and applied forces

to the equilibrium conditions of a mechanical system in terms of mechanical

energy. The use here of Castigliano’s ?rst theorem is for the distinct purpose of

introducing the concept of minimum potential energy without resort to the higher

mathematic principles of the calculus of variations, which is beyond the mathematical level intended for this text.

2.2 LINEAR SPRING AS A FINITE ELEMENT

A linear elastic spring is a mechanical device capable of supporting axial loading

only and constructed such that, over a reasonable operating range (meaning extension or compression beyond undeformed length), the elongation or contraction of the spring is directly proportional to the applied axial load. The constant

of proportionality between deformation and load is referred to as the spring constant, spring rate, or spring stiffness [4], generally denoted as k, and has units

of force per unit length. Formulation of the linear spring as a ?nite element is

accomplished with reference to Figure 2.1a. As an elastic spring supports axial

loading only, we select an element coordinate system (also known as a local coordinate system) as an x axis oriented along the length of the spring, as shown.

The element coordinate system is embedded in the element and chosen, by geometric convenience, for simplicity in describing element behavior. The element

u2

u1

f1

1

Force, f

Hutton: Fundamentals of

Finite Element Analysis

2

k

(a)

f2

k

1

x

Deflection, ? u2 u1

(b)

Figure 2.1

(a) Linear spring element with nodes, nodal displacements, and nodal forces.

(b) Load-de?ection curve.

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

or local coordinate system is contrasted with the global coordinate system. The

global coordinate system is that system in which the behavior of a complete

structure is to be described. By complete structure is meant the assembly of

many ?nite elements (at this point, several springs) for which we desire to compute response to loading conditions. In this chapter, we deal with cases in which

the local and global coordinate systems are essentially the same except for translation of origin. In two- and three-dimensional cases, however, the distinctions

are quite different and require mathematical recti?cation of element coordinate

systems to a common basis. The common basis is the global coordinate system.

Returning attention to Figure 2.1a, the ends of the spring are the nodes and

the nodal displacements are denoted by u1 and u2 and are shown in the positive

sense. If these nodal displacements are known, the total elongation or contraction

of the spring is known as is the net force in the spring. At this point in our development, we require that forces be applied to the element only at the nodes (distributed forces are accommodated for other element types later), and these are

denoted as f1 and f2 and are also shown in the positive sense.

Assuming that both the nodal displacements are zero when the spring is undeformed, the net spring deformation is given by

= u2 ? u1

(2.1)

and the resultant axial force in the spring is

f = k = k(u 2 ? u 1 )

(2.2)

as is depicted in Figure 2.1b.

For equilibrium, f 1 + f 2 = 0 or f 1 = ? f 2 , and we can rewrite Equation 2.2

in terms of the applied nodal forces as

f 1 = ?k(u 2 ? u 1 )

(2.3a)

f 2 = k(u 2 ? u 1 )

(2.3b)

which can be expressed in matrix form (see Appendix A for a review of matrix

algebra) as

k

?k

u1

f1

=

(2.4)

?k

k

u2

f2

or

[k e ]{u} = { f }

where

[k e ] =

k

?k

?k

k

(2.5)

(2.6)

is de?ned as the element stiffness matrix in the element coordinate system (or

local system), {u} is the column matrix (vector) of nodal displacements, and {f}

is the column matrix (vector) of element nodal forces. (In subsequent chapters,

21

Hutton: Fundamentals of

Finite Element Analysis

22

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

the matrix notation is used extensively. A general matrix is designated by

brackets [ ] and a column matrix (vector) by braces { }.)

Equation 2.6 shows that the element stiffness matrix for the linear spring

element is a 2 ? 2 matrix. This corresponds to the fact that the element exhibits

two nodal displacements (or degrees of freedom) and that the two displacements

are not independent (that is, the body is continuous and elastic). Furthermore, the

matrix is symmetric. A symmetric matrix has off-diagonal terms such that k i j =

kji. Symmetry of the stiffness matrix is indicative of the fact that the body is linearly elastic and each displacement is related to the other by the same physical

phenomenon. For example, if a force F (positive, tensile) is applied at node 2

with node 1 held ?xed, the relative displacement of the two nodes is the same as

if the force is applied symmetrically (negative, tensile) at node 1 with node 2

?xed. (Counterexamples to symmetry are seen in heat transfer and ?uid ?ow

analyses in Chapters 7 and 8.) As will be seen as more complicated structural

elements are developed, this is a general result: An element exhibiting N degrees

of freedom has a corresponding N ? N, symmetric stiffness matrix.

Next consider solution of the system of equations represented by Equation 2.4. In general, the nodal forces are prescribed and the objective is to solve

for the unknown nodal displacements. Formally, the solution is represented by

u1

f1

?1

= [k e ]

(2.7)

u2

f2

where [k e ]?1 is the inverse of the element stiffness matrix. However, this inverse

matrix does not exist, since the determinant of the element stiffness matrix is

identically zero. Therefore, the element stiffness matrix is singular, and this also

proves to be a general result in most cases. The physical signi?cance of the

singular nature of the element stiffness matrix is found by reexamination of

Figure 2.1a, which shows that no displacement constraint whatever has been imposed on motion of the spring element; that is, the spring is not connected to any

physical object that would prevent or limit motion of either node. With no constraint, it is not possible to solve for the nodal displacements individually.

Instead, only the difference in nodal displacements can be determined, as this

difference represents the elongation or contraction of the spring element owing

to elastic effects. As discussed in more detail in the general formulation of interpolation functions (Chapter 6) and structural dynamics (Chapter 10), a properly

formulated ?nite element must allow for constant value of the ?eld variable. In

the example at hand, this means rigid body motion. We can see the rigid body

motion capability in terms of a single spring (element) and in the context of several connected elements. For a single, unconstrained element, if arbitrary forces

are applied at each node, the spring not only deforms axially but also undergoes

acceleration according to Newton’s second law. Hence, there exists not only

deformation but overall motion. If, in a connected system of spring elements, the

overall system response is such that nodes 1 and 2 of a particular element displace the same amount, there is no elastic deformation of the spring and therefore

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

no elastic force in the spring. This physical situation must be included in the

element formulation. The capability is indicated mathematically by singularity

of the element stiffness matrix. As the stiffness matrix is formulated on the basis

of deformation of the element, we cannot expect to compute nodal displacements

if there is no deformation of the element.

Equation 2.7 indicates the mathematical operation of inverting the stiffness

matrix to obtain solutions. In the context of an individual element, the singular

nature of an element stiffness matrix precludes this operation, as the inverse of a

singular matrix does not exist. As is illustrated profusely in the remainder of the

text, the general solution of a ?nite element problem, in a global, as opposed to

element, context, involves the solution of equations of the form of Equation 2.5. For

realistic ?nite element models, which are of huge dimension in terms of the matrix

order (N ? N) involved, computing the inverse of the stiffness matrix is a very inef?cient, time-consuming operation, which should not be undertaken except for the

very simplest of systems. Other, more-ef?cient solution techniques are available,

and these are discussed subsequently. (Many of the end-of-chapter problems

included in this text are of small order and can be ef?ciently solved via matrix inversion using “spreadsheet” software functions or software such as MATLAB.)

2.2.1 System Assembly in Global Coordinates

Derivation of the element stiffness matrix for a spring element was based on

equilibrium conditions. The same procedure can be applied to a connected system of spring elements by writing the equilibrium equation for each node. However, rather than drawing free-body diagrams of each node and formally writing

the equilibrium equations, the nodal equilibrium equations can be obtained more

ef?ciently by considering the effect of each element separately and adding the

element force contribution to each nodal equation. The process is described as

“assembly,” as we take individual stiffness components and “put them together”

to obtain the system equations. To illustrate, via a simple example, the assembly

of element characteristics into global (or system) equations, we next consider the

system of two linear spring elements connected as shown in Figure 2.2.

For generality, it is assumed that the springs have different spring constants

k1 and k2. The nodes are numbered 1, 2, and 3 as shown, with the springs sharing

node 2 as the physical connection. Note that these are global node numbers. The

global nodal displacements are identi?ed as U1, U2, and U3, where the upper case

is used to indicate that the quantities represented are global or system displacements as opposed to individual element displacements. Similarly, applied nodal

U2

U1

1

F1

1

k1

U3

2

2

F2

k2

3

F3

Figure 2.2 System of two springs with node numbers,

element numbers, nodal displacements, and nodal forces.

23

Hutton: Fundamentals of

Finite Element Analysis

24

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

u1(1)

u(1)

2

u1(2)

1

f 1(1)

1

u(2)

2

2

2

k1

f (1)

2

f (2)

2

2

3

k2

(a)

f (2)

3

(b)

f 1(1)

F1

(2)

f (1)

2 f2

1

2

f (2)

3

F2

(d)

(c)

F3

(e)

Figure 2.3 Free-body diagrams of elements and nodes for the

two-element system of Figure 2.2.

forces are F1, F2, and F3. Assuming the system of two spring elements to be

in equilibrium, we examine free-body diagrams of the springs individually (Figure 2.3a and 2.3b) and express the equilibrium conditions for each spring, using

Equation 2.4, as

(1)

(1)

u1

f 1

k1

?k 1

=

(2.8a)

(1)

(1)

?k 1

k1

u2

f 2

(2)

(2)

u1

f 2

k2

?k 2

=

(2.8b)

(2)

(2)

?k 2

k2

u

f

2

3

where the superscript is element number.

To begin “assembling” the equilibrium equations describing the behavior

of the system of two springs, the displacement compatibility conditions, which

relate element displacements to system displacements, are written as

u

(1)

1

= U1

u

(1)

2

= U2

u

(2)

1

= U2

u

(2)

2

= U3

(2.9)

The compatibility conditions state the physical fact that the springs are connected at node 2, remain connected at node 2 after deformation, and hence, must

have the same nodal displacement at node 2. Thus, element-to-element displacement continuity is enforced at nodal connections. Substituting Equations 2.9 into

Equations 2.8, we obtain

(1)

f 1

k1

?k 1

U1

=

(2.10a)

(1)

?k 1

k1

U2

f 2

and

(2)

f 2

k2

?k 2

U2

=

(2.10b)

(2)

?k 2

k2

U3

f 3

Here, we use the notation f

node i.

( j)

i

to represent the force exerted on element j at

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

Equation 2.10 is the equilibrium equations for each spring element expressed

in terms of the speci?ed global displacements. In this form, the equations clearly

show that the elements are physically connected at node 2 and have the same displacement U2 at that node. These equations are not yet amenable to direct combination, as the displacement vectors are not the same. We expand both matrix

equations to 3 ? 3 as follows (while formally expressing the facts that element 1

is not connected to node 3 and element 2 is not connected to node 1):

?

?

? f (1) ?

U1

k1 ?k1 0

? 1 ?

U2 = f (1)

?k1 k1 0

(2.11)

?

?

? 2 ?

0

0

0

0

0

?

?

? 0 ?

0

0

0

0

?

?

(2)

0 k2 ?k2

U2 = f 2

(2.12)

?

?

? (2) ?

0 ?k2 k2

U3

f

3

The addition of Equations 2.11 and 2.12 yields

k1

?k1

0

?k1

k1 + k2

?k2

0

?k2

k2

U1

U2

U3

=

?

?

?

?

?

f (1)

1

(2)

f (1)

2 + f 2

f (2)

3

?

?

?

?

?

(2.13)

Next, we refer to the free-body diagrams of each of the three nodes depicted in

Figure 2.3c, 2.3d, and 2.3e. The equilibrium conditions for nodes 1, 2, and 3

show that

f

(1)

1

= F1

f

(1)

2

+ f

(2)

2

= F2

f

(2)

3

= F3

respectively. Substituting into Equation 2.13, we obtain the ?nal result:

?k 1

0

k1

U1

F1

?k 1 k 1 + k 2 ?k 2

U 2 = F2

0

?k 2

k2

U3

F3

(2.14)

(2.15)

which is of the form [K ]{U} = {F}, similar to Equation 2.5. However, Equation 2.15 represents the equations governing the system composed of two connected spring elements. By direct consideration of the equilibrium conditions,

we obtain the system stiffness matrix [K ] (note use of upper case) as

k1

?k 1

0

[K ] = ?k 1 k 1 + k 2 ?k 2

(2.16)

0

?k 2

k2

Note that the system stiffness matrix is (1) symmetric, as is the case with all linear systems referred to orthogonal coordinate systems; (2) singular, since no

constraints are applied to prevent rigid body motion of the system; and (3) the

system matrix is simply a superposition of the individual element stiffness

matrices with proper assignment of element nodal displacements and associated

stiffness coef?cients to system nodal displacements. The superposition procedure is formalized in the context of frame structures in the following paragraphs.

25

Hutton: Fundamentals of

Finite Element Analysis

26

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

EXAMPLE 2.1

Consider the two element system depicted in Figure 2.2 given that

Node 1 is attached to a ?xed support, yielding the displacement constraint U1 = 0.

k1 = 50 lb./in., k2 = 75 lb./in., F2 = F3 = 75 lb.

for these conditions determine nodal displacements U2 and U3.

? Solution

Substituting the speci?ed values into Equation 2.15 yields

?? ? ? ?

50 ?50

0

? 0 ? ? F1 ?

? ?50 125 ?75 ? U2 = 75

? ? ? ?

75

U3

0

?75 75

?

and we note that, owing to the constraint of zero displacement at node 1, nodal force F1

becomes an unknown reaction force. Formally, the ?rst algebraic equation represented in

this matrix equation becomes

?50U 2 = F1

and this is known as a constraint equation, as it represents the equilibrium condition

of a node at which the displacement is constrained. The second and third equations

become

125

?75

?75

75

U2

U3

=

75

75

which can be solved to obtain U2 = 3 in. and U3 = 4 in. Note that the matrix equations

governing the unknown displacements are obtained by simply striking out the ?rst row

and column of the 3 ? 3 matrix system, since the constrained displacement is zero.

Hence, the constraint does not affect the values of the active displacements (we use the

term active to refer to displacements that are unknown and must be computed). Substitution of the calculated values of U2 and U3 into the constraint equation yields the value

F1 = ?150 lb., which value is clearly in equilibrium with the applied nodal forces of

75 lb. each. We also illustrate element equilibrium by writing the equations for each

element as

50

?50

75

?75

(1)

f 1

?150

0

=

lb.

=

(1)

150

3

f 2

(2)

f 2

?75

?75

3

=

lb.

=

(2)

75

75

4

f 3

?50

50

for element 1

for element 2

Example 2.1 illustrates the general procedure for solution of ?nite element models: Formulate the system equilibrium equations, apply the speci?ed constraint

conditions, solve the reduced set of equations for the “active” displacements, and

substitute the computed displacements into the constraint equations to obtain the

unknown reactions. While not directly applicable for the spring element, for

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

27

more general ?nite element formulations, the computed displacements are also

substituted into the strain relations to obtain element strains, and the strains are,

in turn, substituted into the applicable stress-strain equations to obtain element

stress values.

EXAMPLE 2.2

Figure 2.4a depicts a system of three linearly elastic springs supporting three equal

weights W suspended in a vertical plane. Treating the springs as ?nite elements, determine the vertical displacement of each weight.

? Solution

To treat this as a ?nite element problem, we assign node and element numbers as shown

in Figure 2.4b and ignore, for the moment, that displacement U1 is known to be zero by

the ?xed support constraint. Per Equation 2.6, the stiffness matrix of each element is

(preprocessing)

k

(1)

=

k (2) =

k (3) =

3k

3k

?3k

?3k

3k

2k

?2k

?2k

2k

?k

k

k

?k

1

U1

3k

1

2k

2

k

3

W

2

U2

2k

W

3

U3

k

W

4

U4

(a)

(b)

Figure 2.4 Example 2.2: elastic

spring supporting weights.

Hutton: Fundamentals of

Finite Element Analysis

28

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

The element-to-global displacement relations are

u

(1)

1

= U1

u

(1)

2

=u

(2)

1

= U2

u

(2)

2

=u

(3)

1

= U3

u

(3)

2

= U4

Proceeding as in the previous example, we then write the individual element equations as

?? ? ?

?

U

3k ?3k 0 0 ?

? 1?

? ?

?

? ?3k 3k 0 0 ?? U2 ? ?

?

?

=

? 0

U ? ?

0

0 0 ??

?

? 3?

? ?

?

?

U4

0

0

0 0

?? ? ?

?

?

U

0

0

0

0 ?

? 1?

? ?

?

? 0 2k ?2k 0 ?? U2 ? ?

?

?

? 0 ?2k 2k 0 ?? U3 ? = ?

?

? ?

? ?

?

?

U4

0

0

0

0

?

?? ? ?

?

U

0 0 0

0 ?

? 1?

? ?

?

?0 0 0

?? U2 ? ?

0

?

?

=

? 0 0 k ?k ?? U3 ? ?

?

? ?

? ?

?

?

U4

0 0 ?k k

?

?

f (1)

1 ?

?

?

?

f (1)

2

?

0 ?

?

?

0

?

0 ?

(2) ?

?

f 1 ?

?

f (2)

2 ?

?

?

0

?

0 ?

?

?

0 ?

f (3)

1 ?

?

?

(3) ?

f 2

(1)

(2)

(3)

Adding Equations 1–3, we obtain

?? ? ? ?

U

F

3 ?3 0

0 ?

?

? ?

? 1?

? 1?

? ?3 5 ?2 0 ?? U2 ? ? W ?

?

k?

=

? 0 ?2 3 ?1 ?? U3 ? ? W ?

?

?

? ?

? ?

? ?

U4

W

0

0 ?1 1

?

(4)

where we utilize the fact that the sum of the element forces at each node must equal the

applied force at that node and, at node 1, the force is an unknown reaction.

Applying the displacement constraint U1 = 0 (this is also preprocessing), we obtain

?3kU 2 = F1

(5)

as the constraint equation and the matrix equation

?? ? ? ?

5 ?2 0 ? U2 ? ? W ?

k ? ?2 3 ?1 ? U3 = W

? ? ? ?

W

U4

0 ?1 1

?

(6)

for the active displacements. Again note that Equation 6 is obtained by eliminating the

constraint equation from 4 corresponding to the prescribed zero displacement.

Simultaneous solution (the solution step) of the algebraic equations represented by

Equation 6 yields the displacements as

U2 =

W

k

U3 =

2W

k

and Equation 5 gives the reaction force as

F1 = ?3W

(This is postprocessing.)

U4 =

3W

k

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.2 Linear Spring as a Finite Element

29

Note that the solution is exactly that which would be obtained by the usual statics

equations. Also note the general procedure as follows:

Formulate the individual element stiffness matrices.

Write the element to global displacement relations.

Assemble the global equilibrium equation in matrix form.

Reduce the matrix equations according to speci?ed constraints.

Solve the system of equations for the unknown nodal displacements (primary

variables).

Solve for the reaction forces (secondary variable) by back-substitution.

EXAMPLE 2.3

Figure 2.5 depicts a system of three linear spring elements connected as shown. The node

and element numbers are as indicated. Node 1 is ?xed to prevent motion, and node 3 is

given a speci?ed displacement as shown. Forces F2 = ? F and F4 = 2F are applied at

nodes 2 and 4. Determine the displacement of each node and the force required at node 3

for the speci?ed conditions.

? Solution

This example includes a nonhomogeneous boundary condition. In previous examples, the

boundary conditions were represented by zero displacements. In this example, we have

both a zero (homogeneous) and a speci?ed nonzero (nonhomogeneous) displacement

condition. The algebraic treatment must be different as follows. The system equilibrium

equations are expressed in matrix form (Problem 2.6) as

?

?k

4k

?3k

0

k

? ?k

?

? 0

0

0

?3k

5k

?2k

?

?? ? ? ? ?

0

U1 ?

F1 ?

F1 ?

?

?

?

?

?

?

?

?

?

? ? ? ? ?

?

0 ?

? U2 = F2 = ?F

?

?2k ?

U ? ?

F ? ?

F ?

?

? 3?

? ?

? ?

?

? 3?

? 3 ?

2k

U4

F4

2F

Substituting the speci?ed conditions U 1 = 0 and U 3 = results in the system of

equations

?

k

? ?k

?

? 0

0

?k

4k

?3k

0

k

?

?? ? ?

0

0 ?

F1 ?

?

?

?

?

?

?

? ? ?

?

0 ?

? U2 = ?F

?

? F3 ?

?2k ?

? ?

?

?

?

? ?

?

?

2k

2F

U4

F2 F

1

1

0

?3k

5k

?2k

3

3

2

3k

2

4

F4 2F

2k

?

Figure 2.5 Example 2.3: Three-element system with speci?ed

nonzero displacement at node 3.

Hutton: Fundamentals of

Finite Element Analysis

30

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Since U 1 = 0 , we remove the ?rst row and column to obtain

?

4k

? ?3k

0

?3k

5k

?2k

?

?? ? ?

0

? U2 ? ? ?F ?

=

F

?2k ?

? ? ? 3 ?

2F

U4

2k

as the system of equations governing displacements U2 and U4 and the unknown nodal

force F3. This last set of equations clearly shows that we cannot simply strike out the row

and column corresponding to the nonzero speci?ed displacement because it appears in

the equations governing the active displacements. To illustrate a general procedure, we

rewrite the last matrix equation as

?

5k

? ?3k

?2k

?3k

4k

0

?

?? ? ?

?2k ? ? ? F3 ?

0 ? U2 = ?F

?

? ? ?

U4

2F

2k

Next, we formally partition the stiffness matrix and write

?

5k

? ?3k

?2k

?3k

4k

0

?? ?

?2k ? ?

{}

{F }

[K ] [K U ]

=

0 ? U2 =

? ?

{FU }

[K U ] [K UU ] {U }

U4

2k

with

[K ] = [5k]

[K U ] = [?3k

?2k]

?3k

[K U ] = [K U ] T =

?2k

4k 0

[K UU ] =

0 2k

{} = {}

U2

{U } =

U4

{F } = {F3 }

?F

{FU } =

2F

From the second “row” of the partitioned matrix equations, we have

[K U ]{} + [K UU ]{U } = {FU }

and this can be solved for the unknown displacements to obtain

{U } = [K UU ]?1 ({F } ? [K U ]{})

provided that [K UU ]?1 exists. Since the constraints have been applied correctly, this

inverse does exist and is given by

?

[K UU ]?1

1

? 4k

=?

?

0

?

0 ?

?

1 ?

2k

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

Substituting, we obtain the unknown displacements as

?

?

?

3 ?

F

?

?

?

0 ?

?

?

+

?

?

4k

4 ?

? ?F + 3k

=

?

1 ? 2F + 2k

?

?

F

?

?

?

+ ?

?

?

2k

k

?

1

?

U2

? 4k

=?

{U } =

U4

?

0

The required force at node 3 is obtained by substitution of the displacement into the upper

partition to obtain

5

3

F3 = ? F + k

4

4

Finally, the reaction force at node 1 is

F

3

? k

4

4

F1 = ?kU 2 =

As a check on the results, we substitute the computed and prescribed displacements into

the individual element equations to insure that equilibrium is satis?ed.

Element 1

?k

k

k

?k

0

U2

=

?kU2

kU2

=

?

?

? f (1) ?

1

? f (1) ?

2

which shows that the nodal forces on element 1 are equal and opposite as required for

equilibrium.

Element 2

3k

?3k

?3k

3k

U2

U3

?

?

3

F

?3k ? ? + ?

4k

4

?

3k ?

?

?

3

3F

?

(2)

?

?

? k ?

?

??

f 2

4k

4

=

=

?

?

?

f (2)

?

? 3F + 3 k ?

3

4k

4

3k

=

?3k

which also veri?es equilibrium.

Element 3

2k

?2k

?2k

2k

U3

U4

=

2k

?2k

?2k

2k

Therefore element 3 is in equilibrium as well.

F

+

k

=

?2F

2F

=

f

f

(3)

3

(3)

4

2.3 ELASTIC BAR, SPAR/LINK/TRUSS ELEMENT

While the linear elastic spring serves to introduce the concept of the stiffness matrix, the usefulness of such an element in ?nite element analysis is rather limited.

Certainly, springs are used in machinery in many cases and the availability of a

?nite element representation of a linear spring is quite useful in such cases. The

31

Hutton: Fundamentals of

Finite Element Analysis

32

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

spring element is also often used to represent the elastic nature of supports for

more complicated systems. A more generally applicable, yet similar, element is

an elastic bar subjected to axial forces only. This element, which we simply call

a bar element, is particularly useful in the analysis of both two- and threedimensional frame or truss structures. Formulation of the ?nite element characteristics of an elastic bar element is based on the following assumptions:

1.

2.

3.

4.

The bar is geometrically straight.

The material obeys Hooke’s law.

Forces are applied only at the ends of the bar.

The bar supports axial loading only; bending, torsion, and shear are not

transmitted to the element via the nature of its connections to other

elements.

The last assumption, while quite restrictive, is not impractical; this condition is

satis?ed if the bar is connected to other structural members via pins (2-D) or balland-socket joints (3-D). Assumptions 1 and 4, in combination, show that this is

inherently a one-dimensional element, meaning that the elastic displacement of

any point along the bar can be expressed in terms of a single independent variable. As will be seen, however, the bar element can be used in modeling both

two- and three-dimensional structures. The reader will recognize this element

as the familiar two-force member of elementary statics, meaning, for equilibrium, the forces exerted on the ends of the element must be colinear, equal in

magnitude, and opposite in sense.

Figure 2.6 depicts an elastic bar of length L to which is af?xed a uniaxial

coordinate system x with its origin arbitrarily placed at the left end. This is the

element coordinate system or reference frame. Denoting axial displacement at

any position along the length of the bar as u(x), we de?ne nodes 1 and 2 at each

end as shown and introduce the nodal displacements u 1 = u(x = 0) and

u 2 = u(x = L ) . Thus, we have the continuous ?eld variable u(x), which is to be

expressed (approximately) in terms of two nodal variables u 1 and u 2 . To accomplish this discretization, we assume the existence of interpolation functions

N 1 (x ) and N 2 (x ) (also known as shape or blending functions) such that

u(x ) = N 1 (x )u 1 + N 2 (x )u 2

(2.17)

u1

u2

1

2

x

u(x)

L

Figure 2.6 A bar (or truss) element with element

coordinate system and nodal displacement

notation.

x

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

(It must be emphasized that, although an equality is indicated by Equation 2.17,

the relation, for ?nite elements in general, is an approximation. For the bar element, the relation, in fact, is exact.) To determine the interpolation functions, we

require that the boundary values of u(x ) (the nodal displacements) be identically

satis?ed by the discretization such that

u(x = 0) = u 1

u(x = L ) = u 2

(2.18)

Equations 2.17 and 2.18 lead to the following boundary (nodal) conditions:

N 1 (0) = 1

N 2 (0) = 0

(2.19)

N1( L ) = 0

N2( L ) = 1

(2.20)

which must be satis?ed by the interpolation functions. It is required that the displacement expression, Equation 2.17, satisfy the end (nodal) conditions identically, since the nodes will be the connection points between elements and the

displacement continuity conditions are enforced at those connections. As we

have two conditions that must be satis?ed by each of two one-dimensional functions, the simplest forms for the interpolation functions are polynomial forms:

N 1 (x ) = a0 + a1 x

(2.21)

N 2 (x ) = b0 + b1 x

(2.22)

where the polynomial coef?cients are to be determined via satisfaction of the

boundary (nodal) conditions. We note here that any number of mathematical

forms of the interpolation functions could be assumed while satisfying the

required conditions. The reasons for the linear form is explained in detail in

Chapter 6.

Application of conditions represented by Equation 2.19 yields a0 = 1,

b0 = 0 while Equation 2.20 results in a1 = ?(1/L ) and b1 = x /L . Therefore,

the interpolation functions are

N 1 (x ) = 1 ? x /L

(2.23)

N 2 (x ) = x /L

(2.24)

and the continuous displacement function is represented by the discretization

u(x ) = (1 ? x /L )u 1 + (x /L )u 2

(2.25)

As will be found most convenient subsequently, Equation 2.25 can be expressed

in matrix form as

u

u(x ) = [N 1 (x ) N 2 (x )] 1 = [N ] {u}

(2.26)

u2

where [N ] is the row matrix of interpolation functions and {u} is the column

matrix (vector) of nodal displacements.

Having expressed the displacement ?eld in terms of the nodal variables, the

task remains to determine the relation between the nodal displacements and

applied forces to obtain the stiffness matrix for the bar element. Recall from

33

Hutton: Fundamentals of

Finite Element Analysis

34

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

elementary strength of materials that the de?ection of an elastic bar of length L

and uniform cross-sectional area A when subjected to axial load P is given by

PL

=

(2.27)

AE

where E is the modulus of elasticity of the material. Using Equation 2.27, we

obtain the equivalent spring constant of an elastic bar as

P

AE

k=

=

(2.28)

L

and could, by analogy with the linear elastic spring, immediately write the stiffness matrix as Equation 2.6. While the result is exactly correct, we take a more

general approach to illustrate the procedures to be used with more complicated

element formulations.

Ultimately, we wish to compute the nodal displacements given some loading

condition on the element. To obtain the necessary equilibrium equations relating

the displacements to applied forces, we proceed from displacement to strain,

strain to stress, and stress to loading, as follows. In uniaxial loading, as in the bar

element, we need consider only the normal strain component, de?ned as

du

?x =

(2.29)

dx

which, when applied to Equation 2.25, gives

u2 ? u1

?x =

(2.30)

L

which shows that the spar element is a constant strain element. This is in accord

with strength of materials theory: The element has constant cross-sectional area

and is subjected to constant forces at the end points, so the strain does not vary

along the length. The axial stress, by Hooke’s law, is then

u2 ? u1

x = E ? x = E

(2.31)

L

and the associated axial force is

AE

P = x A =

(u 2 ? u 1 )

(2.32)

L

Taking care to observe the correct algebraic sign convention, Equation 2.32 is

now used to relate the applied nodal forces f 1 and f 2 to the nodal displacements

u 1 and u 2 . Observing that, if Equation 2.32 has a positive sign, the element is in

tension and nodal force f 2 must be in the positive coordinate direction while

nodal force f 1 must be equal and opposite for equilibrium; therefore,

AE

f1 = ?

(u 2 ? u 1 )

(2.33)

L

AE

f2 =

(u 2 ? u 1 )

(2.34)

L

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

Equations 2.33 and 2.34 are expressed in matrix form as

AE 1 ?1

u1

f1

=

?1

1

u

f2

L

2

35

(2.35)

Comparison of Equation 2.35 to Equation 2.4 shows that the stiffness matrix for

the bar element is given by

AE 1 ?1

[k e ] =

(2.36)

L ?1 1

As is the case with the linear spring, we observe that the element stiffness matrix

for the bar element is symmetric, singular, and of order 2 ? 2 in correspondence

with two nodal displacements or degrees of freedom. It must be emphasized that

the stiffness matrix given by Equation 2.36 is expressed in the element coordinate system, which in this case is one-dimensional. Application of this element

formulation to analysis of two- and three-dimensional structures is considered in

the next chapter.

EXAMPLE 2.4

Figure 2.7a depicts a tapered elastic bar subjected to an applied tensile load P at one end

and attached to a ?xed support at the other end. The cross-sectional area varies linearly

from A 0 at the ?xed support at x = 0 to A 0 /2 at x = L . Calculate the displacement of the

end of the bar (a) by modeling the bar as a single element having cross-sectional area

equal to the area of the actual bar at its midpoint along the length, (b) using two bar

elements of equal length and similarly evaluating the area at the midpoint of each, and

(c) using integration to obtain the exact solution.

? Solution

(a) For a single element, the cross-sectional area is 3A 0 /4 and the element “spring

constant” is

k=

AE

3A 0 E

=

L

4L

and the element equations are

3A 0 E

4L

1

?1

?1

?1

U1

U2

=

F1

P

The element and nodal displacements are as shown in Figure 2.7b. Applying the

constraint condition U 1 = 0 , we ?nd

U2 =

4PL

PL

= 1.333

3A 0 E

A0 E

as the displacement at x = L .

(b) Two elements of equal length L /2 with associated nodal displacements are

depicted in Figure 2.7c. For element 1, A 1 = 7A 0 /8 so

k1 =

7A 0 E

A1E

7A 0 E

=

=

L1

8( L /2)

4L

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

1

x

7

A

8 o

L

(

A(x) Ao 1

x

2L

A

)

3

A

4 o

5

A

8 o

2

P

(b)

(a)

(c)

1.4

Exact

One element

Two elements

1.2

P

Ao

x

?(x)

0.8

?Y

1.0

0.6

0.4

0.2

0

0

0.1

0.2

0.3

0.4

P

0.5

x

L

0.6

0.7

0.8

0.9

1.0

(e)

(d)

2

Exact

One element

Two elements

1.8

(

?

P

?o ?o Ao

)

1.6

1.4

1.2

1

0.8

0.6

0.4

0.2

0

0

0.1

0.2

0.3

0.4

0.5

x

L

0.6

0.7

0.8

0.9

1.0

(f)

Figure 2.7

(a) Tapered axial bar, (b) one-element model, (c) two-element model, (d) free-body diagram

for an exact solution, (e) displacement solutions, (f) stress solutions.

36

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.3 Elastic Bar, Spar/Link/Truss Element

while for element 2, we have

5A 0

8

A1 =

A2 E

5A 0 E

5A 0 E

=

=

L2

8( L /2)

4L

and k 2 =

Since no load is applied at the center of the bar, the equilibrium equations for the

system of two elements is

?

?? ? ? ?

0 ? U1 ? ? F1 ?

?k2 ? U2 =

0

? ? ? ?

P

k2

U3

?k1

k1 + k2

?k2

k1

? ?k1

0

Applying the constraint condition U 1 = 0 results in

k1 + k2

?k 2

?k 2

k2

U2

U3

=

0

P

Adding the two equations gives

P

4PL

=

k1

7A 0 E

U2 =

and substituting this result into the ?rst equation results in

k1 + k2

48PL

PL

=

= 1.371

k2

35A 0 E

A0 E

U3 =

(c) To obtain the exact solution, we refer to Figure 2.7d, which is a free-body diagram of

a section of the bar between an arbitrary position x and the end x = L . For equilibrium,

x A = P

x

A = A(x ) = A 0 1 ?

2L

and since

the axial stress variation along the length of the bar is described by

x =

A0

P

x

1?

2L

Therefore, the axial strain is

?x =

x

=

E

E A0

P

x

1?

2L

Since the bar is ?xed at x = 0 , the displacement at x = L is given by

L

=

0

=

P

? x dx =

EA 0

L

0

dx

x

1?

2L

L

2 PL

2 PL

2 PL

PL

[?ln(2L ? x )]0 =

[ln(2L ) ? ln L ] =

ln 2 = 1.386

E A0

E A0

E A0

A0 E

37

Hutton: Fundamentals of

Finite Element Analysis

38

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Comparison of the results of parts b and c reveals that the two element solution

exhibits an error of only about 1 percent in comparison to the exact solution from

strength of materials theory. Figure 2.7e shows the displacement variation along the

length for the three solutions. It is extremely important to note, however, that the

computed axial stress for the ?nite element solutions varies signi?cantly from that of

the exact solution. The axial stress for the two-element solution is shown in Figure 2.7f, along with the calculated stress from the exact solution. Note particularly

the discontinuity of calculated stress values for the two elements at the connecting

node. This is typical of the derived, or secondary, variables, such as stress and strain,

as computed in the ?nite element method. As more and more smaller elements are

used in the model, the values of such discontinuities decrease, indicating solution

convergence. In structural analyses, the ?nite element user is most often more interested in stresses than displacements, hence it is essential that convergence of the

secondary variables be monitored.

2.4 STRAIN ENERGY, CASTIGLIANO’S

FIRST THEOREM

When external forces are applied to a body, the mechanical work done by those

forces is converted, in general, into a combination of kinetic and potential energies. In the case of an elastic body constrained to prevent motion, all the work

is stored in the body as elastic potential energy, which is also commonly

referred to as strain energy. Here, strain energy is denoted U e and mechanical

work W. From elementary statics, the mechanical work performed by a force F

as its point of application moves along a path from position 1 to position 2 is

de?ned as

2

r

W = F · d

(2.37)

1

where

d

r = dx i + dy j + dz k

(2.38)

is a differential vector along the path of motion. In Cartesian coordinates, work

is given by

x2

y2

z2

W = Fx dx + Fy dy + Fz dz

x1

y1

(2.39)

z1

where Fx , Fy , and Fz are the Cartesian components of the force vector.

For linearly elastic deformations, de?ection is directly proportional to applied force as, for example, depicted in Figure 2.8 for a linear spring. The slope

of the force-de?ection line is the spring constant such that F = k. Therefore,

the work required to deform such a spring by an arbitrary amount 0 from its

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

Force, F

2.4 Strain Energy, Castigliano’s First Theorem

k

1

Deflection, ?

Figure 2.8 Force-de?ection

relation for a linear elastic

spring.

free length is

0

0

1

W = F d = k d = k02 = U e

2

0

(2.40)

0

and we observe that the work and resulting elastic potential energy are quadratic

functions of displacement and have the units of force-length. This is a general

result for linearly elastic systems, as will be seen in many examples throughout

this text.

Utilizing Equation 2.28, the strain energy for an axially loaded elastic bar

?xed at one end can immediately be written as

Ue =

1 2

1 AE 2

k =

2

2 L

(2.41)

However, for a more general purpose, this result is converted to a different form

(applicable to a bar element only) as follows:

2

1 2

1 AE PL

1 P

P

1

U e = k =

=

AL = ?V

(2.42)

2

2 L

AE

2 A

AE

2

where V is the total volume of deformed material and the quantity 12 ? is strain

energy per unit volume, also known as strain energy density. In Equation 2.42,

stress and strain values are those corresponding to the ?nal value of applied

force. The factor 12 arises from the linear relation between stress and strain as the

load is applied from zero to the ?nal value P. In general, for uniaxial loading, the

strain energy per unit volume u e is de?ned by

?

u e = d?

(2.43)

0

which is extended to more general states of stress in subsequent chapters. We note

that Equation 2.43 represents the area under the elastic stress-strain diagram.

39

Hutton: Fundamentals of

Finite Element Analysis

40

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Presently, we will use the work-strain energy relation to obtain the governing equations for the bar element using the following theorem.

Castigliano’s First Theorem [1]

For an elastic system in equilibrium, the partial derivative of total strain energy

with respect to de?ection at a point is equal to the applied force in the direction

of the de?ection at that point.

Consider an elastic body subjected to N forces Fj for which the total strain

energy is expressed as

j

N

Ue = W =

Fj dj

(2.44)

j=1

0

where j is the de?ection at the point of application of force Fj in the direction of

the line of action of the force. If all points of load application are ?xed except

one, say, i, and that point is made to de?ect an in?nitesimal amount i by an

incremental in?nitesimal force Fi , the change in strain energy is

i

U e = W = Fi i + Fi di

(2.45)

0

where it is assumed that the original force Fi is constant during the in?nitesimal

change. The integral term in Equation 2.45 involves the product of in?nitesimal

quantities and can be neglected to obtain

U e

= Fi

i

(2.46)

which in the limit as i approaches zero becomes

?U

= Fi

? i

(2.47)

The ?rst theorem of Castigliano is a powerful tool for ?nite element formulation, as is now illustrated for the bar element. Combining Equations 2.30, 2.31,

and 2.43, total strain energy for the bar element is given by

2

1

1

u2 ? u1

U e = x ? x V = E

AL

(2.48)

2

2

L

Applying Castigliano’s theorem with respect to each displacement yields

?U e

AE

=

(u 1 ? u 2 ) = f 1

? u1

L

(2.49)

?U e

AE

=

(u 2 ? u 1 ) = f 2

? u2

L

(2.50)

which are observed to be identical to Equations 2.33 and 2.34.

Hutton: Fundamentals of

Finite Element Analysis

2. Stiffness Matrices,

Spring and Bar Elements

Text

© The McGraw?Hill

Companies, 2004

2.4 Strain Energy, Castigliano’s First Theorem

41

The ?rst theorem of Castigliano is also applicable to rotational displacements. In the case of rotation, the partial derivative of strain energy with respect

to a rotational displacement is equal to the moment/torque applied at the point of

concern in the sense of the rotation. The following example illustrates the application in terms of a simple torsional member.

EXAMPLE 2.5

A solid circular shaft of radius R and length L is subjected to constant torque T. The shaft

is ?xed at one end, as shown in Figure 2.9. Formulate the elastic strain energy in terms of

the angle of twist at x = L and show that Castigliano’s ?rst theorem gives the correct

expression for the applied torque.

? Solution

From strength of materials theory, the shear stress at any cross section along the length of

the member is given by

=

Tr

J

where r is radial distance from the axis of the member and J is polar moment of inertia of

the cross section. For elastic behavior, we have

Tr

=

G

JG

=

where G is the shear modulus of the material, and the strain energy is then

1

Ue =

2

1

dV =

2

V

=

?

?

L

Tr

Tr

?

dA?dx

J

JG

A

0

T2

2J 2 G

L

r 2 dA dx =

T2L

2JG

0 A

where we have used the de?nition of the polar moment of inertia

J =

r2 dA

A

R

L

T

Figure 2.9 Example 2.5:

Circular cylinder subjected to

torsion.

Hutton: Fundamentals of

Finite Element Analysis

42

2. Stiffness Matrices,

Spring and Bar Elements

CHAPTER 2

Text

© The McGraw?Hill

Companies, 2004

Stiffness Matrices, Spring and Bar Elements

Again invoking the strength of materials results, the angle of twist at the end of the

### Случайные PDF

Файл |

27.pdf |

ДЗ 1.pdf |

3j_semestr_var_20_1.pdf |

Тема 1.pdf |

Potapcev1.pdf |