Finite Element Analysis of Structures through Unified Formulation

By Erasmo Carrera, Maria Cinefra, Marco Petrolo and Enrico Zappino

The finite element method (FEM) is a computational tool widely used to design and analyse  complex structures. Currently, there are a number of different approaches to analysis using the FEM that vary according to the type of structure being analysed: beams and plates may use 1D or 2D approaches, shells and solids 2D or 3D approaches, and methods that work for one structure are typically not optimized to work for another.

Finite Element Analysis of Structures Through Unified Formulation deals with the FEM used for the analysis of the mechanics of structures in the case of linear elasticity. The novelty of this book is that the finite elements (FEs) are formulated on the basis of a class of theories of structures known as the Carrera Unified Formulation (CUF). It formulates 1D, 2D and 3D FEs on the basis of the same 'fundamental nucleus' that comes from geometrical relations and Hooke's law, and presents both 1D and 2D refined FEs that only have displacement variables as in 3D elements. It also covers 1D and 2D FEs that make use of 'real' physical surfaces rather than ’artificial’ mathematical surfaces which are difficult to interface in CAD/CAE software.

Key features:

• Covers how the refined formulation can be easily and conveniently used to analyse laminated structures, such as sandwich and composite structures, and to deal with multifield problems
• Shows the performance of different FE models through the 'best theory diagram' which allows different models to be compared in terms of accuracy and computational cost
• Introduces an axiomatic/asymptotic approach that reduces the computational cost of the structural analysis without affecting the accuracy
• Introduces an innovative 'component-wise' approach to deal with complex structures
• Accompanied by a website hosting the dedicated software package MUL2 (www.mul2.com)

Finite Element Analysis of Structures Through Unified Formulation is a valuable reference for researchers and practitioners, and is also a useful source of information for graduate students in civil, mechanical and aerospace engineering.

• Language: English
• Publisher: Wiley
• Release date: Jul 29, 2014
• ISBN: 9781118536650

Book preview

1

Introduction

1.1 What is in this Book

This book is devoted to the FE analysis of structures referring to linear elastic materials and small displacement assumptions. Attention is mainly focused on displacement formulations. This book is intended for two kinds of readers:

Those who are not familiar with FEs and would like to learn about them in a unified formulation framework.

Those who are familiar with FEs, but would like to learn how they are formulated with a unified formulation to overcome their limitations based on classical theories of structures.

Compared with other books on the subject, the present book offers the following novel features:

It formulates 1D, 2D and 3D FEs on the basis of the same ‘fundamental nucleus’ that comes from geometrical relations and Hooke's law; the only difference between 3D elements and 1D and 2D elements is that cross-section and through-the-thickness integrals are introduced, respectively. The differential operators are the same in all three cases.

It formulates refined 1D and 2D theories through formulae that remain invariant with respect to the variable order of the expansions used for the unknown displacement variables over the beam cross-section and plate/shell thickness.

It shows that an appropriate refinement of 1D FEs can make them suitable for analysing shell structures.

It presents both 1D and 2D refined FEs that only have displacement variables as in 3D elements. This is obtained using Lagrange polynomials to expand the displacement fields over the beam cross-section and plate/shell thickness.

It presents 1D and 2D FEs that make use of ‘real’ physical surfaces rather than ‘artificial’ mathematical surfaces. Classical 1D and 2D elements in fact need lines and reference surfaces, both of which are artificial entities, to build mathematical models of a given structure. This means that the models from computer-aided design (CAD) tools have to be modified to obtain lines and reference surfaces for plate/shell elements. The use of real surface properties can facilitate the direct construction of FE mathematical models from CAD. On the other hand, modifications of the structure from FE analysis can easily be transferred to CAD, since physical surfaces are used.

It shows how the described refined formulation can be easily and conveniently used to analyse laminated structures, such as sandwich and composite structures.

It shows how refined elements are essential in the case of multifield loadings, e.g. thermal, electrical and magnetic loadings.

It underlines that the use of refined FEs, unlike most available FEM codes, requires one to overcome the constraints of having only 6 DOFs for each node (3 displacements + 3 rotations, as in the Newtonian mechanics of a rigid body).

It introduces an axiomatic/asymptotic approach that reduces the computational cost of the structural analysis without affecting the accuracy.

