[SciPy-user] verify whether a matrix is positive definite or not
Anne Archibald
peridot.faceted@gmail....
Thu Jun 28 09:54:50 CDT 2007
On 28/06/07, Nils Wagner <nwagner@iam.uni-stuttgart.de> wrote:
> Hi all,
>
> I have a parameter-dependent matrix
>
> B(x) = B_0 + x B_1, 0 \le x \le 1
>
> where B_0 and B_1 are symmetric. How can I determine critical values x*
> (if any) such that B(x*) is not positive definite ?
>
>
> from scipy import *
>
> def B(x):
>
> return array(([[11.,8.],[8.,7.]])) - x*array(([[20.,1.],[1.,26]]))
>
> X = linspace(0,1,100)
>
> for x in X:
> print x
> L=linalg.cholesky(B(x),lower=1)
>
> I mean it would be nice if cholesky could return info=1 if the matrix is
> not spd.
> The current behaviour is
>
> Traceback (most recent call last):
> File "test_spd.py", line 11, in ?
> L=linalg.cholesky(B(x),lower=1)
> File "/usr/lib64/python2.4/site-packages/scipy/linalg/decomp.py", line
> 552, in cholesky
> if info>0: raise LinAlgError, "matrix not positive definite"
> numpy.linalg.linalg.LinAlgError: matrix not positive definite
>
> Helpful suggestions would be appreciated.
Well, you can always use try/except to catch the LinAlgError. It's
remotely possible that cholesky might fail to converge and throw a
different LinAlgError which you would want to re-raise.
You can also look at the eigenvalues - the matrix is positive definite
if and only if they're all positive. So making a function that takes a
parameter and returns the least eigenvalue should give you a
relatively smooth function to do root-finding on. With symmetric
matrices, eigenvalue finding ought to be fairly reliable.
For this particular case, note that the positive-definite matrices
form a cone, that is, the sum of two or any positive multiple of a
positive-definite matrix is also positive definite. In particular this
means it's convex. So if you're tracing the line between two
endpoints, as you are here (the endpoints are B_0 and B_0+B_1), you
can check the endpoints and know that the matrix is positive definite
between them if they're both positive definite. If one is positive
definite and the other isn't, then clearly there's some point in
between where they stop being positive definite. If neither is
positive definite, then it's possible that some positive definite
matrix lies between them (good luck finding it; you could try
numerically maximizing the least eigenvalue).
Anne M. Archibald
More information about the SciPy-user
mailing list