6
$\begingroup$

To compute the eigenvector corresponding to dominant eigenvalue of a symmetric matrix $A\in\mathbb{R}^{n\times n}$, one used Power Iteration, i.e., given some random initialization, $u_1\in\mathbb{R}^n$, one iteratively computes $$u_1\leftarrow Au_1,$$after which a normalization is applied to $u_1$. Now, suppose that eigenvectors $u_1, u_2$ are computed in advance, and one wants to compute the eigenvector $u_3$ associated with the third dominant eigenvalue.

In case the initial $u_3$ is orthogonal to both $u_1$ and $u_2$, can it be shown that the series $$u_3\leftarrow Au_3$$ converges in the direction of the eigenvector of $A$ corresponding to third dominant eigenvalue.

Note that the question is motivated by an observation that most implementations of the Power Iteration that computes higher eigenvectors provides Gram-Schmidt (GS) orthonormalization in each iteration, i.e., after each $u_3\leftarrow Au_3$, orthogonalization is applied wrt $v_1, v_2$. In case the orthonogonality of $u_3$ is imposed with its initialization (ie. from the start), is the GS necessary after each matrix-vector multiplication?

$\endgroup$
11
$\begingroup$

In exact arithmetic you shouldn't need to reorthogonalize regularly, but practically you do. Your u1 and u2 are close to (but not exactly) the true eigenvectors, so your initial deflation almost (but not entirely) removed the true eigenvectors from u3. The tiny components you left behind will be amplified by repeated multiplication by A, you will need to continue damping it out by gram-schmidt subtraction.

$\endgroup$
  • $\begingroup$ Good explanation. So, if the eigenvectors $u_1, u_2$ are true, and given in advance, the GS will not be necessary after each matrix-vector multiplication $Au_3$, in case $u_3$ is orthogonal to $u_1, u_2$. How can it be shown that the series $Au_3, A^2u_3, A^3u_3...$ is also orthogonal to $u_1, u_2$? $\endgroup$ – usero Oct 3 '12 at 13:08
  • 1
    $\begingroup$ To show that power series remains orthogonal, expand your initial u3 in terms of (all of) the eigenvectors of A. It will be deficient in u1 and u2, but will still have components along all the rest. Applying A repeatedly to that expansion will make it easy to see that if u3 is initially deficient in u1 and u2, it will remain so forever. But you still have the issue of inexact arithmetic, which will add some random error every time you apply A. So even though u3 is initially deficient in u1 and u2, multiplying A*u3 in inexact arithmetic will yield some miniscule but nonzero u1/u2 components. $\endgroup$ – rchilton1980 Oct 3 '12 at 13:25
  • $\begingroup$ Just one more thing: in case initial $u_3$ is orthogonal to $u_1, u_2$, could one perform multiplication series for $k$ steps, ie. $Au_3, A^2u_3, A^3u_3, ..., A^ku_3$, and then perform GS just at the and, ie. only after $u_3\leftarrow A^ku_3$. I suppose the 'error' would accumulate during the process, but will be removed after $k$ steps. However, I'm not sure: by repeated multiplication one might 'regain' component of e.g, $v_2$ and converge in different direction. Your opinion? $\endgroup$ – usero Oct 3 '12 at 17:34
  • $\begingroup$ I think that should work. But as you push closer to practical implementation, I also think the Lanczos or Arnoldi algorithms will be a better choice. If you have matlab, it tidily wraps up ARPACK's implementation of the Arnoldi method. (If not, ARPACK has a pretty liberal license and is not too bad to figure out). I'm not aware of a comparable OSS Lanczos implementation - Arnoldi/ARPACK is less efficient for your specific real+symmetric problem than Lanczos would be, but it will still be good enough. $\endgroup$ – rchilton1980 Oct 4 '12 at 13:24
  • $\begingroup$ to avoid opening another short question: given a general non-symmetric real matrix $A$, would GS orthonormalization lead to subsequent eigenvectors (eg, we know $u_1$, and want $u_2$)? In such case, the eigenvectors would not be orthonormal. On the other hand, I'm not sure if multiplication $Au_2$ would re-introduce component of $u1$; in case GS is performed, once more orthongonal vectors are obtained, but the eigenvectors $u1,u2$ might not be $\endgroup$ – usero Oct 10 '12 at 9:44
