Advanced Numerical Analysis Chapter 6
6.4.2 Uniform meshes in 2D
We can compute fpx, yq “ ´∆uex explicitly. Moreover, uex “ 0 on BΩ. Then, we consider the homogeneous
Dirichlet problem:
´ ∆u “ f, Ω; (6.30a)
u “ 0, BΩ. (6.30b)
It is immediate to see that uex satisfies the boundary condition. Therefore, by the construction of f, we see
that the exact solution of the above problem is uex.
Exercise 6.8. Show that u R H2
pΩq. (You do not have to show that u P H1` 2
3 ´ǫ
pΩq for any small 0 ă ǫ ă
2
3
.)
Exercise 6.9. Show that ∆pr
α sinpαθqq “ 0.
Exercise 6.10. Compute fpx, yq explicitly.
The 1st step of the FEM is to formulate Problem (6.30) in a weak form. Let V “ H1
0
pΩq. Then the weak
formulation is to find u P V such that
apu, vq “ pf, vq @ v P V. (6.31)
The 2nd step is to define the triangulation Th of Ω. We begin with a uniform triangulation with n and mesh
size h “
1
n
, which is composed of cubes of size h ˆ h. The smartest way of defining the vertex points of the
cubes which are in Ω may be as follows:
pxj , ykq “ hpj ´ n, k ´ nq,
$
&
%
j “ 0, ¨ ¨ ¨ , 2n ; k “ 0, ¨ ¨ ¨ , n;
j “ 0, ¨ ¨ ¨ , n ; k “ n ` 1, ¨ ¨ ¨ , 2n.
(6.32)
Let use denote Rjk “ pxj´1, xj q ˆ pyk´1, ykq.
Exercise 6.11. If g P C
0
pRjkq, for numerical approximation of the integral ş
Rjk
gpx, yq dx use the 2ˆ2 Gauss
quadrature rule. Give the explicit formula for this 2 ˆ 2 Gauss quadrature rule in terms of the coordinates
of Rjk “ pxj´1, xj q ˆ pyk´1, ykq and g.
The 3rd step is to define the finite element space Vh, associated with Th, which approximates V. Notice
that we are solving the Dirichlet boundary value problem, and therefore the basis functions associated with
Th can be given for every interior vertex points pxj , ykq. Denote by Vi
h
and Vb
h
the set of all interior and
boundary vertices of Th.
Let us use the Park-Sheen P1-nonconforming quadrilateral finite element, which consists of piecewise
©Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr) 143
Advanced Numerical Analysis Chapter 6
linear polynomial on each cube that is continuous at the midpoint of interior edges. The benefit of the
Park-Sheen element is that the gradients of the basis function are a constant vector on each cube, so
that the integration on each element is trivial. Explicitly, let Vh be the the Park-Sheen P1-nonconforming
quadrilateral finite element space:
Vh “ tv P L
2
pΩq | v|Q P P1pQq @ Q P Th; rrvssepmq “ 0 @ midpoint m of interior edge e in Th;
vpmq “ 0 @ midpoint m of boundary edge in Thu,
where rrvssepmq denotes the jump of v across the interface e of Ql and Qr at m.
A set of basis functions for Vh consists of those functions φjk associated with the interior vertices pxj , ykq P
Ω defined as follows. First, notice that there are four rectangles which share an interior vertex pxj , ykq in
Th. Let us call them as
RI
jk “ pxj , xj`1q ˆ pyk, yk`1q, RII
jk “ pxj´1, xj q ˆ pyk, yk`1q, RIII
jk “ pxj´1, xj q ˆ pyk´1, ykq, RIV
jk “
pxj , xj`1q ˆ pyk´1, ykq.
The closure of these four rectangles is the support of φjk. That is, φjk is defined as zero on ΩzYι“I,II,III,IV Rι
jk,
and it is the defined such that it is piecewise linear and it has value 1 at the four midpoints of interior edges
of these four rectangles and value 0 at the eight midpoints of exterior edges of these four rectangles. Notice
that the four midpoints of interior edges of these four rectangles are pxj` 1
2
, ykq pxj , yk` 1
2
q, pxj´ 1
2
, ykq, and
pxj , yk´ 1
2
q. The eight midpoints of exterior edges of these four rectangles are
pxj`1, yk` 1
2
q, pxj` 1
2
, yk`1q on BRI
jk;
pxj´ 1
2
, yk`1q, pxj´1, yk` 1
2
q on BRII
jk;
pxj´1, yk´ 1
2
q, pxj´ 1
2
, yk´1q on BRIII
jk ;
pxj` 1
2
, yk´1q, pxj`1, yk´ 1
2
q on BRIV
jk .
Exercise 6.12. Define the basis function φjkpx, yq on each Rι
jk, ι “ I, II, III, IV. Also, write the formula
for the gradient of the basis function, ∇φjkpx, yq on each Rι
jk, ι “ I, II, III, IV.
By finishing defining the FEM basis functions, we have
Vh “ Spantφjk | pxj , ykq P Vi
h
u,
recalling that Vi
h
denotes the set of all interior vertices of Th.
The 4th step is to set up the linear system using the basis functions φjk in Vh. We will abuse the notation
as follows:
jk P Vi
h ðñ pxj , ykq P Vi
h
.
We will integrate on Ω by integrating over each Q P Th and summing up them. That is,
ahpu, vq :“
ÿ
QPTh
p∇u, ∇vqQ “
ÿ
QPTh
ż
Q
∇u ¨ ∇v dx .
The the finite element solution is to find uh “
ř
jkPVi
h
Ujkφjk such that
ahpuh, vq “ pf, vq @ v P Vh.
©Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr) 144
Advanced Numerical Analysis Chapter 6
This is equivalent to finding the coefficient vector pUjkqjkPVi
h
such that
ÿ
jkPVi
h
ahpφjk, φlmqUjk “ pf, φlmq @ lm P Vi
h
. (6.33)
Exercise 6.13. For n “ 10, 20, 40, 80, 160 with h “
1
n
, use the CGM (COnjugate Gradient Method) to solve
(6.33) for pUjkqjkPVi
h
. You do not have to report the vectors pUjkqjkPVi
h
. Then, compute the error
Eh “
gffe
ÿ
jkPVi
h
ż
Rjk
|uhpx, yq ´ uexpx, yq|2 dx (6.34)
Then, compute the error reduction ratios
log2
E2h
Eh
, h “
1
10
,
1
20
,
1
40
,
1
80
,
1
160
.
Here, to appoximate ş
Rjk
gpx, yq dx use the 2 ˆ 2 Gauss quadrature rule.
Exercise 6.14. In the above exercise, we used the Park-Sheen P1-nonconforming FEM. But this time, repeat
the same exercise by using the conforming Q1-element (or bilinear element). For this, the space Vh in (6.33)
should be replaced by
Vh “ tv P C
0
pΩq | v|Q P Q1pQq @ Q P Th; v “ 0 on BΩu.
Here, and in what follows, for domain K Ă R
2
, Q1pKq “ Spant1, x, y, xyu “ ta0 ` a1x ` a2y ` a3xy |
a0, a1, a2, a3 P Ru. Then, in this Q1-conforming element case, define the basis function φjkpx, yq P Vh on
each Rι
jk, ι “ I, II, III, IV so that
φjk|Rι
jk
P Q1pR
ι
jkq @ ι “ I, II, III, IV
and
φjkpxl
, ymq “ δjlδkm @ lm P Vi
h
.
Of course, you need to have the gradient ∇φjk, too.
6.4.3 Graded meshes in 2D
Continuing from the previous subsection, now we proceed to work on graded meshes.
As in the 1D example, we attempt to introduce graded meshes on the corner singular elements.
This procedure is a recursive procedure.
They are three elements of size h ˆ h from Th that share their vertices with p0, 0q :
Rnn “ pxn´1, xnq ˆ pyn´1, ynq “: G
p0q
c
,
Rn`1,n “ pxn, xn`1q ˆ pyn´1, ynq “: G
p0q
r
,
Rn,n`1 “ pxn´1, xnq ˆ pyn, yn`1q “: G
p0q
t
.
©Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr) 145
Advanced Numerical Analysis Chapter 6
Let us denote the above three rectangles by G
p0q
c , Gp0q
r , Gp0q
t and their centers by g
p1q
c , g
p1q
r , g
p1q
t
, respectively.
(“c,r,t” mean center, right, top.) Also the superscript pkq denotes the k-th level of mesh grading. Using the
centers g
p1q
ι , ι “ c, r, t, we subdive each of the rectangles G
p0q
ι , ι “ c, r, t, into two quadrilaterals Q
p1q
ι,1
, Qp1q
ι,2
,
and a cube G
p1q
ι of size h{2. The cubes G
p1q
ι are supposed to have common vertex p0, 0q and edges are parallel
to the x- and y-axes.
In this procedure, there are two new vertices, denoted by g
p1q
h “ hp´1
2
, 0q and g
p1q
v “ hp0, ´
1
2
q (“h,v”
represent horizontal and vertical).
By this, we have 5 more interior vertices g
p1q
ι , ι “ c, r, t, h, v and three original rectangles of size h are
subdivided to 6 quadrilaterals and 3 cubes G
p1q
ι of size h
2
.
We repeat this procedure for G
p1q
ι to get 3 cubes G
p2q
ι of size h
2
2 with 5 more vertices g
p2q
ι , ι “ c, r, t, h, v.
Repeating this procedure for s steps, we have 5s new graded vertices g
pkq
ι , ι “ c, r, t, h, v for k “ 1, 2, ¨ ¨ ¨ , s
in addition to the uniform vertices pxj , ykq
1
s.
For each of these graded vertices, we can creat new basis functions ψ
pkq
ι , ι “ c, r, t, h, v for k “ 1, 2, ¨ ¨ ¨ , s.
Hence we will have modified finite element space associated with the graded meshes:
V
psq
h “ Spantφjk,pj, kq P Vi
h
; ψ
pkq
ι
, ι “ c, r, t, h, v for k “ 1, 2, ¨ ¨ ¨ , su.
Here, observe that the five basis functions φn´1,n`1, φn´1,n, φn´1,n´1, φn,n´1, φn`1,n´1 are modified from
those associated with the uniform mesh due to the modification of the origianl uniform mesh, and they have
new supporting quadrilaterals.
The finite element solution is then to find uh P V
psq
h
of the form
uh “
ÿ
jkPVi
h
Ujkφjk `
ÿs
k“1
ÿ
ι“c,r,t,h,v
Gkιψ
pkq
ι
.
which satisfies
ÿ
jkPVi
h
ahpφjk, φlmqUjk `
ÿs
k“1
ÿ
ι“c,r,t,h,v
ahpψ
pkq
ι
, φlmqGkι “ pf, φlmq @ lm P Vi
h
, (6.35a)
ÿ
jkPVi
h
ahpφjk, ψppq
q
qUjk `
ÿs
k“1
ÿ
ι“c,r,t,h,v
ahpψ
pkq
ι
, ψppq
q
qGkι “ pf, ψppq
q
q (6.35b)
@ p “ 1, ¨ ¨ ¨ , s, q “ c, r, t, h, v.
In (6.35), observe those terms without common support are zeros. So the above system is sparse.
Exercise 6.15. You can choose to use either the Park-Sheen P1-nonconforming element or the conforming
Q1 element. Let the grading step number to be s “ 5.
1. Give the explicit formula of the basis functions and their gradients φn´1,n`1, φn´1,n, φn´1,n´1, φn,n´1, φn`1,n´1.
Also give the explicit formula of the basis functions and their gradients ψ
pkq
ι , ι “ c, r, t, h, v for k “
1, 2, ¨ ¨ ¨ , su.
©Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr) 146
Advanced Numerical Analysis Chapter 6
2. Use the Conjugate Gradient Method to solve (6.35) in order to find the coefficients the coefficients
pUjkqjkPVi
h
and pGkιqι“c,r,t,h,v;k“1,2,¨¨¨ ,s.
Finally compute the error reduction ratios