Structured Doubling Algorithms for Solving g-Palindromic Quadratic - - PowerPoint PPT Presentation

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Structured Doubling Algorithms for Solving g-Palindromic Quadratic Eigenvalue Problems

Eric King-wah Chu

School of Mathematical Sciences Building 28, Monash University VIC 3800, Australia eric.chu@sci.monash.edu.au This is a joint work with T.M. Huang and W.W. Lin

RANMEP2008

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Existing, Alternative Approaches S + S−1, H2 Transformations Structure-Preserving Doubling Algorithms Palindromic Factorization g-SDA Conclusions

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Abstract

The T-palindromic quadratic eigenvalue problem (λ2B + λC + A)x = 0 with A, B, C ∈ Cn×n, B⊤ = A and C ⊤ = C, governs the vibration behaviour of trains. One way to solve the problem is to apply a structure-preserving doubling algorithm (SDA) to the nonlinear matrix equation (NME) X + BX −1A = C and “square-root” the matrix quadratic involved. In this talk, we generalize the SDA for the solution of (odd and even) T- and H-palindromic quadratic eigenvalue problems in a unified fashion.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T-Palindromic Quadratic Eigenvalue Problems

Consider the palindromic quadratic eigenvalue problem (QEP) P(λ)x ≡ (λ2A⊤

1 + λA0 + A1)x = 0

(1) where A1, A0 ∈ Cn×n, and A⊤

0 = A0.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T-Palindromic Quadratic Eigenvalue Problems

Consider the palindromic quadratic eigenvalue problem (QEP) P(λ)x ≡ (λ2A⊤

1 + λA0 + A1)x = 0

(1) where A1, A0 ∈ Cn×n, and A⊤

0 = A0.

Symplectic structure: eigenvalues in pairs λ and 1/λ.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T-Palindromic Quadratic Eigenvalue Problems

Consider the palindromic quadratic eigenvalue problem (QEP) P(λ)x ≡ (λ2A⊤

1 + λA0 + A1)x = 0

(1) where A1, A0 ∈ Cn×n, and A⊤

0 = A0.

Symplectic structure: eigenvalues in pairs λ and 1/λ. RevPm(λ) ≡ λmPm(1/λ) = Pm(λ)⊤

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T-Palindromic Quadratic Eigenvalue Problems

Consider the palindromic quadratic eigenvalue problem (QEP) P(λ)x ≡ (λ2A⊤

1 + λA0 + A1)x = 0

(1) where A1, A0 ∈ Cn×n, and A⊤

0 = A0.

Symplectic structure: eigenvalues in pairs λ and 1/λ. RevPm(λ) ≡ λmPm(1/λ) = Pm(λ)⊤ For QEPs, P(λ) = λ2A⊤

1 + λA0 + A1,

we have RevP(λ) = A⊤

1 + λA0 + λ2A1.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T-Palindromic Quadratic Eigenvalue Problems

Consider the palindromic quadratic eigenvalue problem (QEP) P(λ)x ≡ (λ2A⊤

1 + λA0 + A1)x = 0

(1) where A1, A0 ∈ Cn×n, and A⊤

0 = A0.

Symplectic structure: eigenvalues in pairs λ and 1/λ. RevPm(λ) ≡ λmPm(1/λ) = Pm(λ)⊤ For QEPs, P(λ) = λ2A⊤

1 + λA0 + A1,

we have RevP(λ) = A⊤

1 + λA0 + λ2A1.

Important to preserve the symplectic structure!

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Standard Linearizations

1. I −A1 −A0 x λx

• = λ

I A⊤

1

x λx

• 2.

−A1 I x λx

• = λ

A0 A⊤

1

I x λx

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Standard Linearizations

1. I −A1 −A0 x λx

• = λ

I A⊤

1

x λx

• 2.

−A1 I x λx

• = λ

A0 A⊤

1

I x λx

• ◮ QZ cannot preserve the symplectic structure.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Standard Linearizations

1. I −A1 −A0 x λx

• = λ

I A⊤

1

x λx

• 2.

−A1 I x λx

• = λ

A0 A⊤

1

I x λx

