Take the 2-minute tour ×
MathOverflow is a question and answer site for professional mathematicians. It's 100% free, no registration required.

I'm wondering about which bounds one can put on the condition number of a $n\times n$ square matrix which is obtained from another $n\times n$ square matrix by orthogonalizing the first $m < n$ columns. Obviously, the condition number will become 1 if $m = n$ (and hence minimal). I am wondering whether the condition number is always improving, even if only a few columns are orthogonalized. If this is not the case, I would also be happy to know whether one can put a bound on how much the condition number can become worse.

To be more specific, I have a $n \times n$ matrix \[A = (A_0\, |\, A_1)\] with $A_0 = Q\,R$ with $Q$ a $n \times m$ unitary matrix and $R$ an $m\times m$ upper triangular matrix, $m\leq n$. The columns of $A_0$ can be assumed to be normalized to 1, $||A(:,i)||_p = 1$ for $i = 1 \dots n$. (If it helps, one can also assume the columns of $A_1$ to be normalized to unity).

I now define \[A' = (Q\,|\,A_1).\] The condition number in the $p$-norm is defined as $\kappa_p(A) = ||A||_p \,||A^{-1}||_p$.

Is it possible to give an upper bound on $\kappa_p(A')$ compared to $\kappa_p(A)$? (for example, does $\kappa_p(A') \leq \kappa_p(A)$ hold? [probably not, see below])

Background

The problem arises in the numerical solution of a quantum mechanical scattering problem where we have to solve an ordinary linear system. The columns of $A_0$ are given by the eigenvectors of a different problem that describes modes that are coupled out and in of the system; the matrix $A_1$ corresponds to the scatterer and is in principle unrelated to $A_0$. (There is a lot of additional structure in the problem, but that is probably too special for mathoverflow). It turns out that sometimes the eigenvectors forming $A_0$ are numerically almost linearly dependent (they are eigenvectors of a general eigenproblem, not Hermitian).

We were able to reformulate the physical question such that instead of putting in the eigenvectors explicitely in $A_0$ we could use an orthogonal basis spanning the space of eigenvectors, which we can easily find using the Schur form of the eigenproblem. This allows us to write $A_0 = Q R$ as described above (and we actually do not need to worry about accuracy in this step).

Now, if the whole matrix $A$ is ill-conditioned because of nearly linearly dependen vectors in $A_0$, going to $A'$ will certainly improve the condition number drastically. I am wondering if this step however can also do harm, if we always apply it, regardless of the condition of $A_0$.

Attempts of mine so far

Numerical tests on random, ill-conditioned matrices (constructed via SVD) have shown that if I orthogonalize a few columns, the condition number does not get much better, but also not much worse. However, I saw examples where it did get somewhat worse (10-20% larger condition number), hence $\kappa_p(A') \leq \kappa_p(A)$ is probably not true.

With $A = A' \begin{pmatrix}R &0\\\\0&1\end{pmatrix}$, and since $\kappa_2(\begin{pmatrix}R &0\\\\0&1\end{pmatrix})=\kappa_2(R)$ (the largest singular value of R is $||R||_2 \geq \max_i ||R(:,i)||_2=1$) I can find the bounds $\kappa_2(A)/\kappa_2(R) \leq \kappa_2(A') \leq \kappa_2(A) \kappa_2(R)$. This tells me that if $A_0$ is well-conditioned, the condition number of $A'$ is similar to $A$, but the bounds are way to loose if $A_0$ has nearly linearly dependent columns (where the orthogonalization helps most. Not surprising, as I wasn't able to use that fact that $Q$ is orthogonal.)

share|improve this question
    
I'm a bit confused about your numerical tests. Haven't you found that $\kappa_p(A') \le \kappa_p(A)$ is definitely not true? Or are you simply leaving open the possibility that the numerics are highly inaccurate? –  Mark Meckes Jul 6 '13 at 10:53
    
I chanced across this and would suggest you give a try at asking over at scicomp: scicomp.stackexchange.com ; techniques of this nature are very common in scientific computation. My intuition on the matter would say that you need to know far more about the matrix before being able to quantify a reduction in condition number (I suspect normality of the matrix would play an important role in the determination of a concrete condition number reduction) –  Reid Atcheson Jul 8 '13 at 5:13
add comment

2 Answers

Though you can obtain specific bounds based on the entries of the matrix, there is no reason to expect any nice behavior in the condition number, since condition number is a global property, i.e., depends on the entire matrix, whereas the QR algorithm (or any other standard decomposition algorithms for that matter) are typically local, i.e., they choose a row or column and then extract the row or column updating the rest. For instance, consider the following matrix $$\begin{bmatrix}10^{100} & 0\\ 0 & 10^{100}\end{bmatrix}$$The condition number of this matrix, based on two norm, is $1$. After one step of QR, we get the matrix $$\begin{bmatrix}1 & 0\\ 0 & 10^{100}\end{bmatrix}$$whose condition number, based on two norm, is $10^{100}$.

share|improve this answer
add comment

This method cannot work because the condition number $cond$ (here associated to the $||.||_2$) reflects only the maximum trouble, in the following sense: for the sake of simplicity, assume that the eigenvalues $(\lambda_i)_i$ of $A$ are s.t. $\max|\lambda_i|/\inf|\lambda_i|$ is small. Thus, if $cond(A)$ is great, then there are at least $2$ eigenvectors that are ill conditionned. For instance let $D=diag(-1,1,-2,2)$ with $cond(D)=2$ and $P=diag(1,1,U)$ where $U=\begin{pmatrix}1&1\\1+\epsilon&1\end{pmatrix}$ with $\epsilon=10^{-5}$. Then $cond(p)=4*10^5$ and $cond(P^{-1}DP)=1.6*10^{11}$. Moreover, if $Q=diag(U,U)$ (2 troubles), then $cond(Q^{-1}DQ)=3.2*10^{11}$ that has the same order than the first one !

Conclusion. You can work as follows: find two eigenvectors of $A$ which are close and keep away them from the other (by a change of basis, turn them in orthogonal vectors). Then repeat the previous operation with the obtained matrix.

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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