Introduction to the conforming and nonconforming finite element methods Dongwoo Sheen

Graduate School of Informatics Kyoto University Kyoto 608–8501, Japan

Email: [email protected] http://www.nasc.snu.ac.kr July 14, 2015: 13:14

Finite element methods

Section 0.0

c

2012–2015 Dr. Dongwoo Sheen All rights reserved.

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

Page 2 of 65

Contents 1 Introduction to conforming and nonconforming finite elements 1.1 The model problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.1 Weak formulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.2 Unique Solvability . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.3 Elliptic regularity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1.4 Equivalence of three formulations . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Finite Element Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.1 Th : triangulation of Ω . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.2 Construction of FE space . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.3 Compute FEM solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.4 Convergence study: Galerkin method . . . . . . . . . . . . . . . . . . . . . 1.2.5 Finite elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Construction of basis functions and calculation of matrix components . . . . . . . 1.3.1 The 2–dimensional P1 triangular conforming element (the Courant element) 1.3.2 Calculation of matrix components . . . . . . . . . . . . . . . . . . . . . . . 1.3.3 Implementation of P1 conforming triangular element on a uniform grid . . . 1.4 Higher order finite elements in Rd . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.1 Pm (conforming) simplicial element . . . . . . . . . . . . . . . . . . . . . . . 1.4.2 Qm (conforming) d-linear element . . . . . . . . . . . . . . . . . . . . . . . 1.4.3 Unisolvency and optimality . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4.4 Triangulations and reference elements . . . . . . . . . . . . . . . . . . . . . P 1.5 Assembly of matrix Akj = K∈Th ah (φj , φk ) . . . . . . . . . . . . . . . . . . . . .

5 . . . . . . 5 . . . . . . 5 . . . . . . 6 . . . . . . 7 . . . . . . 7 . . . . . . 9 . . . . . . 9 . . . . . . 10 . . . . . . 11 . . . . . . 12 . . . . . . 13 . . . . . . 14 and the Crouzeix–Raviart P . . . . . . 15 . . . . . . 16 . . . . . . 18 . . . . . . 18 . . . . . . 18 . . . . . . 19 . . . . . . 19 . . . . . . 23

2 NC (Nonconforming) finite elements 25 2.1 P1 triangular NC element (Crouzeix-Raviart, 1973) and the rotated Q1 -rectangular NC elements 26 2.2 The P1 -NC quadrilateral and hexahedral elements . . . . . . . . . . . . . . . . . . . . . . . . 28 2.2.1 Error analysis for linear nonconforming Galerkin method . . . . . . . . . . . . . . . . 29 2.2.2 Implementation of P1 –NC quadrilateral element . . . . . . . . . . . . . . . . . . . . . 30 2.3 Hermite–type finite elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

3 Preliminaries from functional analysis 3.1 Duality . . . . . . . . . . . . . . . . . 3.2 Annihilators . . . . . . . . . . . . . . . 3.3 Big Theorems . . . . . . . . . . . . . . 3.4 Compact operators . . . . . . . . . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

. . . .

39 40 42 42 43

4 Sobolev Spaces 45 4.1 Distributions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 4.2 Sobolev Spaces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 4.2.1 Schwartz space S(Rd ) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 3

Finite element methods

4.3 4.4

4.2.2 Fourier transformation . . . . . . . . . . . . . . . . . . . . . . . 4.2.3 Prolongation (Extension) . . . . . . . . . . . . . . . . . . . . . 4.2.4 The definition of traces and trace theorems . . . . . . . . . . . 4.2.5 Application of Trace Theorem . . . . . . . . . . . . . . . . . . . The embedding theorems and fractional Sobolev spaces . . . . . . . . Application to elliptic problems and finite element subspace of H 1 (Ω)

Section 0.0 . . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

. . . . . .

50 52 53 56 57 60

5 Appendix 61 5.1 Brief review on mathematical analysis and linear algebra . . . . . . . . . . . . . . . . . . . . . 61

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

4

Chapter 1

Introduction to conforming and nonconforming finite elements 1.1

The model problem

Throughout the lecture F denotes the real number field R or the complex number field C. Let Ω ⊂ Rd be a bounded open region with smooth or polygonal boundary ∂Ω. Consider the following boundary value problem: given f ∈ H −1 (Ω) find u ∈ H01 (Ω) such that L u := −∇·A(x)∇u + b(x)u

=

f,

Ω,

(1.1a)

u =

0,

∂Ω.

(1.1b)

• Uniform ellipticity assumption: A(·) is a uniformly positive–definite d × d matrix–valued function defined on Ω: i.e., hA(x)ξ, ξi

b(x)

1.1.1

≥ c0 |ξ|2 ≥ 0

a.e.

∀ξ ∈ F d xΩ.

a.e.

x ∈ Ω;

(1.2a) (1.2b)

Weak formulation

Multiply (1.1a) by an arbitrary v ∈ H01 (Ω) and use the Divergence Theorem, keeping in mind (1.1b). Then we have the following weak problem: to find u ∈ H01 (Ω) such that (A∇u, ∇v) + (bu, v) = hf, viH −1 ,H 1 0

∀v ∈ H01 (Ω).

(1.3)

Here, and in what follows, standard notations from functional analysis are used. • The duality paring h·, ·iV ′ ,V , where V and V ′ are a normed linear space V and its dual space V ′ , respectively. The corresponding norms are denoted by k · kV and k · kV ′ . The duality is defined as hf, viV ′ ,V = f (v) for all v ∈ V such that f : V 7→ F is a continuous linear functional on V. The dual norm is given by hf, viV ′ ,V , kf kV ′ = sup kvkV 06=v∈V for all f ∈ V ′ .

• The space of square integrable functions L2 (Ω) on Ω is denoted by Z 2 L (Ω) = {v : Ω → F : |v(x)|2 dx < ∞}, Ω

5

Finite element methods

Section 1.1

with the L2 (Ω)–inner product: (u, v) = 2

and the L (Ω)–norm: kvk =

Z

∀u, v ∈ L2 (Ω)

u(x)v(x) dx Ω

sZ

∀v ∈ L2 (Ω).

|v(x)|2 dx

Ω

• The Sobolev spaces, H m (Ω) = {v ∈ L2 (Ω) : kvkm < ∞}, are Hilbert spaces with the inner product X  ∂αu ∂αv  . , (u, v)m,Ω = ∂xα ∂xα |α|≤m

Pd Here, α = (α1 , · · · , αd ) ∈ Zd+ is a multi–index with its length |α| = j=1 αj , which is used to denote the partial derivatives ∂αv ∂α ∂α = 1 · · · d v. α ∂x ∂x1 ∂xd The corresponding norms and seminorms are designated by v u X q ∂αv u k α k2 , kvkm,Ω = (v, v)m,Ω and |v|m,Ω = t ∂x |α|=m

respectively.

• In particular, for m = 1, H 1 (Ω) = {v ∈ L2 (Ω) :

∂v ∈ L2 (Ω), j = 1, · · · , N }, ∂xj

is equipped with the inner product, norm, and seminorm: (u, v)1,Ω =

Z

uv dx + Ω

k=1 1

respectively.

d Z X

kuk1,Ω = [(u, u)1,Ω ] 2 ,

and

Ω

∂u ∂v dx, ∂xk ∂xk 1

|u|1,Ω = [(∇u, ∇u)] 2 ,

• H0m (Ω) denotes the completion of C0∞ (Ω) in the norm of H m (Ω) • Often, H −1 (Ω) will mean the dual space of H01 (Ω) with the dual norm kf k−1,Ω =

1.1.2

hf, viH −1 (Ω),H 1 (Ω)

sup 06=v∈H01 (Ω)

kvk1,Ω

.

Unique Solvability

Set V = H01 (Ω) and define a bilinear form a(·, ·) : V × V 7→ F and a linear functional ℓ : V 7→ F by a(u, v) = (A∇u, ∇v) + (bu, v),

and

ℓ(v) = hf, viH −1 (Ω),H 1 (Ω)

∀ u, v ∈ V.

Then, (1.3) is rewritten as to find u ∈ V such that a(u, v) = ℓ(v)

∀v ∈ V.

(1.4)

We now state the following frequently used lemma: c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

6

Finite element methods

Section 1.1

Lemma 1.1 (The Lax-Milgram Lemma). Assume the following 1. V is a Hilbert space; 2. a(·, ·) is continuous on V , i.e. there exists C > 0 such that |a(u, v)| ≤ CkukV kvkV

∀u, v ∈ V ;

3. a(·, ·) is V –coercive on V , i.e. there exists α > 0 such that |a(u, u)| ≥ αkuk2V

∀u ∈ V ;

4. ℓ(·) is continuous on V , i.e. there exists M > 0 such that |ℓ(v)| ≤ M kvkV

∀v ∈ V.

Then, there exists a unique solution u ∈ V to (1.4) such that kukV ≤

M kℓkV ′ . α

Exercise 1.1. Prove the above Lax–Milgram Lemma.

1.1.3

Elliptic regularity

