2
votes

In Matlab I have a real and symmetric n x n matrix A, where n > 6000. Even though A is positive definite it is close to singular. A goes from being positive definite to singular to indefinite for a particular variable which is changed. I am to determine when A becomes singular. I don't trust the determinants so I am looking at the eigenvalues, but I don't have the memory (or time) to calculate all n eigenvalues, and I am only interested in the smallest - and in particular when it changes sign from positive to negative. I've tried

D = eigs(A,1,'smallestabs')

by which I lose the sign of the eigenvalue, and by

D = eigs(A,1,'smallestreal')

Matlab cannot get the lowest eigenvalue to converge. Then I've tried defining a shift value like

for i = 1:10 
   if i == 1
      D(i) = eigs(A,1,0) 
   else
      D(i) = eigs(A,1,D(i-1))
   end
end

where i look in the range of the last lowest eigenvalue. However, the eigenvalues seem to behave oddly, and I am not sure if I actually find the true lowest one.

So, any ideas on how to

  • without doubt find the smallest eigenvalue with 'eigs', or
  • by another way determine when A becomes singular (when changing a variable in A)

are highly appreciated!

1
Shouldn't you use the condition number, rather than the smallest eigenvalue, to assess whether the matrix is close to singular? See condLuis Mendo
Well, the point is that I know it is close to singular, but need to establish a criterion for when it shifts from positive definite over singular to indefinite. Can I do that by the condition number?DavidH
Sure. It’s the ratio between largest and smallest eigenvalues. So it’s better than just observing the smallest eigenvalueLuis Mendo
Well, as the matrix becomes singular, the smallest eigenvalue approaches zero, and the condition number explodes. Since the condition number explodes, the numerical niceness of the matrix turns into outright numerical rudeness, and eigenvalue methods may suffer depending on the spectrum of A. Not knowing any more about A limits the advice that I think can be given. But I agree with @Luis, that cond or rcond may be more robust paths.TroyHaskin
I'm not quite sure about the condition number and 'cond' command. Even though I'm in a range where 'eigs' gives me negative eigenvalues, I only get a positive condition number...?DavidH

1 Answers

1
votes

Solution

I seem to have solved my particular problem. Matlabs command chol have the possibility to return a value p which is zero if the matrix is positive definite. Thus, performing

[~,p] = chol(A)

in my case determines the transition from positive definite to not positive definite (meaning first singular then indefinite), and is also computationally very efficient. In the documentation for chol it is also preferred over eigs to check for positive definiteness. However, there seem to be some confusion about the result if the matrix is only positive semi-definite, so be carefull if this is the case.


Alternative solutions

I've come across several possible solutions which I'd like to state:

Determinant: For a positive definite matrix the determinant is positive. However, for an indefinite matrix it may be negative - this could indicate the transition. Though, generally determinants for large nearly singular matrices are not recommended.

Eigenvalues: For a positive definite matrix the real part of all eigenvalues are positive. If at least one eigenvalue is zero the matrix is singular, and if one becomes negative and the rest is positive it is indefinite. Detecting the shift in sign for the lowest eigenvalue indicates the point the matrix becomes singular. In matlab the lowest eigenvalue may be found by

D = eigs(A,1,'smallestreal')

However, in my case Matlab coudn't perform this. Alternatively you can try searching around zero:

D = eigs(A,1,0)

This however only finds the eigenvalue closest to zero. Even if you make a loop like indicated in my original question above, you are not guaranteed to actually find the lowest. And the precision of the eigenvalues for a nearly singular matrix seems to be low in some cases.

Condition number: Matlabs cond returns the condition number of the matrix by performing

C = cond(A)

which states the ratio of the largest eigenvalue to the lowest. A shift in sign in the condition number thereby states the transition. This, however, didn't work for me, as I only got positive condition numbers even though I had negative eigenvalues. But maybe it will work in other cases.