• ◮ QZ cannot preserve the symplectic structure.

No significant figures in numerical solutions.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Linearizations

Hilliges/Mehl/Mehrmann 2004 2006 λZ + Z ⊤ with Z = A⊤

1

A0 − A1 A⊤

1

A⊤

1

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Linearizations

Hilliges/Mehl/Mehrmann 2004 2006 λZ + Z ⊤ with Z = A⊤

1

A0 − A1 A⊤

1

A⊤

1

• ◮ Preserve symplecticity to some extent.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Linearizations

Hilliges/Mehl/Mehrmann 2004 2006 λZ + Z ⊤ with Z = A⊤

1

A0 − A1 A⊤

1

A⊤

1

• ◮ Preserve symplecticity to some extent.

◮ This linearization, coupled with standard software, lead to

better but still low accuracy.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Linearizations

Hilliges/Mehl/Mehrmann 2004 2006 λZ + Z ⊤ with Z = A⊤

1

A0 − A1 A⊤

1

A⊤

1

• ◮ Preserve symplecticity to some extent.

◮ This linearization, coupled with standard software, lead to

better but still low accuracy.

◮ Pencil 2n × 2n, O((2n)3) flops.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Anti-Triangular Schur Form & QR-Like Method

Mackey/Mackey/Mehl/Mehrmann 2007 Schr¨

• der 2007

λZ + Z ⊤ − → λ     × × × × × × × × × ×     +     × × × × × × × × × ×    

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Anti-Triangular Schur Form & QR-Like Method

Mackey/Mackey/Mehl/Mehrmann 2007 Schr¨

• der 2007

λZ + Z ⊤ − → λ     × × × × × × × × × ×     +     × × × × × × × × × ×    

◮ Congruence transformations using unitary matrices.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Anti-Triangular Schur Form & QR-Like Method

Mackey/Mackey/Mehl/Mehrmann 2007 Schr¨

• der 2007

λZ + Z ⊤ − → λ     × × × × × × × × × ×     +     × × × × × × × × × ×    

◮ Congruence transformations using unitary matrices. ◮ Jacobi method

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Anti-Triangular Schur Form & QR-Like Method

Mackey/Mackey/Mehl/Mehrmann 2007 Schr¨

• der 2007

λZ + Z ⊤ − → λ     × × × × × × × × × ×     +     × × × × × × × × × ×    

◮ Congruence transformations using unitary matrices. ◮ Jacobi method ◮ QR-like method

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Other Palindromic Linearizations

From the companion linearization λL − M ≡ λ

• I

A0 A⊤

1

• −
• I

−A1

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Other Palindromic Linearizations

From the companion linearization λL − M ≡ λ

• I

A0 A⊤

1

• −
• I

−A1

• One step of doubling produces

τ L − M = τ A⊤

1

A0 A⊤

1

• +

A1 A0 A1

• with τ = λ2 and the same eigenvector (x⊤, λx⊤)⊤.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Scaling the first row-blocks by −A−1

1

and the second column-blocks by A−⊤

1

, then swap the roles of L and M, we have a SSF form τ I G A⊤

• −
• A

−H I

• with A ≡ −A−1

1 A⊤ 1 , H ≡ −A0 = H⊤, G ≡ A−1 1 A0A−⊤ 1

= G ⊤.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Scaling the first row-blocks by −A−1

1

and the second column-blocks by A−⊤

1

, then swap the roles of L and M, we have a SSF form τ I G A⊤

• −
• A

−H I

• with A ≡ −A−1

1 A⊤ 1 , H ≡ −A0 = H⊤, G ≡ A−1 1 A0A−⊤ 1

= G ⊤. The SDA can then be applied.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Similarly, another doubling step produces (τZ + Z ⊤)y = 0 , τ = λ4 with Z =

• −A⊤

1 A−1 1 A0

−A⊤

1 A−1 1 A⊤ 1

A⊤

1 − A0A−1 1 A0

−A0A−1

1 A⊤ 1

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

S + S−1 Transformation

Lin/Qian 2007 The matrices M = A1 −A0 −I

• ,

L = I A⊤

1

