Добавил:
Upload Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:
Lecture Notes on Solving Large Scale Eigenvalue Problems.pdf
Скачиваний:
49
Добавлен:
22.03.2016
Размер:
2 Мб
Скачать
☆

2

CHAPTER 1. INTRODUCTION

1.2Example 1: The vibrating string

1.2.1Problem setting

Let us consider a string as displayed in Fig. 1.1. The string is fixed at both ends, at x = 0

u

 

 

u(x,t)

 

x

0

L

Figure 1.1: A vibrating string fixed at both ends.

and x = L. The x-axis coincides with the string’s equilibrium position. The displacement of the rest position at x, 0 < x < L, and time t is denoted by u(x, t).

We will assume that the spatial derivatives of u are not very large:

 

∂u

 

 

 

 

 

 

 

∂x

 

is small.

 

 

 

 

This assumption entails that we may neglect terms of higher order.

Let v(x, t) be the velocity of the string at position x and at time t. Then the kinetic energy of a string section ds of mass dm = ρ ds is given by

 

 

 

 

 

 

 

1

 

 

 

1

 

 

 

∂u

2

(1.1)

 

 

dT =

 

dm v2 =

 

 

ρ ds

 

 

.

2

2

∂t

From Fig. 1.2 we see that ds2 = dx2 + ∂u∂x 2 dx2 and thus

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

ds

= s1 +

∂u

 

2

1

 

 

∂u

2

 

 

 

 

 

 

 

 

= 1 +

 

 

 

 

 

+ higher order terms.

 

dx

∂x

 

2

∂x

Plugging this into (1.1) and omitting also the second order term (leaving just the number 1)

gives

dT = ρ dx ∂u 2 . 2 ∂t

The kinetic energy of the whole string is obtained by integrating over its length,

T =

Z0

dT (x) = 2

Z0

ρ(x)

∂t

dx

 

L

1

L

 

∂u

2

The potential energy of the string has two components

1.2. EXAMPLE 1: THE VIBRATING STRING

3

ds

dx

Figure 1.2: A vibrating string, local picture.

1. the stretching times the exerted strain τ,

 

 

 

 

 

 

τ Z0

ds − τ Z0

dx = τ

Z0

s

 

 

 

 

− 1 dx

1 + ∂x

2

L

L

 

 

L

1

 

 

∂u

 

 

 

 

 

L

∂u

 

2

 

 

 

= τ Z0

 

 

 

 

 

 

+

higher order terms! dx

 

 

2

∂x

2. exterior forces of density f,

 

 

− Z0L fudx.

 

 

 

 

 

 

 

 

 

Summing up, the potential energy of the string becomes

 

 

 

 

L

 

τ

∂u

 

2

 

 

 

 

(1.2)

V = Z0

 

 

 

 

 

 

− fu! dx.

2

 

∂x

To consider the motion (vibration) of the string in a certain time interval t1 ≤ t ≤ t2 we form the integral

(1.3)

I(u) =

Zt1

(T − V ) dt =

2

Zt1

Z0

"ρ(x)

∂t

− τ

∂x

+ 2fu# dx dt

 

 

t2

 

1

t2

 

L

∂u

2

∂u

2

Here functions u(x, t) are admitted that are di erentiable with respect to x and t and satisfy the boundary conditions (BC) that correspond to the fixing,

(1.4)

u(0, t) = u(L, t) = 0,

t1 ≤ t ≤ t2,

as well as given initial conditions and end conditions,

(1.5)

u(x, t1) = u1

(x),

0 < x < L.

u(x, t2) = u2

(x),

 

 

4

CHAPTER 1. INTRODUCTION

According to the principle of Hamilton a mechanical system with kinetic energy T and potential energy V behaves in a time interval t1 ≤ t ≤ t2 for given initial and end positions

such that

Z t2

I = L dt, L = T − V,

t1

is minimized.

Let u(x, t) be such that I(u) ≤ I(w) for all w, that satisfy the initial, end, and boundary conditions. Let w = u + ε v with

v(0, t) = v(L, t) = 0, v(x, t1) = v(x, t2) = 0.

v is called a variation. We now consider I(u + ε v) as a function of ε. Then we have the equivalence

 

 

 

 

 

 

 

 

dI

 

 

 

 

 

 

 

I(u) minimal

 

 

dε (u) = 0 for all admitted v.

 

Plugging u + ε v into eq. (1.3) we obtain

 

 

 

− τ ∂(u∂x

+ 2f(u + ε v)# dx dt

I(u + ε v) = 2 Z

Z "ρ(x)

 

 

∂t

 

 

 

1

 

t2

L

 

∂(u + ε v) 2

 

 

+ ε v) 2

(1.6)

 

 

t1

0

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Z ρ(x) ∂t ∂t − τ ∂x ∂x + 2fv dx dt + O(ε2).

= I(u) + ε Z

 

 

 

 

t2

L

 

 

 

∂u ∂v

 

 

