xi,j as shown in Fig. 3.1 and represent the solution by its values at the grid
points ui,j . Also the differential equation could be discretized using so-called
finite difference schemes, e.g. for the Laplace operator r2 on a square grid,
£
r2u
¤
i,j '
1
h2
³
ui+1,j + ui−1,j + ui,j+1 + ui,j−1 − 4ui,j
´
, (3.3)
where h is the grid point spacing. Thus we would obtain a set of algebraic
equations to solve for the unknowns ui,j .
44 CHAPTER 5. NUMERICAL MODEL
i+1
i
i−1
j
j+1
j−1
fj
Figure 5.1: Sketch of a triangular computational mesh for the nite element method (FEM). The domain
is divided into triangular regions and a basis function j (in text n) is associated to each nodal point.
The basis function is depicted as a piecewise linear function, but in general, it can have a higher-order
spatial dependence. Figure adapted from Ref. 45.
Figure 3.1: The finite difference method.
A computational grid showing a grid point
xi,j (•) and neighbors (±).
Figure 3.2: The finite element method.
Division of the domain into a mesh of
finite elements; showing a piecewise linear
function 'j(x) with compact support.
is a forcing term. We now expand g(r) on a set of basis functions n,
X
In the finite element method the discretization takes a different route
g(r) =
cn n(r); (5.2)
in that one represents the solution as a linear combination of some basis
functions '` defined on the whole domain
where cn are the expansion coecients. Dening the defect d(r),
u(x) =
X
`
n
u` '`(x). (3.4)
d(r) = Lfg(r)g �� F(r); (5.3)
the expansion of g(r) Eq. (5.2) will provide an approximative solution to the strong form
equation (5.1), if the inner product of each basis function and the defect d(r) vanishes
Z
The discretization is then introduced by choosing only a finite set of basis
functions '`, rather than a complete and therefore infinite set.
*
m(r);
"
L
X
��F(r)
The question is now what kind of functions to choose for the basis. A well
m(r) d(r)
dr =
m(r); d(r)
=
cn n(r)
known choice is Fourier expansion using plane waves with the characteristics
that they are infinitely
smooth and mutually orthogonal, and n
that they
extend across the entire domain.
The particular choice of the finite element method is somewhat opposite
as one chooses basis functions '` that are local, in the sense that each '` have
#+
= 0; for all m:
(5.4)
This is a formalized way of saying that the values of the coecients cn is chosen such
as to obtain the best approximative solution to the governing equation (5.1). A solution
fullling Eq. (5.4) is called a weak solution, whereas a solution fullling Eq. (5.1) is called
a strong solution.
If the dierential operator L is linear, Eq. (5.4) becomes
X
n
cn
m(r);L n(r)
=
m(r); F(r)
; for all m; (5.5)
which can be written as a matrix equation
Kc = f; (5.6)
whereK is refereed to as the stiness matrix and its elements are given by Km;n =
m;L n
,
c is a
vector
with the expansion coecient cn, and f is the force vector with elements
fm =
m; F
. In order to nd a weak solution to the governing equation (5.1), we thus
need to solve the matrix equation (5.6) for the unknown expansion coecients in c.
In the following we consider an inhomogeneous PDE on the general divergence form
r J + F = 0; (5.7)