It illustrates the performance of different FE models through the ‘best theory diagram’ (BTD), which allows different models to be compared in terms of accuracy and computational cost.

This book is devoted mostly to the development of FEs and not to procedures that can be used to solve FE problems since these procedures have already been covered in many other books on the FEM (Bathe 1996; Hughes 2000; Oñate 2009; Reddy 2005; Zienkiewicz et al. 2005).

1.2 The Finite Element Method

From a historical point of view, the FEM was introduced when the following two developments became available:

Technological: The development of computers capable of conducting mathematical operations very quickly.

Mathematical: The development of mathematical methods that can be used to solve differential equations in ‘approximated’ or ‘weak’ forms (both ordinary differential equations, ODEs, and partial differential equations, PDEs). The most significant of these methods is probably the weight residual method (WRM).

These two developments are common to all computational mechanics, and are discussed briefly in the following subsections.

1.2.1 Approximation of the Domain

‘Technological’ developments mean that automatic procedures are needed to solve the problem automatically in weak form. The problem can therefore be ‘discretized’ into a number of finite domains, or elements, in which mathematical tools work properly. The problem can then be ‘assembled’ and finally ‘solved’ using a computer.

In the case of the theory of structures, the problem is related to a 3D, 2D or 1D domain, which can be quite complex in terms of geometry, geometrical boundary conditions (conditions on generalized displacements in a given set of points) and mechanical boundary conditions (loadings). Geometrical domain examples are given in Figures 1.1, 1.2 and 1.3, where typical 3D, 2D and 1D domains are shown.

Figure 1.1 Example of 3D structure

Figure 1.2 Example of 2D structure

Figure 1.3 Example of 1D structure (helicopter blade)

The key idea of the FEM is to ‘discretize’ complex domains into simpler ones, as in the following:

Let us consider the 3D domain V in Figure 1.1 which can be discretized into a finite number of regular 3D solids (brick elements), as in Figure 1.4. It should be noted that the original domain can be violated slightly to correspond to its boundary surfaces. This cannot be avoided in the discretization process.

Let us consider a 2D domain, Ω, which coincides with the reference surface of curved or flat panels, see Figure 1.2. The Ω domain can be discretized into a number of regular quadrilateral or triangular figures, or a mixture of these, as can be seen in Figure 1.5.

Let us consider a 1D domain, l, which coincides with the reference line (axis) of a ‘beam’ structure, see Figure 1.3. The l domain can be discretized into a finite number of regular lines, as in Figure 1.6.

Figure 1.4 Example of a 3D FEM model

Figure 1.5 Example of 2D FEM model

Figure 1.6 Example of 1D FEM model

The above discretizations are common to all typical structures. One well-known example is the truss structure, shown in Figure 1.7. This kind of structure is significant from a historical point of view, since it represents one of the first applications of the FEM to the analysis of structures.

Figure 1.7 Example of truss structure

1.2.2 The Numerical Approximation

The solution of a structural problem consists of finding a solution to a given set of governing equations, which are defined in the V, Ω or l domain. The given set of governing equations of the unknown functions Φ

(1.1) numbered Display Equation

can appear, from a mathematical point of view, as one of the following cases:

O is a differential operator, for instance

(1.2) numbered Display Equation

O is an integral operator

(1.3) numbered Display Equation

O is an algebraic operator

(1.4) numbered Display Equation

O is given by any combination of the above.

The FEM solves the above equations, in a weak sense, at the element (or subdomain) level. The unknown function is usually assumed to be a combination of a finite number of Φ in a given set of points, which are ‘the nodes’ of the element. In the 3D case

(1.5) numbered Display Equation

The application of the WRM, for the governing equations (namely Equations (1.1–1.4)) of the element, leads to the following system of algebraic equations:

(1.6) numbered Display Equation

which is written for each element, E, and which is the governing equation for the element in the FEM sense. The governing FEM equation of the whole structure can be found, or ‘assembled’, by imposing compatibility/equilibrium conditions on the values of Φi at the nodes:

(1.7) numbered Display Equation

Equation (1.6), the equilibrium equation, represents the result of a ‘mathematical’ problem. Writing and solving Equations (1.7) for a structure is a matter of automatic calculation. New ‘technology’ arising from the introduction of computers allows the problem to be solved easily.

