Rational Krylov Methods for Solving Nonlinear Eigenvalue Problems - - PowerPoint PPT Presentation

rational krylov methods for solving nonlinear eigenvalue
SMART_READER_LITE
LIVE PREVIEW

Rational Krylov Methods for Solving Nonlinear Eigenvalue Problems - - PowerPoint PPT Presentation

Rational Krylov Methods for Solving Nonlinear Eigenvalue Problems Roel Van Beeumen rvanbeeumen@lbl.gov Computational Research Division Lawrence Berkeley National Laboratory BASCD 2016 Stanford December 3, 2016 Quadratic eigenvalue


slide-1
SLIDE 1

Rational Krylov Methods for Solving Nonlinear Eigenvalue Problems

Roel Van Beeumen

rvanbeeumen@lbl.gov

Computational Research Division Lawrence Berkeley National Laboratory

BASCD 2016 Stanford – December 3, 2016

slide-2
SLIDE 2

Quadratic eigenvalue problem

Vibration analysis in structural analysis gives rise to (λ2M + λC + K)x = 0 where λ is an eigenvalue x is an eigenvector M is the mass matrix C is the damping matrix K is the stiffness matrix

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 1

slide-3
SLIDE 3

Motivation: Nonlinear damping

Clamped beam:

  • λ2M + λC + K
  • x = 0

|C|

|λ| Clamped sandwich beam:

  • λ2M + λC(λ) + K
  • x = 0

|λ|

|C(λ)| for λ on the imaginary axis

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 2

slide-4
SLIDE 4

Motivation: Active damping

Active damping in cars:

input

  • utput

System Controller

Delay eigenvalue problem

  • λ2M + λC + K + e−λτE
  • x = 0
  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 3

slide-5
SLIDE 5

Motivation: Nonlinear boundary conditions

Cavity design of a linear accelerator ∇ × 1 µ∇ × E

  • − λεE = 0

Maxwell’s equations + nonlinear waveguide boundary conditions:  K − λM +

k

  • j=1

i

  • λ − κ2

c,jWj

  x = 0 where λ = ω2/c2 κc,j are the cutoff values

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 4

slide-6
SLIDE 6

Outline

1 Motivation 2 Solving Nonlinear Eigenvalue Problems

Approximation Linearization pencils Solving linear eigenvalue problem

3 Numerical Experiment

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 5

slide-7
SLIDE 7

Nonlinear eigenvalue problem (NLEP)

NLEP

The nonlinear eigenvalue problem: A(λ)x = 0 where λ ∈ Ω ⊆ C: eigenvalue x ∈ Cn\{0}: eigenvector A : Ω → Cn×n: matrix-valued function

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 6

slide-8
SLIDE 8

Nonlinear eigenvalue problem (NLEP)

NLEP

The nonlinear eigenvalue problem: A(λ)x = 0 where λ ∈ Ω ⊆ C: eigenvalue x ∈ Cn\{0}: eigenvector A : Ω → Cn×n: matrix-valued function Note that the NLEP is nonlinear in eigenvalue λ, linear in eigenvector x.

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 6

slide-9
SLIDE 9

Solving NLEPs

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution

1 approximation via interpolation 2 linearization 3 solving linear eigenvalue problem

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 7

slide-10
SLIDE 10

Approximation

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 1: Polynomial interpolation A(λ) ≈ Pd(λ) = A0 + A1λ + A2λ2 + · · · + Adλd

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 8

slide-11
SLIDE 11

Approximation

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 1: Polynomial interpolation A(λ) ≈ Pd(λ) = A0 + A1λ + A2λ2 + · · · + Adλd

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 8

slide-12
SLIDE 12

Approximation

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 1: Polynomial interpolation A(λ) ≈ Pd(λ) = A0 + A1λ + A2λ2 + · · · + Adλd

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 8

slide-13
SLIDE 13

Approximation

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 1: Polynomial interpolation A(λ) ≈ Pd(λ) = A0 + A1λ + A2λ2 + · · · + Adλd

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 8

slide-14
SLIDE 14

Approximation

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 1: Polynomial interpolation A(λ) ≈ Pd(λ) = A0 + A1λ + A2λ2 + · · · + Adλd

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 8

slide-15
SLIDE 15

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-16
SLIDE 16

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-17
SLIDE 17

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + A2n2(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-18
SLIDE 18

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A3n3(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-19
SLIDE 19

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A4n4(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-20
SLIDE 20

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A5n5(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-21
SLIDE 21

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A6n6(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-22
SLIDE 22

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A7n7(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-23
SLIDE 23

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A8n8(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-24
SLIDE 24

Approximation: Newton interpolation

Dynamic polynomial interpolation (Newton) A(λ) ≈ A0n0(λ) + A1n1(λ) + · · · + A9n9(λ)

1 2 3 4 −0.5 0.5 1

λ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 9

slide-25
SLIDE 25

Approximation: Polynomial versus Rational

Scalar nonlinear function: A(λ) = 0.2 √ λ − 0.6 sin(2λ) = 0 with target set: Σ = [0.01, 4]

0.5 1 1.5 2 2.5 3 3.5 4 0.5

λ A(λ)

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 10

slide-26
SLIDE 26

Approximation: Polynomial versus Rational

Scalar nonlinear function: A(λ) = 0.2 √ λ − 0.6 sin(2λ) = 0 with target set: Σ = [0.01, 4]

0.5 1 1.5 2 2.5 3 3.5 4 0.5

Leja points

λ A(λ)

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 10

slide-27
SLIDE 27

Approximation: Polynomial versus Rational

Scalar nonlinear function: A(λ) = 0.2 √ λ − 0.6 sin(2λ) = 0 with target set: Σ = [0.01, 4]

0.5 1 1.5 2 2.5 3 3.5 4 0.5

Leja points

λ A(λ)

0.5 1 1.5 2 2.5 3 3.5 4 0.5

Leja–Bagby points

λ A(λ)

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 10

slide-28
SLIDE 28

Approximation: Polynomial versus Rational

Scalar nonlinear function: A(λ) = 0.2 √ λ − 0.6 sin(2λ) = 0

20 40 60 80 100 10-15 10-10 10-5 100 number of interpolation nodes

interpolation error

  • pol. Leja
  • rat. Leja–Bagby
  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 11

slide-29
SLIDE 29

Approximation: Polynomial versus Rational

Scalar nonlinear function: A(λ) = 0.2 √ λ − 0.6 sin(2λ) = 0

20 40 60 80 100 10-15 10-10 10-5 100 number of interpolation nodes

interpolation error

  • pol. Leja
  • rat. Leja–Bagby

20 40 60 80 100 10-15 10-10 10-5 100 iteration

convergence of eigenvalues

Newton Rational Krylov Fully Rational Krylov

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 11

slide-30
SLIDE 30

Linearization pencils

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 2: Linearization Pd(λ)x = 0 ⇓ L(λ)x = (A − λB)x = 0

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 12

slide-31
SLIDE 31

Linearization pencils

Linearization: idea Second order ODE M ¨ q + C ˙ q + Kq = 0 System of first order ODEs

  • ˙

q1 = q2 M ˙ q2 = −Cq2 − Kq1

I M ˙ q1 ˙ q2

  • =

I −K −C q1 q2

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 13

slide-32
SLIDE 32

Linearization pencils

Linearization: idea Second order ODE M ¨ q + C ˙ q + Kq = 0 System of first order ODEs

  • ˙

q1 = q2 M ˙ q2 = −Cq2 − Kq1

I M ˙ q1 ˙ q2

  • =

I −K −C q1 q2

  • Quadratic eigenvalue problem
  • Mλ2 + Cλ + K
  • x = 0

Linear eigenvalue problem

λ I M x1 x2

  • =

I −K −C x1 x2

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 13

slide-33
SLIDE 33

Linearization pencils

Step 2: Companion linearization

PEP

Pd(λ)x = 0 ⇒

GEP

L(λ)x = (A − λB)x = 0 Pd(λ)x =

  • A0 + A1λ + A2λ2 + · · · + Adλd

x = 0

       A0 A1 A2 · · · Ad−1 I I ... I       

  • A

       x1 x2 x3 . . . xd       

x

= λ

       · · · −Ad I ... ... I I       

  • B

       x1 x2 x3 . . . xd       

x

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 14

slide-34
SLIDE 34

Linearization pencils: Uniform representation

L(λ) =

  • A0

A1 · · · Ad−1 M ⊗ In

  • A

−λ

  • B0

B1 · · · Bd−1 N ⊗ In

  • B

Monomial basis [Mackey, Mackey, Mehl, Mehrmann 2006] Chebyshev basis [Effenberger, Kressner 2012] Lagrange basis [VB, Michiels, Meerbergen 2015] Newton/Hermite basis [Amiraslani, Corless, Lancaster 2009] Rational monomial basis [Nakatsukasa, Tisseur 2014] Rational Newton basis [G¨

uttel, VB, Meerbergen, Michiels 2014]

Spectral discretization [Jarlebring, Meerbergen, Michiels 2010] . . .

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 15

slide-35
SLIDE 35

Solving linear eigenvalue problem

NLEP

A(λ)x = 0 ⇓

PEP

Pd(λ)x = 0 ⇓

GEP

L(λ)x = 0 ⇓ Solution Step 3: Solve generalized linear eigenvalue problem L(λ)x = (A − λB)x = 0 by the rational Krylov method.

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 16

slide-36
SLIDE 36

Compact Rational Krylov (CORK) framework

L(λ) = A −λ B full exploitation of structure compact Krylov subspace reduction in memory cost reduction in computation cost

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 17

slide-37
SLIDE 37

Compact Rational Krylov (CORK) framework

CORK

Generic class of numerical methods with a lot flexibility: polynomial or rational approximation in different bases interpolation nodes (and poles) shifts in the rational Krylov method dynamic / static / hybrid Compact Rational Krylov variants:

1 Newton Rational Krylov (NRK) method [VB, Meerbergen, Michiels 2013] 2 Fully Rational Krylov (FRK) method [G¨

uttel, VB, Meerbergen, Michiels 2014]

3 Infinite Arnoldi (InfA) method [Jarlebring, Michiels, Meerbergen 2012] 4 . . .

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 18

slide-38
SLIDE 38

Gun problem

Maxwell’s equations + nonlinear waveguide boundary conditions: A(λ)x =

  • K − λM + i
  • λ − α2

1W1 + i

  • λ − α2

2W2

  • x = 0,

where K, M, W1, W2 ∈ R9956×9956 λ = ω2/c2

ω the wavenumber c the speed of light

α1 and α2 are the cutoff values = ‘Gun’ problem from the NLEVP collection

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 19

slide-39
SLIDE 39

Gun problem

50 100 150 200 250 300 350 50 100 150 Re √ λ Im √ λ target set Σ

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 20

slide-40
SLIDE 40

Gun problem

50 100 150 200 250 300 350 50 100 150 Re √ λ Im √ λ target set Σ nodes σ poles ξ

A(λ) approximated by a rational Newton polynomial of degree d = 35

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 20

slide-41
SLIDE 41

Gun problem

50 100 150 200 250 300 350 50 100 150 Re √ λ Im √ λ target set Σ nodes σ poles ξ shifts

A(λ) approximated by a rational Newton polynomial of degree d = 35

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 20

slide-42
SLIDE 42

Gun problem

50 100 150 200 250 300 350 50 100 150 Re √ λ Im √ λ target set Σ nodes σ poles ξ shifts eigenvalues λ

A(λ) approximated by a rational Newton polynomial of degree d = 35

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 20

slide-43
SLIDE 43

Compact Rational Krylov (CORK) framework

20 40 60 80 100 20 40 60 80 100 120 140 iteration memory (MB)

n = 10.000 d = 35 k = 50

FRK CORK

O((k · d) · n) ⇒ O((k + d) · n)

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 21

slide-44
SLIDE 44

Conclusions

How to solve NLEP’s

A(λ)x = 0

1 Approximation via interpolation 2 Linearization 3 Solve linear eigenvalue problem

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 22

slide-45
SLIDE 45

Compact Rational Krylov (CORK) framework

CORK

linearization pencils Krylov methods compact subspace reduction memory cost reduction

  • rthog. cost

system solves

  • riginal dim.
  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 23

slide-46
SLIDE 46

References

[1] VB, Meerbergen, and Michiels. Compact rational Krylov methods for nonlinear eigenvalue problems SIAM Journal on Matrix Analysis and Applications, 36(2), 820–838, 2015 [2] G¨ uttel, VB, Meerbergen, and Michiels. NLEIGS: A class of fully rational Krylov methods for nonlinear eigenvalue problems SIAM Journal on Scientific Computing, 36 (6), A2842–A2864, 2014 [3] VB, Meerbergen, and Michiels. A rational Krylov method based on Hermite interpolation for nonlinear eigenvalue problems SIAM Journal on Scientific Computing, 35 (1), A327-A350, 2013

CORK & NLEIGS Matlab toolboxes available on my homepage:

http://www.roelvanbeeumen.be

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 24

slide-47
SLIDE 47

References

[1] VB, Meerbergen, and Michiels. Compact rational Krylov methods for nonlinear eigenvalue problems SIAM Journal on Matrix Analysis and Applications, 36(2), 820–838, 2015 [2] G¨ uttel, VB, Meerbergen, and Michiels. NLEIGS: A class of fully rational Krylov methods for nonlinear eigenvalue problems SIAM Journal on Scientific Computing, 36 (6), A2842–A2864, 2014 [3] VB, Meerbergen, and Michiels. A rational Krylov method based on Hermite interpolation for nonlinear eigenvalue problems SIAM Journal on Scientific Computing, 35 (1), A327-A350, 2013

CORK & NLEIGS Matlab toolboxes available on my homepage:

http://www.roelvanbeeumen.be

Thank you for your attention !

  • R. Van Beeumen (Berkeley Lab)

Rational Krylov methods for NLEPs Stanford – December 3, 2016 24