I am trying to understand a limitation in a routine in the interval arithmetic software Intlab. From matrices starting from a given size (in my particular problems), a linear system appearing in an interative process simply fills the whole memory in Matlab.

The problem is as follows:

  • take $A$, $B$ two matrices (symmetric positive definite in my case) and $\tilde \lambda, \tilde x$ an approximate eigenpair: $A\tilde x \approx \tilde \lambda B \tilde x$ (found using eigs in Matlab).

  • in intlab the routine verifyeig for finding an interval enclosure for the eigenvalue and the eigenvector has as first step the following linear system (following the paper "Verification methods: Rigorous results using floating-point arithmetic" by S. Rump)

-- $R = A-\tilde \lambda B$

-- $v$ is the index where $\tilde x$ has the largest magnitude

-- The column $v$ in $R$ is replaced by $-B\tilde x$.

-- The first iteration in verifyeig solves the linear system

$$ Ry = (A\tilde x-\tilde \lambda B\tilde x)$$

While everything works fine for small matrices using just the backslash, for larger matrices the resolution of the system fails and fills the whole memory.

Having recently observed that in some particular cases, simple implementations are way more efficient than the backslash (exploiting the structure), I wonder if there is a way to improve the resolution of this system.

It seems a particular system, since its matrix $R$ is obtained from a sparse matrix $A-\lambda B$ (almost singular, since $\lambda$ is almost an eigenvalue) plus a column matrix $-Bx$ where $x$ is an "almost eigenvector".

Do you have suggestions for algorithms that might work better than default ones for this kind of problems? Thank you in advance.

Seems like throwing some appropriate sledgehammer on the problem will do the trick. tfqmr in Matlab seems work fine. I will also take a look at some preconditioning ideas.

Nevertheless, if you have some more refined idea that can be tried I am still interested.

  • $\begingroup$ You can try generic iterative solvers. I would have suggested conjugate gradients but you mention you substitute the $v$-th column which can make $R$ non-symmetric. If you can rewrite it as an SPD you can use CG. Otherwise try using some methods that work for non-symmetric matrices, e.g. the biconjugate gradient stabilized method. $\endgroup$
    – lightxbulb
    Commented Jun 10 at 8:34


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.