1.3 Calculation of the Area of a Surface with a Complex Geometry via the FEM

Let us consider a surface Ω with a complex geometry (boundary Γ), such as that in Figure 1.8. The aim here is to calculate the area A of surface Ω. The complexity of the geometry of the surface does not allow a closed-form formula to be used to compute the area; that is, to express the area in terms of the geometrical parameters Ω and Γ:

(1.8) numbered Display Equation

A ‘numerical’ solution can be obtained by ‘discretizing’ the surface into simpler ‘elements’ or ‘subdomains’, such as quadrilaterals or triangles, whose area can easily be computed, see Figure 1.9. The number of elements, Ne, is finite. The words ‘finite’ and ‘elements’ give the name to the FEM.

Figure 1.8 Surface Ω

Figure 1.9 Coarse mesh

Figure 1.10 Refined mesh

The area, A, of a generic element i is denoted by Ai. If Ne is the number of elements, the unknown area is

(1.9) numbered Display Equation

where is not the exact value of the area, but only an approximation. As usual, the approximation process leads to an error

(1.10) numbered Display Equation

where E is the grey area in Figure 1.9, which can be positive or negative. The error arises because of the approximation of the surface boundary, Γ, where a continuous curved line is approximated by the sum of the straight lines. However, such an error can be reduced by introducing smaller triangles or quadrilaterals, see Figure 1.10, which allow a better simulation of the boundary to be made.

Figure 1.11 An axially loaded bar

This simple problem clearly illustrates the nature of ‘discetrization’ problems, and the need to conduct several calculations of the area for high Ne. When Ne is increased, an automatic tool, e.g. ‘a computer’, is needed to compute the area of a generic complex surface.

This simple example is not sufficient to show the mathematical difficulties involved in solving a problem at the element level. The area of a triangle does not in fact introduce any approximation, and can be computed exactly using well-known formulae.

Unfortunately, this is not the case for FEs of structures. This problem requires the application of a mathematical approximation process to solve the PDEs that usually govern elastic problems related to beams, plates, shells and solids.

1.4 Elasticity of a Bar

The simplest structural element is a bar. Let us consider a bar loaded by an axial loading, q(y), as shown in Figure 1.11:

it has a 1D behaviour, which means that the problem variables can be expressed in terms of only one coordinate, in this case the y coordinate;

it can only carry loadings applied along its axis.

The only stress acting on the generic cross-section (σyy = σ) is assumed constant over the section itself. The stress resultant can be defined as

(1.11) numbered Display Equation

If an infinitesimal portion of the bar, with length dy and cross-section A (see Figure 1.12), is considered, the equilibrium along the bar axis leads to

(1.12) numbered Display Equation

which means that the variation in the axial forces (dN/dy) is balanced by the applied axial loading (q). It is clear that the solution of the previous equation, i.e. the distribution of N along the axis, is subordinate to the form of the applied loading q. If q is the only load applied on the bar, the following cases can be of particular interest:

if q = 0,

(1.13) numbered Display Equation

where N0 is the value of N at y = 0.

if q is constant,

(1.14) numbered Display Equation

in more complex cases, N can vary linearly or parabolically along the axis, and can be expressed as

(1.15) numbered Display Equation

where L is the length of the bar.

Hooke's law and the strain displacement relation are now introduced:

(1.16) numbered Display Equation

where E is Young's modulus and A is the area of the cross-section of the bar. E and A are considered constant. The variation of N(y), in terms of axial displacements, is

(1.17) numbered Display Equation

The equilibrium equation can therefore be written in terms of displacements

(1.18) numbered Display Equation

The derivative of the displacement at each point of the bar can be evaluated by integrating Equation (1.18),

(1.19) numbered Display Equation

and, by integrating once more, the displacement becomes

(1.20)

numbered Display Equation

where s and r vary between 0 and L. The q = 0 case leads to a linear displacement field along the bar axis.

In the simplest case, q = 0 and (duy/dy)|y = 0 = 0, and if only a concentrated load is considered, one obtains

(1.21) numbered Display Equation

where c is a constant that can be computed by imposing the given boundary condition. The linear form of uy is consistent with the constant value of the axial stress resultant N.

1.5 Stiffness Matrix of a Single Bar