∂u ∂v

 

 

 

 

 

 

t1

0

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Thus,

= Zt1

Z0

ρ ∂t2

 

 

− τ ∂x2 + 2 f v dx dt = 0

 

 

 

∂ε

 

 

 

 

 

∂I

t2

L

 

 

∂2u

 

 

 

∂2u

 

 

for all admissible v. Therefore, the bracketed expression must vanish,

 

 

 

 

 

 

∂2u

 

 

 

 

∂2u

 

 

 

 

 

(1.7)

 

 

 

 

− ρ

 

 

 

+ τ

 

 

 

= 2 f.

 

 

 

 

 

 

∂t2

 

∂x2

 

 

 

This last di erential equation is named Euler-Lagrange equation.

If the force is proportional to the displacement u(x, t) (like, e.g., in Hooke’s law) then we get a di erential equation of the form

−ρ(x)∂2u + ∂ p(x)∂u + q(x)u(x, t) = 0. (1.8) ∂t2 ∂x ∂x

u(0, t) = u(1, t) = 0

which is a special case of the Euler-Lagrange equation (1.7). Here, ρ(x) plays the role of a mass density, p(x) of a locally varying elasticity module. We do not specify initial and end conditions for the moment. Note that there are no external forces present in (1.7).

From physics we know that ρ(x) > 0 and p(x) > 0 for all x. These properties are of importance also from a mathematical view point! For simplicity, we assume that ρ(x) = 1.

X∞
akwk(x) = u0(x),
k=0
X∞ p
λk bkwk(x) = u1(x),
k=0

1.2. EXAMPLE 1: THE VIBRATING STRING

5

1.2.2The method of separation of variables

For the solution u in (1.8) we make the ansatz

(1.9) u(x, t) = v(t)w(x).

Here, v is a function that depends only on the time t, while w depends only on the spacial variable x. With this ansatz (1.8) becomes

(1.10) v′′(t)w(x) − v(t)(p(x)w′(x))′ − q(x)v(t)w(x) = 0.

Now we separate the variables depending on t from those depending on x,

v′′(t)

=

 

1

(p(x)w′(x))′ + q(x).

v(t)

w(x)

 

 

This equation holds for any t and x. We can vary t and x independently of each other without changing the value on each side of the equation. Therefore, each side of the equation must be equal to a constant value. We denote this value by −λ. Thus, from the left side we obtain the equation

(1.11)

− v′′(t) = λv(t).

√

 

√

 

 

 

This equation has the well-known solution v(t) = a · cos(

λt) + b · sin( λt) where λ > 0

is assumed. The right side of (1.10) gives a so-called Sturm-Liouville problem

 

 

 

 

(1.12)

−(p(x)w′(x))′ + q(x)w(x) = λw(x),

w(0) = w(1) = 0.

 

A value λ for which (1.12) has a non-trivial (i.e. nonzero) solution w is called an eigenvalue; w is a corresponding eigenfunction. It is known that all eigenvalues of (1.12) are

positive. By means of our ansatz (1.9) we get

h √ √ i u(x, t) = w(x) a · cos( λt) + b · sin( λt)

as a solution of (1.8). It is known that (1.12) has infinitely many real positive eigenvalues

0 < λ1 ≤ λ2

≤ · · · , (λkk→∞∞). (1.12) has a non-zero solution, say wk(x), only for these

 

−→

 

 

 

 

 

 

particular values λk. Therefore, the general solution of (1.8) has the form

(1.13)

u(x, t) = k=0 wk(x) hak · cos(

 

λk t) + bk · sin(

λk t)i .

 

∞

 

 

 

 

 

 

 

X

p

 

 

p

 

 

The coe cients ak and bk require that

are determined by initial and end conditions. We could, e.g.,

u(x, 0) =

∂u∂t (x, 0) =

where u0 and u1 are given functions. It is known that the wk form an orthogonal basis in the space of square integrable functions L2(0, 1). Therefore, it is not di cult to compute the coe cients ak and bk.

In concluding, we see that the di cult problem to solve is the eigenvalue problem (1.12). Knowing the eigenvalues and eigenfunctions the general solution of the time-dependent problem (1.8) is easy to form.

Eq. (1.12) can be solved analytically only in very special situation, e.g., if all coe cients are constants. In general a numerical method is needed to solve the Sturm-Liouville problem (1.12).

6

CHAPTER 1. INTRODUCTION

1.3Numerical methods for solving 1-dimensional problems

In this section we consider three methods to solve the Sturm-Liouville problem.

1.3.1Finite di erences

We approximate w(x) by its values at the discrete points xi = ih, h = 1/(n + 1), i = 1, . . . , n.

 

 

 

 

 

 

 

x

 

 

 

 

 

 

 

0

xi−1 xi

x

i+1

L

 

Figure 1.3: Grid points in the interval (0, L).

At point xi we approximate the derivatives by finite di erences. We proceed as

follows. First we write

 

 

 

 

 

 

 

 

