<<

. 47
( 67 .)



>>

tions w. It therefore should also hold for w = g( x). For this particular choice of
w, the integral expression yields
Numerical solution of one-dimensional diffusion equation
236

b
g2 ( x) dx = 0 ’ g( x) = 0 on a ¤ x ¤ b. (14.7)
a
This follows immediately from the observation that the square of a function g( x)
is always greater than or equal to zero for any value of x, i.e. g2 ( x) ≥ 0, such that
the integral of g2 ( x) over the domain a ¤ x ¤ b must be greater than or equal to
zero, i.e.
b
g2 ( x) dx ≥ 0, (14.8)
a
and can only be equal to zero if g( x) = 0 for all a ¤ x ¤ b.
Effectively, the method of weighted residuals transforms the requirement that a
function, say g( x), must be equal to zero on a given domain at an in¬nite number
of points into a single evaluation of the integral, that must be equal to zero.
Using the method of weighted residuals, the differential equation, Eq. (14.1), is
transformed into an integral equation:
b d du
+f dx = 0,
w c (14.9)
dx dx
a

which should hold for all weighting functions w( x). The term between the square
brackets contains second order derivatives d2 u/dx2 of the function u. As has been
outlined in the introduction, approximate solutions of u are sought by de¬ning an
interpolation of u on the domain of interest and transforming the integral equation
into a discrete set of linear equations. De¬ning interpolation functions that are
both second-order differentiable and still integrable is far from straightforward,
in particular in the multi-dimensional case on arbitrarily shaped domains. Fortu-
nately, the second-order derivatives can be removed by means of an integration by
parts:
b b b
dw du
du
’ dx + wf dx = 0.
c (14.10)
wc
dx dx dx
a a
a
This introduces the boundary terms:
b
du du du
= ’w( a) c + w( b) c
wc . (14.11)
dx dx dx
a x=a x=b
At the boundary either u is prescribed (i.e. the essential boundary condition at
x = a) or the derivative c du/dx is prescribed (i.e. the natural boundary condition
at x = b). Along the boundary where u is prescribed the corresponding ¬‚ux, say
pu , with
du
pu = c , (14.12)
dx x=a
237 14.4 Polynomial interpolation


dx x=b = P is known. For
is unknown, while along the other boundary the ¬‚ux c du
the time being the combination of boundary terms in Eq. (14.11) is written as
b
du
B= w c = ’w( a) pu + w( b) P, (14.13)
dx a
realizing that pu is considered as yet unknown. The term B is introduced as an
abbreviation for the boundary contribution.
Consequently, the following integral equation, known as the weak form,
results:
b b
dw du
dx = wf dx + B.
c (14.14)
dx dx
a a
Eq. (14.14) is the point of departure for the Finite Element Method, and is called
the weak form, because the differentiability requirements imposed on u have been
reduced: in the original differential equation, Eq. (14.1), the second-order deriva-
tive d2 u/dx2 appears, while in the weak form, Eq. (14.14), only the ¬rst-order
derivative du/dx has to be dealt with.



14.4 Polynomial interpolation

Suppose that at a ¬nite number of points xi in the domain , the function values
ui = u( xi ) are known, then a polynomial approximation, denoted by uh , of u( x)
can be constructed. The polynomial approximation uh of degree n ’ 1 can be
constructed by
uh ( x) = a0 + a1 x + a2 x2 + · · · + an’1 xn’1 , (14.15)
if u is known at n points. The coef¬cients ai can be identi¬ed uniquely and
expressed in terms of ui , by solving the set of equations:
⎡ ¤⎡ ¤⎡ ¤
2 · · · xn’1
1 x1 x1 u1
a0
1
⎢ 1 x x2 · · · xn’1 ⎥ ⎢ a ⎥ ⎢ u ⎥
⎢ ⎥⎢ 1⎥ ⎢ 2⎥
2
⎢ 2 · · · xn’1 ⎥ ⎢ a ⎥⎢ ⎥
2 2
⎢ 1 x 3 x3 ⎥⎢ 2 ⎥ = ⎢ u3 ⎥ . (14.16)
⎢ ⎥⎢ ⎥⎢ ⎥
3
⎢. . . ⎥⎢ . ⎥ ⎢ . ⎥
.. . . . ¦⎣ . ¦ ⎣ . ¦
. .
⎣. . . .
. . .
···
2 n’1 an’1 un
1 xn xn xn
An example is given in Fig. 14.2, where, in the domain xi’1 ¤ x ¤ xi+2 , a third-
order polynomial (dashed curve) is used to approximate a given function (solid
curve) based on the function values ui’1 , ui , ui+1 and ui+2 .
Clearly, the coef¬cients ai are linearly dependent on the values ui , therefore the
polynomial may be rewritten in terms of ui by
Numerical solution of one-dimensional diffusion equation
238


