I'm trying to do an in place Cholesky factorization, but either that feature is not actually implemented in scipy or there is something I am not understanding. I am posting here in case it is the latter. Here is a simplified example of what I'm doing:
import numpy
import scipy.linalg
numpy.random.seed(0)
X = numpy.random.normal(size=(10,4))
V = numpy.dot(X.transpose(),X)
R = V.copy()
scipy.linalg.cholesky(R,False,overwrite_a=True)
print V
print R
I think what should happen is for R to be overwritten with an upper triangular matrix. However, when I run this code, my V and R come out identical (cholesky is not overwriting R). Am I misunderstanding the purpose of overwrite_a, or making some other mistake? If this is just a bug in scipy, is there a workaround or another way I can do an in place Cholesky factorization in Python?
Try it again, with
>>> scipy.__version__
'0.11.0'
>>> np.__version__
'1.6.2'
it works perfectly:
X = np.random.normal(size=(10,4))
V = np.dot(X.transpose(),X)
print V
R = V.copy()
print R
C = scipy.linalg.cholesky(R,False,overwrite_a=True)
print V
print R
Output is:
[[ 11.22274635 5.10611692 0.70263037 3.14603088] # V before
[ 5.10611692 8.94518939 -3.17865941 1.64689675]
[ 0.70263037 -3.17865941 7.35385131 -2.23948391]
[ 3.14603088 1.64689675 -2.23948391 8.25112653]]
[[ 11.22274635 5.10611692 0.70263037 3.14603088] # R before
[ 5.10611692 8.94518939 -3.17865941 1.64689675]
[ 0.70263037 -3.17865941 7.35385131 -2.23948391]
[ 3.14603088 1.64689675 -2.23948391 8.25112653]]
[[ 11.22274635 5.10611692 0.70263037 3.14603088] # V after
[ 5.10611692 8.94518939 -3.17865941 1.64689675]
[ 0.70263037 -3.17865941 7.35385131 -2.23948391]
[ 3.14603088 1.64689675 -2.23948391 8.25112653]]
[[ 3.35003677 1.52419728 0.20973811 0.93910339] # R after
[ 0. 2.57332704 -1.35946252 0.08375069]
[ 0. 0. 2.33703292 -0.99382158]
[ 0. 0. 0. 2.52478036]]
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With