CHAPTER 1
INTRODUCTION
One of the most important advances in applied mathematics in the 20th century has been
the development of the Finite Element Method as a general mathematical tool for
obtaining approximate solutions to boundary-value problems. The theory of finite
elements draws on almost every branch of mathematics and can be considered as one of
the richest and most diverse bodies of the current mathematical knowledge.
1. 1 Mathematical Modeling of Physical Systems
Due to the complexity of physical systems, some approximation must be made in the
process of turning physical reality into a mathematical model. It is important to decide at
what points in the modeling process these approximations are made. This, in turn,
determines what type of analytical or computational scheme is required in the solution
process. Let us consider a diagram of the two common branches of the general modeling-
solution process given in Figure 1:
For many real world problems the second approach is in fact the only possibility. For
instance suppose that the aim is to find the thermo-mechanical stresses in an air-cooled
turbine blade depicted in Figure 2.
The complex three-dimensional geometry of the blade along with the combined thermal
and mechanical loadings makes the analysis of the blade a formidable task. Nevertheless,
many powerful commercial finite element packages are available that can be
implemented to perform this task with relative ease.
1
1.2 FEM Analysis Process
A model-based simulation process using FEM consists of a sequence of steps. This
sequence takes two basic configurations depending on the environment in which FEM is
used. These are referred to as the Mathematical FEM and the Physical FEM.
Mathematical FEM
The centerpiece in the process steps of the Mathematical FEM is the mathematical model
which is often an ordinary or partial differential equation in space and time. Using the
methods provided by the Variational Calculus, a discrete finite element model is
generated from of the mathematical model. The resulting FEM equations are processed
by an equation solver, which provides a discrete solution. In this process we may also
think of an ideal physical system, which may be regarded as a realization of the
mathematical model. For example, if the mathematical model is the Poisson’s equation,
realizations may be a heat conduction problem. In Mathematical FEM this step is
unnecessary and indeed FEM discretizations may be constructed without any reference
to physics.
The concept of error arises when the discrete solution is substituted in the mathematical
and discrete models. This replacement is generically called verification. The solution
error is the amount by which the discrete solution fails to satisfy the discrete equations.
This error is relatively unimportant when using computers. More relevant is the
discretization error, which is the amount by which the discrete solution fails to satisfy the
mathematical model.
2
Figure 1.1 Comparison of Analytical and Computational Model
Physical FEM
The processes of idealization and discretization are carried out concurrently to produce
the discrete model. Indeed FEM discretizations may be constructed and adjusted without
reference to mathematical models, simply from experimental measurements. The concept
of error arises in the physical FEM in two ways, known as verification and validation,
respectively. Verification is the same as in the Mathematical FEM: the discrete solution
is replaced into the discrete model to get the solution error. As noted above, this error is
not generally important. Validation tries to compare the discrete solution against
observation by computing the simulation error, which combines modeling and solution
errors. Since the latter is typically insignificant, the simulation error in practice can be
identified with the modeling error. Comparing the discrete solution with the ideal
physical system would in principle quantify the modeling errors.
Figure 1.2 Finite Element Discretization
3
Figure 1.3 Finite Element Solution
Figure 1.4 Deviation of the solution from the mathematical model
4
Figure 1.5 Deviation of the solution from physical system
1.3 Application of of FEM
Figure 1.6 Application Examples of FEM
5
FEM PROCEDURE
6
7
CHAPTER 2
FEM FORMULATION
2.1 Direct Formulation
2.2 Weighted Residual Methods
Example with a single governing equation with only one independent variable
f [T ( x)] 0 in
T is the function sought , function of x only
is the domain of the region governed by f
Boundary conditions
g1 [T ( x)] 0 in 1
g 2 [T ( x)] 0 in 2
1 and 2 are parts of the boundary of
Approximation of the solution with a T ' function:
n
T ' T ' ( x; a1 , a 2 , , a n ) ai N i ( x)
i 1
which has one or more unknown( but constant) parameters a1 , a 2 , a n satisfies exactly
the boundary conditions. No surprise if the approximation does not satisfy the equation
exactly! We will get a residual error:
f [T ' ( x; a1 , a 2 , , a n )] R ( x; a1 , a 2 , , a n )
The method of weighted residuals requires that the parameters a1 , a 2 , a n be
determined satisfying:
W ( x)R( x; a
i 1 , a 2 , , a n ) d 0
where the functions Wi ( x ) are the n arbitrary weighting functions
The conditions of the weighting functions is generally left to a personal judgment
The most popular weighted-residual methods are:
8
1) Point collection
2) Subdomain collection
3) Least squares
4) Galerkin
Point Collection
The weighted functions Wi ( x) are ( x xi ) and defined such that
b
a
( x xi ) dx 1 for x xi
b
a
( x xi )dx 0 for x xi
Substitution of this choice of wi (x) gives:
( x x ) R( x; a , a
i 1 2 , a n ) d 0 for i 1,2, , n
which is evaluated at n collection points x1 , x 2 , , x n results n algebraic equations in n
unknowns
R ( x1 ; a1 , a 2 , , a n ) 0
R ( x 2 ; a1 , a 2 , , a n ) 0
R ( x n ; a1 , a 2 , , a n ) 0
Subdomain Collection
The weighted functions Wi (x ) are:
1 for x in 1
W1 ( x )
0 for x not in 1
1 for x in 2
W2 ( x )
0 for x not in 2
Substitution of this choice of wi (x) gives the following n integral equations
R( x; a , a
1
1 2 , a n ) d 1 0
R( x; a , a
2
1 2 , a n ) d 2 0
R( x; a , a
n
1 2 , a n ) d n 0
Least Squares
The method of least squares requires that the integral I of the square of the residual R be
minimized. That is:
I [ R( x; a
1 , a 2 , , a n )] 2 d be a minimum, or equivalently
9
I
ai
ai [ R( x; a , a
1 2 , , a n )]2 d a
i
[ R ( x; a1 , a 2 , , a n )]2 d
Carrying out the differentiations and simplifying, we have:
R R R
R a
i
d R a
2
d R a
n
d 0
which means that the weighting functions are:
R
R a d 0
i
i 1,2, , n
Galerkin
The weighting functions wi (x) are N i (x )
Therefore for the Galerkin method of the weighted residuals we have:
N
i ( x) R ( x; a1 , a 2 , , a n ) d 0 i 1,2, , n
Remember that our approximation of the T ( x ) function is:
n
T T ' ( x; a1 , a 2 , , a n ) ai N i ( x)
i 1
10
2.2 Galerkin’s Method
Let us consider a steady state continuous physical system described the following system
of PDEs:
F(u) + f = 0 on domain
G(u) + f = 0 on bounder
Example: Poisson’s equation
2T 2T q ' ' '
0 on
x 2 y 2 k
Boundary condition
a) Dirchlet BC (Natural BC)
T=Ts on 1
b) Neuman BC (Essential BC)
T
h(T T ) on 2
n
T
c) q on 3
n
- domain - boundary
Figure 2.1 The domain of the boundary
The residual is defined as follows:
R(u) = F(u) + f
The residual vanishes when the solution is substituted.
The weighted residual method consists in finding functions u that satisfy the
following integral equation
W R(u )dR W F (u )
i
i f d 0 i 1...n
11
Where Wi is weighting function and u is the solution that satisfy the boundary condition.
Example
Integral of Poisson’s Equation
2u 2u
Wi x , y 2 2 f d 0
x y
u – Must be twice differentiable and it should satisfy all the boundary conditions on
u and f
Integration by Parts
Gradient or green Theorem
F dxdy n Fds
F F
(i
x
j
y
) dxdy
(ni n y j ) Fds
Divergence Theorem
.Gdxdy n.Gds
G x G y
(
x
y
)dxdy
( n x G x n y G y ) ds
G
w( 2G ) dxdy w.Gdxdy w
n
ds
G w
w
x
dxdy
x
Gdxdy n x wGds
G w
w
x
dxdy
y
Gdxdy
n x wGds
In one dimension:
x2 x2
du dw
x W dx dx x dx udx Wu
x2
x1
1 1
x2 x2
d 2u dW du du
x dx 2 dx dx dx dx W dx
W
x1
1
In two dimension:
12
u w
W
x
dxdy
x
udxdy Wudy
w
= udxdy Wulds
x
l cos
u W
W y
dxdy
y
udxdy Wudx
W
udxdy Wumds
y
m sin
u 2
W u u
W x 2 dxdy A x x dxdy x lds
2u 2u 2W 2W u W
x 2 y 2 x 2 y 2
W u dxdy W
n
u ds
n
Weak Integral Form
A given integral form may be transformed to obtain a so-called weak form through
integration by parts. By this process, the order of the highest derivative can be reduced.
Boundary conditions other than u can also be specified. However, the integration by parts
introduce derivates of the weighing function W. Thus, the continuity condition of W are
more severe
Example
Weak integral of Poisson’s equation
2u 2u N i u N i u u u
i x 2 y 2 f dxdy x x y y N i f d N i n ds N i n ds 0
N
f
.
Ni must be twice differentiable
N u N i u
i N i f d N i f U ds 0
x x y y
Ni must be twice differentiable
N nno N j N i nno N j nno
i u j u j N f
i d N f
i N j u j ds 0 i 1...nno
x j 1 x y j 1 y j 1
K U F
13
N i N j N i N j
K i , j d N i N j ds 0
i 1...nno j 1...nno
x x y y
Fi N i f d f N i ds 0 i 1...nno
Where:
K is stiffness Matrix
F is load vector
U is nodal values of function of interest
14
Variational Formulation
A functional is linear
u u
(u ,
x
) (a u a
1 2
x
) d
A functional is quadratic if all terms are of second order, for example
u 2
(a
1 (
x
) a 2 u 2 ) d
For a purely quadratic functional
u
u u
1
u D d
2
x x
D is a symmetric matrix independent of u
u
(u )
u
u D d
x x
is positive definite if D is positive definite matrix that is all the value of D are
positive
Example
Find the variation of the following one dimensional functional
du x2 1 du
(u ,
dx
) x1
( ( ) 2 u f )dx
2 dx
Its variation
x2 1 du 2
x1
( ) u f ) dx
2 dx
Using properties of
x2 du du
x1
( ( )
dx dx
u f ) dx
x2 d du
x1
(
dx
(u )
dx
u f )dx
The second variation of can be obtained
x2 du 2
2 ( )
x1
( (
dx
)) dx
x2 d (u ) 2
( ) dx 0
x1 dx
d
w (u )
dx
(u ) 0
x1 dw du
0 x2
(
dx dx
w f )dx
15
R 0
R u L(u ) f d 0
L is a linear operator
f and f are independent of u
These conditions are sufficient for a functional to exist
Example
Formulate functional for poission’s equation
2u 2u
F (u ) f 2 2 f 0
x x
The corresponding integral form is obtained previously
w u w u
R (
x x
y y
wf ) d w( u f ) d 0
Choosing w u
( u ) u ( u ) u
R ( x x y y u f )d u ( u f )d 0
Defining the functional
u u 1 u 1 u 1
(u , ,
x y
)
( ( ) 2 ( ) 2 uf ) d ( u 2 uf ) d
2 x 2 y 2
R 0
A solution u for R 0 also renders the functional stationary 0 . At this condition the
functional is either a minimum or a maximum.
To illustrate the process let us consider now a specific example.
Suppose we specify the problem by requiring the stationary of a functional
1 T 2 1 T 2 .
k k Q T d q d (9.72)
x
2 2 y v rq
In which k and Q depend only on position and T such that T 0 on
where and q are bounding the domain .
We now perform the variation. This can be written following rules of
differentiation as
T T T T
k
x x
k Qv T d (q T ) d
y y
rq
As
16
T
T
x x
We can integrate by parts and, noting that T 0 on , and obtain
T T
T k k Qv d
x x y y
T
T k q d 0
rq
n
We immediately observe that the Euler equations are
T T
T k k Qv in
x y y y
T
T k q=0 on q
n
If T is so prescribed that T T on and T 0 on that boundary, then the problem
is precisely the one we have already discussed and the functional specifies the two-
dimensional heat conduction problem in an alternative way.
In this case we have ‘guessed’ the functional but the reader will observe that the variation
operation could have been carried out for any functional specified and corresponding
Euler equations could have been established.
Let us continue the problem to obtain an approximate solution of the linear heat
conduction problem. Taking, as usual,
T T N i ai (9.76)
We substitute this approximation into the expression for the functional and obtain
2 2
1 N i 1 N i
k a i d k ai d
2
x 2 y
Q v N i a i d q N i a i d
rq
On differentiation with respect to atypical parameter a j we have
17
N i N j N i N j
a j
k
x xa i d k
y ai y d
N d
Q v j qN j d
rq
and a system of equations for solution of the problem is
K a f
with
N i N j N i N j
K ij K ji
k
x x
d k
y y
d
f i N j Q d N j q d
rq
18