8
$\begingroup$

In theory, yes. In practice, rounding errors will usually result in (initially slow) convergence to $u_1$.

At essentially the same cost one can run the Lanczos algorithm, which will have much faster convergence, and produce the three dominant eigenvalues unless two of these eigenvalues are essentially the same. For Lanczos, selective reorthogonalization is sufficient to get good results.

$\endgroup$
5
$\begingroup$

Yes. If $A\in\mathbb{R}^{n\times n}$ is symmetric, there exists an orthogonal basis of eigenvectors $v_1,\dots,v_n$ (sorted by decreasing magnitude) of $A$. If you expand your starting vector $u_3$ with respect to this basis, you get $$ u_3 = \alpha_1 v_1 + \alpha_2 v_2 + \dots + \alpha_n v_n,$$ where $\alpha_1 = \alpha_2 =0$ because $u_3$ is orthogonal to $v_1$ and $v_2$. If $u_3$ is not orthogonal to $v_3$, you can now apply the usual argument for the convergence of the power method (with $\lambda_3$ and $v_3$ in place of $\lambda_1$ and $v_1$) to show convergence to $v_3$.

This is true in exact arithmetic only; in floating point arithmetic, rounding errors will lead to loss of orthogonality, and orthogonalizing $u^k := A^k u_3$ against $v_1$ and $v_2$ is required (although it is not necessary to apply this after every single iteration).

$\endgroup$
  • $\begingroup$ Thanks. However, the problem occurs if the initialization to $u_3$ is also orthogonal to true eigenvector $u_3$. One more point: in case of a random initialization for the top 3 eigenvectors, and the initialization $u_3$ being linearly dependent on $u_1, u_2$, the power method will result in $0_n$ in place of $u_3$? In such case, one should check the final result for zero vectors. (or?) $\endgroup$ – usero Oct 3 '12 at 13:12
  • 1
    $\begingroup$ Of course, the starting vector needs to contain some contribution from the dominant eigenvector (i.e., have nonzero projection on its span), since multiplication by $A$ can only amplify, but never generate them. If you start with a linear combination $u$ of the first two eigenvectors $v_1$ and $v_2$ and orthogonalize $u$ against $v_1$ and $v_2$, you of course get the zero vector, and the power method will not do anything. $\endgroup$ – Christian Clason Oct 3 '12 at 13:35
  • $\begingroup$ By the way, it's highly improbable for a random vector to be linear combination of $v_1$ and $v_2$ only - this would be equivalent to a random point in $\mathbb{R}^n$ lying in the $(x_1,x_2)$-plane. (You can see this for $n=3$ by drawing a picture, and it will only get more unlikely as $n$ becomes larger.) $\endgroup$ – Christian Clason Oct 3 '12 at 13:38
  • $\begingroup$ A question related to the last paragraph of your answer: the main purpose of repeated orthonormalization is to prevent 're-gain' of components $\alpha_1v_1$, $\alpha_2v2$; in case $\alpha_1\neq 0$ appears after a certain application of $A$, the convergence is in direction of $v_1$ (or?). This would mean that the orthogonalization needs to be frequent; however, the case of only initial orthogonalization of $u_3$ to $u_1, u_2$, and one more orthogonalization deferred to the end (to the final vector) may result in zero vector (as the final might have converged to $u_1$ or $u_2$). Correct? $\endgroup$ – usero Oct 4 '12 at 9:31
  • $\begingroup$ That is correct, the orthogonalization needs to be frequent enough to avoid $\alpha_3$ becoming small enough to be masked by rounding errors. If you need to orthogonalize against many vectors (and multiplication with $A$ is fast), doing this only every second or third iteration is already a significant speedup. (Of course, as Arnold Neumaier pointed out, deflation might not be the best option for computing the eigenvector.) Since there are only two vectors in your case, I'd just orthogonalize after every iteration. $\endgroup$ – Christian Clason Oct 4 '12 at 13:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.