# Condition Number dependent algorithms for matrix operations

Using the Conjugate gradient method we can solve a linear system $$Ax=b$$, where $$A\in\mathbb R^{n\times n}$$ in time $$O(n^2 \sqrt{\kappa})$$, where $$\kappa=\frac{\sigma_\mathrm{max}(A)}{\sigma_\mathrm{min}(A)}$$ is the condition number. For well conditioned matrices this is better than inverting $$A$$ in $$O(n^\omega)$$ time (as well as various matrix factorization algorithms.)

My intuition is that well conditioned matrices are in some sense trivial, and this should be a benefit, not just for linear solving, but also other matrix related problems. In particular I'm interested in whether there are condition number dependent algorithms (like $$O(n^{2} f(\kappa))$$ time) for any of the following problems:

For the first, of course we can solve $$n$$ linear systems, but that's $$n^3\sqrt{\kappa}$$ time.

The Lyapunov matrix equation (or Sylvester Equation) $$A\Sigma + \Sigma A^T = R$$ can be solved by the alternating-direction implicit (ADI) method. But I'm not sure how fast that actually is. In this paper (Theorem 6 and 7) they seem to suggest a condition number related bound, but unless I misunderstand something, they need roughly $$n-\varepsilon/\kappa$$ iterations (each $$n^2$$ time). Not really better than $$n^3$$.

For matrix determinant it seems from this mathoverflow answer that $$n^\omega$$ is not even possible. So maybe it's unrealistic to beat $$n^3$$ here, even though some structured matrices can be solved faster.

• BTW, answer here gives a way to bound convergence of CG linear solver when $A$ matrix eigenvalues are contained in two clusters of width $\epsilon$. High condition number is a special case of this analysis with eigenvalues contained in just one cluster of small width Nov 3 at 13:35

This was meant as a comment, but too long.

If we fix matrix norm to be 1, we can ask - which algorithms terminate faster when smallest eigenvalue is far from zero?

We can rule out algebraic algorithms like in Ryan Williams answer -- those perform the same number of steps regardless of numeric values. Instead, only consider iterative algorithms which continue until certain numeric error criterion is met.

My feeling is that many such algorithms terminate faster for well-conditioned matrices.

Here's a sketch of how to show this for linear equation solving. The trick is to reformulate gradient descent as operating on the space of eigenvalues rather than dimensions.

The the problem of solving $$Ax=b$$ using gradient descent with learning rate $$1$$ on least squares objective. Since gradient descent is rotationally symmetric we can take $$H=A^T A$$ to be diagonal WLOG, with diagonal entries $$\{h_1,\ldots,h_d\}$$. Suppose we care about relative loss rather than absolute, hence can take $$\|A\|=1$$. If the initial loss is bounded, we have a bound on $$(x^*_i-x^0_i)^2$$ for solution $$x^*$$ and initial guess $$x^0$$. For simplicity, assume this quantity is exactly $$1$$.

After $$t$$ steps, our least squares loss is the following

$$f(t)=\sum_{i=1}^d h_i (1-h_i)^{2t}$$

We can write this sum down as integral over $$\mathrm{d}h$$, the empirical density of $$h_i$$

$$f(t)=\int \mathrm{d}h\ h (1-h)^{2t}$$

You can view this as a convolution with a kernel that is moving to the left as $$t$$ increases. Fixing range of the spectrum fixes the range of $$t$$ for which the loss is changing. For instance, visualize this for $$h=\{1/200,1/20,1/2\}$$ Since we have norm $$1$$, upper bound on condition gives a lower bound on smallest eigenvalue $$\lambda_d$$. The "gradient descent kernel" is centered at $$\frac{1}{2t}$$, hence one would expect for the error to decay exponentially with $$t$$ after $$t>\frac{1}{2\lambda_d}$$ steps, hence needing $$O\left(\frac{\log \epsilon }{2\lambda_d}\right)$$ steps for accurate approximation

Regarding Lyapunov's equation $$AX+XA'=B$$ with Hermitian $$A$$ associated with a matrix of eigenvectors $$U$$ and column vector of eigenvalues $$s$$, you can write solution as follows (proof):

$$X=U \left( \frac{U' BU}{s + s'} \right) U'$$

If our application doesn't care about choice of basis, we would just need to compute the central part

$$\frac{U' BU}{s + s'}$$

Now the question would be -- is this cheap to approximate for well-conditioned $$A$$?