Jacobi rotation - Alchetron, The Free Social Encyclopedia

In numerical linear algebra, a Jacobi rotation is a rotation, Q_kℓ, of a 2-dimensional linear subspace of an n-dimensional inner product space, chosen to zero a symmetric pair of off-diagonal entries of an n×n real symmetric matrix, A, when applied as a similarity transformation:

A ↦ Q k ℓ T A Q k ℓ = A ′ . [ ∗ ⋯ ∗ ⋱ a k k ⋯ a k ℓ ⋮ ⋮ ⋱ ⋮ ⋮ a ℓ k ⋯ a ℓ ℓ ⋱ ∗ ⋯ ∗ ] → [ ∗ ⋯ ∗ ⋱ a k k ′ ⋯ 0 ⋮ ⋮ ⋱ ⋮ ⋮ 0 ⋯ a ℓ ℓ ′ ⋱ ∗ ⋯ ∗ ] .

It is the core operation in the Jacobi eigenvalue algorithm, which is numerically stable and well-suited to implementation on parallel processors.

Only rows k and ℓ and columns k and ℓ of A will be affected, and that A′ will remain symmetric. Also, an explicit matrix for Q_kℓ is rarely computed; instead, auxiliary values are computed and A is updated in an efficient and numerically stable way. However, for reference, we may write the matrix as

Q k ℓ = [ 1 ⋱ 0 c ⋯ s ⋮ ⋱ ⋮ − s ⋯ c 0 ⋱ 1 ] .

That is, Q_kℓ is an identity matrix except for four entries, two on the diagonal (q_kk and q_ℓℓ, both equal to c) and two symmetrically placed off the diagonal (q_kℓ and q_ℓk, equal to s and −s, respectively). Here c = cos ϑ and s = sin ϑ for some angle ϑ; but to apply the rotation, the angle itself is not required. Using Kronecker delta notation, the matrix entries can be written

q i j = δ i j + ( δ i k δ j k + δ i ℓ δ j ℓ ) ( c − 1 ) + ( δ i k δ j ℓ − δ i ℓ δ j k ) s .

Suppose h is an index other than k or ℓ (which must themselves be distinct). Then the similarity update produces, algebraically,

a h k ′ = a k h ′ = c a h k − s a h ℓ a h ℓ ′ = a ℓ h ′ = c a h ℓ + s a h k a k ℓ ′ = a ℓ k ′ = ( c 2 − s 2 ) a k ℓ + s c ( a k k − a ℓ ℓ ) = 0 a k k ′ = c 2 a k k + s 2 a ℓ ℓ − 2 s c a k ℓ a ℓ ℓ ′ = s 2 a k k + c 2 a ℓ ℓ + 2 s c a k ℓ .

Numerically stable computation

To determine the quantities needed for the update, we must solve the off-diagonal equation for zero (Golub & Van Loan 1996, §8.4). This implies that

c 2 − s 2 s c = a ℓ ℓ − a k k a k ℓ .

Set β to half of this quantity,

β = a ℓ ℓ − a k k 2 a k ℓ .

If a_kℓ is zero we can stop without performing an update, thus we never divide by zero. Let t be tan ϑ. Then with a few trigonometric identities we reduce the equation to

t 2 + 2 β t − 1 = 0.

For stability we choose the solution

t = sgn ⁡ ( β ) | β | + β 2 + 1 .

From this we may obtain c and s as

c = 1 t 2 + 1 s = c t

Although we now could use the algebraic update equations given previously, it may be preferable to rewrite them. Let

ρ = s 1 + c ,

so that ρ = tan(ϑ/2). Then the revised update equations are

a h k ′ = a k h ′ = a h k − s ( a h ℓ + ρ a h k ) a h ℓ ′ = a ℓ h ′ = a h ℓ + s ( a h k − ρ a h ℓ ) a k ℓ ′ = a ℓ k ′ = 0 a k k ′ = a k k − t a k ℓ a ℓ ℓ ′ = a ℓ ℓ + t a k ℓ

As previously remarked, we need never explicitly compute the rotation angle ϑ. In fact, we can reproduce the symmetric update determined by Q_kℓ by retaining only the three values k, ℓ, and t, with t set to zero for a null rotation.

References

Jacobi rotation Wikipedia

(Text) CC BY-SA