• define the ⊤-symplectic pencil M − λL, with

MJ M⊤ = LJ L⊤, J = I −I

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

S + S−1 Transformation

Lin/Qian 2007 The matrices M = A1 −A0 −I

• ,

L = I A⊤

1

• define the ⊤-symplectic pencil M − λL, with

MJ M⊤ = LJ L⊤, J = I −I

• After the S + S−1 transformation:
• M − λ

L ≡ (MJ L⊤ + LJ M⊤) − λLJ L⊤

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

S + S−1 Transformation, cont.

◮ A matrix H ∈ C2n×2n is T-skew-Hamiltonian, if

(HJ )⊤ = −HJ .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

S + S−1 Transformation, cont.

◮ A matrix H ∈ C2n×2n is T-skew-Hamiltonian, if

(HJ )⊤ = −HJ .

◮

• M − λ

L = A1 − A⊤

1

A0 −A0 A1 − A⊤

1

• − λ

−A1 A⊤

1

• =
• A0

A⊤

1 − A1

A1 − A⊤

1

A0

• − λ

−A1 −A⊤

1

• J

≡ (K − λN)J

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

S + S−1 Transformation, cont.

◮ A matrix H ∈ C2n×2n is T-skew-Hamiltonian, if

(HJ )⊤ = −HJ .

◮

• M − λ

L = A1 − A⊤

1

A0 −A0 A1 − A⊤

1

• − λ

−A1 A⊤

1

• =
• A0

A⊤

1 − A1

A1 − A⊤

1

A0

• − λ

−A1 −A⊤

1

• J

≡ (K − λN)J

◮ Both K and N are ⊤-skew-Hamiltonian.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

⊤-skew-Hamiltonian pencil K and N: K N             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate the strictly lower triangular part of N(1 : 4, 1 : 4): K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate K(6, 1): K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × ⊗ × × × × × × ⊗ × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate K(7, 1): K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × ⊗ × × × × × × × ⊗ × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate K(8, 1) using Gs1 transformation: K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × × × × × × × × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate K(4, 1): K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × × × × × × ⊗ × × × × × × × × × × ⊗ × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Annihilate K(3, 1): K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × ⊗ × × × × × × × × × × × × ⊗ × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Generalized Givens and Householder Transformations

Repeat . . . K ← QKZ N ← QNZ             × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×                         × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×            

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Structure-Preserving Doubling Algorithm for DAREs

Chu/Fan/Lin/Wang 2004 Discrete linear system: xj+1 = Axj + Buj Optimal control: min

uj J = 1

2

∞

• j=1
• x⊤

j Hxj + u⊤ j Ruj

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Structure-Preserving Doubling Algorithm for DAREs

Chu/Fan/Lin/Wang 2004 Discrete linear system: xj+1 = Axj + Buj Optimal control: min

uj J = 1

2

∞

• j=1
• x⊤

j Hxj + u⊤ j Ruj

• ◮ Discrete Algebraic Riccati Equation (DARE):

X = A⊤X(I + GX)−1A + H

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for DAREs, cont.

◮ Solution of DARE via the EVP: (λL − M)z = 0, with

M =

• A

−H I

• ,

L = I G A⊤

• Matrix pair (M, L) is in Standard Symplectic Form (SSF).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for DAREs, cont.

◮ Solution of DARE via the EVP: (λL − M)z = 0, with

M =

• A

−H I

• ,

L = I G A⊤

• Matrix pair (M, L) is in Standard Symplectic Form (SSF).

◮ (M, L) is symplectic:

MJ M⊤ = LJ L⊤, J = I −I

• Eric Chu

Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for DAREs, cont.

◮ Solution of DARE via the EVP: (λL − M)z = 0, with

M =

• A

−H I

• ,

L = I G A⊤

• Matrix pair (M, L) is in Standard Symplectic Form (SSF).

◮ (M, L) is symplectic:

MJ M⊤ = LJ L⊤, J = I −I

• ◮ Symplectic in SSF is stronger.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Stable Solution

