Divide-and-conquer eigenvalue algorithm: Difference between revisions

Revision as of 19:03, 4 May 2011

Latest revision as of 12:10, 24 June 2024

(17 intermediate revisions by 12 users not shown)

Line 1:

{{Short description|Algorithm on Hermitian matrices}}

{{Multiple issues|

{{No footnotes|date=May 2024}}

{{More citations needed|date=May 2024}}

}}

'''Divide-and-conquer eigenvalue algorithms''' are a class of [[eigenvalue algorithm]]s for [[Hermitian matrix|Hermitian]] or [[real number|real]] [[Symmetric matrix|symmetric matrices]] that have recently (circa 1990s) become competitive in terms of [[Numerical stability|stability]] and [[Computational complexity theory|efficiency]] with more traditional algorithms such as the [[QR algorithm]]. The basic concept behind these algorithms is the [[Divide and conquer algorithm|divide-and-conquer]] approach from [[computer science]]. An [[eigenvalue]] problem is divided into two problems of roughly half the size, each of these are solved [[Recursion|recursively]], and the eigenvalues of the original problem are computed from the results of these smaller problems.

HereThis wearticle presentcovers the simplestbasic versionidea of a divide-and-conquerthe algorithm, similar to the oneas originally proposed by Cuppen in 1981., which Manyis detailsnot thatnumerically lie outside the scope of this article will be omitted; however,stable without consideringadditional these details, the algorithm is not fully stablerefinements.

==Background==

As with most eigenvalue algorithms for Hermitian matrices, divide-and-conquer begins with a reduction to [[Tridiagonal matrix|tridiagonal]] form. For an <math>m \times m</math> matrix, the standard method for this, via [[Householder reflection]]s, takes <math>\frac{4}{3}m^{3}</math> [[flops]]floating point operations, or <math>\frac{8}{3}m^{3}</math> if [[eigenvector]]s are needed as well. There are other algorithms, such as the [[Arnoldi iteration]], which may do better for certain classes of matrices; we will not consider this further here.

In certain cases, it is possible to ''deflate'' an eigenvalue problem into smaller problems. Consider a [[block diagonal matrix]]

Line 10 ⟶ 16:

The eigenvalues and eigenvectors of <math>T</math> are simply those of <math>T_{1}</math> and <math>T_{2}</math>, and it will almost always be faster to solve these two smaller problems than to solve the original problem all at once. This technique can be used to improve the efficiency of many eigenvalue algorithms, but it has special significance to divide-and-conquer.

For the rest of this article, we will assume the input to the divide-and-conquer algorithm is an <math>m \times m</math> real symmetric tridiagonal matrix <math>T</math>. Although theThe algorithm can be modified for Hermitian matrices, we do not give the details here.

==Divide==

Line 18 ⟶ 24:

:[[Image:Almost block diagonal.png]]

The size of submatrix <math>T_{1}</math> we will call <math>n \times n</math>, and then <math>T_{2}</math> is <math>(m - n) \times (m - n)</math>. Note that the remark about <math>T</math> beingis almost block diagonal is true regardless of how <math>n</math> is chosen (i.e., there are many ways to so decompose the matrix). However, it makes sense, from anFor efficiency standpoint,we totypically choose <math>n \approx m/2</math>.

We write <math>T</math> as a block diagonal matrix, plus a [[Rank (linear algebra)|rank-1]] correction:

Line 35 ⟶ 41:

:<math>T = \begin{bmatrix} Q_{1} & \\ & Q_{2} \end{bmatrix} \left( \begin{bmatrix} D_{1} & \\ & D_{2} \end{bmatrix} + \beta z z^{T} \right) \begin{bmatrix} Q_{1}^{T} & \\ & Q_{2}^{T} \end{bmatrix}</math>

The remaining task has been reduced to finding the eigenvalues of a diagonal matrix plus a rank-one correction. Before showing how to do this, let us simplify the notation. We are looking for the eigenvalues of the matrix <math>D + w w^{T}</math>, where <math>D</math> is diagonal with distinct entries and <math>w</math> is any vector with nonzero entries. In this case <math>w = \sqrt{|\beta|}\cdot z</math>.

