Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

Miniproject – Numerical Linear Algebra – 3NMNLA(4)

Background

I Consider the following Neumann-Laplacian eigenvalue problem:

NLEVP :

where

∂2 v ∂2 v

∂x2 ∂y2 ,

n ≥ gradv(x, y) = (nl, n2 ) ≥ ╱ , 、 = nl + n2 ,

with n being the unit outward normal vector to the boundary Γ of the bounded open domain Ω . The boundary condition is known as a natural or Neumann boundary condition, hence the terminology and acronym used for the eigenvalue problem.  The eigenvalues are also sometimes denoted by µQ , to indicate dependence on the

domain Ω .

It is easy to see that for any domain Ω, u = 0 is an eigenvalue with corresponding eigenfunction v(x, y) = 1. It can be shown that the eigenvalues are real, non-negative, with finite multiplicity, while the corresponding eigenfunctions can be taken to be real.  They form a countable positive sequence (uk)keN, usually labeled in increasing order:

0 = ul  < u2 k u3 k ≥ ≥ ≥ k uk k ≥ ≥ ≥ .

In several applications, the eigenvalue of interest is the first non-zero one, i.e., u2 .  Our aim is to approximate this value using a single vector iteration.

II The finite element method is used to discretise LEVP as the generalised eigenvalue problem

GEVP :    Av = λMv,

where A, M è Rn×n  are sparse symmetric matrices, with A positive semidefinite and M positive-definite.  In particular, it is found that GEVP has an eigenpair (0 , 1), which is the discrete counterpart of the eigenpair (ul, v(x, y))  =  (0, 1) highlighted above.   It can be shown that the eigenvalues of GEVP are real with real corresponding eigenvectors (see Q23, Examples 4). We will use a similar labeling for the eigenvalues, with the smallest being the first:

0 = λl  < λ2 k ≥ ≥ ≥ k λn .

We wish to devise an efficient and reliable procedure for computing approximations to u2  for any domain Ω, by computing λ2  using a sequence of eigenvalue problems GEVP of increasing size.  One possible approach is to modify GEVP so that λ2 becomes the smallest eigenvalue (with λl = 0 removed as an eigenvalue). This ensures that the Inverse Iteration method is defined and can therefore be used for estimating λ2  and therefore u2 .

III Since the Inverse Iteration method requires the solution of a large sparse linear system at each step, we will

also aim to incorporate a preconditioned iterative method in our procedure.   The modified GEVP referred to above yields a symmetric indefinite matrix.  Therefore, we need to use a dedicated iterative method, such as the Minimum Residual Method (MINRES) or the Symmetric LQ method (SYMMLQ). These methods are implemented in Matlab – type help minres,  help  symmlq, for details on usage.  Both have the requirement that the preconditioner is a symmetric and positive-definite matrix.

1.  (a) Consider the rectangular domain Ω = (0, a) ! (0, b), for some a, b > 0.

(i) Check that NLEVP has eigensolutions (eigenpairs) (rk,i, wk,i) given by

π 2 k2 π 2l2 kπx lπy

rk,i  =   a2     +   b2    ,    wk,i(x, y) = cos   a   cos   b  ,        k, l è N ↓ ↓0{ .

(ii) Deduce that u2 = π 2 min ,a-2 , b-2 、.

(iii) Show that u2  has multiplicity one if any only if a b.

(b) Consider the following modified generalised eigenvalue problem

GEVP’ : ┌ATe   e┐ u = g ┌M α┐ u,

where e è Rn / ↓0{ , α > 0 are quantities to be chosen suitably.

Let (λi , vi) denote the eigenpairs of GEVP.

(i) Let e satisfy eT vi = 0 for i = 2, 3, . . . , n. Show that for this choice of e and for any α , ╱λi , ┌0(v)i ┐、

are eigenpairs of GEVP’ for i = 2, 3, . . . , n.

(ii) Let e satisfy eT vl = β 0. Show that the remaining two eigenpairs for GEVP’ have eigenvectors

of the form u = ┌c(v)l ┐ for values of c which you should identify.  Deduce that the corresponding

eigenvalues are real and non-zero for any α > 0.

(iii) By considering the orthogonality relations in Q23, Examples 4, suggest an expression for e that

fulfils the requirements in parts (i) and (ii). Given this choice of e, write down a value of α that ensures the smallest eigenvalue of GEVP’ in absolute value is λ2 .

(c) Design a sparse, symmetric and positive definite preconditioner P for the matrix on the left of GEVP’; P should be an approximation to this matrix:  this may be obtained via splittings, factorisations, or generally modifications in the structure of the matrix and/or its factors.

2. We wish to use the symmetric version of the Inverse Iteration method coupled with suitable preconditioning for estimating u2 . We apply the method to GEVP’ with α, e chosen as in Q1(b)(iii).

(a) Modify the file spowergensym .m to provide an implementation of the symmetric Inverse Iteration method

for the smmetric generalised eigenvalue problem.  At every step, the linear system arising should be solved using a suitable preconditioned iterative method with the preconditioner introduced in Q1(c). The input, output lists and stopping criterion should be left unchanged. The output eigenvector should satisfy the orthogonality relations in Q23, Examples 4. You should name your file ipowergen .m.

(b) The matrices provided in the MAT file (download the file and type load  GEVP) are pairs (Ai, Mi), i =

1, 2, 3, 4  corresponding  to  uniform  discretisations  of NLEVP  for  the  case Ω  =  (0 , 2) ! (0, 1).   Use ipowergen to investigate the performance of the Inverse Iteration method equipped with your pre- conditioner. Write down the number of iterations k required by the Inverse Iteration method to satisfy the stopping criterion. Write down also the average number of preconditioned iterations required per step of the Inverse Iteration method.  Note: you should use a suitable stopping criterion with your iterative solver (which may have to be related to the stopping criterion for the Inverse Iteration method).

(c) Comment on the results from part (b). In particular, discuss the dependence on the problem size of the performance of (i) the preconditioner employed and (ii) the Inverse Iteration method.

(d) Compute the errors ei = │u2 { λi} │, where λi}  is the approximation resulting from employing the pair of matrices (Ai, Mi). Assuming {λi} } are the first terms in a sequence of approximations, estimate the rate of convergence to u2 .