g(xi+ 21 ) − g(xi− 12 )

 

 

 

 

 

 

 

 

 

 

 

 

 

 

d

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

g(xi) ≈

 

 

 

 

 

 

.

 

 

 

 

 

 

 

 

 

 

 

 

 

dx

 

 

 

 

h

 

 

 

 

 

 

For g = pdw

we get

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

dx

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

g(x

i+

1 ) = p(x

i+

1 )

w(xi+1) − w(xi)

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2

 

 

2

 

h

 

 

 

 

 

 

 

and finally, for i = 1, . . . , n,

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

d

dw

1

 

 

 

 

 

w(xi+1)

w(xi)

 

 

w(xi)

w(xi

 

1)

 

−

 

p

 

(xi)

≈ −

 

p(xi+ 21 )

 

 

 

 

−

 

− p(xi− 21 )

 

−

−

 

 

dx

dx

h

 

 

 

 

h

 

 

h

= h2 hp(xi−

21 )wi−1 + (p(xi− 12 ) + p(xi+ 21 ))wi − p(xi+ 21 )wi+1i .

1

 

 

Note that at the interval endpoints w0 = wn+1 = 0.

We can collect all equations in a matrix equation,

 

p(x1 ) + p(x3 )

2

2

 

+ q(x1)

 

h2

 

 

3

 

 

 

 

p(x2 )

 

 

2

 

 

− h

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

or, briefly,

 

p(x3 )

− h2

 

2

 

 

p(x3 ) + p(x5 )

2

2

+ q(x2)

h2

 

p(x5 )

−

2

 

 

h2

 

p(x5 )

− h22

...

 

 

w2

 

 

w2

 

 

 

 

 

 

w1

 

 

 

w1

 

 

 

 

 

.

 

.

 

 

 

 

.

 

 

 

.

 

 

..

 

 

.

 

.

 

 

 

 

w3

 

= λ

 

w3

 

,

 

.

 

 

 

 

 

 

 

wn

 

wn

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(1.14)

Aw = λw.

By construction, A is symmetric and tridiagonal. One can show that it is positive definite as well. Note that this matrix has just a few nonzeros: out of the n2 elements of A only 3n − 2 are nonzero. This is a first example of a sparse matrix.

1.3.2The finite element method

We write (1.12) in the form

1.3. NUMERICAL METHODS FOR SOLVING 1-DIMENSIONAL PROBLEMS

7

Find a twice di erentiable function w with w(0) = w(1) = 0 such that

Z 1

−(p(x)w′(x))′ + q(x)w(x) − λw(x) φ(x)dx = 0

0

for all smooth functions φ that satisfy φ(0) = φ(1) = 0.

To relax the requirements on w we integrate by parts and get the new so-called weak form of the problem:

Find a di erentiable function w with w(0) = w(1) = 0 such that

Z 1

(1.15) p(x)w(x)′φ′(x) + q(x)w(x)φ(x) − λw(x)φ(x) dx = 0

0

for all di erentiable functions φ that satisfy φ(0) = φ(1) = 0.

Remark: Requiring di erentiability is too strong and does not lead to a mathematically suitable formulation. In particular, the test functions that will be used below are not di erentiable in the classical sense. It is more appropriate to require w and φ to be weakly di erentiable. In terms of Sobolev spaces: w, φ H01([0, 1]). An introduction to Sobolev spaces is, however, beyond the scope of these notes.

 

Ψi

 

 

 

 

1

 

 

 

 

 

 

 

 

 

 

 

 

 

x

 

 

 

 

 

 

 

 

0

 

xi−1 xi

x

i+1

L

 

Figure 1.4: A basis function of the finite element space: a hat function. We now write w as the linear combination

 

 

 

 

 

 

n

 

 

 

 

 

(1.16)

 

 

 

 

w(x) =

ξi Ψi(x),

 

 

 

 

 

 

 

 

 

=1

 

 

 

 

 

 

 

 

 

 

Xi

 

 

 

 

where

 

 

 

 

 

| +

 

 

 

 

 

(1.17)

Ψ

(x) =

1

− |

x − xi

= max

0, 1

−

|x − xi|

,

 

i

 

 

h

{

 

h

}

is the function that is linear in each interval (xi, xi+1) and satisfies

 

 

 

 

 

 

 

 

1,

i = k,

 

 

 

 

 

Ψi(xk) = δik := 0,

i 6= k.

 

 

An example of such a basis function, a so-called hat function, is given in Fig. 1.4.

We now replace w in (1.15) by the linear combination (1.16), and replace testing ‘against all φ’ by testing against all Ψj . In this way (1.15) becomes

0

1

−p(x)(

n

ξi Ψi′ (x))Ψj′ (x) + (q(x) − λ)

n

ξi Ψi(x)Ψj (x)! dx, for all j,

Z

 

 

X

 

X

 

 

 

i=1

 

i=1

 

Соседние файлы в предмете [НЕСОРТИРОВАННОЕ]