- •1 Introduction
- •1.1 What makes eigenvalues interesting?
- •1.2 Example 1: The vibrating string
- •1.2.1 Problem setting
- •1.2.2 The method of separation of variables
- •1.3.3 Global functions
- •1.3.4 A numerical comparison
- •1.4 Example 2: The heat equation
- •1.5 Example 3: The wave equation
- •1.6 The 2D Laplace eigenvalue problem
- •1.6.3 A numerical example
- •1.7 Cavity resonances in particle accelerators
- •1.8 Spectral clustering
- •1.8.1 The graph Laplacian
- •1.8.2 Spectral clustering
- •1.8.3 Normalized graph Laplacians
- •1.9 Other sources of eigenvalue problems
- •Bibliography
- •2 Basics
- •2.1 Notation
- •2.2 Statement of the problem
- •2.3 Similarity transformations
- •2.4 Schur decomposition
- •2.5 The real Schur decomposition
- •2.6 Normal matrices
- •2.7 Hermitian matrices
- •2.8 Cholesky factorization
- •2.9 The singular value decomposition (SVD)
- •2.10 Projections
- •2.11 Angles between vectors and subspaces
- •Bibliography
- •3 The QR Algorithm
- •3.1 The basic QR algorithm
- •3.1.1 Numerical experiments
- •3.2 The Hessenberg QR algorithm
- •3.2.1 A numerical experiment
- •3.2.2 Complexity
- •3.3 The Householder reduction to Hessenberg form
- •3.3.2 Reduction to Hessenberg form
- •3.4 Improving the convergence of the QR algorithm
- •3.4.1 A numerical example
- •3.4.2 QR algorithm with shifts
- •3.4.3 A numerical example
- •3.5 The double shift QR algorithm
- •3.5.1 A numerical example
- •3.5.2 The complexity
- •3.6 The symmetric tridiagonal QR algorithm
- •3.6.1 Reduction to tridiagonal form
- •3.6.2 The tridiagonal QR algorithm
- •3.7 Research
- •3.8 Summary
- •Bibliography
- •4.1 The divide and conquer idea
- •4.2 Partitioning the tridiagonal matrix
- •4.3 Solving the small systems
- •4.4 Deflation
- •4.4.1 Numerical examples
- •4.6 Solving the secular equation
- •4.7 A first algorithm
- •4.7.1 A numerical example
- •4.8 The algorithm of Gu and Eisenstat
- •4.8.1 A numerical example [continued]
- •Bibliography
- •5 LAPACK and the BLAS
- •5.1 LAPACK
- •5.2 BLAS
- •5.2.1 Typical performance numbers for the BLAS
- •5.3 Blocking
- •5.4 LAPACK solvers for the symmetric eigenproblems
- •5.6 An example of a LAPACK routines
- •Bibliography
- •6 Vector iteration (power method)
- •6.1 Simple vector iteration
- •6.2 Convergence analysis
- •6.3 A numerical example
- •6.4 The symmetric case
- •6.5 Inverse vector iteration
- •6.6 The generalized eigenvalue problem
- •6.7 Computing higher eigenvalues
- •6.8 Rayleigh quotient iteration
- •6.8.1 A numerical example
- •Bibliography
- •7 Simultaneous vector or subspace iterations
- •7.1 Basic subspace iteration
- •7.2 Convergence of basic subspace iteration
- •7.3 Accelerating subspace iteration
- •7.4 Relation between subspace iteration and QR algorithm
- •7.5 Addendum
- •Bibliography
- •8 Krylov subspaces
- •8.1 Introduction
- •8.3 Polynomial representation of Krylov subspaces
- •8.4 Error bounds of Saad
- •Bibliography
- •9 Arnoldi and Lanczos algorithms
- •9.2 Arnoldi algorithm with explicit restarts
- •9.3 The Lanczos basis
- •9.4 The Lanczos process as an iterative method
- •9.6 Partial reorthogonalization
- •9.5 An error analysis of the unmodified Lanczos algorithm
- •9.7 Block Lanczos
- •9.8 External selective reorthogonalization
- •Bibliography
- •10 Restarting Arnoldi and Lanczos algorithms
- •10.2 Implicit restart
- •10.3 Convergence criterion
- •10.4 The generalized eigenvalue problem
- •10.5 A numerical example
- •10.6 Another numerical example
- •10.7 The Lanczos algorithm with thick restarts
- •10.8 Krylov–Schur algorithm
- •10.9 The rational Krylov space method
- •Bibliography
- •11 The Jacobi-Davidson Method
- •11.1 The Davidson algorithm
- •11.2 The Jacobi orthogonal component correction
- •11.2.1 Restarts
- •11.2.2 The computation of several eigenvalues
- •11.2.3 Spectral shifts
- •11.3 The generalized Hermitian eigenvalue problem
- •11.4 A numerical example
- •11.6 Harmonic Ritz values and vectors
- •11.7 Refined Ritz vectors
- •11.8 The generalized Schur decomposition
- •11.9.1 Restart
- •11.9.3 Algorithm
- •Bibliography
- •12 Rayleigh quotient and trace minimization
- •12.1 Introduction
- •12.2 The method of steepest descent
- •12.3 The conjugate gradient algorithm
- •12.4 Locally optimal PCG (LOPCG)
- •12.5 The block Rayleigh quotient minimization algorithm (BRQMIN)
- •12.7 A numerical example
- •12.8 Trace minimization
- •Bibliography
1.4. EXAMPLE 2: THE HEAT EQUATION |
9 |
1.3.4A numerical comparison
We consider the above 1-dimensional eigenvalue problem
(1.20) |
− |
((1 + x)w′(x))′ + w(x) = λw(x), w(0) = w(1) = 0, |
|
|
and solve it with the finite di erence and finite element methods as well as with the global functions method. The results are given in Table 1.1.
Clearly the global function method is the most powerful of them all. With 80 basis functions the eigenvalues all come right. The convergence rate is exponential.
With the finite di erence and finite element methods the eigenvalues exhibit quadratic convergence rates. If the mesh width h is reduced by a factor of q = 2, the error in the eigenvalues is reduced by the factor q2 = 4.
1.4Example 2: The heat equation
The instationary temperature distribution u(x, t) in an insulated container satisfies the equations
∂u(x, t) |
− u(x, t) = 0, |
x Ω, t > 0, |
|||
|
|
||||
|
∂t |
||||
(1.21) |
|
∂u(x, t) |
= 0, |
x ∂Ω, t > 0, |
|
|
|
|
∂n |
||
u(x, 0) = u0(x), x Ω.
Here Ω is a 3-dimensional domain2 with boundary ∂Ω. u0(x), x = (x1, x2, x3)T R3, is a given bounded, su ciently smooth function.
(1.22) u =
X ∂2u
∂xi2
is called the Laplace operator and ∂n∂u denotes the derivative of u in direction of the outer normal vector n. To solve the heat equation the method of separation of variables is employed. We write u in the form
(1.23) |
|
u(x, t) = v(t)w(x). |
|
|||
If a constant λ can be found such that |
|
|||||
|
|
|
|
|
|
|
(1.24) |
w(x) + λw(x) = 0, w(x) 6= 0, x in Ω, |
|||||
|
∂w(x, t) |
|
|
|
||
|
|
= 0, |
x on ∂Ω, |
|||
|
|
|
||||
|
|
∂n |
|
|||
then the product u = vw is a solution of (1.21) if and only if |
||||||
(1.25) |
|
|
dv(t) |
+ λv(t) = 0, |
|
|
|
|
|
||||
|
|
|
dt |
|
||
the solution of which has the form a·exp(−λt). By separating variables, the problem (1.21) is divided in two subproblems that are hopefully easier to solve. A value λ, for which (1.24)
2In the sequel we understand a domain to be bounded and simply connected.
10 |
CHAPTER 1. INTRODUCTION |
Finite di erence method
k |
λk(n = 10) |
|
λk(n = 20) |
λk(n = 40) |
λk(n = 80) |
1 |
15.245 |
|
15.312 |
15.331 |
15.336 |
2 |
56.918 |
|
58.048 |
58.367 |
58.451 |
3 |
122.489 |
|
128.181 |
129.804 |
130.236 |
4 |
206.419 |
|
224.091 |
229.211 |
230.580 |
5 |
301.499 |
|
343.555 |
355.986 |
359.327 |
6 |
399.367 |
|
483.791 |
509.358 |
516.276 |
7 |
492.026 |
|
641.501 |
688.398 |
701.185 |
8 |
578.707 |
|
812.933 |
892.016 |
913.767 |
9 |
672.960 |
|
993.925 |
1118.969 |
1153.691 |
10 |
794.370 |
|
1179.947 |
1367.869 |
1420.585 |
|
|
|
|
|
|
|
|
|
|
||
|
|
Finite element method |
|
||
|
|
|
|
|
|
|
|
|
|
|
|
k |
λk(n = 10) |
|
λk(n = 20) |
λk(n = 40) |
λk(n = 80) |
1 |
15.447 |
|
15.367 |
15.345 |
15.340 |
2 |
60.140 |
|
58.932 |
58.599 |
58.511 |
3 |
138.788 |
|
132.657 |
130.979 |
130.537 |
4 |
257.814 |
|
238.236 |
232.923 |
231.531 |
5 |
426.223 |
|
378.080 |
365.047 |
361.648 |
6 |
654.377 |
|
555.340 |
528.148 |
521.091 |
7 |
949.544 |
|
773.918 |
723.207 |
710.105 |
8 |
1305.720 |
|
1038.433 |
951.392 |
928.983 |
9 |
1702.024 |
|
1354.106 |
1214.066 |
1178.064 |
10 |
2180.159 |
|
1726.473 |
1512.784 |
1457.733 |
|
|
|
|
||
|
|
Global function method |
|
||
|
|
|
|
|
|
|
|
|
|
|
|
k |
λk(n = 10) |
|
λk(n = 20) |
λk(n = 40) |
λk(n = 80) |
1 |
15.338 |
|
15.338 |
15.338 |
15.338 |
2 |
58.482 |
|
58.480 |
58.480 |
58.480 |
3 |
130.389 |
|
130.386 |
130.386 |
130.386 |
4 |
231.065 |
|
231.054 |
231.053 |
231.053 |
5 |
360.511 |
|
360.484 |
360.483 |
360.483 |
6 |
518.804 |
|
518.676 |
518.674 |
518.674 |
7 |
706.134 |
|
705.631 |
705.628 |
705.628 |
8 |
924.960 |
|
921.351 |
921.344 |
921.344 |
9 |
1186.674 |
|
1165.832 |
1165.823 |
1165.822 |
10 |
1577.340 |
|
1439.083 |
1439.063 |
1439.063 |
|
|
|
|
|
|
|
|
|
|
|
|
Table 1.1: Numerical solutions of problem (1.20)
1.5. EXAMPLE 3: THE WAVE EQUATION |
11 |
has a nontrivial (i.e. a nonzero) solution is called an eigenvalue; w then is called a corresponding eigenfunction.
If λn is an eigenvalue of problem (1.24) with corresponding eigenfunction wn, then
e−λntwn(x)
is a solution of the first two equations in (1.21). It is known that equation (1.24) has infinitely many real eigenvalues 0 ≤ λ1 ≤ λ2 ≤ · · · , (λnt→∞−→ ∞). Multiple eigenvalues are counted according to their multiplicity. An arbitrary bounded piecewise continuous function can be represented as a linear combination of the eigenfunctions w1, w2, . . .. Therefore, the solution of (1.21) can be written in the form
|
∞ |
|
(1.26) |
u(x, t) = |
cne−λntwn(x), |
|
n=1 |
|
|
X |
|
where the coe cients cn are determined such that |
||
|
|
∞ |
(1.27) |
u0(x) = |
cnwn(x). |
|
|
n=1 |
|
|
X |
The smallest eigenvalue of (1.24) is λ1 = 0 with w1 = 1 and λ2 > 0. Therefore we see from (1.26) that
(1.28) |
u(x, t) |
−→ |
c1. |
t→∞
Thus, in the limit (i.e., as t goes to infinity), the temperature will be constant in the whole container. The convergence rate towards this equilibrium is determined by the smallest positive eigenvalue λ2 of (1.24):
∞ |
∞ |
|
X |
X |
|
ku(x, t) − c1k = k |
cne−λntwn(x)k ≤ |e−λnt|kcnwn(x)k |
|
n=2 |
n=2 |
∞ |
|
|
|
|
≤ e−λ2t |
X |
|
kcnwn(x)k ≤ e−λ2tku0(x)k. |
|
n=2
Here we have assumed that the value of the constant function w1(x) is set to unity.
1.5Example 3: The wave equation
The air pressure u(x, t) in a volume with acoustically “hard” walls satisfies the equations
(1.29) |
∂2u(x, t) |
− u(x, t) = 0, |
|
x Ω, t > 0, |
|||
∂t2 |
|
||||||
(1.30) |
|
|
∂u(x, t) |
= 0, |
|
x ∂Ω, t > 0, |
|
|
|
|
|
|
|||
|
|
|
∂n |
|
|||
(1.31) |
|
|
|
u(x, 0) = u0 |
(x), |
x Ω, |
|
(1.32) |
|
∂u(x, 0) |
= u1 |
(x), |
x Ω. |
||
|
|
|
|
||||
|
|
|
∂t |
||||
Sound propagates with speed − u, along the (negative) gradient from high to low pressure.
12 |
CHAPTER 1. INTRODUCTION |
To solve the wave equation we proceed as with the heat equation in section 1.4: separation of u according to (1.23) leads again to equation (1.24) but now together with
(1.33) |
d2v(t) |
+ λv(t) = 0. |
|
|
dt2 |
||
We know this equation from the analysis of the vibrating string, see (1.11). From there we know that the general solution of the wave equation has the form
(1.13) |
u(x, t) = k=0 wk(x) hak · cos( |
λk t) + bk · sin( |
|
λk t)i . |
|||
|
∞ |
|
|
|
|
|
|
|
X |
p |
|
|
p |
|
|
where the wk, k = 1, 2, . . ., are the eigenfunctions of the eigenvalue problem (1.24). The coe cients ak and bk are determined by (1.31) and (1.32).
If a harmonic oscillation is forced on the system, an inhomogeneous problem
(1.34) |
∂2u(x, t) |
− u(x, t) = f(x, t), |
∂t2 |
is obtained. The boundary and initial conditions are taken from (1.29)–(1.32). This problem can be solved by setting
|
|
|
|
|
∞ |
v˜n(t)wn(x), |
||||
|
u(x, t) := |
|||||||||
(1.35) |
|
|
|
|
n=1 |
|
|
|
|
|
|
|
|
|
X |
|
|
|
|
|
|
|
|
|
|
|
∞ |
|
|
|
|
|
|
|
|
|
|
X |
φn(t)wn(x). |
||||
|
f(x, t) := |
|||||||||
|
|
|
|
|
n=1 |
|
|
|
|
|
With this approach, v˜n has to satisfy equation |
||||||||||
|
|
|
d2v˜n |
|
|
|
|
|
||
(1.36) |
|
|
|
|
+ λnv˜n = φn(t). |
|||||
|
|
dt2 |
||||||||
If φn(t) = a sin ωt, then the solution becomes |
|
|
|
|
|
|||||
|
|
|
|
|
|
|
|
1 |
|
|
|
|
|
|
|
|
|
||||
(1.37) |
v˜n = An cos pλnt + Bn sin pλnt + |
|
a sin ωt. |
|||||||
λn − ω2 |
||||||||||
An and Bn are real constants that are determined by the initial conditions. If ω gets close
√ √
to λ1, then the last term can be very large. In the limit, if ω = λn, v˜n gets the form
|
p |
|
p |
|
|
(1.38) |
v˜n = An cos |
λnt + Bn sin λnt + at sin ωt. |
|||
In this case, v˜n is not bounded in time anymore. This phenomenon is called resonance. Often resonance is not desirable; it may, e.g., mean the blow up of some structure. In order to prevent resonances eigenvalues have to be known. Possible remedies are changing the domain (the structure).
Remark 1.1. Vibrating membranes satisfy the wave equation, too. In general the boundary conditions are di erent from (1.30). If the membrane (of a drum) is fixed at its boundary, the condition
(1.39) |
u(x, t) = 0 |
is imposed. This boundary conditions is called Dirichlet boundary conditions. The boundary conditions in (1.21) and (1.30) are called Neumann boundary conditions. Combinations of these two can occur. 