Let us consider the simplest case of an axially loaded bar of length L. The whole bar is the FE under investigation. Points y = 0 and y = L are the nodes of the bar element. See Figure 1.13 for the notation.

Figure 1.12 Equilibrium of an infinitesimal portion of a bar

Figure 1.13 Local reference system and node numeration of a bar

If u1 and u2 are the values of the displacements at the nodes, then

(1.22) numbered Display Equation

where N1(y) and N2(y) are the two linear Lagrange polynomials

(1.23) numbered Display Equation

which satisfy the conditions

(1.24)

numbered Display Equation

N1 and N2 are known as ‘shape functions’ in FE procedures.

The stress resultant, N, is

(1.25) numbered Display Equation

where

(1.26)

numbered Display Equation

Figure 1.14 Relation between stress and resultant

Equation (1.25) represents the relationship between the stress resultant and the strain in the generic section that corresponds to y, see Figure 1.14. If a force is applied at node 1, in the y-direction, or y direction equilibrium leads to

(1.27) numbered Display Equation

therefore

(1.28) numbered Display Equation

while at node 2, one obtains

(1.29) numbered Display Equation

in terms of axial displacement

(1.30) numbered Display Equation

These two equilibrium conditions can be expressed in matrix form:

(1.31) numbered Display Equation

By introducing the displacement vector (where T denotes vector/matrix transposition)

(1.32) numbered Display Equation

and the force vector

(1.33) numbered Display Equation

Equation (1.31) becomes

(1.34) numbered Display Equation

where is

(1.35) numbered Display Equation

This matrix is known as the stiffness matrix in FE procedures. It should be noted that this matrix is symmetric and positive semidefinite.1

Example 1.5.1 Let us consider the bar in Figure 1.15. The bar is axially loaded at node 2 and is clamped at node 1. The displacement vector has two contributions, the free displacements and the constrained displacements , which are

(1.36) numbered Display Equation

and the force vector

(1.37) numbered Display Equation

where is the vector of the forces applied to the free nodes, while is the reaction force in the constrained node, which is denoted by . The problem that has to be solved can be written in the following form:

(1.38) numbered Display Equation

In this case, the problem becomes

(1.39) numbered Display Equation

The displacement vector is obtained by solving the first equation:

(1.40) numbered Display Equation

The reaction force at node 1 can be computed using the second equation:

(1.41) numbered Display Equation

Figure 1.15 Example of an axially loaded bar, physical and FEM model

1.6 Stiffness Matrix of a Bar via the PVD

The most powerful tool that can be used to derive FE equations in both static and dynamic cases and for linear and nonlinear problems is without doubt the principle of virtual work, PVW. In this book, PVW is used for FE problems with only displacement unknowns. In this case, it is referred to as the principle of virtual displacements (PVD). PVD is applied in this book to solids, beams, plates and shells for classical and refined FEM formulations.

Figure 1.16 Linear Lagrange polynomials: triangle similarity

Equations (1.31), which were introduced in Section 1.7, are now obtained using the PVD. For convenience, the displacement at point y is rewritten in the following rather universal notation:

(1.42) numbered Display Equation

where

(1.43) numbered Display Equation

is the matrix of the shape functions and

(1.44) numbered Display Equation

is the vector of the unknown displacements. In this case, N1 and N2 are linear Lagrange polynomials. These can also be derived by considering the conditions

(1.45) numbered Display Equation

(1.46) numbered Display Equation

In order to verify these conditions, the shape functions have to fulfil the following requirements:

(1.47)

numbered Display Equation

Because of the linearity of the functions, it is also possible to derive N(y) in a simple geometrical manner. Figure 1.16 shows the application of the Thales theorem, which can be written as

(1.48) numbered Display Equation

The explicit form of the displacement becomes

(1.49) numbered Display Equation

If matrix notation is introduced, the displacement can be rewritten as

(1.50) numbered Display Equation

and the explicit forms of the shape functions are

(1.51) numbered Display Equation

which coincide with the linear Lagrange polynomials.

The axial strain is

(1.52) numbered Display Equation

If the differential operator is introduced and the displacement, uy, is expressed in terms of nodal displacements, the axial strain becomes

(1.53) numbered Display Equation

where matrix is

(1.54) numbered Display Equation

In explicit form, the strain