October 30, 2017 | Author: Anonymous | Category: N/A
The implicit QZ algorithm for the palindromic eigenvalue problem David S. Watkins watkins ......
The implicit QZ algorithm for the palindromic eigenvalue problem David S. Watkins
[email protected]
Department of Mathematics Washington State University
Luminy, October 2007 – p.
Context
Luminy, October 2007 – p.
Context LQG problem in discrete time (optimal control)
Luminy, October 2007 – p.
Context LQG problem in discrete time (optimal control) ⇒ symplectic eigenvalue problem
Luminy, October 2007 – p.
Context LQG problem in discrete time (optimal control) ⇒ symplectic eigenvalue problem eigenvalue symmetry: λ, λ−1
Luminy, October 2007 – p.
Context LQG problem in discrete time (optimal control) ⇒ symplectic eigenvalue problem eigenvalue symmetry: λ, λ−1 2
1.5
1
0.5
0
−0.5
−1
−1.5
−2 −2
−1.5
−1
−0.5
0
0.5
1
1.5
2
Luminy, October 2007 – p.
Context LQG problem in discrete time (optimal control) ⇒ symplectic eigenvalue problem eigenvalue symmetry: λ, λ−1 2
1.5
1
0.5
0
−0.5
−1
−1.5
−2 −2
−1.5
−1
−0.5
0
0.5
1
1.5
2
want stable invariant subspace Luminy, October 2007 – p.
Palindromic Eigenvalue Problem
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0 generalized eigenvalue problem A − λB
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0 generalized eigenvalue problem A − λB B = AT
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0 generalized eigenvalue problem A − λB B = AT (C − λC T )x = 0
⇔
xT (C T − λC) = 0
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0 generalized eigenvalue problem A − λB B = AT (C − λC T )x = 0
⇔
xT (C T − λC) = 0
eigenvalue symmetry: λ, λ−1
Luminy, October 2007 – p.
Palindromic Eigenvalue Problem Mackey / Mackey / Mehl / Mehrmann (C − λC T )x = 0 generalized eigenvalue problem A − λB B = AT (C − λC T )x = 0
⇔
xT (C T − λC) = 0
eigenvalue symmetry: λ, λ−1 Formulate control problems as palindromic eigenvalue problems.
Luminy, October 2007 – p.
QR-like algorithms
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin)
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin) explicit QR-like algorithm
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin) explicit QR-like algorithm (It works!)
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin) explicit QR-like algorithm (It works!) implicit QR-like algorithm (bulge chase)
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin) explicit QR-like algorithm (It works!) implicit QR-like algorithm (bulge chase) (It works!)
Luminy, October 2007 – p.
QR-like algorithms for palindromic problems C. Schröder (Berlin) explicit QR-like algorithm (It works!) implicit QR-like algorithm (bulge chase) (It works!) Are they equivalent?
Luminy, October 2007 – p.
The Bulge Chase
Luminy, October 2007 – p.
The Bulge Chase compact form needed
Luminy, October 2007 – p.
The Bulge Chase compact form needed
∗ C= ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p.
∗ T C = 0 ∗ 0 ∗ ∗ 0 ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p.
∗ T C = 0 ∗ 0 ∗ ∗ 0 ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
not always attainable!
Luminy, October 2007 – p.
The Bulge Chase
Luminy, October 2007 – p.
The Bulge Chase single-shift case
Luminy, October 2007 – p.
The Bulge Chase single-shift case (for simplicity)
Luminy, October 2007 – p.
The Bulge Chase single-shift case (for simplicity) Pick a shift µ.
Luminy, October 2007 – p.
The Bulge Chase single-shift case (for simplicity) Pick a shift µ. 0 .. . T x = (C − µC )e1 = 0 α β
Luminy, October 2007 – p.
The Bulge Chase single-shift case (for simplicity) Pick a shift µ. 0 .. . T x = (C − µC )e1 = 0 α β
Build a Givens rotation (for example) that annihilates α. Luminy, October 2007 – p.
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p.
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p.
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 1
The Bulge Chase
∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
But what happens when we get to the middle?
Luminy, October 2007 – p. 2
Bulge Chase in the Pencil
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Bulge Chase in the Pencil
T C − λC = ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Bulge Chase in the Pencil
T C − λC = ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Bulge pencils on a collision course!
Luminy, October 2007 – p. 2
Bulge Chase in the Pencil
T C − λC = ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Half a step further:
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Swap the shifts
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Swap the shifts
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Reverse the process
T C − λC = ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 2
Reverse the process
T C − λC = ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 3
Reverse the process
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 3
Bulge Chase is complete
∗ T C − λC = ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 3
This works very well.
Luminy, October 2007 – p. 3
This works very well. Question:
Luminy, October 2007 – p. 3
This works very well. Question: How do we explain it?
Luminy, October 2007 – p. 3
This works very well. Question: How do we explain it? Answer:
Luminy, October 2007 – p. 3
This works very well. Question: How do we explain it? Answer: I don’t have time . . .
Luminy, October 2007 – p. 3
This works very well. Question: How do we explain it? Answer: I don’t have time . . . . . . but I’ll try.
Luminy, October 2007 – p. 3
This works very well. Question: How do we explain it? Answer: I don’t have time . . . . . . but I’ll try. Compare with standard QZ.
Luminy, October 2007 – p. 3
QZ algorithm,
Luminy, October 2007 – p. 3
QZ algorithm, implicit version
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular)
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm p(z) = (z − µ1 ) · · · (z − µm )
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm p(z) = (z − µ1 ) · · · (z − µm ) x = p(AB −1 )e1
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm p(z) = (z − µ1 ) · · · (z − µm ) x = p(AB −1 )e1 Make a bulge,
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm p(z) = (z − µ1 ) · · · (z − µm ) x = p(AB −1 )e1 Make a bulge, then chase it.
Luminy, October 2007 – p. 3
QZ algorithm, implicit version A − λB (Hessenberg, triangular) pick shifts µ1 , . . . µm p(z) = (z − µ1 ) · · · (z − µm ) x = p(AB −1 )e1 Make a bulge, then chase it. ˆ Get Aˆ − λB
Luminy, October 2007 – p. 3
QZ algorithm,
Luminy, October 2007 – p. 3
QZ algorithm, explicit version
Luminy, October 2007 – p. 3
QZ algorithm, explicit version p(AB −1 ) = QR
Luminy, October 2007 – p. 3
QZ algorithm, explicit version p(AB −1 ) = QR
˜ p(B −1 A) = Z R
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B Explicit QZ step is complete.
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B Explicit QZ step is complete. ˆ −1 = Q∗ (AB −1 )Q AˆB
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B Explicit QZ step is complete. ˆ −1 = Q∗ (AB −1 )Q (QR iteration on AB −1 ) AˆB
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B Explicit QZ step is complete. ˆ −1 = Q∗ (AB −1 )Q (QR iteration on AB −1 ) AˆB ˆ −1 Aˆ = Z ∗ (B −1 A)Z B
Luminy, October 2007 – p. 3
QZ algorithm, explicit version ˜ p(AB −1 ) = QR p(B −1 A) = Z R ˆ = Q∗ BZ Aˆ = Q∗ AZ, B Explicit QZ step is complete. ˆ −1 = Q∗ (AB −1 )Q (QR iteration on AB −1 ) AˆB ˆ −1 Aˆ = Z ∗ (B −1 A)Z (QR iteration on B −1 A) B
Luminy, October 2007 – p. 3
Explicit = Implicit?
Luminy, October 2007 – p. 3
Explicit = Implicit? AB −1 and B −1 A are upper Hessenberg.
Luminy, October 2007 – p. 3
Explicit = Implicit? AB −1 and B −1 A are upper Hessenberg. Utilize the Hessenberg form.
Luminy, October 2007 – p. 3
Explicit = Implicit? AB −1 and B −1 A are upper Hessenberg. Utilize the Hessenberg form. implicit-Q theorem, or . . .
Luminy, October 2007 – p. 3
Explicit = Implicit? AB −1 and B −1 A are upper Hessenberg. Utilize the Hessenberg form. implicit-Q theorem, or . . . work directly with the Krylov subspaces.
Luminy, October 2007 – p. 3
Explicit = Implicit? AB −1 and B −1 A are upper Hessenberg. Utilize the Hessenberg form. implicit-Q theorem, or . . . work directly with the Krylov subspaces. time permitting . . .
Luminy, October 2007 – p. 3
Back to the palindromic case:
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity)
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity) Cˆ −T Cˆ = GT (C −T C)G−T
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity) Cˆ −T Cˆ = GT (C −T C)G−T
(similarity)
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity) Cˆ −T Cˆ = GT (C −T C)G−T ˜ Need p(CC −T ) = GR
(similarity)
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity) Cˆ −T Cˆ = GT (C −T C)G−T (similarity) ˜ and p(C −T C) = G−T R Need p(CC −T ) = GR
Luminy, October 2007 – p. 3
Back to the palindromic case: Bulge chase gives Cˆ = G−1 CG−T Cˆ Cˆ −T = G−1 (CC −T )G (similarity) Cˆ −T Cˆ = GT (C −T C)G−T (similarity) ˜ and p(C −T C) = G−T R Need p(CC −T ) = GR or something like that.
Luminy, October 2007 – p. 3
What we actually get:
Luminy, October 2007 – p. 3
What we actually get: r(CC −T ) = GL
Luminy, October 2007 – p. 3
What we actually get: r(CC −T ) = GL m Y z − µk where r(z) = µk z − 1 k=1
Luminy, October 2007 – p. 3
What we actually get: r(CC −T ) = GL m Y z − µk where r(z) = µk z − 1 k=1
L is lower triangular.
Luminy, October 2007 – p. 3
What we actually get: r(CC −T ) = GL m Y z − µk where r(z) = µk z − 1 k=1
L is lower triangular. r(C −T C) = G−T R
Luminy, October 2007 – p. 3
What we actually get: r(CC −T ) = GL m Y z − µk where r(z) = µk z − 1 k=1
L is lower triangular. r(C −T C) = G−T R r(CC −T )−T = r(C −T C),
R = L−T
Luminy, October 2007 – p. 3
Our Approach:
Luminy, October 2007 – p. 3
Our Approach: want r(CC −T ) = GL
Luminy, October 2007 – p. 3
Our Approach: want r(CC −T ) = GL Define L = G−1 r(CC −T )
Luminy, October 2007 – p. 3
Our Approach: want r(CC −T ) = GL Define L = G−1 r(CC −T ) Then show that L is lower triangular.
Luminy, October 2007 – p. 3
Our Approach: want r(CC −T ) = GL Define L = G−1 r(CC −T ) Then show that L is lower triangular. This is not entirely straightforward.
Luminy, October 2007 – p. 3
Our Approach: want r(CC −T ) = GL Define L = G−1 r(CC −T ) Then show that L is lower triangular. This is not entirely straightforward. Utilize the Hessenberg form.
Luminy, October 2007 – p. 3
Recall that
∗ C= ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗ ∗
0 ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 4
This implies
Luminy, October 2007 – p. 4
This implies
CC −T
=
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 4
This implies
CC −T
=
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
partly lower Hessenberg (n/2 − 1 columns)
Luminy, October 2007 – p. 4
and
Luminy, October 2007 – p. 4
and
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ −T C C= ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
Luminy, October 2007 – p. 4
and
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗ −T C C= ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
∗ ∗ ∗ ∗ ∗ ∗ ∗ ∗
partly upper Hessenberg (n/2 − 2 columns)
Luminy, October 2007 – p. 4
Use both.
Luminy, October 2007 – p. 4
Use both. L11 L12 L= L21 L22
Luminy, October 2007 – p. 4
Use both. L11 L12 L= L21 L22 Using CC −T ,
Luminy, October 2007 – p. 4
Use both. L11 L12 L= L21 L22 Using CC −T , deduce L12 = 0 . . .
Luminy, October 2007 – p. 4
Use both. L11 L12 L= L21 L22 Using CC −T , deduce L12 = 0 . . . and L22 is lower triangular.
Luminy, October 2007 – p. 4
L=
L11 L21 L22
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular R11 R12 −T R=L = R22
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular R11 R12 −T R=L = R22 Using C −T C,
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular R11 R12 −T R=L = R22 Using C −T C, deduce that R11 is upper triangular.
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular R11 R12 −T R=L = R22 Using C −T C, deduce that R11 is upper triangular. Hence L11 is lower triangular.
Luminy, October 2007 – p. 4
L=
L11 L21 L22
just need L11 lower triangular R11 R12 −T R=L = R22 Using C −T C, deduce that R11 is upper triangular. Hence L11 is lower triangular. done!
Luminy, October 2007 – p. 4
Conclusion:
Luminy, October 2007 – p. 4
Conclusion: Schröder’s palindromic bulge-chasing algorithm
Luminy, October 2007 – p. 4
Conclusion: Schröder’s palindromic bulge-chasing algorithm effects a combination QL and QR iteration
Luminy, October 2007 – p. 4
Conclusion: Schröder’s palindromic bulge-chasing algorithm effects a combination QL and QR iteration driven by a rational function m Y z − µk r(z) = µk z − 1 k=1
Luminy, October 2007 – p. 4
Conclusion: Schröder’s palindromic bulge-chasing algorithm effects a combination QL and QR iteration driven by a rational function m Y z − µk r(z) = µk z − 1 k=1
with shifts µ1 , . . . , µm in the numerator and shifts µ−1 1 , . . . , µ−1 m in the denominator.
Luminy, October 2007 – p. 4
Conclusion: Schröder’s palindromic bulge-chasing algorithm effects a combination QL and QR iteration driven by a rational function m Y z − µk r(z) = µk z − 1 k=1
with shifts µ1 , . . . , µm in the numerator and shifts µ−1 1 , . . . , µ−1 m in the denominator.
That’s all, folks!
Luminy, October 2007 – p. 4