u




x
xi “1 xi xi +1 xi +2

Figure 14.2
Solid line: u( x), dashed line: polynomial approximation of u( x).


n
uh ( x) = Ni ( x) ui , (14.17)
i=1

where the functions Ni ( x) are polynomial expressions of order n ’ 1 in terms of
the coordinate x. These functions Ni ( x) are called shape functions because they
de¬ne the shape of the interpolation of uh , for instance linear, quadratic etc. To
illustrate this, consider a ¬rst-order polynomial on the domain [ x1 , x2 ], with u( x)
known at x1 and x2 . In that case
x ’ x1 x ’ x1
uh ( x) = 1 ’ u1 + u2 , (14.18)
x2 ’ x1 x2 ’ x1
implying that
x ’ x1 x ’ x1
N1 = 1 ’ N2 =
, . (14.19)
x2 ’ x1 x2 ’ x1
Rather than approximating u( x) with a single polynomial of a certain degree over
the entire domain , the domain may also be divided into a number of non-
overlapping subdomains, say e . Within each subdomain e a local polynomial
approximation of u may be constructed. A typical example is a piecewise lin-
ear approximation within each subdomain e , see Fig. 14.3. Consider one of the
subdomains e = [ xi , xi+1 ]. Then within e the function u( x) is approximated by

uh ( x) | = N1 ( x) ui + N2 ( x) ui+1 , (14.20)
e

with, in conformity with Eq. (14.19):
x ’ xi x ’ xi
N1 = 1 ’ N2 = . (14.21)
xi+1 ’ xi xi+1 ’ xi
More generally, if within a subdomain e an n-th order polynomial approximation
of u is applied, the subdomain should have n + 1 points at which the function u is
239 14.5 Galerkin approximation

u ui




ui +1

x
xi “1 xi xi +1 xi +2
©e

Figure 14.3
Solid line: u( x), dashed line: piecewise linear approximation of u( x).



known. For instance, to use a second-order (quadratic) polynomial, the subdomain
should cover at least three consecutive points, for example e =[ xi , xi+2 ], such
that
uh ( x) | = N1 ( x) ui + N2 ( x) ui+1 + N3 ( x) ui+2 . (14.22)
e

The subdomain e within which a certain polynomial approximation is used is
referred to as an element. The points at which the values of u are de¬ned are
called nodes.
The shape functions Ni may not be chosen arbitrarily. The most stringent
requirement is that uh must be interpolated continuously over the total domain
(the ¬rst-order derivative of uh ( x) should exist). Suppose that the nodes xj within
an element are numbered j = 1, . . . , n and that the associated shape functions Ni
are numbered i = 1, . . . , n. In that case, for consistency, the shape functions must
be chosen such that
Ni ( xj ) = δij , (14.23)
with
δij = 0 if i = j
= 1 if i = j. (14.24)



14.5 Galerkin approximation

To transform the weak form into a linear set of equations in order to derive an
approximate solution the following steps are taken.

Step 1. Element division As shown in the previous section the domain may
be split into a number of subdomains e , elements, and within each element a
Numerical solution of one-dimensional diffusion equation
240

polynomial interpolation can be made of the function u. An example of such an
element division is given in Fig. 14.4. This distribution of elements is called a
mesh. Then, the integration over the domain can be performed by summing up
the integrals over each element. Consequently, Eq. (14.14) yields:
Nel Nel
dw du
dx = wf dx + B,
c (14.25)
dx dx
e e
e=1 e=1

where Nel denotes the number of elements.



Step 2. Interpolation Suppose that the domain has been divided into three
linear elements, as depicted in Fig. 14.4. Then the nodal values ui may be collected
in an array u:

⎡ ¤
u1
⎢ ⎥
⎢ u2 ⎥
u=⎢ ⎥. (14.26)
⎣ u3 ¦


u4

The unknowns associated with each of the elements are collected in the arrays
e
ue , such that



u1 u2 u3
u1 = u2 = u3 =
, , . (14.27)
∼ ∼ ∼
u2 u3 u4

So, it is important to realize that each particular element array ue contains a subset

of the total, or global, array u. Within each element array ue a local numbering
∼ ∼
may be used, such that for this particular example with linear elements:


<<

. 47
( 67 .)



>>