If the domain Ω and the coefficients A(x) and b(x) are smooth enough and f ∈ L2 (Ω), then the above solution u ∈ H01 (Ω) to (1.4) belongs to H 2 (Ω) ∩ H01 (Ω). u) for all u ∈ H01 (Ω) and c ∈ L∞ (div; Ω). Here, Exercise 1.2. Show that (c · ∇u, u) = − 21 ((∇·c)u, ∞ ∞ d ∞ L (div; Ω) = f ∈ [L (Ω)] | ∇·f ∈ [L (Ω)] .

Exercise 1.3. Consider the problem:

− ∇·A(x)∇u + c(x) · ∇u + b(x)u = u =

f, 0,

Ω, ∂Ω,

(1.5a) (1.5b)

where A and b satisfy the above assumptions (1.2) and c(x) fulfills the following assumption: 1 ∇·c(x) ≤ b(x) 2

a.e. ∈ Ω.

Then, using the above Lax–Milgram Lemma, show that the above problem (1.5) has a unique solution in H01 (Ω).

1.1.4

Equivalence of three formulations

Introduce the quadratic functional E : V → R given by E(v) =

1 1 a(v, v) − hf, viV ′ ,V = a(v, v) − ℓ(v), 2 2

(1.6)

where E is called a complementary potential energy. Consider the minimization problem: min E(v),

(1.7)

v∈V

or equivalently, find u ∈ V such that The following lemma is useful.

E(u) ≤ E(v)

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

∀v ∈ V. 7

Finite element methods

Section 1.1

Lemma 1.2. If w is continuous on Ω and Z w v dx = 0 Ω

∀v ∈ H01 (Ω).

Then w(x) = 0 for all x ∈ Ω. Proof. Exercise. Hint: use ǫ − δ definition of continuous functions and C0∞ (Ω) ⊂ H01 (Ω), where C0∞ (Ω) denotes the space of infinitely differentiable functions with compact support in Ω. The support supp(v) of v means the set {x ∈ Ω | v(x) 6= 0}. Since Ω is a subset of Rd , the compactness means supp(v) is bounded.  Remark 1.1. Recall elliptic regularity for higher dimensional elliptic problems: Let u be the solution of a(u, v) = (f, v), v ∈ V , where a(v, w) =

 d  X ∂v ∂w aij + (c v, w) , ∂xi ∂xj i,j=1

V = {v ∈ H 1 (Ω) : v = 0

on ∂Ω}

If Ω is a convex polygonal domain and f ∈ L2 (Ω), then the variational solution u ∈ H 2 (Ω)∩H01 (Ω), and kuk2,Ω ≤ C kf k0,Ω

[Gri85]

(1.8)

If ∂Ω is smooth, and f ∈ H k (Ω), then the variational solution u ∈ H k+2 (Ω) ∩ H01 (Ω), and kukk+2,Ω ≤ C(k, Ω) kf kk,Ω

for k = 0, 1, · · ·

(1.9)

Theorem 1.1. Let u, v ∈ H 1 (Ω) and A ∈ L∞ (Ω). Moreover, assume that ∇·(A∇u) ∈ L2 (Ω). Then the normal trace ν · (A∇u) ∈ H −1/2 (∂Ω), and the following integration by parts formula holds: (A∇u, ∇v) = −(∇·(A∇u), v) + hν · (A∇u), viH −1/2 (∂Ω),H 1/2 (∂Ω)

∀v ∈ H 1 (Ω).

(1.10)

Theorem 1.2. The three problems (1.1), (1.4) and (1.7) are equivalent under a sufficient regularity condition such as ∇·(A∇u) ∈ L2 (Ω) and f ∈ L2 (Ω). Proof. Step 1. (1.1)⇒(1.4). This was already shown above. Step 2. (1.4)⇒(1.1). Assume that v ∈ V satisfies (1.4). Generally, this doesn’t hold. But, if ∇·A∇u ∈ L2 (Ω), one can use the integration by parts Theorem 1.1. Thus, Z ∀v ∈ V. (L u − f ) v dx = 0, Ω

By the above Lemma 1.2, L u − f = 0 in V ′ . Step 3. (1.4)⇒(1.7). Suppose u is a solution of (1.4). Then, for an arbitrary v ∈ V , we have   1 1 a(v, v) − ℓ(v) − a(u, u) − ℓ(u) E(v) − E(u) = 2 2 1 = (a(u, u) − 2a(u, v) + a(v, v)) 2 1 = a(u − v, u − v) ≥ 0, 2 where we used a(u, u − v) = ℓ(u − v), since u solves (1.4). c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

8

Finite element methods

Section 1.2

Step 4. (1.7)⇒(1.4). Suppose u is a solution of (1.7). Then for any ǫ ∈ R, E(u) ≤ E(u + ǫv) ⇐⇒ g(ǫ). Since g(ǫ) has a minimum at ǫ = 0, g ′ (0) = 0. Therefore,   d 1 a(u + ǫv, u + ǫv) − hf, u + εviV ′ ,V g ′ (0) = dǫ 2 ǫ=0 = a(u, v) − hf, viV ′ V = 0.  Theorem 1.3. A solution to (1.4) is unique if it exists, under the stated assumptions. Proof. Suppose u1 and u2 ∈ V solve (1.4). Then a(u1 , v) a(u2 , v) By subtracting, a(u1 − u2 , v) = 0 Choose v = u1 − u2 ∈ V . Then

= =

hf, viV ′ ,V , hf, viV ′ ,V

(1.11) (1.12)

for all v ∈ V.

for all v ∈ V . a(u1 − u2 , u1 − u2 ) = 0,

(1.13)

from which by the V –coercivity of a, ku1 − u2 kV ≤ Ca(u1 − u2 , u1 − u2 ) = 0 Hence, u1 (x) − u2 (x) ⇐⇒ constant for all x ∈ Ω. Since u1 (0) − u2 (0) = 0, u1 (x) − u2 (x) = 0, for all x ∈ Ω. 

1.2

Finite Element Methods

The finite element methods to solve (1.3) usually consist of the following procedure. 1. (Th )h : Triangulation of the domain Ω 2. (Vh )h : Construction of FE space: Vh → V as h → 0 3. FE solution: Find uh ∈ Vh for each h 4. Error Analysis: Convergence of uh → u as h → 0 5. Numerical Simulation

1.2.1

Th : triangulation of Ω Decompose Ω into a union of finite number of elementary geometry Standard 2−dimensional FEM mesh−grid

Standard 2−dimensional FEM mesh−grid

y

y

(1,1) (0,1)

(1,1) (0,1)

x

(0,0)

(1,0)

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

x

(0,0)

(1,0)

9

Finite element methods

Section 1.2 Example of Conforming Two−Dimensional Mesh

Example of Two-Dimensional Mesh

1.2.2

Construction of FE space Based on Th construct finite dimensional subspace Vh

Example of Conforming Two−Dimensional Mesh

Piecewise linear basis function at (x , y ) i j

(x (x ,y ) i j+1

,y ) i+1 j+1

5 4 6 (x

,y ) i-1 j

(x ,y ) i j

(x

,y ) i+1 j

3

1 2

(x ,y ) i j-1 (x

,y ) i-1 j-1

• Let Vh ⊆ V = H01 (Ω) given by Vh = {v ∈ C 0 (Ω) : v|K is piecewise linear for all K ∈ Th } for triangular decomposition Th or Vh = {v ∈ C 0 (Ω) : v|K is piecewise bilinear for all K ∈ Th } for quadrilateral decomposition Th • In general, the function spaces are based on piecewise polynomials such that v|K ∈ Pℓ (K) for some ℓ ≥ 0. • Let ϕj , j = 1 · · · , N, be the basis functions for Vh Vh

=

Span{ϕj , j = 1, · · · , N }

=

{v =

N X j=1

cj ϕj : cj ∈ F , j = 1, · · · , N }

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

10

Finite element methods

Section 1.2

• Let ψj , j = 1 · · · , M, be the basis functions for Wh

1.2.3

{w =

=

Wh

M X j=1

cj ψj : cj ∈ F , j = 1, · · · , M }

Compute FEM solution

1. Find uh ∈ Vh such that a(uh , v) = ℓ(v) 2. Find uh =

PN

j=1

∀v ∈ Wh .

αj ϕj such that N X

a(ϕj , v)αj = ℓ(v)

j=1

3. Find uh =

PN

j=1

(1.14)

∀v ∈ Wh .

αj ϕj such that N X

a(ϕj , ϕk )αj = ℓ(ϕk )

j=1

4. Denote b = (ℓ(ϕ1 ), · · · , ℓ(ϕM ))t

∀ϕk ∈ Wh .

α = (α1 , · · · , αN )t ; Akj = a(ϕj , ϕk ), k = 1, · · · , M, j = 1, · · · , N

5. Solve the M × N linear system: Aα = b.

(1.15)

6. If M > N, α minimizes kAα′ − bk • Vh : solution space,

Wh : test function space

• Galerkin method (Ritz method) if Vh =Wh in (2) • Petrov-Galerkin method: if Vh 6=Wh in (2) • If either Vh or Wh is not a subspace of H01 (Ω), Nonconforming FEM. Linear solvers 1. Direct method: Gaussian elimination (LU-decomposition) method, Cholesky decomposition method, QR-decomposition method 2. Iterative method: Gauss-Jacobi and Gauss-Seidel methods, SOR, Richardson iteration method, Gradient method, Conjugate Gradient Method, BiConjugate Gradient Method, Domain Decomposition Method, MultiGrid solver,.... c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

11

Finite element methods

1.2.4

Section 1.2

Convergence study: Galerkin method Error equation

Subtracting (1.3) from (1.14), the error equation holds: a(uh − u, ϕ) = 0 ∀ϕ ∈ Vh .

(1.16)

Lemma 1.3 (Ce´a’s lemma). kuh − uk1 ≤

C inf kv − uk1 . α v∈Vh

(1.17)

In particular, v = Ph u, “the Vh -interpolation of u”, leads to kuh − uk1 ≤

C ku − Ph uk1 . α

(1.18)

Proof. From (1.16), coercivity, and boundedness; for any v ∈ Vh , αkuh − uk21 ≤ a(uh − u, uh − u) = a(uh − v + (v − u), uh − u) = a(v − u, uh − u)

≤ Ckv − uk1 ku − uh k1 .

Therefore, by dividing by u − uh k1 , one gets (1.17). In particular, v = Ph u, “the Vh -interpolation of u”, leads to (1.18)  By the interpolation property ku − Ph uk1 ≤ C1 h|u|2 , kuh − uk1 ≤ C1

C h|u|2 . α

Hence, FE error : ellipticity * approximation properties in the FE space

Duality argument for L2 (Ω) error estimate The aim is to show kuh − uk ≤ Ch2 |u|2 , under elliptic regularity conditions. In order to obtain an L2 -error estimate, the following inclusions and identifications are used: H01 (Ω) ֒→ L2 (Ω) ⇐⇒ L2 (Ω)′ ֒→ H −1 (Ω). If φ ∈ L2 (Ω), obviously φ ∈ H −1 (Ω). Then, for all v ∈ H01 (Ω), one has hφ, viH −1 (Ω),H 1 (Ω) = (φ, v)L2 (Ω) . 0

We now state and prove the following lemma: c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

12

Finite element methods

Section 1.2

Lemma 1.4 (Aubin-Nitsche Lemma). Let H be Hilbert space with the norm k · kH and V be a subspace of H which is a Hilbert space with the norm k·kV . Suppose V ֒→ H is continuous. Then for any subspace Vh ⊂ V fulfilling a(uh − u, v) = 0

v ∈ Vh ,

(1.19)

the following estimation is valid: kuh − ukH ≤ Ckuh − ukV sup

g∈H

inf v∈Vh kv − ϕg kV , kgkH

where ϕg ∈ V is the unique solution of a(w, ϕg ) = hg, wiV ′ ,V

∀w ∈ V,

(1.20)

where the bilinear and liner forms satisfy the assumptions in the Lax-Milgram Lemma. Proof. Let g ∈ H ⇐⇒ H ′ ֒→ V ′ . In (1.20), the choice of w = uh − u leads to for any v ∈ Vh , by (1.19) | hg, uh − uiV ′ ,V |

= |a(uh − u, ϕg )| = |a(uh − u, ϕg − v)| ≤ C kuh − ukV kϕg − vkV .

Now, exploiting the Gelfand triplet: (g, uh − u)H = hg, uh − uiV ′ ,V for g ∈ H, one has kuh − ukH = sup

g∈H

inf v∈Vh kϕg − vkV (g, uh − u)H ≤ C kuh − ukV sup . kgkH kgkH g∈H

This proves the lemma.  Corollary 1.1. Let Th be shape-regular and u ∈ V ⇐⇒ H01 (Ω) be the solution of a(u, v) = hf, viV ′ ,V for every v ∈ H01 (Ω). Then, kuh − ukL2 (Ω) ≤ Ch kuh − ukH 1 (Ω) . Moreover, assume that u ∈ H 2 (Ω) ∩ H01 (Ω) satisfies a(u, v) = hf, viV ′ ,V for every v ∈ H 2 (Ω), then kuh − uk0 ≤ Ch2 |u|2 . Proof. Set H = L2 (Ω) and V = H01 (Ω) in the Aubin-Nitsche Lemma, with k · kH = k·k0 and k·kV = k·k1 respectively. Notice that H01 (Ω) ֒→ L2 (Ω) is continuous. Thus by the Aubin-Nitsche Lemma and elliptic regularity,   inf v∈Vh kϕg − vk1 kuh − uk0 ≤ C kuh − uk1 sup kgk0 g∈L2 (Ω) ≤

Ch kuh − uk1 ≤ Ch2 |u|2 .



1.2.5

Finite elements Finite element:(K, PK , ΣK )

1. K: “element”–geometric object, triangle, quadrilateral, simplex, hexahedron, .... 2. PK : “local finite element space”–vector space of polynomials defined on K c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

13

Finite element methods

Section 1.3

3. ΣK : “DOFs”–(Degrees of freedom) to determine an element p ∈ PK uniquely • The Courant element or P1 triangular element and the bilinear element or the Q1 element.

K=T

V1 02 0 1

V1

PK = Span{1, x1, x2}

e3 e1

e2

V2

T

e1

R

K=R

ΣK = {φ(Vj ), j = 1, 2, 3}

PK = Q 1

e3

e2

0 1

0 V1 3

1 0 0 1 V1

V3

ΣK = {φ(Vj ), j = 1, · · · , 4}

e4

V4

• P2 triangular element and Q2 -rectangular element

K=T

V1 02 0 1

e3

PK = P2

T 1 0 0 V1 3

O

1 0 0m2 1

e2

1 0 0 1 V1

V3

e1

R

m1(1, 0)

K=R PK = Q 2

ΣK = {φ(Vj ), φ(mj ), j = 1, 2, e3}3

e1

V1

m2(0, 1) m3(−1, 0)

m3 1 0 0 1

1 0 0m1 1

e2

V2

e4

m4(0, −1)

ΣK = {φ(Vj ), φ(j ), V4 j = 1, · · · , 4; φ(O)}

Here, and in what follows, the following notations for polynomial spaces are used: Pk Qk

= =

αd 1 Span{xα 1 · · · xd | α1 + · · · + αd ≤ k, αj ≥ 0},

1 Span{xα 1

d · · · xα d

| 0 ≤ αj ≤ k}.

(1.21) (1.22)

1.3

Construction of basis functions and calculation of matrix components

1.3.1

The 2–dimensional P1 triangular conforming element (the Courant element) and the Crouzeix–Raviart P1 nonconforming element

Consider a triangle K with vertices Vj = (pj , qj ), j = 1, 2, 3. Let Pj = (xj , yj ), j = 1, 2, 3, be three points on the triangle K. We consider the two cases: 1. Pj = Vj , j = 1, 2, 3 : the P1 triangular conforming element 2. Pj = Mj =

Vj−1 +Vj ,j 2

= 1, 2, 3 : the Crouzeix–Raviart P1 nonconforming element

The basis function φ(j, x, y) associated with Pj is constructed as follows: φ(j, x, y) =

(x − xk )(yk − yl ) − (xk − yl )(y − yk ) , (xj − xk )(yk − yl ) − (xk − xl )(yj − yk )

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

{j, k, l} = {1, 2, 3}.

(1.23) 14

Finite element methods

Section 1.3

Then obviously, φ(j, Pk ) = δjk , where δjk is the Kronecker delta. The gradient of φ(j, x, y) is given as follows:   1 yk − yl ∇φ(j, x, y) = , {j, k, l} = {1, 2, 3}. (1.24) (xj − xk )(yk − yl ) − (xk − yl )(yj − yk ) −(xk − yl ) The area of K is given by the cross product of vectors:  1 1 1 1 |K| = |(V1 − V3 ) × (V2 − V3 )| = |V1 − V3 | |V2 − V3 || sin θ| = det p1 2 2 2 q1

where θ is the angle between the two vectors.

1.3.2

1 p2 q2

 1 p3  , q3

(1.25)

Calculation of matrix components

Now let us compute the matrix component Akj = a(φj , φk ) in (1.15). Z Akj = a(φj , φk ) = (A∇φj ) · ∇φk + b(x)φj φk dx Ω X Z = (A(x)∇φj ) · ∇φk + b(x)φj φk dx K∈Th

X

=:

K

aK (φj , φk ) =:

X

AK kj .

K∈Th

K∈Th

Exercise 1.4. (Due the start of class of July 3) Submit the program code with input data, if there is any, and analysis. 2 Consider the elliptic equation (1.1) with discontinuous coefficient. The domain is taken as Ω = (0, 1)
1 3 2 with Ω− = 4 , 4 and the coefficients A(x) = aχΩ− (x)I + χΩ\Ω− (x)I and b(x) = 1 for all x ∈ Ω. Here, χ is the characteristic function and I denotes the identity matrix. Define ( x − x2 , x ∈ [0, 14 ] ∪ [ 34 , 1],   φ(x) = 3 1 1 1 2 x ∈ [ 41 , 43 ]. 16 + a 16 − (x − 2 ) ,

1 Then set f (x, y) = −∇·(A(x, y)∇(φ(x)φ(y)) + φ(x)φ(y) Try with a = 100 . Divide Ω into Nx × Ny uniform rectangles and then divide each rectangle into two triangles with the diagonal from the top left to the bottom right. Try with Nx × Ny = 50 × 50, 100 × 100, 200 × 200. In each case, calculate the L2 (Ω)–error

sZ

Ω

|uh (x, y) − u(x, y)|2 dxdy,

h = 1/50, 1/100, 1/200.

The linear system can be solved by using the multi-grid method or the Conjugate Gradient Method. Let u(x, y) = U (r, θ). Then we have ∂U ∂r 1 ∂U r ∂θ



=

 ∂u 

=

 cos θ sin θ





cos θ − sin θ

sin θ cos θ

  ∂u  ∂x ∂u ∂y

(1.26)

Hence, ∂x ∂u ∂y

− sin θ cos θ

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)



∂U ∂r 1 ∂U r ∂θ



.

(1.27) 15

Finite element methods

Section 1.3

Thus, ∆u =



∂  ∂x ∂ ∂y

·

 ∂u  ∂x ∂u ∂y

=



cos θ sin θ

− sin θ cos θ



∂ ∂r 1 ∂ r ∂θ

  cos θ · sin θ

− sin θ cos θ



∂U ∂r 1 ∂U r ∂θ



(1.28)

Or,  ∂u  ∂x ∂u ∂y

  − sinr θ ∂U cos θ ∂U ∂r ∂θ = cos θ ∂U sin θ ∂U ∂r + r ∂θ

(1.29)

Consequently, ∆u =



     ∂ ∂U ∂ ∂U sin θ ∂ sin θ ∂U cos θ ∂ cos θ ∂U cos θ + sin θ sin θ (1.30) − − + + ∂r r ∂θ ∂r r ∂θ ∂r r ∂θ ∂r r ∂θ   2 2 2 ∂U 1 ∂ U 1 ∂U 1 ∂ U 1 ∂ ∂ U r + 2 + + 2 = . (1.31) ∂r2 r ∂r r ∂θ2 r ∂r ∂r r ∂θ2

= =

1.3.3

cos θ

Implementation of P1 conforming triangular element on a uniform grid

(xj−1 , yk+1 )

(xj , yk+1 )

We begin with the following simplest example:

2 Kjk

3 Kjk (xj−1 , yk )

1 Kjk

(xj , yk ) (xj+1 , yk )

4 Kjk

=

f,

Ω = (0, 1)2 ,

(1.32a)

u =

0,

∂Ω.

(1.32b)

− ∆u

Let Th be the triangulation of Ω into the 2N × 2N uniform triangles with vertices Vjk = (xj , yk ) = h(j, k), j, k = 0, 1, · · · , N, with h = N1 .

6 Kjk

5 Kjk

(xj , yk−1 ) 1 2 6 In each triangle Kjk , Kjk , · · · , Kjk which has Vjk as a vertex, we have

 1 − h (x + y − (j + k + 1)h)     − h1 (y − (k + 1)h)    1 h (x − (j − 1)h) ϕjk (x, y) = 1  h (x + y − (j + k − 1)h)   1    h (y − (k − 1)h)  − h1 (x − (j + 1)h)

on on on on on on

1 Kjk , 2 Kjk , 3 Kjk , 4 Kjk , 5 Kjk , 6 Kjk .

PN −1 PN −1 Let uh (x, y) = j=1 k=1 αjk ϕjk (x, y) be an element of Vh . Find αjk for j = 1, · · · , N − 1 and k = 1, · · · , N − 1 such that N −1 X

αjk

j,k=1

Z

Ω

∇ϕjk · ∇ϕlm dxdy =

for all l = 1, · · · , N − 1, m = 1, · · · , N − 1. Let Z Alm,jk = ∇ϕjk · ∇ϕlm dxdy Ω

If we write α (and b) as the following order:

Z

f ϕlm dxdy Ω

and blm =

Z

f ϕlm dxdy. Ω

(α11 , α21 , · · · , xN −1,1 , α12 , α22 , · · · , αN −1,2 , · · · , α1,N −1 , α2,N −1 , · · · , αN −1,N −1 )t , c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

16

Finite element methods

Section 1.3

then the matrix A satisfying Aα = b will be of the form:  A11,11 A11,21 ··· A11,N −1,1 ···  A21,11 A · · · A ··· 21,21 21,N −1,1   . .. .. .. . .. ..  . . .   AN −1,1,11 A · · · A · · · N −1,1,21 N −1,1,N −1,1   .. .. .. .. ..  . . . . .   A1,N −1,11 A1,N −1,21 ··· A1,N −1,N −1,1 ···   .. .. .. . .. .  . . . . . AN −1,N −1,11 AN −1,N −1,21 · · · AN −1,N −1,N −1,1 · · ·

A11,1,N −1 A21,1,N −1 .. .

··· ··· .. .

A11,N −1,N −1 A21,N −1,N −1 .. .

AN −1,1,1,N −1 .. .

··· .. . ··· .. .

AN −1,1,N −1,N −1 .. .

···

AN −1,N −1,N −1,N −1

A1,N −1,1,N −1 .. . AN −1,N −1,1,N −1

A1,N −1,N −1,N −1 .. .



       .      

In the above matrix, the nonzero components may be found only if |j −l|+|k −m| ≤ 1 or |j −l| = |k −m| = 1. To compute these nonzero terms, notice that ∇ϕjk in each region:

∇ϕjk =

1. Then Ajk,jk =

P6

R

j j=1 Kjk

              

1 h (−1, 1) 1 h (0, −1) 1 h (1, 0) 1 h (1, −1) 1 h (0, 1) 1 h (−1, 0)

on on on on on on

1 Kjk , 2 Kjk , 3 Kjk , 4 Kjk , 5 Kjk , 6 Kjk .

|∇ϕjk |2 dxdy = [( h1 )2 + ( h1 )2 + 2( h1 )2 + ( h1 )2 + ( h1 )2 + 2( h1 )2 ] ×

h2 2

= 4.

3 2. For (l, m) = (j − 1, k) [or (j + 1, k)], we can compute Alm,jk by calculating integrals only on Kjk and 4 1 6 Kjk [or Kjk and Kjk ]. Therefore, Alm,jk = −( h1 )2 ·

h2 2

× 2 = −1.

3. Similarly, for (l, m) = (j, k − 1) [or (j, k + 1)], we can compute Alm,jk by calculating integrals only on 2 5 1 2 4 Kjk and Kjk [or Kjk and Kjk ]. Therefore, Alm,jk = −( h1 )2 · h2 × 2 = −1. 4. Finally, for (l, m) = (j − 1, k − 1) [or (j + 1, k + 1)], we can compute Alm,jk by calculating integral only 3 5 6 2 on Kjk and Kjk [or Kjk and Kjk ]. In these cases, the gradients of ϕjk and ϕlm are perpendicular, so Alm,jk = 0. Combining the above calculations, we  B −I · · ·  −I B . . . A=  . .. ..  .. . . O O ···

can see that  O  O  ..  where . B

and I is the (N − 1) × (N − 1) identity matrix. Similarly, blm can be calculated by Z

ϕlm f dxdy = Ω

XZ jk



−1 · · ·  .. −1 4 .  B= . . . .. ..  .. 0 0 ···

ϕlm (x, y)f (x, y) dxdy. = Kjk

4

6 Z X j=1

j Klm

0



 0  ..  , . 4

ϕlm (x, y)f (x, y) dxdy.

The jk-th row of Aα = b is given by 4αjk − (αj,k+1 + αj,k−1 + αj−1,k + αj+1,k ) = c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

Z

ϕjk f dxdy. Ω

17

Finite element methods

Section 1.4

For the constant force term, say, f (x, y) = 1 for all (x, y) ∈ Ω, one can calculate bjk = as the sum of volumes of six tetrahedrons: bjk = 6 ·

R

Ω

f (x, y)ϕjk (x, y) dxdy

h2 = h2 . 6

Notice that this is identical to the linear system where finite difference method applied, by dividing both sides by h2 .. In particular, if f = 0, then αjk =

1 (αj,k+1 + αj,k−1 + αj−1,k + αj+1,k ), 4

which is called the Discrete Mean Value Property. Generally, if u is harmonic (i.e. ∆u = 0), then u has the Mean Value Property: u(x0 ) =

1 2πr

Z

u(s)ds = S(x0 ;r)

1 πr2

Z

u(x)dx B(x0 ;r)

where S(x0 ; r) and B(x0 ; r) denote circle and disk of radius r centered at x0 .

1.4 1.4.1

Higher order finite elements in Rd Pm (conforming) simplicial element

Let K be a non-degenerate simplex with vertices V1 , · · · , Vd , Vd+1 with identification V0 = Vd+1 . Denote by Pd+1 J the set of indices α = (α1 , · · · , αd+1 ) in Zd+1 with j=1 αj = m Then consider the set of points + P=

X αj Vj . m

α∈J

• K = d − simplex • PK = Pm (K) with dim(PK ) = • ΣP K = {φ(V ) : V ∈ P}

1.4.2



d+m m



as the degrees of freedom

Qm (conforming) d-linear element

Similar as simplicial elements. Let K be the d–cube (−1, 1)d . Then consider the set of points P=



2 (α1 , · · · , αd ) − (1, · · · , 1) ∈ Rd | αj ∈ {0, 1, · · · , m} ∀j m



• K = d − cube = (−1, 1)d • PK = Qm (K) with dim(PK ) = (m + 1)d • ΣP K = {φ(V ) : V ∈ P}

as the degrees of freedom.

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

18

Finite element methods

1.4.3

Section 1.4

Unisolvency and optimality

It is necessary to show that the DOFs determine the function in PK . Also in order to construct optimal finite elements, one usually require that Pm (K) ⊂ PK . In this case the Bramble–Hilbert lemma shows the element is optimal. We state the following results without proof at the moment. Theorem 1.4. All the finite elements defined above have the “unisolvency” property. That is, if the DOFs of a given function f ∈ PK are zero, then the function f ⇐⇒ 0 in K. Also we have the continuity property as follows: Theorem 1.5. All the finite elements defined above have the “continuity” property. That is, if the DOFs restricted on a face f of K of a given function f ∈ PK are zero, then the function f ⇐⇒ 0 on the face f .

1.4.4

Triangulations and reference elements

Let (Th )0 0 such that for every h, hK ≤ κ for all K ∈ Th , ρK the family of triangulations is called “shape-regular”. If, in addition, there exists c > 0 such that hK ≥ ch

∀K ∈ Th

∀h,

the family of triangulations is called “quasi-uniform”. We now introduce the following reference elements: b be the reference simplicial element in Rd with coordinates Vbj , j = 0, 1, · ∗ d 1. d-simplicial element: Let K such that Vb0 = 0 and Vbj − Vb0 = b ej be the j th unit vector for j = 1, · · · , d. Then, for each triangle b → K such that FK (b or a simplex K ∈ Th , there is a surjective one–to–one affine map FK : K x) = b + bK . BK x In R2 and R2 the elements are called triangular and simplicial elements, respectively.

b be the d–cube (−1, 1)d . Then, for any d-rectangle K ∈ Th , there is a 2. d-rectangle-type element: Let K b → K such that FK (b surjective one–to–one d-linear map FK : K x) = BK (b x), where each component b → K is a d-linear map. of BK : K In R2 and R2 the elements are called quadrilateral and hexahedral elements, respectively. Notice that in two dimension the faces of a quadrilateral are straight lines segments, but in three dimension the faces of a hexahedron may not be flat.

b the reference finite element is denoted by the triple: On the reference element K b b ). b Pb b , Σ (K, K K

b b or ΣK , are evaluations of function φb ∈ Pb b Remark 1.2. The DOFs(degrees of freedom), denoted by Σ K K b or φ ∈ PK . Hence the set ΣKb or ΣK which denotes the DOFs is a subset of dual space of PbKb or PK , respectively. That is, b b ) = dim(Pb b ) or ΣK ⊂ (PK )′ with #(ΣK ) = dim(PK ) b b ⊂ (Pb b )′ with #(Σ Σ K K K K

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

19

Finite element methods

Section 1.4

If the DOFs consist of the evaluation of function at points or the integrations of function on a subset of the element only, the finite element is called “Lagrange–type finite element”. If they contain any evaluation of derivative of function or any integration of derivative, which is not included in the Lagrange–type DOFs, the finite element is called “Hermite–type finite element”. All the finite elements introduced above are of Lagrange–type. There are a number of Hermite–type finite elements. We will introduce an important Hermite–type finite element, so–called the 21 DOFs “Argyris element” and the 6 DOFs “Morley nonconforming element” Using the reference finite element, one can define the finite element (K, PK , ΣK ) on each element K ∈ Th as follows: b • K = FK (K) −1 • PK = {b p ◦ FK | pb ∈ PbKb }

b K ), • ΣK = {f ∈ (PK )′ | Σ(φ) = Σ(F

b ∈Σ bb Σ K

∀φ ∈ PK }

−1 Notice that in ΣK there are at least three types of DOFs to distinguish as follows: with φ = φb ◦ FK , which b is equivalent to φ = φ ◦ FK , b = φ(b b xj ) ⇐⇒ Σ(φ) = (φ ◦ FK )(b b φ) 1. Point–value type: Σ( xj ) ∂bx R R −1 b x) db b = φ(b b φ) dσ(x) σ (b x) ⇐⇒ Σ(φ) = e (φ ◦ FK )(x) ∂x 2. Integral type: Σ( e b Pd b ∂xℓ K) b = ∂Σ b φ) xj ) ⇐⇒ Σ(φ) = ℓ=1 ∂(φ◦F (b xj ) ∂b 3. Derivative type: Σ( ∂b xk (b ∂xℓ xk

b Indeed, they are Let {φbj | j = 1, · · · , J} be the standard basis functions on the reference element K. constructed as follows: denote b j ∈ (Pb b )′ , b b = {Σ Σ K K

b j , find φbj ∈ Pb b such that Associated to each Σ K b k (φbj ) = δjk , Σ

j = 1, . . . , dim(PbKb )}.

δjk being the Kronecker delta function.

(1.33)

Hence, the standard finite element basis functions are given by φbj , j = 1, · · · , dim(PbKb ).

(1.34)

The basis functions on each K ∈ Th is defined by the pull–back

−1 {φK,j (x) = (φbj ◦ FK )(x) | j = 1, · · · , J}.

Then the global finite element space Vh can be defined by Vh = {vh ∈ C 0 (Ω) | vh |K =

J X

αK,j φK,j

j=1

∀K ∈ Th }.

(1.35)

Notice that in order to meet the criterion that vh ∈ C 0 (Ω), one needs to have αK,j = αK ′ ,j ′

∀j ′ such that xj (∈ K) = xj ′ (∈ K ′ ).

Definition 1.1. Consider a finite element (K, PK , ΣK ) and the global finite element space Vh . Let ℓ be the maximum order of derivatives in the definition of the DOFs ΣK . For each K ∈ Th , define the local PK –interpolation operator IK : C ℓ (K) → PK such that Σℓ (IK v − v) = 0

∀Σℓ ∈ ΣK

∀v ∈ C ℓ (K).

(1.36)

Then the global Vh –interpolation operator Ih : Ω → Vh is defined piecewisely such that Ih v|K (x) = IK v c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

∀K ∈ Th . 20

Finite element methods

Section 1.4

Ih v and IK v are called the global Vh – and local PK –interpolants of v, respectively. Let S ⊂ Rd be a bounded set with diam (S) = h. Bramble and Hilbert use the following normalized norms on W m,p (S) and W m,∞ (S), which are equivalent to the usual Sobolev norms: |||u|||m,p,S

=

|||u|||m,∞,S

=

m X

h

k=1 m X k

k=0

kp−d

|u|pk,p,S

! p1

h |u|k,∞,S .

(1.37) (1.38)

Theorem 1.6 (Bramble–Hilbert I.). Let Q = W k,p (S)/Pk−1 be a quotient space. Then the semi–norm h−d/p |u|k,p,S is a norm on Q equivalent to the quotient norm k[u]kQ = inf h−d/p kvkk,p,S , v∈[u]

[u] denotes the equivalent class.

Furthermore, there exists C > 0 independent of h and u such that d

d

hk− p kukk,p,S ≤ k[u]kQ ≤ Chk− p kukk,p,S

∀u ∈ W k,p (S).

Theorem 1.7 (Bramble–Hilbert II.). Let S be as in the above theorem. Suppose that L : W k,p (S) → R be a bounded linear functional such that 1. |L(u)| ≤ C|||u|||k,p,S

∀u ∈ W k,p (S) for C > 0 independent of h and u,

2. Pk−1 (S) ⊂ Ker(L). Then there exists C > 0 such that |Lv| ≤ C1 hk−d/p |u|k,p,S

∀u ∈ W k,p (S) for C1 > 0 independent of h and u.

Theorem 1.8 (Bramble–Hilbert III.). Let S be as in the above theorem. Suppose that L : W k,p (S) → R be a bounded linear functional such that 1. |L(u)| ≤ C|||u|||j,∞,S

∀u ∈ C j (S) for C > 0 independent of h and u,

2. Pk−1 (S) ⊂ Ker(L). Then there exists C1 > 0 such that, for all |Lv| ≤ C1 hk−d/p |u|k,p,S

1 p

−

k−j d

< 0, such that W k,p (S) ֒→ C j (S) and

∀u ∈ W k,p (S) for C1 > 0 independent of h and u.

Theorem 1.9 (Bramble–Hilbert IV.). Let S be as in the above theorem. Suppose u ∈ C s , s = [s]+α ∈ (0, 1], where [s] denotes the Gauss integer which is the largest integer ≤ s. Suppose that L : C 0 (S) → R be a bounded linear functional such that 1. |L(u)| ≤ C|u|0,∞,S

∀u ∈ C 0 (S) for C > 0 independent of h and u,

2. Pk−1 (S) ⊂ Ker(L). Then there exists C1 > 0 such that, for all |Lv| ≤ C1 hs sup

x,y∈S

X ∂ α u(x) − ∂ α u(y)| |x − y|α

|α|=[s]

1 p

−

k−j d

< 0, such that W k,p (S) ֒→ C j (S) and

0≤s 0 independent of h and u. 21

Finite element methods

Section 1.5

Applications of Bramble–Hilbert Theorems are as follows: In Theorem 1.7, suppose K ∈ Th be of size diam (K) = h. Then, choose k = 2 and PK = P1 (K). Consider the difference beween the following standard local P1 –interpolant and a function : #(vertices of K)

Lu(x) = IK u(x) − u(x) =

X

u(xj )φK,j (x) − u(x)

j=1

Then, verify that |Lu(x)| = |IK u(x) − u(x)| ≤ C|||u|||2,p,K .

Thus from Theorem 1.7 it follows that

d

|Lu(x)| ≤ C1 h2− p |u|2,p,K . Hence, kIK u − uk0,K =

Z

p

K

|Lu(x)| dx

 p1

d

1

≤ C1 h2− p |u|2,p,K |K| p = C1 h2 |u|2,p,K .

∂u In the meanwhile, for the energy–norm estimate, set wj = ∂x choose k = 1 and PK = P0 (K). Consider j the difference between the P0 –interpolation (i.e., the mean–average operator on K) and a function: Z 1 Lwj (x) = wj (x) dx − wj (x). |K| K

Then, verify that Thus from Theorem 1.7 it follows that

|Lwj (x)| ≤ C|||wj |||1,p,K . d

|Lwj (x)| ≤ C1 h1− p |wj |1,p,K . Hence, |IK u − u|1,p,K

  p1 d Z X 1 d = |Lwj (x)|p dx ≤ C1 h1− p |u|2,p,K |K| p = C1 h1 |u|2,p,K . j=1

K

Definition 1.2. For m ≥ 1, the broken norms (or mesh-dependent norms) are defined as follows: ! 21 X 2 ∀v ∈ Vh . kvkm,T kvkm,h :=

(1.39)

T ∈Th

Summarizing the above, we have Theorem 1.10. Let m = 2. Suppose that the finite element spaces PK contain P1 , as all the finite elements introduced above. Then, kIh u − uk0,h kIh u − uk1,h

≤ Ch2 |u|2,Ω , ≤ Ch|u|2,Ω ,

for all u ∈ H 2 (Ω). Invoking the Ce´ a lemma, we have the following error estimates between the finite element solution and the exact equation: kuh − uk1,h ≤ C kIh u − uk1,h

≤

Ch|u|2,Ω .

Theorem 1.11. Let m ≥ 2, {Th }h be a family of shape-regular triangulations of Ω, and Vh is a finite b ⊂ Pb b . Then, there exists a constant C > 0 such that element space such that Pm−1 (K) K ku − Ih ukk,h ≤ Chm−k |u|m,Ω

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

∀u ∈ H m (Ω) for 0 ≤ k ≤ min{1, m}.

22

Finite element methods

1.5

Assembly of matrix Akj =

P

Section 1.5 K∈Th

ah (φj , φk )

For this topics see Long Chen’s implementation of iFEM. Long Chen’s lecture notes “http://www.math.uci.edu/ chenlong/226/Ch3FEMCode.pdf”

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

23

Finite element methods

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

Section 1.5

24

Chapter 2

NC (Nonconforming) finite elements All the finite element spaces Vh introduced earlier are subspaces of H 1 (Ω) ∩ C 0 (Ω). This means that any vh ∈ Vh are continuous across the interfaces Γℓm := ∂Kℓ ∩ ∂Km for all elements Kℓ and Km in Th . Denote by [[·]]Γℓm the jump across the interface Γℓm defined by [[u]]Γℓm = γΓℓm (u|Kℓ ) − γΓℓm (u|Km ),

(2.1)

where u|Kℓ denotes the restriction of u ∈ L2 (Ω) to Kℓ , and γΓℓm uKℓ the trace of uKℓ onto Γℓm . Recall that Kℓ ’s are open sets. Conforming finite elements require that [[u]]Γℓm (x) = 0

∀x ∈ Γℓm .

(2.2)

Breaking this continuity restriction would substantially widen the realm of finite element spaces. However, there should be certain rules to relate functions defined on neighboring elements Kℓ and Km . These are essentially the degrees of freedom. Let us classify DOFs (the degrees of freedom) in two categories: 1. Interface DOFs: the DOFs that transfer the information on the element to neighboring element 2. Interior DOFs: the DOFs that determine information only on the element Observe that conforming finite elements have interior DOFs such that (2.2) holds. Remark 2.1. DOFs are also classified into two kinds as follows: 1. Lagrange type of finite elements: DOFs consist of point values of a function (usual Pm−1 or Qm−1 finite elements) 2. Hermitian type of finite elements: DOFs consist of point values of the derivatives of the function upto certain order (Morley element, for instance) Nonconforming finite elements may have two types of interface DOFs 1. Gauss points DOFs: ΣP K = {φ(gj ) : gj are Gauss points on the faces of K} R 2. Orthogonality DOFs: ΣP K = { f φq dσ : q ∈ Pm−1 (f ), f is a face of K}

which relax (2.2) as follows

1. [[u]]Γℓm (gj ) = 0 ∀Gauss points on the faces Γℓm R 2. Γℓm [[u]]Γℓm q dσ = 0 ∀q ∈ Pm−1 (f ), f is a face of K}

respectively.

25

Finite element methods

2.1

Section 2.1

P1 triangular NC element (Crouzeix-Raviart, 1973) and the rotated Q1 -rectangular NC elements

1. P1 triangular NC element (Crouzeix-Raviart, 1973) • K = triangle

• PK = Span{1, x, y}

• DOF: {φ(mj ), mj midpoints of edges, j = 1, 2, 3}

K=T

V2 e3 m1

1 0 0 1

1 m3 0 0 1

T

e1

PK = Span{1, x, y}

1 0 0 1

V3

ΣK = {φ(mj ), j = 1, 2, 3}

e2 V1

m2

Exercise 2.1. (Due to July 10) Let K be the triangle with vertices (0, 0), (1, 0), (0, 1), (0, 0). Find three basis functions φj (x, y), j = 1, · · · , 3, for PK such that φj (mk ) = δjk , where δjk is the Kronecker delta function. 2. The rotated Q1 -rectangular NC elements (Han (1984), Rannacher–Turek (1992), Chen (1992), Douglas– Santos–Sheen–Ye (1999) e2

V2

V1 m2

m4 e4

ΣM K ΣIK (a) H. Han (1984):

K=R

m1

e3

V3

e1

R

m3

PK = Span{1, x, y, θℓ(x) − θℓ(y)} V4

= {φ(mj ), j = 1, · · · , 4} R = { ej φ dσ, j = 1, · · · , 4}

• K = R = [−1, 1]2 • PK = Span{1, x, y, x2 − 5/3x4 , y 2 − 5/3y 4 } R • DOF: ΣM K = {φ(mj ), mj midpoints of edges, j = 1, 2, 3, 4} plus K φ

(b) Rannacher-Turek (1992, rotated Q1 element, also Z. Chen): • K = [−1, 1]2 • PK = Span{1, x, y, x2 − y 2 } c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

26

Finite element methods

Section 2.1

• DOF1 ΣP K = {φ(m R j ), mj midpoints of edges, j = 1, 2, 3, 4}; DOF2: ΣIK = { ej φdσ, ej four edges, j = 1, 2, 3, 4};

(c) DSSY element(Douglas-Santos-Sheen-Ye, 1999): Cai-Douglas-Ye, 1999 stable Stokes element: e2

V2

V1 m2

K=R

e1

R

m3

m1

e3 m4 V3

e4

ΣM K ΣIK

PK = Span{1, x, y, θℓ(x) − θℓ(y)} V4

= {φ(mj ), j = 1, · · · , 4} R = { ej φ dσ, j = 1, · · · , 4}

• K = [−1, 1]2 • PK = Span{1, x, y, θℓ (x) − θℓ (y)}, ℓ = 0, 1, 2 R • DOF1 = DOF2; |e1j | ej φ dσ = φ(mj )

 2  t , t2 − 5/3t4 , θℓ (t) =  2 t − 25/6t4 + 7/2t6 ,

ℓ = 0; ℓ = 1; ℓ = 2.

Exercise 2.2. (Due to July 10) Let K = (−1, 1)2 . Find four basis functions φj (x, y), j = 1, · · · , 4, for PK (DSSY element) such that φj (mk ) = δjk , where δjk is the Kronecker delta function. (d) For truly quadrilaterals, (Cai-Douglas-Santos-S.-Ye, CALCOLO, 2000): • K = [−1, 1]2 • PK = Span{1, x, y, θℓ (x) − θℓ (y), xy}, ℓ = 0, 1, 2 R • DOF: {φ(mj ), mj midpoints of edges, j = 1, 2, 3, 4} and Rb φ(x, y)xy dxdy;

3. P2 triangular NC element (Fortin-Soulie, 1983; Lee-Sheen 2005) gj ’s denotes Gauss points on the face

K=T

V2 1 0 0 1 0 1

1 0 0 1 g1

e1 g2 00 11 00 11 g3 0 1 0 0 1 V3 1

PK = P2

g6

e3 T

ΣK 1 0 0 1 0 1

g5 g4 e2

1 0 0 1 0 1

= {φ(gj ), j = 1, · · · , 6;

R

K φ dx}

V1

4. P2 rectangular NC element (or “rotated Q2 -rectangular NC element) (Lee-Sheen 2005) gj ’s denotes Gauss points on the face c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

27

Finite element methods

Section 2.2

e2

V2 g4

V1

g3

g5

R

g2 e1

O

e3

g1

g6 g8

g7 e4

V3

K=R

PK = P2 ⊕ Span{x2y, xy2} R ΣK = {φ(gj ), j = 1, · · · , 8; K φ dx}

V4

The P1 -NC quadrilateral and hexahedral elements

2.2

(C. Park, Thesis 2002, Park-Sheen SIAM J. Numer. Anal. 2003) • K = Q: a quadrilateral V2

• PK = S (Q) = Span{1, x, y}

e2 m2

• ΣK = {φMj , j = 1, · · · , 4}

V1

• Basis function

m3 m1

φj (m) =

e3 e1

V3

e4

V4

m4



1 2, 0,

if m = Mj , Mj+1 , if m = Mj+2 , Mj+3

• S (Q) = Span{φ1 , · · · , φ4 } and dim(S (Q)) = 3

Lemma 2.1. If u ∈ S(Q), then u(M1 ) + u(M3 ) = u(M2 ) + u(M4 ). Conversely, if uj is a given value at Mj , for j = 1, · · · , 4, satisfy u1 + u3 = u2 + u4 , then there exists a unique function u ∈ S(Q) such that u(Mj ) = uj , j = 1, · · · , 4. The P1 -NC hexahedral element v6

v7

• K = hexahedron

1 c3 0 1 0

v1

v2

• PK = S (Q) = Span{1, x, y, z} 1 0 1 0

c 1 0 1 5 0

v7 1 0 0 1

c1

v3

1 0 0 c6 1 0 1

c2 v8

c41 0

• ΣK = {φ(Mj ), j = 1, · · · , 6} • S (Q) = Span{φ1 , · · · , φ8 } and dim(S (Q)) = 4

0 1

v4

Lemma 2.2. If u ∈ S(Q), then u(M1 ) + u(M6 ) = u(M2 ) + u(M5 ) = u(M3 ) + u(M4 ). Conversely, if uj is a given value at Mj , for 1 ≤ j ≤ 6, satisfying u1 + u6 = u2 + u5 = u3 + u4 , then there exists a unique function u ∈ S(Q) such that u(Mj ) = uj , 1 ≤ j ≤ 6. With these elements we can prove optimal convergence. Theorem 2.1. The above element is unisolvent. Optimal error estimates holds for elliptic problems: c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

28

Finite element methods

2.2.1

Section 2.2

Error analysis for linear nonconforming Galerkin method

Define the global nonconforming finite element space N C h0 on Th by N C h0

=

{v ∈ L2 (Ω) | v|K ∈ PK

h[[v]]e , qie = 0,

∀K ∈ Th ,

∀q ∈ Pm−1 (e), hve , qie = 0, ∀e ∈ ∂} R where [[f ]]e means the jump of f across the interface e and h[[v]]e , qie = e [[v]]e q dσ. Define the bilinear form ah : [N C h0 + H 1 (Ω)]2 → R by X aK (u, v), where aK (u, v) denotes the integrations restricted to the domain K ah (u, v) = K∈Th

Define the broken energy norm: for w ∈ H 1 (Ω) + N C h0 p |w|1,h = ah (w, w)

Recall that the weak solution u ∈ H01 (Ω) satisfies a(u, v) = hf, vi

∀v ∈ H01 (O).

(2.3)

The nonconforming Galerkin method is to find uh ∈ N C h0 such that ah (uh , vh ) = hf, vh i

∀vh ∈ N C h0 .

(2.4)

Lemma 2.3 (The second Strang Lemma). If u ∈ H01 (Ω) and uh ∈ N C0h are the solutions of (2.3) and (2.4) respectively, then ! |ah (u, w) − (f, w)| (2.5) ah (u − uh , u − uh ) ≤ C inf |u − v|1,h + sup |w|1,h v∈N C0h w∈N C0h Remark 2.2. In the right-hand side of (2.5), the second term is called a consistency error term. Notice that in the consistency error term, ah (u, w) − (f, w) = ah (u, w) − ah (uh , w).

(2.6)

If w ∈ C 0 (Ω), then the consistency error term become zero. Proof. First, by the triangle inequality |u − uh |1,h ≤ |u − zh |1,h + |zh − uh |1,h

∀zh ∈ N C0h .

Next, since ah (uh , uh − zh ) = (f, uh − zh ), one has |uh − zh |21,h

= = = ≤

ah (uh − zh , uh − zh )

ah (u − zh , uh − zh ) + ah (uh − u, uh − zh ) ah (u − zh , uh − zh ) + (f, uh − zh ) − ah (u, uh − zh )

|u − zh |1,h |uh − zh |1,h + (f, uh − zh ) − ah (u, uh − zh ).

Then |uh − zh |1,h ≤ |u − zh |1,h +

(f, uh − zh ) − ah (u, uh − zh ) , |uh − zh |1,h

which combined with the triangle inequality gives |u − uh |1,h ≤ 2|u − zh |1,h +

(f, uh − zh ) − ah (u, uh − zh ) . |uh − zh |1,h

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

29

Finite element methods

Section 2.2

If we choose zh ∈ N C h such that |u − zh |1,h =

inf

v∈N C0h

|u − v|1,h ,

then |u − uh |1,h ≤ 2 inf

v∈N C0h

|u − v|1,h + sup

w∈N C0h

|(f, w) − ah (u, w)| . |w|1,h



For the nonconforming finite elements introduced above PK including Pm−1 (K), m = 2, 3, by the Bramble-Hilbert Theorem 1.7, kIh u − ukm,Kb ≤ C|u|m,Kb

b ∀u ∈ H m (K),

from which the first part in the second Strang lemma is estimated as follows: kIh u − ukk,h ≤ Chm−k |u|m,Ω ,

2.2.2

0 ≤ k ≤ m.

(2.7)

Implementation of P1 –NC quadrilateral element

Consider the special case: Ω = (0, 1) is decomposed into the uniform N × N uniform rectangles with coordinates (xj , yk ) = h(j, k), j, k = 0, 1, · · · , N, where h = N1 . The basis functions φjk (x, y) contain four rectangles whose vertices contain (xj , yk ). The basis function values are shown in the following, but in a shifted region (−h, h)2 . To find components of the mass matrix, using symmetry and translation invariance, one has the following: • (l, m) = (j, k) 4 (φjk , φjk ) = 2 4h

Z

h 0

Z

h 0



3 −x − y + h 2

2

dx dy =

20 2 h , 48

• (l, m) = (j + 1, k) 2 (φjk , φj+1,k ) = 2 4h

Z

h 0

Z

h 0



−x − y + 32 h 2h



(x − y +

h 18 2 ) dx dy = h , 2 48



−x − y + 23 h 2h



(x + y −

h −11 2 ) dx dy = h . 2 48

and • (l, m) = (j + 1, k + 1) (φjk , φj+1,k+1 ) =

1 4h2

Z

h 0

Z

h 0

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

30

Finite element methods

Section 2.2

( h2 , h)

(0, h)

(−h, h)

(h, h)

−11

1 Kjk

2 Kjk

−x−y+ 23 h 2h

x−y+ 23 h 2h

(−h, 0)

x+y+ 23 h 2h

−11 (h, h2 )

18

h2 48

(h, 0)

( h2 , 0)

( h2 , h) 1 Kjk

2 Kjk

(h, h2 )

(0, 0)

18

20 ( h2 , 0)

18

−x+y+ 32 h 2h

3 Kjk

3 Kjk

4 Kjk

(0, −h)

(−h, −h)

(h, −h)

−11

φ(x, y)

4 Kjk

18 (φjk , φℓm)

−11

For the components of mass matrix, using symmetry and translation invariance,

• (l, m) = (j, k) 4 (∇φjk , ∇φjk ) = 2 4h

Z

Z

h 0



h 0

   −1 −1 dx dy = 2, · −1 −1

• (l, m) = (j + 1, k)

(∇φjk , ∇φj+1,k ) =

2 4h2

Z

h 0

Z

h 0



   −1 1 · dx dy = 0, −1 −1

and

• (l, m) = (j + 1, k + 1)

(∇φjk , ∇φj+1,k+1 ) =

1 4h2

Z

h 0

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

Z

h 0



   1 −1 1 · dx dy = − . −1 1 2 31

Finite element methods

Section 2.3

(h, h)

( h2 , h)

(0, h)

(−h, h)

1 Kjk   −1 1 2h −1

2 jk  K 1 1 2h −1

(−h, 0)

(0, 0)

  1 1

1 2h

3 Kjk

4 ( h2 , 0)

3 Kjk

(h, −h)

−1

∇ φ(x, y)

2.3

(h, h2 )

  −1 1

(0, −h)

(−h, −h)

(0, h2 )

2

4 Kjk

−1

1 Kjk

0

(h, 0) 1

( h2 , h)

0 2 Kjk

(h, h2 )

( h2 , 0) 1 2h

−1

0

4 Kjk

0 (∇ φjk , ∇ φℓm)

−1

Hermite–type finite elements

In this section we introduce Hermite–type finite elements for solving the fourth–order partial differential equations (the biharmonic equation) ∆2 u(x) u(x) ∂u ∂ν

= =

f (x), 0,

x ∈ Ω, x ∈ ∂Ω,

=

0,

x ∈ ∂Ω.

(2.8)

∂u ∂v v − u dσ ∂ν ∂ν

(2.9)

Green’s second identity: Z

Ω

∆u v dx −

Z

u ∆v dx = Ω

Z

∂Ω

By Green’s second identity, we can have Z Z Z (∆2 u, v) = (f, v) −→ ((∆2 u)v dx − ∆u ∆v dx =

∂v ∂(∆u) v − ∆u dσ ∂ν ∂ν Ω Ω ∂Ω H 2 (Ω) ∂v 2 ∞ ∀v ∈ H0 (Ω) = C0 (Ω) (i.e. v|∂Ω = 0, = 0) ∂ν ∂Ω

The weak problem is given as finding u ∈ H02 (Ω) such that (∆u, ∆v) = (f, v)

∀v ∈ H02 (Ω)

(2.10)

Finite element method is to find Vh ⊂ H02 (Ω) and then to find uh ∈ Vh such that (∆uh , ∆v) = (f, v)

∀v ∈ Vh

(2.11)

In particular, we want to find a finite element space Vh ⊂ C 1 (Ω) ∩ H02 (Ω).

Argyris element (C 1 -quintic triangular element or 21-DOFs triangular element)

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

32

Finite element methods

Section 2.3

• T = triangle • PT = P5 (T )

 ∂p ∂2p ∂2p ∂2 ∂p (ai ), (ai ), (a ), (a ), • ΣT = {∂ p(ai ), |α| ≤ 2, i = 1, 2, 3 : ∂ν p(aij . 1 ≤ i ≤ j ≤ 3} = p(ai ), i i ∂x1 ∂x2 ∂x21 ∂x1 ∂x2 ∂x  ∂p (aij ), 1 ≤ i ≤ j ≤ 3 i = 1, 2, 3 ; ∂ν α

• dim(P5 (T )) = 21.

Argyris triangle

Figure 2.3.1: Argyris element Let vh ∈ Vh , where Vh is constructed using Argyris element. We want to show that v h |T 1 − v h | T 2 = 0

along T¯1 ∩ T¯2 .

Let t denotes the variable along T¯1 ∩ T¯2 so that v(t) is polynomial of degree ≤ 5 fulfilling v(t1 ) = v(t2 ) = v ′ (t1 ) = v ′ (t2 ) = v ′′ (t1 ) = v ′′ (t2 ) = 0 =⇒ v(t) has factors (t − t1 )3 and (t − t2 )3 =⇒ v = 0. Of course, we have v ′ (t) = 0. Next, to show that the the normal directional derivative of v equals zero, let ! ∂vh ∂vh w(t) = − (t), ∂ν T1 ∂ν T2

where ν denotes the normal vector from T1 to T2 . Note that w(t) is a polynomial of degree ≤ 4 and w(t1 ) = w(t2 ) = w(t12 ) = 0, where t12 represents the midpoint between a1 and a2 . ! ∂v ∂v ∂ h − (tj ) w′ (tj ) = ∂t ∂ν ∂ν T1

=

T2

∇(( ∇vh |T1 · ν) − ( ∇vh |T2 · ν)) · τ = 0,

j = 1, 2.

Thus w(t) = 0 along T¯1 ∩ T¯2 . And hence vh ∈ C 1 (T 1 ∪ T 2 ). This proves that Vh ⊂ C 1 (Ω) ∩ H02 (Ω). Quadratic nonconforming elements on rectangles with Heejeong Lee 1. Triangular case c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

33

Finite element methods

Section 2.3

K=T

V2 m1

PK = P2 m3

T

V3

∂φ ΣK = {φ(Vj), ∂ν (mj),j = 1,2,3}

m2

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

V1

34

Finite element methods

Section 2.3

V2

K=T

1 0 0 1 1 0 0 0 g1 1 g1 6

PK = P2

e1 g2 00 11 00 11 g3 0 1 0 1 V3 10

T

e3 R

g5 1010ΣK = {φ(gj ),j = 1,··· ,6; K φdx} 0 1

g4 e2

1 0 0 1 0 1

V1

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

35

Finite element methods

Section 2.3

Incomplete biquadratic element

a2

a1

R a3

a4

• Morley element

• Fortin-Soulie element 2. Rectangular case • Incomplete biquadratic element

• Reduced incomplete biquadratic element (with Heejeong Lee) Exercise 2.3 (Due July 24 (everything before the coding), and July 31 (coding)). Choose any nonconforming finite element as you like. Then consider the same PDE as in the previous programming exercise. 1. Use N × N rectangular meshes with N = 50, 100, 200. 2. Find an explicit formula for the standard local basis functions on a reference domain. Using the reference basis functions, find an explicit formula for the global basis functions φ′jk , where φjk is associated with the point (xj , yk ) = h(j, k), h = N1 . 3. Give an explicit description how you would compute ah (φjk , φlm ) and (f, φjk ). ′ 4. Formulate a linear system Ax = b, where x denotes the unknown coefficients αjk s. (Submit up to this by July 24.)

5. Use any numerical linear algebra to solve Ax = b. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

36

Finite element methods

Section 2.3

Quadratic nonconforming rectangle

g4

V2

g3

V1 g2

g5

R g6 V3

g1 g7

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

g8

V4

37

Finite element methods

Section 2.3

6. Compare the numerical solutions with the numerical solutions obtained by conforming finite elements. 7. Submit the codes and analysis and comments by July 31.

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

38

Chapter 3

Preliminaries from functional analysis Let Ω be a nonempty open set in Rd , and set • for 1 ≤ p < ∞, Lp (Ω) = {f : Ω → C : kf k0,p,Ω < ∞}, where kf k0,p,Ω :=

R

Ω

|f |p dx

• C k (Ω) = {f : Ω → C : all the derivatives of f up to order k are continuous},

(1/p)

,

• C0k (Ω) = {f ∈ C κ (Ω) : supp(f ) is compact in Ω}, • C ∞ (Ω) = {f : Ω → C : f is infinitely differentiable}, • D(Ω) = C0∞ (Ω) = {f ∈ C ∞ (Ω) : supp(f ) is compact in Ω}, where • supp(f ) = {x ∈ Ω : f (x) 6= 0}. • Note that ∂Ω = Ω \ Ω • A subset K of Ω is compact iff every open cover {Oα }α∈A of K containsSa finite subcover, that is, there m exists a finite collection {Oα1 , · · · , Oαm }, α1 , · · · , αm ∈ A such that i=1 Oαi ⊃ K. In case U ⊂ Ω and U is compact, we denote U ⋐ Ω. Notations C m (Ω) D(Ω)

= =

Rd+

=

Rd−

=

{f : Ω → C : there exist an open set O containing Ω and f˜ ∈ C m (O) such that f = f˜|Ω }, {f : Ω → C : there exist an open set O containing and f˜ ∈ D(O) such that f = f˜|Ω },

{x = (x′ , xd ) ∈ Rd : xd > 0, x′ ∈ Rd−1 },

{x = (x′ , xd ) ∈ Rd : xd < 0, x′ ∈ Rd−1 }.

Example 3.1. Notice that the spaces of polynomials, trigonometric functions, and exponential functions are dense in C ∞ (Ω) under suitable norms. The following function    1 exp , if |x − a| < r, |x − a|2 − r2 ϕ(x) =  0, otherwise. belongs to C0∞ (Rd ).

39

Finite element methods

Section 3.1

Definition 3.1. A function k·k : X 7→ R+ is called a norm if kx + yk ≤ kxk + kyk

(i)

∀x, y ∈ X

kαxk = |α| kxk kxk > 0

(ii) (iii)

∀x, ∈ X ∀x, x 6= 0.

∀α ∈ R

Definition 3.2. A function | · | : X 7→ R+ is called a semi-norm if (i) and (ii) hold in Definition 3.1. Theorem 3.1. [Rud91a, p. 30] A topological vector space X is normable if and only if the origin has a convex bounded neighborhood. For an open subset Ω of Rd , p

L (Ω)

:=

{f : Ω 7→ R | kf k0,p,Ω :=

L∞ (Ω)

:=

Z

Z

Ω

|f (x)|p dx < ∞} p

Ω

|f (x)| dx

 p1

{f : Ω 7→ R | ess sup |f | < ∞} Ω

kf k0,∞,Ω := ess sup |f | Ω

p

Then L (Ω) is Banach space for p > 0 for the norm k·k0,p,Ω . (i.e., a vector space eqipped with a norm with which the space is complete.) If p = 2, then L2 (Ω) is Hilbert space with the scalar product (·, ·)0,2,Ω . (i.e., a vector space with an inner product (scalar product) which generates a norm with which the space is complete.) Z q (f, g)0,2,Ω := f (x) g(x) dx, kf k0,2,Ω = (f, f )0,2,Ω Ω

Definition 3.3. For Banach spaces V, W , denote by L (V, W ) the set of all bounded linear mappings from V into W . If W = V , we use the abbreviation L (V ) = L (V, V ).

3.1

Duality

Definition 3.4. Let V be a vector space with inner product (·, ·) and let A be a bounded linear operator on V . A∗ is called an adjoint of A, if for all u, v ∈ V , (Au, v) = (u, A∗ v)

or

(u, Av) = (A∗ u, v).

If A = A∗ , then A is called self-adjoint. Definition 3.5. Let V be a normed linear space. The dual space V ′ is the set of all bounded linear functionals v ′ : V → C, which forms a complete normed linear space, (a Banach space) with the norm given by kv ′ kV ′

:= =

|v ′ (v)| v∈V,v6=0 kvkV sup

sup

v∈V,v6=0,kvkV ≤1

=

sup v∈V,kvkV =1

(3.1a) |v ′ (v)|

|v ′ (v)|

(3.1b) (3.1c)

The duality paring h·, ·iV ′ ,V between V ′ and V will be denoted by hv ′ , viV ′ ,V := v ′ (v) c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

∀v ∈ V ∀v ′ ∈ V ′ . 40

Finite element methods

Section 3.2

Definition 3.6. [Yos95, p.120] A sequence (vn ) in a normed linear space V is said to converge weakly to v∞ ∈ V if limn→∞ hv ′ , vn iV ′ ,V exists and it is equal to hv ′ , v∞ iV ′ ,V for every v ′ ∈ V ′ . Theorem 3.2. [Yos95, p.121] A sequence (vn ) in a normed linear space V is said to converge weakly to v∞ ∈ V iff 1. supn≥1 kvn kV < ∞ and 2. limn→∞ hv ′ , vn iV ′ ,V = hv ′ , v∞ iV ′ ,V for all v ′ ∈ D′ for any strongly dense subset D′ of V ′ . Definition 3.7. [Yos95, p.125] A sequence (vn′ ) in the dual space of V ′ of a normed linear space V is said ′ ′ , viV ′ ,V for every v ∈ V. ∈ V ′ if limn→∞ hvn′ , viV ′ ,V exists and it is equal to hv∞ to converge weakly∗ to v∞ Theorem 3.3. [Yos95, p.121] A sequence (vn′ ) in the dual space of V ′ of a normed linear space V converge ′ weakly∗ to v∞ ∈ V ′ if iff 1. supn≥1 kvn′ kV ′ < ∞ and ′ 2. limn→∞ hvn′ , viV ′ ,V = hv∞ , viV ′ ,V for all v on D for any strongly dense subset D of V.

Theorem 3.4. If V is a normed linear space, the closed unit ball B ′ of V ′ is weak∗ -compact. The weak∗ topology is the weakest topology that makes all linear functional v ∗ → hv ∗ , viV ′ ,V continuous. Theorem 3.5. If V is a normed linear space, V ′ is a complete normed linear space, (a Banach space) with the dual norm k · kV ′ Theorem 3.6. Let V be a normed linear space. The map V ′ → C : v ′ 7→ hv ′ , viV ′ ,V is a bounded linear functional on V ′ with the norm kvkV . Let V be a Banach space, and V ′′ the dual space of V ′ , which consists of the set of bounded linear functional on V ′ . Thus, every v ∈ V can be regarded as an element in v ′′ ∈ V ′′ . Morever, due to the above theorem, the embedding I : V → V ′′ defined by I (v) = v ′′ is an isometric isomorphism and I (V ) is a closed subspace of V ′′ . In case I (V ) = V ′′ , V is called a reflexive Banach space. Definition 3.8. Let V, W be Banach spaces. Let T : D(T ) ⊂ V → W be an unbounded operator with the domain D(T ) dense in V . An unbounded adjoint operator T ∗ : D(T ∗ ) ⊂ W ′ → V ′ is defined as follows: First, the domain D(T ∗ ) is defined by follows: D(T ∗ )

=

{w′ ∈ W ′ : there exists a c such that | hw′ , T viW ′ ,W | ≤ ckvkV

∀v ∈ D(T )}.

(3.2)

It is clear that D(T ∗ ) is a subspace of W ′ . Since the map v → hw′ , T viW ′ ,W : D(T ) → C is bounded and D(T ) is dense in V , by Hahn-Banach theorem there exists a unique linear map f : V → C such that |f (v)| ≤ ckvkV ∀v ∈ V. Denote such a map (depending on w′ ) by T ∗ w′ = f, Then the following is fulfilled: hw′ , T viW ′ ,W = hT ∗ w′ , viV ′ ,V

∀v ∈ D(T ), ∀w′ ∈ D(T ∗ ).

In particular if V and W are identical Hilbert space and T is a bounded linear map in V , then (T ∗ x, y) = (x, T y)

for all x, y ∈ V.

Theorem 3.7. If V and W are normed linear spaces. Then for each T ∈ L (V, W ), kT ∗ kL (W ′ ,V ′ ) = kT kL (V,W ) . c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

41

Finite element methods

3.2

Section 3.4

Annihilators

Definition 3.9. Suppose M ⊂ V is a subspace of a Banach space V , and N ⊂ V ′ is a subspace of V ′ . The annihilatos M ⊥ and N ⊥ are defined as follows: M⊥ ⊥

N

= {v ′ ∈ V ′ : hv ′ , viV ′ ,V = 0 = {v ∈ V : hv ′ , viV ′ ,V = 0

∀v ∈ M },

∀v ′ ∈ N }.

Theorem 3.8. Suppose M ⊂ V is a subspace of a Banach space V , and N ⊂ V ′ is a subspace of V ′ . Then 1.

⊥

(M ⊥ ) = M

k·kV

. ( ⊥ N )⊥ = N

weak∗ V ′

.

2. In addidtion, if M is closed, then M ′ = V ′ /M ⊥ ; 3. For T ∈ L (V, W ),

N (T ∗ ) = R(T )⊥ ;

(V /M )′ = M ⊥ .

N (T ) =

⊥

R(T ∗ ).

4. For T ∈ L (V, W ), N (T ∗ ) is weak∗ -closed in W ′ 5. For T ∈ L (V, W ), R(T ) is dense in W if and only if T ∗ is one-to-one. 6. For T ∈ L (V, W ), R(T ∗ ) is weak∗ -dense in V ′ if and only if T is one-to-one.

3.3

Big Theorems

Theorem 3.9 (Bohnenblust-Sobczyk Extension Theorem). Let V be a linear space with a seminorm | · |V . Let M be a linear subspace of V and f is a linear functional defined on M with |f (v)| ≤ |v|V for all v ∈ M. Then there exists a linearly extended fuctional fe defined on V such that |fe(v)| ≤ |v|V for all v ∈ V.

Theorem 3.10 (Hahn-Banach Extension Theorem in Normed Linear Spaces). Let V be a normed linear space with norm k · kV . Let M be a linear subspace of V and f is a continuous linear functional defined on M with |f (v)| ≤ kvkV for all v ∈ M. Then there exists a continuous linearly extended fuctional fe defined on V such that kfekV ′ = kf kV ′ .

Theorem 3.11 (Banach Open Mapping Theorem in Banach Spaces). Let V and W be Banach spaces. Let T ∈ L (V, W ) with T (V ) = W . Then T maps every open set of V onto an open set of W. Theorem 3.12 (Closed graph theorem). Let V, W be Banach spaces. Then if a map T : V → W is continuous, then the graph G = {(v, T v) : v ∈ V } is closed in V × W. Conversely, if a map T : V → W is linear and the graph G = {(v, T v) : v ∈ V } is closed in V × W, then the map T : V → W is continuous. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

42

Finite element methods

3.4

Section 3.4

Compact operators

Definition 3.10. For Banach spaces V, W , a linear map T : V → W is said to be compact if T (B(0; 1)) is compact in W, where B(0; 1) denotes the open unit ball in V. Notice the following facts: 1. A compact operator T : V → W is evidently bounded. 2. T : V → W is compact if and only if every bounded sequence {v n } in V contains a subsequence {v nk } such that {T v nk } converges to some w in W. Definition 3.11. For a Banach space V , suppose that T ∈ L (V ). The spectrum σ(T ) is the set of all λ ∈ C such that T − λI is not invertible. Any λ ∈ σ(T ) should then satisfy at least one of the following statement: 1. R(T − λI) 6= V ; 2. N (T − λI) 6= {0}; in this case λ is called an eigenvalue of T with dim(N (T − λI)) ≥ 1, and any v ∈ N (T − λI) is called an eigenvector of T associated with the eigenvalue λ. The following theorems are standard. For instance see Rudin [Rud91b]. Theorem 3.13. Let V, W be Banach spaces. 1. If T ∈ L (V, W ) and R(T ) is a finite dimensional subspace of W , then T is compact. 2. If T ∈ L (V, W ) is compact and R(T ) is a closed subspace of W , then R(T ) is a finite dimensional subspace of W. 3. If T ∈ L (V ) is compact, then the null space N (T − λI)) is of finite dimension for all nonzero λ. 4. If T ∈ L (V ) is compact in an infinite dimensional Banach space V , then 0 ∈ σ(T ). 5. If S, T ∈ L (V ) and T is compact, the composition S ◦ T is compact. 6. If T ∈ L (V, W ), then T is compact if and only if T ∗ is compact. [Schauder’s Theorem] 7. If T ∈ L (V ) is compact, then T − λI has closed range for all nonzero λ. 8. Suppose that T ∈ L (V ) is compact, and r > 0. Let E be a set of eigenvalues λ of T such that |λ| > r. Then E is a finite set, and R(T − λI) 6= V for each λ ∈ E. 9. Suppose that T ∈ L (V ) is compact. Then for all nonzero λ, dim(N (T − λI)) = dim(N (T ∗ − λI)) = dim(V /R(T − λI)) = dim(V ′ /R(T ∗ − λI)) < ∞. 10. Suppose that T ∈ L (V ) is compact. If λ ∈ σ(T ) is nonzero, then λ is an eigenvalue of T and T ∗ . 11. Suppose that T ∈ L (V ) is compact. Then σ(T ) is compact, at most countable, and has possibly at most one limit point 0. Each nonzero eigenvalue has finite multiplicity. Theorem 3.14 (Fredholm Alternative). Let T ∈ L (V ) be a compact operator in a normed linear space V . Then either 1. the homogeneous equation x − Tx = 0 has a nontrivial solution or c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

43

Finite element methods

Section 3.4

2. for each y ∈ V , the nonhomogeneous equation x − Tx = y is uniquely solvable. In this case, (I − T )−1 is also bounded. Proof. See [GT83].  Theorem 3.15 (Fredholm Alternative in Hilbert spaces). Let T ∈ L (V ) be a compact operator in a (real) Hilbert space. Then there exists a countable set Λ ⊂ R with no limit point except possibly λ = 0, such that 1. for each nonzero λ ∈ / Λ the equations λx − T x = y,

λx − T ∗ x = y

(3.3)

hava a unique solution x ∈ V for each y ∈ V . Moreover (λI − T )−1 and (λI − T ∗ )−1 are bounded. 2. for each λ ∈ Λ,

λx − T x = y is solvable if and only if y ⊥ N (λI − T ∗ ),

and λx − T ∗ x = y is solvable if and only if y ⊥ N (λI − T );

moreover, N (λI − T ) and N (λI − T ∗ ) are of finite dimension. Proof. See [GT83]. 

Definition 3.12. For more general linear operator T : V toV with D(T ) ⊂ V, set Tλ = λI − T,

λ ∈ C.

The resolvent set of T is defined by ρ(T ) = {λ ∈ C | R(Tλ ) is dense in V and (Tλ )−1 exists and is continuous} For λ ∈ ρ(T ), the resolvent of T is defined by (Tλ )−1 and denoted by R(λ : T ). The spectrum of T is the complement of ρ(T ) in C, i.e. σ(T ) = C \ ρ(T ). There are three kinds of spectra which are disjoint: Point spectrum Pσ (T ) = {λ ∈ C | Tλ has no inverse}: = the set of eigenvalues of T ; Continuous spectrum Cσ (T ) = {λ ∈ C | Tλ has discontinuous inverse Tλ−1 , D(Tλ−1 ) is dense in V }; Residual spectrum Rσ (T ) = {λ ∈ C | Tλ has an inverse Tλ−1 , D(Tλ−1 ) is not dense in V }. Theorem 3.16. [Yos95, p.209] If T : V → V is a closed linear operator, for any λ ∈ ρ(T ) the resolvent R(λ; T ) is an everywhere defined continuous linear operator, i.e. R(λ; T ) ∈ L (T ).

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

44

Chapter 4

Sobolev Spaces Multi-index and partial derivatives Let α = (α1 , · · · , αd ), αi ≥ 0 integers, be a multi-index, and set |α| = α1 + · · · + αd . Denote by ∂ α ϕ the partial derivative  α1  αd  ∂ ∂ |α| ϕ ∂ α D ϕ= . ··· ϕ = α1 d ∂x1 ∂xd ∂x1 · · · ∂xα d

4.1

Distributions

Definition of Distribution Definition 4.1. (φj )j ⊂ C0∞ (Ω) is said to converge in the sense of the space D(Ω) to the function φ ∈ C0∞ (Ω) if (i) there exists a compact subset K of Ω such that supp(φj ) ⊂ K

∀j;

(ii) for all multi-index α, ∂ α φj → ∂ α φ uniformly on K, (i.e. ∀ ǫ > 0, ∃ an integer N (ǫ) such that n ≥ N (ǫ) −→ |∂ α φm (x) − ∂ α φ(x)| < ǫ C0∞ (Ω)

∀x ∈ K.)

From now on, designate by D(Ω) the space equipped with the topology structure (complete locally convex topological vector space) given by (i) and (ii). Definition 4.2. The dual space D ′ (Ω) of D(Ω) is called the space of distributions. The dual space D ′ (Ω) is obviously the linear space (or vector space) of all continuous linear mapping T : D(Ω) → C. Thus if S, T ∈ D ′ (Ω), c ∈ C, then (T + S)(φ) (cT )(φ)

= =

T (φ) + S(φ), c(T (φ)),

T (φ + cψ)

=

T (φ) + cT (ψ)

for all φ, ψ ∈ D(Ω). T ∈ D ′ (Ω) if and only if lim T (φj ) = T (φj ) whenever D(Ω) − lim φj = φ in D(Ω).

j→∞

j→∞

′

T ∈ D (Ω) if and only if for any compact subset K ⊂ Ω there exist CK > 0 and NK > 0 such that T (φ) ≤ CK sup |∂ α φ| for all |α| ≤ NK , whenever φ ∈ D(Ω) with supp(φ) ⊂ K. K

The topology structure on D ′ (Ω) is given by the following convergence criterion: Tj → T in D ′ (Ω) ⇐⇒ Tj (φ) → T (φ) in C for all φ ∈ D(Ω). That is, the convergence of Tj → T in D ′ (Ω) is weak* convergence. 45

Finite element methods

Section 4.1

Definition 4.3. A function u defined almost everywhere on Ω is said to be locally integrable on Ω if u ∈ L1 (K) for every measurable compact subset K of Ω. In this case, we denote by u ∈ L1loc (Ω).

Example 4.1. f (x) = |x|, x ∈ (−∞, ∞) f is not differentiable at 0, but f ′ is defined almost everywhere on (−∞, ∞) and f ′ ∈ L1loc (−∞, ∞). Examples of distribution

Example 4.2. Regular distributions : For all u ∈ L1loc (Ω), there is an associated distribution Tu ∈ D ′ (Ω). Indeed, let u ∈ L1loc (Ω), and define the map Tu : D(Ω) → C by Z u φ dx, ∀φ ∈ D(Ω). Tu (φ) = Ω

Then for c ∈ C, φ, ψ ∈ D(Ω),

Tu (φ + cψ)

= =

Z

ZΩ Ω

=

u(φ + cψ)dx Z uφdx + c uψdx Ω

Tu (φ) + cTu (ψ).

Therefore Tu : D(Ω) → C is a linear map. To show that Tu is continuous, suppose that φj → φ in D(Ω) as j → ∞. Then there exists a K ⋐ Ω such that supp(φj ) ⊂ K, and supp(φ) ⊂ K. Thus, Z |Tu (φj ) − Tu (φ)| ≤ sup |φj (x) − φ(x)| |u(x)|dx → 0 x∈K

K

as n → ∞, since |φj − φ| → 0 uniformly. This shows that Tu is a distribution.

Example 4.3. Dirac delta distribution. For a ∈ Ω, define a linear form δa : D(Ω) → R by δa (φ) = φ(a)

∀φ ∈ D(Ω).

Exercise 4.1. Check that δa is a distribution. Example 4.4. Heaviside function Define H : R → R by ( 1, x≥0 H(x) = 0, x < 0. Then H ∈ L1loc (R). Therefore, by Example 4.2 there is an element TH ∈ D ′ (R). Indeed it is given by Z Z ∞ φ(x)dx. TH (φ) = H(x)φ(x)dx = R

0

Differentiation of distribution Definition 4.4. Given a distribution T ∈ D ′ (Ω), its weak or distributional derivative (or derivative in the ∂T ∈ D ′ (Ω) is defined by distribution sense) ∂x j ∂T ∂φ (φ) = −T ( ) ∂xj ∂xj

∀φ ∈ D(Ω).

In general, for T ∈ D ′ (Ω), ∂ α T ∈ D ′ (Ω) is given by

∂ α T (φ) = (−1)|α| (∂ α φ)

α

′

′

∀φ ∈ D(Ω).

(4.1)

The partial differentiation ∂ : D (Ω) → D (Ω) defined by (4.1) is continuous in the following sense: if Tj → T in D ′ (Ω), then ∂ α Tj → ∂ α T in D ′ (Ω).

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

46

Finite element methods

Section 4.2

Example 4.5. u(x) = 21 |x|, x ∈ (−∞, ∞) (i) Tu′ : D(R) 7→ R : let ϕ ∈ D(R) be arbitrary. Then set K = supp(ϕ), which is compact. "Z # Z Z 1 1 ′ ′ ′ ′ Tu (ϕ) = − |x| ϕ (x) dx = x ϕ (x) dx − x ϕ (x) dx 2 (−∞,0)∩K K 2 (0,∞)∩K # " Z Z 1 ϕ(x) dx ϕ(x) dx + − = 2 (0,∞)∩K (−∞,0)∩K Z 1 = (H(x) − ) ϕ(x) dx. 2 K Consequently, Tu′ = TH− 12 in D ′ (R) ′ (ii) TH : D(R) 7→ R : ∀ϕ ∈ D(R), ′ TH (ϕ)

(−1) TH (ϕ′ ) = −

= =

ϕ(0) = δ0 (ϕ)

Z

∞

ϕ′ (x) dx 0

where δ0 is Dirac delta distribution(function). Example 4.6. If f : R → R has a continuous and bounded derivative in R \ {x1 , · · · , xm }, with possible jumps at xk , k = 1, · · · , m. Let Jk = f (xk+ ) − f (xk− ), k = 1, · · · , m. Then m

X d Tf = T f J k δxk , df + dx dx k=1

f df where dx is regarded as an L1loc function defined everywhere by the derivative of f except at xk ’s. Indeed, setting x0 = −∞, xm+1 = ∞, for all φ ∈ D(R), Z m+1 m X Z xk df X dφ d Jk φ(xk ) Tf (φ) = − (x)φ(x)dx + f (x) dx = dx dx dx R k=1 xk−1 k=1 Z f m X df Jk φ(xk ). (x)φ(x)dx + = R dx k=1

2

L1loc (Ω),

Let v ∈ L (Ω). Then since v ∈ we identify it with its corresponding distribution Tv ∈ D ′ (Ω) defined in Example 4.2. From now on, all the derivatives are understood in the distributional sense.

4.2

Sobolev Spaces

Definition 4.5. For a nonnegative integer m and p ∈ [1, ∞] the Sobolev norms are defined by   p1 X kukm,p,Ω :=  kDα ukpp,Ω  , 1 ≤ p < ∞,

(4.2a)

0≤|α|≤m

kukm,∞,Ω

:=

max kDα uk∞,Ω .

0≤|α|≤m

(4.2b)

With the above Sobolev norms, there are two approaches to define Sobolev spaces of order (m, p) as follows: H m,p (Ω) W

m,p

(Ω)

W0m,p (Ω)

:= := :=

{u ∈ C m (Ω) : kukm,p < ∞} p

α

p

k·km,p

,

(4.3a)

{u ∈ L (Ω) : D u ∈ L (Ω) for 0 ≤ |α| ≤ m},

(4.3b)

C0∞ (Ω)

(4.3c)

W m,p (Ω)

,

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

47

Finite element methods

Section 4.2

Theorem 4.1. In general, H m,p (Ω) ⊂ W m,p (Ω).

The following salient theorem [MS64] removes pains to distinguish the above two definitions. Theorem 4.2 (H=W, Meyers-Serrin (1964)). For p ∈ [1, ∞),

H m,p (Ω) = W m,p (Ω).

Exercise 4.2. Let Ω = (−1, 1) and u(x) = |x|. Then show that u ∈ W 1,∞ (Ω) but u 6∈ H 1,∞ (Ω).

Due to Theorem 4.2, we will follow the definition of the original Sobolev spaces (due to Sobolev himself) W m,p (Ω). In particular if p = 2, we will use the notation H m (Ω) for W m,2 (Ω) so that H m (Ω) = W m,2 (Ω). Theorem 4.3. W m,p (Ω) is a Banach space for p ∈ [1, ∞] with the Sobolev norm k·km,p,Ω .

Theorem 4.4. W m,p (Ω) is a separable Banach space for p ∈ [1, ∞) and a reflexive Banach space for p ∈ (1, ∞) Corollary 4.1. H m (Ω) is a separable Hilbert space for all nonnegative integer m with the inner product X Z (u, v)m,Ω = Dα uDα v dx 0≤|α|≤m

Ω

Recall that Definition 4.6. A topological vector space V is called separable if there exists a countable dense subset S of V. Thus Theorem 4.4 implies that there exists a countable subset (fn )n of W m,p (Ω) such that for any f ∈ W m,p (Ω), there exists (αn )n ∈ C such that

∞

X

αnn = 0.

f −

n=1

m,p,Ω

Denote by | · |m,p,Ω the Sobolev seminorm defined by  p1  X |f |m,p,Ω =  k∂ α f k0,p,Ω  . |α|=m

Notice that

|f |m,p,Ω = 0

if

f (x) = polynomial of order m − 1 in Ω.

The spaces H 1 (Ω) and W 1,p (Ω) Example 4.7. The Sobolev space of order (1, 2) on Ω is given by H 1 (Ω) = {v ∈ L2 (Ω) :

∂v ∈ L2 (Ω), j = 1, · · · , d}, ∂xj

equipped with the inner product (u, v)1,Ω = and the norm

Z

uvdx + Ω

d Z X

k=1

Ω

∂u ∂v dx, u, v ∈ H 1 (Ω). ∂xk ∂xk 1

||u||1,Ω = [(u, u)1,Ω ] 2 . c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

48

Finite element methods

Section 4.2

Example 4.8. If m = 1,

kf k1,p,Ω =



kf kp 0,p,Ω

n X

∂f p

+

∂x

0,p,Ω

j=1

 p1



 =

Z

Ω

 p1 p n X ∂f dx , |f |p + ∂x j=1

where the derivatives are understood in the sense of distribution.

Example 4.9. H 1 (Ω) is a Hilbert space for the inner product (·, ·)1,Ω . ∂v

Proof. Suppose {vj } ⊂ H 1 (Ω) is a Cauchy sequence. Then {vj } and { ∂xkj } for k = 1, · · · , d are Cauchy sequences in L2 (Ω). Since L2 (Ω) is complete, there exist v ∈ L2 (Ω) and v (k) ∈ L2 (Ω), 1 ≤ k ≤ d, such that ∂v vj → v in L2 (Ω) and ∂xkj → v (k) in L2 (Ω).

∂v We have to prove that v (k) = ∂x in the sense of distribution on Ω. Since I : L2 (Ω) → D ′ (Ω) defined k by I (φ) = φ for all φ ∈ L2 (Ω) is continuous,

I vj → I v ∂vj I → I v (k) ∂xk Since

∂ ∂xk

in

D ′ (Ω)

in

D ′ (Ω) k = 1, · · · , d.

: D ′ (Ω) → D ′ (Ω) is continuous, I

∂vj ∂I vj ∂I v ∂v = → =I ∂xk ∂xk ∂xk ∂xk

Since the limit is unique in D ′ (Ω), we have I k = 1, · · · , d.

∂v ∂xk

in D ′ (Ω).

= I v (k) ∈ D ′ (Ω); and hence

∂v ∂xk

= v (k) ∈ L2 (Ω) for each

Include a commuting diagram here. Therefore v ∈ H 1 (Ω), and ||vj − v||1,Ω =

"Z

Ω

|vj − v|2 dx +

d Z X

k=1

∂v 2 ∂vj − | dx | ∂x ∂x k k Ω

# 21

→ 0 asj → ∞,

which implies that vj → v in H 1 (Ω) as j → ∞.  Example 4.10. H 1 (Ω) is separable (i.e. there exists a countable dense subset in H 1 (Ω).) For a general domain Ω, D(Ω) is not dense in H 1 (Ω). The closed set Rd \ Ω should satisfy (m, p′ )-polar condition, for example; i.e., a distribution T ∈ D(Ω) having support in Rd \ Ω is the zero distribution. However, we recall Definition 4.5, from which we have Definition 4.7. H01 (Ω) = D(Ω)

||·||1,Ω

Corollary 4.2. D(Rd ) is dense in H 1 (Rd ) i.e., H01 (Rd ) = H 1 (Rd ). Proof.  Corollary 4.3. For 1 ≤ p < ∞, W m,p (Ω) is separable, Theorem 4.5. D(Ω) is dense in W m,p (Ω) if and only if the complement Ωc is (m, p′ )-polar. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

49

Finite element methods

4.2.1

Section 4.2

Schwartz space S(Rd ) S(Rd )

:=

{f ∈ C ∞ (Rd ) | sup |xβ ∂ α f (x)| < ∞ x∈Rd

∀α, β}

{f ∈ C ∞ (Rd ) | sup |(1 + |x|2 )N ∂ α f (x)| < ∞

=

x∈Rd

∀α, N = 0, 1, 2, · · · }. 2

which is the space of rapidly decreasing functions at ∞. Examples include any function in C0∞ (Rd ), e−|x| , and so on. For each polynomial p(x) and multi-index α, N : S(Rd ) → R+ defined by Np (f ) := sup |p(x)∂ α f (x)| x∈Rd

d

forms a seminorm on S(R ). Indeed, due to the rapidly decreasing nature, the seminorms are norms on S(Rd ). Proposition 4.1. [Yos95, p. 146] S(Rd ) is closed with respect to the application of linear partial differential operators with polynomial coefficients. Proposition 4.2. [Yos95, p. 146] C0∞ (Rd ) is a dense subset of S(Rd ) with respect to the topology of S(Rd ).

4.2.2

Fourier transformation

For any f ∈ S(R)∪L1 (R), the Fourier transform fb is defined by fb(ω) := formula Z ∞ 1 fb(ω)eiωt dω. f (t) = 2π −∞

R∞

−∞

f (t)e−iωt dt, with the inversion

For higher dimension, Fourier transforms are defined analogously: for f ∈ S(Rd ) ∪ L1 (Rd ), the fb : Rd 7→ C is defined by Z ∞ Z ∞ b f (x)e−iω·x dx1 · · · dxd ··· F(f )(ω) = f (ω) := −∞

−∞

with the inversion formula F

−1

(fb)(x) :=



1 2π

d Z

∞ −∞

···

Z

∞ −∞

Here, ω = (ω1 , · · · , ωd ) and x = (x1 , · · · , xd ).

fb(ω)eiω·x dω1 · · · dωd .

(4.4)

Proposition 4.3. [Yos95, p. 146] In particular, both Fourier transformation and its inverse transformation are linear and continuous from S(Rd ) to S(Rd ). That is, the Fourier transform F : S → S which maps φ to φb is an isomorphism. Theorem 4.6 (Fourier Integral Theorem). [Yos95, p. 147]

F−1 ◦ F = IS(Rd ) .

One of merits of Fourier transformation is to transfer derivatives in f (x) into algebraic factor in fb(ω). Given f ∈ W m,2 (Ω) d Z d Z d 2 X X ∂f 2 ∂f dω < ∞ ∂xj dx < ∞ ⇐⇒ ∂x j j=1 Ω j=1 R ⇐⇒

⇐⇒ c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

d Z X j=1

Z

R

R

2 ωj2 fb(ω) dω < ∞

2 |ω|2 fb(ω) dω < ∞

50

Finite element methods

Section 4.2

Theorem 4.7 (Parseval formula). kfbk20 = (2π)d kf k20 .

(4.5)

Theorem 4.8 (Plancherel Theorem). The Fourier transform, defined originally on L1 (Rd )∩L2 (Rd ) extends uniquely to a map from L2 (Rd ) to L2 (Rd ) satisfying

for all f, g ∈ L2 (Rd ).

(fb, gb) = (2π)d (f, g),

kfbk20 = (2π)d kf k20 ,

(4.6)

We will extend the Fourier transformation to (temperate) distributions. The space S is endowed with the topology defined by the seminorm |φ|S := sup |xβ ∂ α φ(x)| α,β,x

makes S a Fr´echet space. A space V is called a Fr´echet space if it is a locally convex space whose topology is induced by a complete invariant metric dist(·, ·). Definition 4.8. A temperate distribution S′ is the set of all continuous linear functional u on S. We are now in a position to define the Fourier transformation for temperate distributions. Definition 4.9. For u ∈ S′ , the Fourier transform u b is defined by b u b(φ) = u(φ)

∀φ ∈ S.

(4.7)

We then have the following theorems.

Theorem 4.9. The Fourier transformation F : S′ → S′ is an isomorphism (with the weak topology) with the Fourier inversion of the form (4.4) for all u ∈ S′ . Example 4.11. [Rud91a, p.190] 1. All distributions with compact support are temperate distributions. 2. Lp (Rs ) ⊂ S(Rd ) for p ∈ [1, ∞]. 3. Pk (Rs ) ⊂ S(Rd ), k = 0, 1, · · · . The polynomials are temperate distributions. 4. A Borel measure µ on Rd such that Z

Rd

1 dµ(x) < ∞, (1 + |x|2 )k

for some positive integer k, is a temperate distribution. 5. A measurable function g on Rd such that Z

Rd

|g(x)|p dx < ∞, (1 + |x|2 )k

for some positive real number k, is a temperate distribution. Theorem 4.10. Let u ∈ S′ and v ∈ D′ (Rd ) with supp(v) compact. Then u ∗ v ∈ S′ and the Fourier transform F(u ∗ v) is equal to u bvb. Property 4.1. For u ∈ S′ , the following holds:

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

51

Finite element methods

Section 4.2

d 1. [ b. dxj u = iωj u

d b. 2. x d j u = −i dωj u

For any real number s, the Sobolev space can be defined through Fourier transformation: Z H s (Rd ) := {f : Rd 7→ R | (1 + |ω|2 )s |fb(ω)|2 dω < ∞} Rd

Example 4.12. Every distribution with compact support is in H s (Rd ).

Definition 4.10. Let X = ∪s∈R H s (Rd ). A linear operator T : X → X is said to have order t if TH s (Rd ) : H s (Rd ) → H s−t (Rd ) is a continuous surjective mapping. Example 4.13.

1. Every partial differential operator ∂ α is an operator with order |α|.

2. For t ∈ R, the mapping T : u 7→ v given by t

vb(ω) = (1 + |ω|2 ) 2 u b(ω)

is a linear isometry of H s (Rd ) onto H s−t (Rd ) with order t. The inverse mapping T −1 : u 7→ v given by t vb(ω) = (1 + |ω|2 )− 2 u b(ω)

is a linear isometry of H s (Rd ) onto H s+t (Rd ) with order −t.

4.2.3

Prolongation (Extension)

Theorem 4.11. If v ∈ H01 (Ω), the function ( belongs to H 1 (Rd ).

ve(x) =

Proof. If ϕ ∈ D(Ω), then

belongs to D(Rd ). Moreover,

ϕ(x) e =

(

v(x) 0

x∈Ω x ∈ Rd \ Ω

ϕ(x) 0

x∈Ω x ∈ Rd \ Ω

kϕk e 1,Rd = kϕk1,Ω .

The map D(Ω) → H 1 (Rd ) : ϕ 7→ ϕ e is linear and continuous with respect to the norm induced by k · k1,Ω . Therefore this map extends to the map H01 (Ω) → H 1 (Rd ) : v 7→ ve continuously and linearly. 

Lemma 4.1 (Poincar´e inequality). If Ω is bounded, then there exists a constant C = C(Ω) such that

||v||0,Ω

1/2

d X

∂v 2



≤C

∂xj 

j=1

0,Ω

∀ v ∈ H01 (Ω).

Remark 4.1. v = 1 on Ω does not satisfy Poincar´e inequality. Thus v ∈ H 1 (Ω) is not sufficient. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

52

Finite element methods

Section 4.2

Corollary 4.4 (Poincar´e Inequality). kvk1,2,Ω



d X

∂v 2

 ≤C

∂xj

0,2,Ω

j=1

where C is independent of v, but depends only on Ω.

 12 

∀v ∈ H01 (Ω),

(4.8)

Remark 4.2. This corollary is a general version for Lemma 4.1 since 2

2

2

2

2

2

kvk1,Ω = kvk0,Ω + k∇vk0,Ω ≤ C 2 k∇vk0,Ω + k∇vk0,Ω = (C 2 + 1) k∇vk0,Ω . In H01 (Ω), by Corollary 4.5 k·kH 1 (Ω) ≃ | · |H 1 (Ω) , i.e., kvkH 1 (Ω) ≤ C1 |v|H 1 (Ω) ≤ C2 kvkH 1 (Ω) ,

∀v ∈ H01 (Ω)

Corollary 4.5. If Ω is bounded, the seminorm 

|v|1,Ω := 

d X j=1

k

is equivalent to the Sobolev norm || · ||1,Ω in H01 (Ω). (i.e. Ckvk1 ≤ |v|1 ≤ C ′ kvk1 for some C, C ′ ∈ R)

4.2.4

1/2

∂v 2  k ∂xj 0,Ω

The definition of traces and trace theorems

For 1-dimensional case, let Ω = (a, b). Then if v ∈ H 1 (Ω) then v ∈ C 0 (Ω). In this case, it is not difficult to define v|Γ i.e., the values v(a) and v(b) for all v ∈ H 1 (Ω). However, for d-dimensional case, d ≥ 2, it is not trivial to define boundary values of v ∈ H 1 (Ω). Indeed, for R < 1, let Ω = {x ∈ R2 : 0 < |x| < R} and v(x) = | log |x||p . Then  q p−1 xj ∂v = p log x21 + x22 2 ∂xj x1 + x22

||v||21,Ω

= =

Z

|v|2 + |∇v|2 dx Z Z R | log r|2p rdr + 2p2 π 2π 0


0. Let Ω = Rd+ . Then Γ = {x = (x′ , 0); x′ ∈ Rd−1 }. Lemma 4.2. D(Rd+ ) is dense in H 1 (Rd+ ). Proof.  Definition 4.11. For an integer m ≥ 1, an open set Ω ⊂ Rd is said to be m-regular if its boundary Γ is a manifold of class C m of dimension d − 1, Ω being locally at one side of Γ. Or equivalently, if there exists a are C m maps and finite open cover {Oj }Jj=1 of Ω and invertible ϕj = Oj → B(0; 1) such that ϕj and ϕ−1 j ϕj (Oi ∩ Ω)

ϕj (Oi ∩ Γ)

= B(0; 1) ∩ Rd+ = {y = (y ′ , yd ) ∈ Rd : |y ′ | < 1, yd > 0} = {y = (y ′ , yd ) ∈ Rd : |y ′ | < 1, yd = 0}.

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

53

Finite element methods

Section 4.2

Lemma 4.3. If Ω is 1-regular, then there exists a continuous linear extension operator E : H 1 (Ω) → H 1 (Rd ) such that E v(x) = v(x) ∀x ∈ Ω Proof. (1st step) We begin by considering the case Ω = Rd+ . If v ∈ D(Rd+ ), define E v on Rd by (reflection) ( v(x′ , xd ), if xd ≥ 0 ′ E v(x , xd ) = v(x′ , −xd ), if xd ≤ 0. Then E v is continuous, E v ∈ H 1 (Rd ) and  ∂v  (x′ , xd ),    ∂xj   ∂v ∂ (x′ , xd ), E v(x′ , xd ) = ∂x  ∂xj j    ∂v ′  − (x , −xd ), ∂xj

Then for all v ∈ D(Rd+ ), we have

kE vk1,Rd =

j = 1, · · · , d − 1, if xd > 0, j = d, if xd < 0, j = d.

√ 2 kvk1,Rd . +

Since E can be extended by linearly and continuously to E : H 1 (Rd+ ) → H 1 (Rd ). is dense in H (2nd step) Now assume that Ω is an open 1-regular domain in Rd . We know that there exists a partition of unity {αj }Jj=0 subordinate to the cover {Oj }Jj=1 of the boundary Γ of Ω; i.e., 1

D(Rd+ )

(Rd+ ),

αj ∈ D(Oj ), 1 ≤ j ≤ J;

J X j=0

αj = 1 on Ω, 0 ≤ αj ≤ 1,

and Γ ⊂ ∪Jj=1 Oj .

1

Thus, for v ∈ H (Ω) we can write v=

J X

αj v,

j=0

we can define E (αj v) for all j = 0, 1, · · · , J, and then define E v by Ev =

J X

E (αj v).

j=0

First, we have d E (α0 v) = α g 0 v : the extension of α0 v by 0 to R \ Ω.

Then for j = 1, · · · , J, consider

wj = (αj v) ◦ (ϕ−1 j |B+ ),

B+ = B(0; 1) ∩ Rd+ .

We have wj ∈ H 1 (B+ ) and wj = 0 in a neighborhood of {y ∈ ∂B+ : yd > 0}. One can extend wj by 0 in Rd+ \ B+ . Let w fj denote this extension. Then w fj ∈ H 1 (Rd+ ), which is again d extended to R by reflection (as in the first step) f w fj ∈ H 1 (Rd )

f and supp(w fj ) ⊂ B(0; 1).

^ f f Then let w fj ◦ ϕ ∈ H 1 (Rd ): the extension of w fj ◦ ϕ by 0 in Rd \ Oj . We have Then v 7→

PJ

j=0

^ f E (αj v) = w fj ◦ ϕj ,

1 ≤ j ≤ J.

E (αj v) is a 1-extension. 

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

54

Finite element methods

Section 4.2

Lemma 4.4. If Ω is 1-regular, D(Ω) is dense in H 1 (Ω). Proof. Let v ∈ H 1 (Ω). Then by the Lemma 4.3, there exists a function E v ∈ H 1 (Rd ) such that E v = v on Ω. By Theorem 4.10, H01 (Rd ) = H 1 (Rd ), there exists a sequence {wj }j ⊂ D(Rd ) such that wj → E v in H 1 (Rd ). Let vj be the restriction of wj to Ω. Then (vj )j ⊂ D(Ω)

vj → v in H 1 (Ω).

s.t.

Therefore D(Ω) is dense in H 1 (Ω).  Define kvk0,Γ =

Z

2

Γ

|v(x)| dσ

1/2

.

(4.9)

By using a partition of unity, {αk }K k=0 , we have 2 d−1 ), 1 ≤ k ≤ K}, L2 (Γ) = {v : Γ → R : (αk^ v) ◦ ϕ−1 k (·, 0) ∈ L (R d−1 where (αk^ v) ◦ ϕ−1 \ {y ′ ∈ Rd−1 : |y ′ | < 1}. Then, the map given by k (·, 0) is the extension by 0 to R

"

K X

^ −1 2

v→ 7

(αk v) ◦ ϕk

0,Rd−1

k=0

#1/2

(4.10)

defines a norm which is equivalent to the norm given by k·k0,Γ defined in (4.9). Consider the mapping γ0 : D(Ω) → L2 (Γ) defined by γ0 v = v|G

∀v ∈ D(Ω).

Lemma 4.5. If Ω is 1-regular, then there exists C > 0 such that |γ0 v|0,Γ ≤ C||v||1,Ω

∀ v ∈ D(Ω).

where C is independent of v, but dependent only on Ω. Proof. Let v ∈ D(Ω). By the partition of unity {αk }K k=0 , we put wk = (αk v) ◦ ϕ−1 k , 0 ≤ k ≤ K. Since ||v(·, 0)||0,Rd−1 ≤ ||v||1,Rd+ for all v ∈ D(Rd+ ), fk k1,Rd . kw fk (·, 0)k0,Rd−1 ≤ kw +

Since αk and

ϕ−1 k

are smooth, we can find

kw fk k1,Rd ≤ Ck kvk1,Ω . +

Thus

||w fk (·, 0)||0,Rd−1 ≤ Ck kvk1,Ω .

By the equivalence of the two norms k·k0,Γ given by (4.9) and (4.10), we are done.  Theorem 4.12 (Trace Theorem). Suppose that Ω is 1-regular. Then D(Ω) is dense in H 1 (Ω) and γ0 : v 7→ γ0 v = v|Γ from D(Ω) to L2 (Γ) can be extended linearly and continuously to a map, again denoted by γ0 γ0 : H 1 (Ω) → L2 (Γ). Remark 4.3. The trace theorem holds for domain Ω ⊂ Rd which is Lipschitz domain or piecewise 1-regular. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

55

Finite element methods

4.2.5

Section 4.2

Application of Trace Theorem

Theorem 4.13. Suppose Ω is open bounded in Rd with its boundary Γ being piecewise C 1 . Then H01 (Ω) H01 (Ω)

= =

Ker(γ0 ) i.e., {v ∈ H 1 (Ω) : γ0 v = v|Γ = 0 on Γ}.

Proof. (1) H01 (Ω) ⊂ {v ∈ H 1 (Ω) : v|Γ = 0}. Indeed, let v ∈ H01 (Ω). Then there exists {ϕj } ⊂ D(Ω) such that ϕj → v in H 1 (Ω). Since γ0 is continuous, L2 (Γ) − lim γ0 ϕj = γ0 (H 1 (Ω) − lim vj ) = γ0 v j→∞

j→∞

Since γ0 ϕj = 0 on Γ for all j, we have γ0 v = 0. (2) To show H01 (Ω) ⊃ {v ∈ H 1 (Ω) : v|Γ = 0}. Using the local property and partition of unity, it suffices to prove the case where v ∈ H 1 (Rd+ ) and has a compact support in Rd+ satisfying v(·, 0) = 0. It suffices to show that v ∈ H01 (Rd+ ). Let ve be the extension of v by 0 to Rd \ Rd+ . One can verify that ve ∈ H 1 (Rd ) by showing that g ∂e v ∂v = , ∂xj ∂xj

in the sense of distribution (we already know that implies

∂e v ∂xj

1 ≤ j ≤ d, ∂v ∂xj

∈ L2 (Rd )). In order to show (4.11), we need

(4.11)

∂v 2 d ∈ L2 (Rd+ ) and hence g ∂xj ∈ L (R ). Thus (4.11)

Lemma 4.6. If v, ϕ ∈ H 1 (Rd+ ), we have the following version of Green’s formula Z Z ∂ϕ ∂v dx = − ϕdx, 1 ≤ j ≤ d − 1, v d d ∂x ∂x j j R+ R+ Z Z Z ∂v ∂ϕ dx = − ϕdx − v(x′ , 0)ϕ(x′ , 0)dx′ . v d ∂xd ∂x d−1 d R R Rd + + We begin by obvserving that for ϕ ∈ D(Rd ), Z Z ∂ϕ ∂ϕ ∂e v ∂ϕ ve v (ϕ) = −e v( )=− dx = − dx, 1 ≤ j ≤ d. d ∂xj ∂xj ∂x ∂x d j j R R+

Since v(·, 0) = 0, Lemma 4.6 implies that Z g ∂v ∂v ∂e v (ϕ) = ϕ dx = (ϕ), 1 ≤ j ≤ d. d ∂xj ∂x ∂x j j R+

This proves (4.11), and hence ve ∈ H 1 (Rd ). For all h > 0, consider the translattion operator τh : H 1 (Rd ) → H 1 (Rd ) defined by τh ve(x′ , xd ) = ve(x′ , xd + h).

∂ ∂e v Since ∂x τh ve = τh ∂x , 1 ≤ j ≤ d, τh ve ∈ H 1 (Rd ). Moreover, one can verify that τh ve → ve in H 1 (Rd ) as j j h → 0. We introduce a regularization {ρε } ⊂ D(Rd ). Then ρε ∗ τh ve → τh ve in H 1 (Rd ) as ε → 0. And for sufficiently small ε > 0, ρε ∗ τh ve has a compact support in Rd+ , from which one can construct a sequence of functions in D(Rd+ ) that converges to v in H 1 (Rd+ ). Therefore, v ∈ H01 (Rd+ ).

Proof of Lemma 4.6. It is easy to see that (i) and (ii) hold for v, ϕ ∈ D(Rd+ ). Since D(Rd+ ) is dense in H 1 (Rd+ ) and the trace is continuous (γ0 : H 1 (Rd+ ) → L2 (Rd−1 ). 

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

56

Finite element methods

Section 4.3

Let ν = (ν1 , · · · , νd ) denote the unit outward normal to Γ. Theorem 4.14. Let Ω be open bounded with a piecewise C 1 boundary Γ. Then the map v 7→ (γ0 v, γ1 v) = ∂v |Γ ) from D(Ω) to L2 (Γ) × L2 (Γ) can be extended to a continuous linear map (v|Γ , ∂ν H 2 (Ω) → L2 (Γ) × L2 (Γ). Theorem 4.15. Let Ω be bounded open with piecewise C 1 boundary Γ. Then if u, v ∈ H 1 (Ω), Z

Ω

∂u vdx = − ∂xj

Z

∂v dx + u ∂x j Ω

Z

νj γ0 (uv)dσ, Γ

1 ≤ j ≤ d.

Proof. If u, v ∈ H 1 (Ω), there exist uj , vk ∈ D(Ω) such that uj → u and vk → v, in H 1 (Ω) For all uj , vk ∈ D(Ω), we have Z

Ω

∂uj vk dx = − ∂xl

Z

∂vk dx + uj ∂xl Ω

Z

νl γ0 (uj vk )dσ, Γ

1 ≤ l ≤ d.

By density and continuity, we have Green’s formula. 

4.3

The embedding theorems and fractional Sobolev spaces

Let Ω ∈ Rd . For m ∈ Z+ , define the Sobolev norms 

kukm,p := 

X

0≤|α|≤m

 p1

k∂ α ukpLp (Ω) 

and kukm,∞ :=

max k∂ α ukL∞ (Ω) .

0≤|α|≤m

Then for m ∈ Z+ and 1 ≤ p ≤ ∞, the spaces W m,p (Ω) = {u ∈ Lp (Ω) : kukm,p < ∞} are Banach spaces. The fractional Sobolev spaces are then defined as follows. If s = m + σ, σ ∈ [0, 1), and with m ∈ Z+ , the W m+σ,p (Ω) consists of all functions u ∈ W m,p (Ω) satisfying Z Z Ω

p

Ω

|∂ α u(x) − ∂ α u(y)| dx dy < ∞, |x − y|d+σp

|α| = m.

Example 4.14. Let Γ = (−1, 1) and let u(x) = 0 for −1 < x < 0 and u(x) = 1 for 0 < x < 1. We want to see in which class u belongs to. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

57

Finite element methods

Section 4.3

Then certainly u ∈ L2 (Γ) and

Z Z Γ

2

|∂ α u(x) − ∂ α u(y)| dx dy |x − y|1+2σ Γ Z 1Z 1 α 2 |∂ u(x) − ∂ α u(y)| = dx dy |x − y|1+2σ −1 −1 Z 0Z 1 1 dx dy =2 1+2σ −1 0 (x − y) Z 0Z 1 1 =2 dx dy 1+2σ (x − y) −1 0 Z 0 2 1 1 = − dy 2σ 2σ −1 y (1 − y)2σ Z 1 1 1 0 = − dy σ −1 (−y)2σ (1 − y)2σ Z 1 0 1 1 = − dy σ −1 (−y)2σ (1 − y)2σ Z 1 1 1 1 = − dy 2σ σ 0 (y) (1 + y)2σ 1 (2 − 21−2σ ) < ∞ = σ(1 − 2σ)

1

provided 1 − 2σ > 0. Therefore u ∈ W 0, 2 −ǫ for all ǫ > 0. Theorem 4.16. [Sobolev Embedding Theorem] Suppose Ω ⊂ Rd is an open Lipschitz domain. Then for p ∈ [1, ∞) and integers m ≥ n, set the critical Sobolev index γ = γ(m, n, p) =  n,q  W (Ω) n,q W m,p (Ω) ⊂ Wloc (Ω)   n C (Ω),

1 m−n − . p d

for q = γ1 ∀q ∈ [1, ∞)

if γ > 0, if γ = 0, if γ < 0.

(4.12)

Moreover, if Ω is bounded, then

W m,p (Ω) ⊂ C n (Ω) if γ < 0,

(4.13)

and the following embedding is compact: ′

W m,p (Ω) ⊂ W n,q (Ω)

(4.14)

for all q ′ ∈ [1, γ1 ) if γ > 0, or for all q ′ ∈ [1, ∞) if γ = 0. Thus if Ω ∈ R2 is an open Lipschitz domain, then H 1 (Ω) ⊂ Lqloc (Ω)

∀q ∈ [1, ∞)

since γ =

and in general H 1 (Ω) 6⊂ C 0 (Ω)

since γ =

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

1 1 − = 0. 2 2

1 1 − = 0. 2 2 58

Finite element methods

Section 4.4

but W 1,p (Ω) ⊂ C 0 (Ω)

∀p > 2.

since γ =

1 1 − < 0. p 2

If Ω ∈ R3 is an open Lipschitz domain, then H 1 (Ω) ⊂ L6 (Ω)

since γ =

1 1 2 − = . 2 3 6

and in general H 1 (Ω) 6⊂ C 0 (Ω)

since γ =

1 1 1 − = , 2 3 6

H 2 (Ω) ⊂ C 0 (Ω)

but W 1,p (Ω) ⊂ C 0 (Ω)

∀p > 3 since γ =

since γ =

1 2 − < 0, 2 3

1 1 − < 0. p 3

Theorem 4.17. Suppose that Ω is a bounded open subset of Rd with Lipschitz boundary. Let m1 , m2 , m ∈ Z+ be nonnegative integers and p1 , p2 , p ∈ [1, ∞) be real numbers such that m1 ≥ m,

m2 ≥ m.

Suppose that either 1.

m1 +m2 −m d

≥

1 p1

+

1 p2

− p1 ,

mj −m d

>

1 pj

− p1 ,

j = 1, 2; or

2.

m1 +m2 −m d

>

1 p1

+

1 p2

− p1 ,

mj −m d

≥

1 pj

− p1 ,

j = 1, 2.

Then the multiplication mapping u, v → uv is continuous bilinear map from W m1 ,p1 (Ω) × W m2 ,p2 (Ω) into W m,p (Ω). Compactness results Theorem 4.18. Let Ω be bounded and open in RN with a piecewise C 1 boundary Γ. Then the canonical injection I : H 1 (Ω) → L2 (Ω) is compact (i.e. if (uj )j ⊂ H 1 (Ω) with ||uj ||1,Ω ≤ C for all j for some constant C > 0, then (uj )j is relatively compact in L2 (Ω).) Proof. A classical result for Hilbert spaces: Any bounded set of a Hilbert space V is weakly relatively compact. This means if (fj )j ⊂ V is bounded in V , there exists f ∈ V and a subsequence (fmj )j such that fmj → f weakly in V i.e., (g, fmj ) → (g, f ) ∀ g ∈ V.

For the case where Ω is open and 1-regular, there is a continuous linear extension operator E : H 1 (Ω) → H 1 (RN ) such that E u = u a.e. on Ω ∀ u ∈ H 1 (Ω) By truncating E u, E u has a support contained in a fixed compact set K, K ⊃ Ω.  Remark 4.4. If Ω is unbounded, the above theorem may not hold. For instance, let Ω = RN . Then H 1 (RN ) ֒→ L2 (RN ) is not compact. Remark 4.5. The assumption on the regularity of the domain for the above theorem can be weakened. Theorem 4.19. Let Ω be open, bounded in RN . Then H01 (Ω) ֒→ L2 (Ω) is compact. c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

59

Finite element methods

4.4 Let Ω =

Section 4.4

Application to elliptic problems and finite element subspace of H 1 (Ω) SK

k=1

Ωk such that

1. Ωk is open in RN , ∂Ωk being piecewise C 1 , k = 1, · · · , K. 2. Ωi ∩ Ωj = ∅ if i 6= j. Theorem 4.20. Let u ∈ C 0 (Ω) such that u|Ωi ∈ H 1 (Ωi ), i = 1, · · · , M . Then u ∈ H 1 (Ω) Proof. For each j, define uj ∈ L2 (Ω) by uj |Ωr =

sense of D ′ (Ω). Then

∂u ∂xj

∂ ∂xj (u|Ωk ),

1 ≤ k ≤ K. We want to see uj =

∂u ∂xj

in the

∈ L2 (Ω). To show this, for all ϕ ∈ D(Ω),

uj (ϕ)

=

=

Z

K Z X

K Z X

∂ (u|Ωk )ϕ dx ∂xj Ω k=1 Ωk k=1 Ωk  Z Z K Z X ∂ϕ ∂ϕ u uϕνj dσ = − dx − dx u|Ωk − ∂xj Ω ∂xj Γk Ωk uj ϕ dx =

uj ϕ dx =

k=1

= Therefore, uj =

∂u ∂xj

−u(

∂ϕ ∂u )= (ϕ). ∂xj ∂xj

∈ L2 (Ω), which shows that u ∈ H 1 (Ω). 

Recall that ν = (ν1 , · · · , νd ) denotes the unit outward normal to Ω on its boundary Γ. Theorem 4.21 (Green’s Theorem). Let Ω be a bounded open subset of Rd with a Lipschitz-continuous boundary Γ. 1. If u, v ∈ H 1 (Ω), then

Z

u Ω

∂v dx = − ∂xj

Z

Ω

∂u v dx + ∂xj

Z

γ0 (uv)νj ds.

(4.15)

Γ

2. If u, v ∈ H 1 (Ω) and u ∈ H 2 (Ω), then we have   Z Z Z ∂u v νj dσ. − ∆uvdx = ∇u · ∇vdx − γ0 ∂xj Ω Ω Γ

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

60

Chapter 5

Appendix 5.1

Brief review on mathematical analysis and linear algebra

Unless otherwise stated, by a field F we mean the real number field R or the complex number field C. Definition 5.1. A vector space V over a field F consists of a set of vectors with the two “operations:” (i) the vector addition “+” and (ii) the scalar multiplication such that V is closed under these operations fulfilling the following properties: 1. (Closed under “+”) ∀u, v ∈ V ,

u+v ∈V

2. (Closed under the scalar multiplcation) ∀α ∈ F

∀u ∈ V ,

αu ∈ V

3. (Existence of zero element, or the identity element for +) ∃0 ∈ V such that 0 + u = u = u + 0

∀u ∈ V

4. (Existence of the inverse element for +) ∀u ∈ V , ∃ − u ∈ V such that u + (−u) = 0 = (−u) + u 5. (The identity element of F ) 1u = u 6. (Commutative law) u + v = v + u

∀u ∈ V

∀u, v ∈ V

7. (Associative law) (u + v) + w = u + (v + w) 8. (Associative law) (αβ)u = α(βu)

∀u, v, w ∈ V

∀α, β ∈ F

9. (Distributive law) α(u + v) = αu + αv 10. (Distributive law) (α + β)u = αu + βu

∀u ∈ V

∀α ∈ F

∀u, v ∈ V

∀α, β ∈ F

∀u ∈ V

Definition 5.2. Let V be a vector space over F . A nonempty subset W of V is called a subspace if the following closure rules hold for W : 1. (Closed under “+”) ∀u, v ∈ W ,

u+v ∈W

2. (Closed under the scalar multiplcation) ∀α ∈ F

∀u ∈ W ,

αu ∈ W .

Example 5.1. F is a vector space over itself. Example 5.2. Let d be a positive integer. F d is a vector space over F , where F d = {hv1 , · · · , vd i |vj ∈ F Example 5.3. Let 1 ≤ p ≤ ∞, and Ω ⊂ Bd be an open region. Define (R 1 1 ≤ p < ∞, |f (x)|p dx p Ω kf k0,p,Ω = ess.supx∈Ω |f (x)| p=∞ 61

∀j = 1, · · · , d}

Finite element methods

Section 5.1

and set Lp (a, b) = {f : Ω 7→ R | kf k0,p,Ω < ∞}.

Then Lp (a, b) is a vector space over R. Indeed, the Minkowski inequality kf + gk0,p,Ω ≤ kf k0,p,Ω + kgk0,p,Ω

∀f, g ∈ Lp (a, b).

is the key to verify the closure property under +. Definition 5.3. Let V and W be vector spaces over F . A map T : V 7→ W is called linear if the following holds: 1. (Preserve the “+” structure) ∀u, v ∈ W ,

T (u + v) = T (u) + T (v)

2. (Preserve the scalar multiplcation structure) ∀α ∈ F

∀u ∈ W ,

T (αu) = αT (u).

Notational simplification: If T is a linear map, usually T (u) is abbreviated as T u. Definition 5.4. Let V be a vector space over F . An inner product (·, ·) is a map V × V 7→ F defined by (v, w) ∈ F . Definition 5.5. Let V be a vector space over F . A norm k · kV : V → R is a function such that 1. kukV > 0

∀0 6= v ∈ V ;

2. ku + vkV ≤ kukV + kvkV 3. kαukV = |α| kukV

∀u, v ∈ V ;

∀α ∈ F

∀u ∈ V.

If a vector space V over F is equipped with a norm k · kV is called a “normed linear space.” Definition 5.6. Let V be a vector space over F . A seminorm | · |V : V → R is a function such that 1. |u|V ≥ 0

∀v ∈ V ;

2. |u + v|V ≤ kukV + |v|V 3. |αu|V = |α| |u|V

∀u, v ∈ V ;

∀α ∈ F

∀u ∈ V.

A metric space is more general than a normed linear space. Indeed, Definition 5.7. Let M be a set. A metric d is a real–valued function defined on M × M satisfying the following properties: 1. d(x, y) ≥ 0 2. d(x, y) = 0

∀x, y ∈ M ⇐⇒ x = y

3. d(x, y) = d(y, x)

(non–negativity); ∀x, y ∈ M

∀x, y ∈ M

4. d(x, z) ≤ d(x, y) + d(y, z)

(identity of indiscernibles); (symmetry);

∀x, y, z ∈ M

(triangle inequality).

If a metric d is equipped, the ordered pair (M, d) is called a metric space. A normed linear space is a metric space by defining the metric d(x, y) = ckx − ykV for all x, y ∈ V with any positive constant c > 0. Definition 5.8. A sequence (xj )∞ j=1 in a metric space (M, d) is called a Cauchy sequence if for every ǫ > 0 there exists an integer N = N (ǫ) such that d(xj , xk ) < N (ǫ) c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

∀j, k ≥ N. 62

Finite element methods

Section 5.1

Definition 5.9. A metric space (M, d) is said to be complete if every Cauchy sequence in (M, d) converges to a an element in (M, d). Definition 5.10. A complete normed linear space (V, k · kV ) is called a Banach space. Definition 5.11. A complete normed linear space (V, k · kV ) with a real or complex innter–product (·, ·)V is p called a Hilbert space where the norm is induced from the inner–product, i.e. kvk = (v, v) ∀v ∈ V.

A Hilbert space V is a generalization of the Euclidian space Rd or Cd where V can be a fairly abstract set of elements. The Sobolev spaces H k (Ω) = W k,2 (Ω) are Hilbert spaces while W k,p (Ω) for p ≥ 1 are Banach spaces with the corresponding Sobolev norms. Definition 5.12. Let (H, k · kH ) be a Hilbert space and (V, k · kV ) be a subspace which is dense and the inclusion ι : V → H is continuous, i.e if limj→∞ kvj − vkV = 0 then limj→∞ kvj − vkH = 0. Recall the dual map ι∗ : V ′ → H ′ . Now identifying H ′ with H, the following inclusions and and identifications are called a Gelfand triplet: V ֒→ H ⇐⇒ H ′ ֒→ V ′ . The meaning is as follows. If φ ∈ H ⇐⇒ H ′ , the inclusion implies φ ∈ V ′ . Then, for all v ∈ V, one has hφ, viV ′ ,V = hφ, viH ′ ,H = (φ, v)H .

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

63

Finite element methods

c

Dongwoo Sheen, Ph.D. (http://www.nasc.snu.ac.kr)

Section 5.1

64

Bibliography [Gri85]

P. Grisvard. Elliptic Problems in Nonsmooth Domains. Pitman, Boston, London, 1985.

[GT83]

D. Gilbarg and N. S. Trudinger. Elliptic Partial Differential Equations of Second Order. Springer– Verlag, Berlin, Heidelberg, second edition, 1983.

[MS64]

N. Meyers and J. Serrin. H = W . Proc. Nat. Acad. Sci., USA, 51:1055–1056, 1964.

[Rud91a] W. Rudin. Functional analysis. International Series in Pure and Applied Mathematics. McGrawHill, Singapore, 2nd edition, 1991. [Rud91b] W. Rudin. Functional analysis. International Series in Pure and Applied Mathematics. McGrawHill, Singapore, 2nd edition, 1991. [Yos95]

K. Yosida. Functional Analysis. Springer-Verlag, Berlin, 1995.

65