<<

. 59
( 67 .)



>>



a=0




(a)

a=1




(b)

a = 10




(c)

a = 25




(d)

a = 100




(e)

Figure 16.10
Contour lines of constant u values for a range of values of a. With increasing a the effect of
convection increases.
291 Exercises

cases 10 contour lines are shown ranging from 0.1 to 1 with increments of 0.1.
The effect of an increasing velocity is clearly demonstrated.



Exercises

16.1 The weak form of the two-dimensional convection-diffusion equation is
given by
‚u
+ wv · ∇u + ∇w · ( c∇u)
w d
‚t
= + wn · c∇u = d .
wf d

After discretization, the element inertia matrix is de¬ned as

Me = N NT d ,
∼∼
e

while the element matrix related to the convective part is given by
dN T dN T
Ce = + vy ∼

N vx d.

dx dy
e

Consider the element as depicted in the ¬gure below.

y
4 3


x



1 2


This is a bi-linear element. The element spans the spatial domain ’2 ¤
x ¤ 2 and 2 ¤ y ¤ 2.
(a) Compute the element inertia matrix M e using a 2—2 Gauss integration
rule. It is recommended to use MATLAB for this computation.
(b) Suppose that the location of the integration points coincides with the
nodes of the element. What is the element inertia matrix M e in this
case?
(c) Compute the matrix Ce if vx = 1 and vy = 0, using a 2 — 2 Gauss
integration rule.
(d) Suppose that along the edge of the element located at x = 2 a constant
¬‚ux q = n · c∇u is prescribed. Then compute the column
The three-dimensional convection-diffusion equation
292


Nq d ,

e

for this edge, contributing to the element column f e .

16.2 Consider the mesh depicted in the ¬gure below. The solution vector u con-

tains the nodal solutions in the sequence de¬ned by the node numbers,
hence:

uT = [ u1 u2 · · · u16 ] .




10 12 14 16

3 6 9



5 7 9 15

2 5 8



13
2 4 8

1 4 7



6 11
1 3


Let the element topology array of the third element be given by
top(3,:)=[5 7 12 10 1 1].
(a) Is this topology array unique?
(b) What is the pos array for this element?
(c) De¬ne the dest array.
(d) The solution u is stored in the array sol. How are the nodal solutions

of element 8 extracted from sol using the array pos?
(e) Suppose that the array nodes contains the node numbers of the
left edge, say nodes=[1 2 5 10]. How are the nodal solutions
of these nodes extracted from the solution array sol via the array
dest?
(f) Consider the third element whose element topology array is given by
top(3,:)=[5 7 12 10 1 1]. If the nodal solution array ue for ∼
this element is given by:

uT = [ 1 4 3 7] ,
∼e


compute the solution u for ξ = 1 and · = 3 .
4 4
(i) Suppose the element topology array is given by
top(3,:)=[12 10 5 7 1 1]. What is the corresponding ele-
ment nodal solution array ue ?

293 Exercises

16.3 Consider the square domain 0 ¤ x ¤ 1, 0 ¤ y ¤ 1. Along the line x = 0
the boundary condition u = 0 is imposed, and along the line x = 1 the
boundary condition u = 1 is prescribed. Along the other boundaries the
natural boundary condition n · c∇u = 0 is speci¬ed. Consider the steady
convection diffusion problem on this domain.
(a) Suppose that the convection-diffusion problem represents the tem-
perature equation. What physical meaning does the above natural
boundary condition have?
(b) Adjust the m-¬le demo_cd such that the above problem is solved
using 5 — 5 linear elements. Select the diffusion constant c = 1 and
the convective velocity v = ex (Remark: see the element m-¬le elcd
and the associated m-¬le elcd_a). Specify the structured array mat
according to this.
(c) Extract the solution along the line y = 0. Hint, ¬rst select the nodes
using usercurves containing the nodes along all the curves that
have been de¬ned for the mesh generator crmesh. Then use dest
to extract the relevant solution components.
(d) Compare the solution along y = 0 with the one-dimensional solution,
using the program fem1dcd.
16.4 Consider the domain as sketched in the ¬gure. Along the boundary 1 the
speci¬cation u = 1 is chosen, while along 7 u = 0 is chosen as a boundary
condition. Along all other boundaries n · c∇u = 0 is imposed. Choose the
convective velocity v such that:

y<0: v= 0
y≥0: v = y ex

Furthermore, choose c = 1 throughout the domain. Model this problem
using femlin_cd. Approximately beyond which value of x along the
boundary 6 holds u > 0.8?

3
“6
2

“7 “5
1
y




“3
0
“4 “2
“1
’1

’2
0 2 4 6 8 10
x
The three-dimensional convection-diffusion equation
294

16.5 If the ¬‚ux n · c∇u = P along is prescribed, the integral (element level):
p

NP d ,

e

needs to be computed. It may be transformed to the local coordinate system
ξ ∈[ ’1 1], such that:
‚x
1
NP d = NP dξ ,
‚ξ
∼ ∼
’1
e

where it is assumed that e is oriented in x-direction.
If linear elements are chosen, the shape functions along the boundary are
given by:

2 (1 ’ ξ )
1
N= .
2 (1 + ξ )
1




Use an isoparametric element and let L denote the length of e, then
demonstrate that for constant P:
LP 1
NP d = .

1
2
e
17 Shape functions and numerical
integration

17.1 Introduction

In the previous chapter the shape functions Ni have hardly been discussed in any
detail. The key purpose of this chapter is ¬rst to introduce isoparametric shape
functions, and second to outline numerical integration of the integrals appearing
in the element coef¬cient matrices and element column. Before this can be done
it is useful to understand the minimum requirements to be imposed on the shape

<<

. 59
( 67 .)



>>