IfThe case of a zero entry is simple, since if w_i is zero, (<math>e_i</math>,d_i) is an eigenpair (<math>e_i</math> is in the standard basis) of <math>D + w w^{T}</math> since

<math>(D + w w^{T})e_i = De_i = d_i e_i</math>.

Line 52 ⟶ 58:

This equation is known as the ''secular equation''. The problem has therefore been reduced to finding the roots of the [[rational function]] defined by the left-hand side of this equation.

All general eigenvalue algorithms must be iterative,{{Citation needed|date=April 2024}} and the divide-and-conquer algorithm is no different. Solving the [[nonlinear]] secular equation requires an iterative technique, such as the [[Newton's method|Newton–Raphson method]]. However, each root can be found in [[Big O notation|O]](1) iterations, each of which requires <math>\Theta(m)</math> flops (for an <math>m</math>-degree rational function), making the cost of the iterative part of this algorithm <math>\Theta(m^{2})</math>.

==Analysis==

As is common for divide and conquer algorithms, weW will use the [[Master theorem (analysis of algorithms)|master theorem for divide-and-conquer recurrences]] to analyze the running time. Remember that above we stated we choose <math>n \approx m/2</math>. We can write the [[recurrence relation]]:

:<math>T(m) = 2 \times T\left(\frac{m}{2}\right) + \Theta(m^{2})</math>

In the notation of the Master theorem, <math>a = b = 2</math> and thus <math>\log_{b} a = 1</math>. Clearly, <math>\Theta(m^{2}) = \Omega(m^{1})</math>, so we have

:<math>T(m) = \Theta(m^{2})</math>

Remember that aboveAbove, we pointed out that reducing a Hermitian matrix to tridiagonal form takes <math>\frac{4}{3}m^{3}</math> flops. This dwarfs the running time of the divide-and-conquer part, and at this point it is not clear what advantage the divide-and-conquer algorithm offers over the QR algorithm (which also takes <math>\Theta(m^{2})</math> flops for tridiagonal matrices).

The advantage of divide-and-conquer comes when eigenvectors are needed as well. If this is the case, reduction to tridiagonal form takes <math>\frac{8}{3}m^{3}</math>, but the second part of the algorithm takes <math>\Theta(m^{3})</math> as well. For the QR algorithm with a reasonable target precision, this is <math>\approx 6 m^{3}</math>, whereas for divide-and-conquer it is <math>\approx \frac{4}{3}m^{3}</math>. The reason for this improvement is that in divide-and-conquer, the <math>\Theta(m^{3})</math> part of the algorithm (multiplying <math>Q</math> matrices) is separate from the iteration, whereas in QR, this must occur in every iterative step. Adding the <math>\frac{8}{3}m^{3}</math> flops for the reduction, the total improvement is from <math>\approx 9 m^{3}</math> to <math>\approx 4 m^{3}</math> flops.

Line 69 ⟶ 75:

==Variants and implementation==

The algorithm presented here is the simplest version. In many practical implementations, more complicated rank-1 corrections are used to guarantee stability; some variants even use rank-2 corrections. ^{[citation{{Citation needed]}|date=September 2011}}

There exist specialized root-finding techniques for rational functions that may do better than the Newton-Raphson method in terms of both performance and stability. These can be used to improve the iterative part of the divide-and-conquer algorithm.

The divide-and-conquer algorithm is readily [[Parallel algorithm|parallelized]], and [[linear algebra]] computing packages such as [[LAPACK]] contain high-quality parallel implementations.{{Citation needed|date=September 2023}}

==References==

Line 84 ⟶ 90:

| title = Applied Numerical Linear Algebra

| year = 1997}}.

* {{cite journal |first1=J.J.M. |last1=Cuppen |title=A Divide and Conquer Method for the Symmetric Tridiagonal Eigenproblem |journal=[[Numerische Mathematik]] |volume=36 |pages=177–195 |date=1981 |issue=2 |doi=10.1007/BF01396757 |s2cid=120504744 }}

{{Numerical linear algebra}}

[[Category:Numerical linear algebra]]

[[Category:Divide-and-conquer algorithms]]