With MJ M⊤ = LJ L⊤, consider y⊤M = λy⊤L ⇒ y⊤MJ M⊤ = λy⊤LJ M⊤ ⇔ (y⊤LJ )L⊤ = λ(y⊤LJ )M⊤ so λ ∈ σ(M, L) ⇔ λ−1 ∈ σ(M, L), or σ(M, L) = σ(L, M).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Stable Solution

With MJ M⊤ = LJ L⊤, consider y⊤M = λy⊤L ⇒ y⊤MJ M⊤ = λy⊤LJ M⊤ ⇔ (y⊤LJ )L⊤ = λ(y⊤LJ )M⊤ so λ ∈ σ(M, L) ⇔ λ−1 ∈ σ(M, L), or σ(M, L) = σ(L, M).

◮ Stable deflating subspace gives rise to the stable symmetric

positive definite solution of the DARE.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Stable Solution

With MJ M⊤ = LJ L⊤, consider y⊤M = λy⊤L ⇒ y⊤MJ M⊤ = λy⊤LJ M⊤ ⇔ (y⊤LJ )L⊤ = λ(y⊤LJ )M⊤ so λ ∈ σ(M, L) ⇔ λ−1 ∈ σ(M, L), or σ(M, L) = σ(L, M).

◮ Stable deflating subspace gives rise to the stable symmetric

positive definite solution of the DARE.

◮ For a stable λ, |λ| < 1 and λ2m → 0 fast, as m → ∞. ◮ λ2m, λ−2m grow further apart (relative to the unit circle).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Stable Solution

With MJ M⊤ = LJ L⊤, consider y⊤M = λy⊤L ⇒ y⊤MJ M⊤ = λy⊤LJ M⊤ ⇔ (y⊤LJ )L⊤ = λ(y⊤LJ )M⊤ so λ ∈ σ(M, L) ⇔ λ−1 ∈ σ(M, L), or σ(M, L) = σ(L, M).

◮ Stable deflating subspace gives rise to the stable symmetric

positive definite solution of the DARE.

◮ For a stable λ, |λ| < 1 and λ2m → 0 fast, as m → ∞. ◮ λ2m, λ−2m grow further apart (relative to the unit circle). ◮ Null space of (M−1L)2m (m → ∞) wanted.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

◮ Find

• M,

L

• such that
• M,

L L −M

• = 0.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

◮ Find

• M,

L

• such that
• M,

L L −M

• = 0.

◮ Note that

ML = LM, thus matrix multiplication is “commutative”!

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

◮ Find

• M,

L

• such that
• M,

L L −M

• = 0.

◮ Note that

ML = LM, thus matrix multiplication is “commutative”!

◮ Doubling: (

M, L) ≡ ( MM, LL), with Mx = λLx ⇒

• MMx = λ

MLx ⇒ Mx = λ LMx = λ2 LLx = λ2 Lx

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

◮ Find

• M,

L

• such that
• M,

L L −M

• = 0.

◮ Note that

ML = LM, thus matrix multiplication is “commutative”!

◮ Doubling: (

M, L) ≡ ( MM, LL), with Mx = λLx ⇒

• MMx = λ

MLx ⇒ Mx = λ LMx = λ2 LLx = λ2 Lx

◮ Doubling preserves symplecticity (not in SSF, in general).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

M =

• A

−H I

• , L =

I G A⊤

• .

◮ Find

• M,

L

• such that
• M,

L L −M

• = 0.

◮ Note that

ML = LM, thus matrix multiplication is “commutative”!

◮ Doubling: (

M, L) ≡ ( MM, LL), with Mx = λLx ⇒

• MMx = λ

MLx ⇒ Mx = λ LMx = λ2 LLx = λ2 Lx

◮ Doubling preserves symplecticity (not in SSF, in general). ◮ There are many different possibilities for

• M,

L

• .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Preserving SSF

◮ Using QR: (Benner/Byers 2001) preserves symplecticity, not

in SSF.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Preserving SSF

◮ Using QR: (Benner/Byers 2001) preserves symplecticity, not

in SSF.

◮ Use generalized Gaussian elimination, preserves the SSF.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Preserving SSF

◮ Using QR: (Benner/Byers 2001) preserves symplecticity, not

