Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Cholesky/isHermitian not recognizing valid input #32202

Closed
KevinMDuarte opened this issue May 31, 2019 · 1 comment
Closed

Cholesky/isHermitian not recognizing valid input #32202

KevinMDuarte opened this issue May 31, 2019 · 1 comment

Comments

@KevinMDuarte
Copy link

When trying to compute the cholesky factorization of A'BA , Julia produces a PosDefException: matrix is not Hermitian; I suppose the reason for this error is caused by a numerical roundness issue obtained while computing the product of these matrices.

using LinearAlgebra
A = rand(n,n);
C = rand(n,n);
B = C'*C; #B is symmetric a
AtBA = A'*B*A;
# NOT WORKING : 
cholesky(AtBA) 
# WORKING : 
LinearAlgebra.LAPACK.potrf!('U',AtBA)

I believe this is a similar issue as the one mentioned in issue #28885 where it is said that this "feature" is desired to be kept like this. However, it would be nice to have this inaccuracy caused issue mentioned in the documentation for the cholesky/isposdef/ishermitian calls since newcomers to the Julia language might be confused why their theoretically posdef matrix cannot be factorized using cholesky in Julia unless using other tricks such as using the potrf call from the LAPACK library.

versioninfo()

Julia Version 1.1.1
Commit 55e36cc (2019-05-16 04:10 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin15.6.0)
  CPU: Intel(R) Core(TM) i5-4250U CPU @ 1.30GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, haswell)
@fredrikekre
Copy link
Member

If you know your matrix is analytically Hermitian you should use the Hermitian wrapper, and call cholesky as

using LinearAlgebra
A = rand(n,n)
C = rand(n,n)
B = C'*C
AtBA = A'*B*A
cholesky(Hermitian(AtBA))

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants