Transcript
Nonlinear Eigenvalue Problems: Theory and Applications David Bindel1 Department of Computer Science Cornell University
3 March 2016
1
Joint work with Amanda Hood Arizona
1 / 41
The Computational Science & Engineering Picture
Application
Analysis
Computation
MEMS
Linear algebra
HPC / cloud
Smart grids
Approximation theory
Simulators
Networks
Symmetry + structure
Little languages
Arizona
2 / 41
Why eigenvalues?
Arizona
3 / 41
Why eigenvalues?
The polynomial connection The optimization connection The approximation connection The Fourier/quadrature/special function connection The dynamics connection Arizona
4 / 41
Vibrations are everywhere, and so too are the eigenvalues associated with them. As mathematical models invade more and more disciplines, we can anticipate a demand for eigenvalue calculations in an ever richer variety of contexts. – Beresford Parlett, The Symmetric Eigenvalue Problem Arizona
5 / 41
Why eigenvalues? As a student in a first ODE course: y=eλt v
y 0 − Ay = 0 −−−−−−−−→ (λI − A)v = 0. Me: “How do I compute this?” p(λ) = det(λI − A) = 0.
And then I learned better...
Arizona
6 / 41
Eigenvalues or polynomial roots? The cosy relationship that mathematicians enjoyed with polynomials suffered a severe setback in the early fifties when electronic computers came into general use. Speaking for myself I regard it as the most traumatic experience in my career as a numerical analyst. — J. H. Wilkinson
Representation matters. Arizona
7 / 41
Why nonlinear eigenvalues? y=eλt v
y 0 − Ay = 0 −−−−−−−−→ (λI − A)v = 0 y=eλt v
y 00 + By 0 + Ky = 0 −−−−−−−−→ (λ2 I + λB + K)v = 0 y=eλt v
y 0 − Ay − By(t − 1) = 0 −−−−−−−−→ (λI − A − Be−λ )v = 0 y=eλt v
T (d/dt)y = 0 −−−−−−−−→ T (λ)v = 0 Higher-order ODEs Dynamic element formulations Delay differential equations Boundary integral equation eigenproblems Radiation boundary conditions Arizona
8 / 41
Big and little
v=y 0
y 00 + By 0 + Ky = 0 −−−−−→
0 v B K v + =0 y −I 0 y
B K v (λ I + λB + K)u = 0 −−−−−→ λI + =0 −I 0 u 2
v=λu
Tradeoff: more variables = more linear.
Arizona
9 / 41
My motivation
T (ω)v ≡ K − ω 2 M + G(ω) v = 0
Arizona
10 / 41
My motivation
T (ω)v ≡ K − ω 2 M + G(ω) v = 0
Wanted: Perturbation theory justifying a terrible estimate of G(ω)
Arizona
11 / 41
The big picture
T (λ)v = 0,
v 6= 0.
where T : Ω → Cn×n analytic, Ω ⊂ C simply connected Regularity: det(T ) 6≡ 0
Nonlinear spectrum: Λ(T ) = {z ∈ Ω : T (z) singular}. What do we want? Qualitative information (e.g. no eigenvalues in right half plane) Error bounds on computed/estimated eigenvalues Assurances that we know all the eigenvalues in some region
Arizona
12 / 41
Solver strategies
T (λ)v = 0,
v 6= 0.
A common approach: 1
Approximate T locally (linear, rational, etc)
2
Solve approximate problem.
3
Repeat as needed.
How should we choose solver parameters? What about global behavior? What can we trust? What might we miss?
Arizona
13 / 41
Analyticity to the rescue
T (λ)v = 0,
v 6= 0.
where T : Ω → Cn×n analytic, Ω ⊂ C simply connected Regularity: det(T ) 6≡ 0
Nonlinear spectrum: Λ(T ) = {z ∈ Ω : T (z) singular}. Goal: Use analyticity to compare and to count
Arizona
14 / 41
Winding and Cauchy’s argument principle
Γ
Ω Z 0 1 f (z) WΓ (f ) = dz 2πi Γ f (z) = # zeros − # poles Arizona
15 / 41
Winding, Rouché, and Gohberg-Sigal |f | |f + g| Γ Ω Analytic Winding # Theorem
sΓ f, g : Ω → C R f 0 (z) 1 2πi Γ f (z) dz
T, E : Ω → Cn×n R 1 −1 0 tr 2πi Γ T (z) T (z) dz
Rouché (1862): |g| < |f | on Γ =⇒ same # zeros of f, f + g
Gohberg-Sigal (1971): kT −1 Ek < 1 on Γ =⇒ same # eigs of T, T + E
Arizona
16 / 41
Comparing NEPs
Suppose T, E : Ω → Cn×n analytic
Γ ⊂ Ω a simple closed contour
T (z) + sE(z) nonsingular ∀s ∈ [0, 1], z ∈ Γ Then T and T + E have the same number of eigenvalues inside Γ. Pf: Constant winding number around Γ.
Arizona
17 / 41
Nonlinear pseudospectra Λ (T ) ≡ {z ∈ Ω : kT (z)−1 k > −1 } 20
15
10
5
0
−5
−10
−15
−20 −20
−15
−10
−5
0
Arizona
5
10
15
20
18 / 41
Pseudospectral comparison
U
Ω Ω E analytic, kE(z)k < on Ω . Then Λ(T + E) ∩ Ω ⊂ Λ (T ) ∩ Ω
Also, if U a component of Λ and U¯ ⊂ Ω , then
|Λ(T + E) ∩ U | = |Λ(T ) ∩ U | Arizona
19 / 41
Pseudospectral comparison
U
Ω Ω Most useful when T is linear Even then, can be expensive to compute! What about related tools? Arizona
20 / 41
The Gershgorin picture (linear case) A = D + F,
D = diag(di ),
ρi =
X j
|fij |
ρ2
ρ1
ρ3 d1
d2
Arizona
d3
21 / 41
Gershgorin (+) Write A = D + F , D = diag(d1 , . . . , dn ). Gershgorin disks are: X Gi = z ∈ C : |z − di | ≤ |fij | . j
Useful facts:
Pf:
S Spectrum of A lies in m i=1 Gi S i∈I Gi disjoint from other disks =⇒ contains |I| eigenvalues. S A − zI strictly diagonally dominant outside m i=1 Gi . Eigenvalues of D − sF , 0 ≤ s ≤ 1, are continuous.
Arizona
22 / 41
Nonlinear Gershgorin Write T (z) = D(z) + F (z). Gershgorin regions are X Gi = z ∈ C : |di (z)| ≤ |fij (z)| . j
Useful facts: Spectrum of T lies in
Sm
i=1 Gi
S Bdd connected component of m i=1 Gi strictly in Ω =⇒ same number of eigs of D and T in component =⇒ at least one eig per component of Gi involved Pf: Strict diag dominance test + continuity of eigs
Arizona
23 / 41
Example I: Hadeler T (z) = (ez − 1)B + z 2 A − αI,
A, B ∈ R8×8
20
15
10
5
0
−5
−10
−15
−20 −20
−15
−10
−5
0
Arizona
5
10
15
20
24 / 41
Comparison to simplified problem
Bauer-Fike idea: apply a similarity! T (z) = (ez − 1)B + z 2 A − αI T˜(z) = U T T (z)U = (ez − 1)DB + z 2 I − αE = D(z) − αE
Gi = {z : |βi (ez − 1) + z 2 | < ρi }.
Arizona
25 / 41
Gershgorin regions 20
15
10
5
0
−5
−10
−15
−20 −20
−15
−10
−5
0
Arizona
5
10
15
20
26 / 41
A different comparison Approximate ez − 1 by a Chebyshev interpolant: T (z) = (ez − 1)B + z 2 A − αI T˜(z) = q(z)B + z 2 A − α T (z) = T˜(z) + r(z)B Linearize T˜ and transform both: T˜(z) 7→ DC − zI
T (z) 7→ DC − zI + r(z)E Restrict to Ω = {z : |r(z)| < }: ˆ i = {z : |z − µi | < ρi }, Gi ⊂ G Arizona
ρi =
X j
|eij | 27 / 41
Spectrum of T˜ 10
8
6
4
2
0
−2
−4
−6
−8
−10 −10
−8
−6
−4
−2
0
Arizona
2
4
6
8
10
28 / 41
ˆ i for < 10−10 G 10
8
6
4
2
0
−2
−4
−6
−8
−10 −10
−8
−6
−4
−2
0
Arizona
2
4
6
8
10
29 / 41
ˆ i for = 0.1 G 10
8
6
4
2
0
−2
−4
−6
−8
−10 −10
−8
−6
−4
−2
0
Arizona
2
4
6
8
10
30 / 41
ˆ i for = 1.6 G 10
8
6
4
2
0
−2
−4
−6
−8
−10 −10
−8
−6
−4
−2
0
Arizona
2
4
6
8
10
31 / 41
Example II: Resonance problem
V0
a
b
ψ(0) = 0 d2 − 2 +V −λ ψ =0 on (0, b), dx √ ψ 0 (b) = i λψ(b), Arizona
32 / 41
Reduction via shooting
ψ(0) = 0,
V0
R0a (λ)
a R0a (λ)
ψ(0) ψ(a) = , ψ 0 (0) ψ 0 (a)
ψ(b) ψ(b) Rab (λ) 0 = 0 , ψ (b) ψ (b)
b Rab (λ)
√ ψ 0 (b) = i λψ(b)
Arizona
33 / 41
Reduction via shooting First-order form: du 0 1 ψ(x) = u, where u(x) ≡ 0 . V −λ 0 ψ (x) dx On region (c, d) where V is constant: u(d) = Rcd (λ)u(c),
Rcd (λ) = exp (d − c)
0 1 V −λ 0
Reduce resonance problem to 6D NEP: R0a (λ) −I 0 u(0) 0 Rab (λ) −I u(a) = 0. T (λ)uall ≡ 1 0 0√ 0 u(b) 0 0 0 −i λ 1 Arizona
34 / 41
Expansion via rational approximation Consider the equation
A B u − λI =0 C D v
Partial Gaussian elimination gives the spectral Schur complement A − λI − B(D − λI)−1 C u = 0 Idea: Given T (λ) = A − λI − F (λ),
find a rational approximation F (λ) ≈ B(D − λI)−1 C. Arizona
35 / 41
Expansion via rational approximation u(0) u(a) u(b)
×
Arizona
=0 .. .
36 / 41
Analyzing the expanded system
Tˆ(z) is a Schur complement in K − zM So Λ(Tˆ) is easy to compute.
Or: think T (z) is a Schur complement in K − zM + E(z) Compare Tˆ(z) to T (z) or compare K − zM + E(z) to K − zM
Arizona
37 / 41
Analyzing the expanded system Q: Can we find all eigs in a region not missing anything? Concrete plan ( = 10−8 ) T = shooting system Tˆ = rational approximation Find region D with boundary Γ s.t. D ⊂ Ω (i.e. kT − Tˆk < on D) Γ does not intersect Λ (T )
=⇒ Same eigenvalue counts for T , Tˆ =⇒ Eigs of Tˆ in components of Λ (T ) Converse holds if D ⊂ Ω/2
Can refine eigs of Tˆ in D via Newton. Arizona
38 / 41
Resonance approximation 50
−10 −6
0
40
−8
−8
−6 −4 −2
20
−1
0
30
10
−10
−10
0
0
−20 0
−1
−30
0
50
100
150
200
Arizona
250
300
−6
−10
−8
−6 −4 −2
−50 −50
−8
0
−40
350
400 39 / 41
Resonance approximation 2
−2 −6 −4
0 1
0
0
0
−10
−8
−2
−1
0
−6
−4
−2
−3
−2
−8
−5 −10
−8
−6
−10
−4
−4
−2
0 Arizona
2
4
6
8
10 40 / 41
For more
Localization theorems for nonlinear eigenvalues. David Bindel and Amanda Hood, SIAM Review 57(4), Dec 2015
Arizona
41 / 41