in SSF.

◮ Use generalized Gaussian elimination, preserves the SSF. ◮ SDA iteration:

• A

← A(I + GH)−1A

• G

← G + AG(I + HG)−1A⊤

• H

← H + A⊤(I + HG)−1HA

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Preserving SSF

◮ Using QR: (Benner/Byers 2001) preserves symplecticity, not

in SSF.

◮ Use generalized Gaussian elimination, preserves the SSF. ◮ SDA iteration:

• A

← A(I + GH)−1A

• G

← G + AG(I + HG)−1A⊤

• H

← H + A⊤(I + HG)−1HA

◮ Many nice properties; operation count 15% that of the

QR-based algorithm (operating on the 4n × 2n matrix [L⊤, −M⊤]⊤).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for Palindromic QEPs

◮ SDA can be applied after deflation, etc. (no details here)

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for Palindromic QEPs

◮ SDA can be applied after deflation, etc. (no details here) ◮ Theoretically, may need to invert ill-conditioned matrices;

Cayley transform can be applied but never required numerically.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for Palindromic QEPs

◮ SDA can be applied after deflation, etc. (no details here) ◮ Theoretically, may need to invert ill-conditioned matrices;

Cayley transform can be applied but never required numerically.

◮ Quadratic convergence (linear with eigenvalues on unit circle).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for Palindromic QEPs

◮ SDA can be applied after deflation, etc. (no details here) ◮ Theoretically, may need to invert ill-conditioned matrices;

Cayley transform can be applied but never required numerically.

◮ Quadratic convergence (linear with eigenvalues on unit circle). ◮ Accurate numerical results.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

SDA for Palindromic QEPs

◮ SDA can be applied after deflation, etc. (no details here) ◮ Theoretically, may need to invert ill-conditioned matrices;

Cayley transform can be applied but never required numerically.

◮ Quadratic convergence (linear with eigenvalues on unit circle). ◮ Accurate numerical results. ◮ Twice as expensive as the factorization approach (see below).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Factorization

Chu/Huang/Lin/Wu 2007 The palindromic factorization: P(λ) = (λY − A1)Y −1(λA1 − Y )⊤ with Y = Y ⊤ satisfying the NME-P: Y + A1Y −1A⊤

1 = −A0

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Factorization

Chu/Huang/Lin/Wu 2007 The palindromic factorization: P(λ) = (λY − A1)Y −1(λA1 − Y )⊤ with Y = Y ⊤ satisfying the NME-P: Y + A1Y −1A⊤

1 = −A0 ◮ The SDA2 (Lin/Xu 2006) can then be applied.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Factorization

Chu/Huang/Lin/Wu 2007 The palindromic factorization: P(λ) = (λY − A1)Y −1(λA1 − Y )⊤ with Y = Y ⊤ satisfying the NME-P: Y + A1Y −1A⊤

1 = −A0 ◮ The SDA2 (Lin/Xu 2006) can then be applied. ◮ The palindromic EVP is then “square-rooted”.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Palindromic Factorization

Chu/Huang/Lin/Wu 2007 The palindromic factorization: P(λ) = (λY − A1)Y −1(λA1 − Y )⊤ with Y = Y ⊤ satisfying the NME-P: Y + A1Y −1A⊤

1 = −A0 ◮ The SDA2 (Lin/Xu 2006) can then be applied. ◮ The palindromic EVP is then “square-rooted”. ◮ Cyclic Reduction (CR) can be applied to the NME-P, or the

related quadratic equation: A⊤

1 X 2 + A0X + A1 = 0

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Eigenvalue Distribution

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Convergence (Qk+1 − Qk and Rk)

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

(Relative) Residuals of Approximate Eigenpairs (λ, x)

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Reciprocal Properties of Eigenvalues (|λiλ2n−i+1 − 1|)

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

(∗, ε)-Homomorphism

Definition

A function g : Cn×n → Cn×n is called a (∗, ε)-homomorphism if g(α1Φ1 + α2Φ2) = α∗

1g(Φ1) + α∗ 2g(Φ2)

g(Φ1Φ2) = εg(Φ2)g(Φ1) for all Φ1, Φ2 ∈ Cn×n and α1, α2 ∈ C. Furthermore, g preserves the singularity, i.e., det(Φ) = 0 ⇔ det(g(Φ)) = 0 Here “∗” denotes “H” (Hermition/conjugate transpose) or “T” (transpose) and ε = ±1.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Theorem

Let g be a (∗, ε)-homomorphism. Then it holds (i) g(0) = 0, (ii) g(I) = εI, and (iii) g(Φ−1) = g(Φ)−1.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Theorem

Let g be a (∗, ε)-homomorphism. Then it holds (i) g(0) = 0, (ii) g(I) = εI, and (iii) g(Φ−1) = g(Φ)−1.

Proof.

(i) g(Φ) = g(Φ + 0) = g(Φ) + g(0). Therefore, g(0) = 0. (ii) Let Φ be nonsingular. Then g(Φ) = g(Φ · I) = εg(I)g(Φ). From det(g(Φ)) = 0 follows that g(I) = εI. (iii) g(I) = g(Φ−1 · Φ) = εg(Φ)g(Φ−1). Therefore, from (ii) we get g(Φ−1) = g(Φ)−1.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-Palindromic QEPs

Definition

A quadratic eigenvalue problem (QEP) Q(λ) ≡ (λ2B + λC + A)x = 0, (2) where A, B, C ∈ Cn×n, λ ∈ C, x = 0 ∈ Cn, is called a g-palindromic QEP if there is a ∗-homomorphism g such that g(B) = A, g(C) = C and g(A) = B. Moreover, A and B are called g-related (denoted by A

g

∼ B) and C is g-symmetric.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Symplecticity

Theorem

Let Q(λ) be a g-palindromic quadratic pencil. Then λ ∈ σ(Q(λ)) if and only if 1/λ∗ ∈ σ(Q(λ)).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Symplecticity

Theorem

Let Q(λ) be a g-palindromic quadratic pencil. Then λ ∈ σ(Q(λ)) if and only if 1/λ∗ ∈ σ(Q(λ)).

Proof.

Without loss of generality, assume that λ = 0. Then = det(λ2B + λC + A) = det(g(λ2B + λC + A)) = det

• (λ∗)2g(B) + λ∗g(C) + g(A)
• =

det

• (λ∗)2A + λ∗C + B
• .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-NME

A quadratic pencil can be factorized as λ2B + λC + A = (λB + X)X −1(λX + A) = λ2B + λ(X + BX −1A) + A with X satisfying the nonlinear matrix equation (NME) X + BX −1A = C

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-NME

A quadratic pencil can be factorized as λ2B + λC + A = (λB + X)X −1(λX + A) = λ2B + λ(X + BX −1A) + A with X satisfying the nonlinear matrix equation (NME) X + BX −1A = C If we can find a solution X for the NME structurally, then the palindromic QEP is solved.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-NME

A quadratic pencil can be factorized as λ2B + λC + A = (λB + X)X −1(λX + A) = λ2B + λ(X + BX −1A) + A with X satisfying the nonlinear matrix equation (NME) X + BX −1A = C If we can find a solution X for the NME structurally, then the palindromic QEP is solved. The g-SDA solves g-NMEs, thus g-palindromic QEPs.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-SDA

For a given g-NME, we define M = A C −I

• ,

L = −D I B

• .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-SDA

For a given g-NME, we define M = A C −I

• ,

L = −D I B

• .

We have the g-SDA:

• A

= A(C − D)−1A,

• B = B(C − D)−1B
• C

= C − B(C − D)−1A,

• D = D + A(C − D)−1B

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

g-SDA

For a given g-NME, we define M = A C −I

• ,

L = −D I B

• .

We have the g-SDA:

• A

= A(C − D)−1A,

• B = B(C − D)−1B
• C

= C − B(C − D)−1A,

• D = D + A(C − D)−1B
• A and

B are g-related, and C and D are g-symmetric.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Doubling

Theorem

(i) The pencil M − λ L has the doubling property; i.e., if M U V

• = L

U V

• S,

where U, V ∈ Cn×m and S ∈ Cm×m, then

• M

U V

• =

L U V

• S2.

(ii) The quadratic pencil λ2 B + λ C + A corresponding to M − λ L is still a g-palindromic quadratic pencil.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Proof. (i) From M U V

• = L

U V

• S and M∗L = L∗M, we have
• M

U V

• =

M∗M U V

• = M∗L

U V

• S

= L∗M U V

• S = L∗L

U V

• S2

=

• L

U V

• S2

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

(ii) Using the properties of g, we have g( A) = ε2g(A)(g(C) − g(D))−1g(A) = B(C − D)−1B = B, g( B) = ε2g(B)(g(C) − g(D))−1g(B) = A(C − D)−1A = A, g( C) = g(C) − ε2g(A)(g(C) − g(D))−1g(B) = C − B(C − D)−1A = C, g( D) = g(D) + ε2g(B)(g(C) − g(D))−1g(A) = D + A(C − D)−1B = D.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

(ii) Using the properties of g, we have g( A) = ε2g(A)(g(C) − g(D))−1g(A) = B(C − D)−1B = B, g( B) = ε2g(B)(g(C) − g(D))−1g(B) = A(C − D)−1A = A, g( C) = g(C) − ε2g(A)(g(C) − g(D))−1g(B) = C − B(C − D)−1A = C, g( D) = g(D) + ε2g(B)(g(C) − g(D))−1g(A) = D + A(C − D)−1B = D. Therefore, A

g

∼ B, C and D are g-symmetric and λ2 B + λ C + A is again a g-palindromic quadratic pencil. ✷

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Convergence with Unimodular Eigenvalues

Definition

A solution X of a g-NME is said to possess property (P), if (i) ρ(X −1A) ≤ 1; and (ii) the partial multiplicities of each unimodular eigenvalue of X −1A is half of that of the corresponding unimodular eigenvalue of the associated pair (M, L).

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Theorem

Assume that the g-NME and the dual g-NME have the solutions X and Y with property (P), respectively. Suppose the sequence {Ak, Bk, Ck, Dk} generated by the g-SDA is well-defined. Then (i) Ak = O(2−k) → 0, as k → ∞, (ii) Bk = O(2−k) → 0, as k → ∞, (iii) Ck − X = O(2−k) → 0, as k → ∞, (iv) Dk − Y = O(2−k) → 0, as k → ∞. Furthermore, X and Y are g-symmetric, i.e., g(X) = X and g(Y ) = Y .

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T/H-Palindromic QEPs

If the (∗, ε)-homomorphism is defined by g(Φ) = +Φ∗. Then the g-palindromic QEP becomes (i) T-palindromic QEP (∗ = “T”): (λ2A⊤ + λC + A)x = 0 with C ⊤ = +C. (ii) H-palindromic QEP (∗ = “H”): (λ2AH + λC + A)x = 0 with C H = +C.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

T/H-Anti-Palindromic QEPs

If the (∗, ε)-homomorphism is defined by g(Φ) = −Φ∗. Then the g-palindromic QEP becomes (iii) T-anti-palindromic QEP (∗ = “T”): (λ2A⊤ + λC − A)x = 0 with C ⊤ = −C. (iv) H-anti-palindromic QEP (∗ = “H”): (λ2AH + λC − A)x = 0 with C H = −C.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

The g-SDA for cases (i) and (ii): A0 = A, C0 = C = +C ∗, D0 = 0, Ak+1 = Ak(Ck − Dk)−1Ak, Ck+1 = Ck − A∗

k(Ck − Dk)−1Ak,

Dk+1 = Dk + Ak(Ck − Dk)−1Ak.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

The g-SDA for cases (iii) and (iv): A0 = A, C0 = C = −C ∗, D0 = 0, Ak+1 = Ak(Ck − Dk)−1Ak, Ck+1 = Ck + A∗

k(Ck − Dk)−1Ak,

Dk+1 = Dk − Ak(Ck − Dk)−1A∗

k.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Other g-Palindromic QEPs

For the ∗-(anti-)palindromic QEP: (λ2A∗ + λC ∓ A)x = 0 with C ∗ = ±C,

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

Other g-Palindromic QEPs

For the ∗-(anti-)palindromic QEP: (λ2A∗ + λC ∓ A)x = 0 with C ∗ = ±C, The quadratic pencil can be factorized λA∗ + λC ∓ A = (λA∗ + X)X −1(λX ∓ A) where X satisfies the g-NME X ∓ A∗X −1A = C, C ∗ = ±C.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

If we perform one step of the g-SDA on the above g-NME, then X satisfies X ± A∗X −1 A = C,

• C ∗ = ±

C, where

• A = AC −1A,

C = C ∓ A∗C −1A, D = ±AC −1A∗.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

If we perform one step of the g-SDA on the above g-NME, then X satisfies X ± A∗X −1 A = C,

• C ∗ = ±

C, where

• A = AC −1A,

C = C ∓ A∗C −1A, D = ±AC −1A∗. The g-NME corresponds to the ∗-(anti-)palindromic QEP (λ2 A∗ + λ C ± A)x = 0.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

*-Even/Odd Palindromic QEPs

For the QEP Q(λ)x ≡ (λ2M + λG + K)x = 0, M∗ = ±M, K ∗ = ±K, G ∗ = ∓G,

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

*-Even/Odd Palindromic QEPs

For the QEP Q(λ)x ≡ (λ2M + λG + K)x = 0, M∗ = ±M, K ∗ = ±K, G ∗ = ∓G, It is well-known that Q(λ) has the factorization Q(λ) = (λM + MS + G)(λI − S) if and only if S is a solution of the quadratic matrix equation MS2 + GS + K = 0.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

If λ ∈ σ(Q(λ)), then −λ∗ ∈ σ(Q(λ)). If xi and yi are, respectively, the right and left eigenvectors corresponding to λi of the solvant S, i.e., Sxi = λixi, y∗

i S = λiy∗ i ,

(3) then xi and (λiM + MS + G)−∗yi are eigenvectors corresponding to λi and −λ∗

i , respectively, of the QEP.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

It seems difficult to find the solvant S directly whose eigenvalues are on the right half-plane. Instead, the Cayley transformation S = (I + Y )(I − Y )−1 is used. The solvant S then satisfies εA∗Y 2 + CY + A = 0, (4) where A = M + K + G, C = 2(M − K), ε = ±1. With Y = −X −1A, we have the NME: X + εA∗X −1A = C. The g-SDA can then be applied.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

◮ Palindromic linearization/QZ, Jacobi or QR-like methods for

anti-triangular Schur form.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

◮ Palindromic linearization/QZ, Jacobi or QR-like methods for

anti-triangular Schur form.

◮ S + S−1 method stable and accurate. ◮ SDA theoretically less stable, more efficient, preserving SSF.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

◮ Palindromic linearization/QZ, Jacobi or QR-like methods for

anti-triangular Schur form.

◮ S + S−1 method stable and accurate. ◮ SDA theoretically less stable, more efficient, preserving SSF. ◮ Better finite elements, deflation, correction or refinement,

model reduction, optimization of design parameters, higher-order problems, large-scaled problems, multi-frequency problems; other applications in laser & semi-conductor simulations, optimization, control system design, etc.

◮ Related (⊤ or ∗, odd or even) palindromic problems.

Eric Chu Structured Doubling Algorithms

Existing Approaches S + S−1 Transformation Doubling Factorization g-SDA Conclusions

◮ Palindromic linearization/QZ, Jacobi or QR-like methods for

anti-triangular Schur form.

◮ S + S−1 method stable and accurate. ◮ SDA theoretically less stable, more efficient, preserving SSF. ◮ Better finite elements, deflation, correction or refinement,

model reduction, optimization of design parameters, higher-order problems, large-scaled problems, multi-frequency problems; other applications in laser & semi-conductor simulations, optimization, control system design, etc.

◮ Related (⊤ or ∗, odd or even) palindromic problems. ◮ g-SDA for g-palindromic quadratic eigenvalue problems, via

g-NME and factorization of the matrix quadratic.

Eric Chu Structured Doubling Algorithms