[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
chol(a)
From: |
Salvatore Filippone |
Subject: |
chol(a) |
Date: |
Tue, 15 Nov 2011 11:38:41 +0100 |
Hello there,
While checking exercises for some students, I came across a rather
surprising (for me) behaviour.
Consider the following SPD matrix:
octave:1> a=[12,2,3,4;
> 2,14,5,3;
> 3,5,16,6;
> 4,3,6,16]
a =
12 2 3 4
2 14 5 3
3 5 16 6
4 3 6 16
If I perform the Cholesky factorization on the lower or upper part, I
get the expected results;
octave:2> chol(a)
ans =
3.46410 0.57735 0.86603 1.15470
0.00000 3.69685 1.21725 0.63117
0.00000 0.00000 3.71057 1.14045
0.00000 0.00000 0.00000 3.60107
octave:4> chol(a,'lower')
ans =
3.46410 0.00000 0.00000 0.00000
0.57735 3.69685 0.00000 0.00000
0.86603 1.21725 3.71057 0.00000
1.15470 0.63117 1.14045 3.60107
Now, it is my understanding that internally Octave uses LAPACK; if I
ask for the factorization of an upper triangle in LAPACK, the code
will only look at the upper part of the input.
This is consistent with what I get if I pass triu(a) as in
> chol(triu(a))
ans =
3.46410 0.57735 0.86603 1.15470
0.00000 3.69685 1.21725 0.63117
0.00000 0.00000 3.71057 1.14045
0.00000 0.00000 0.00000 3.60107
In LAPACK, I can do a perfectly symmetric choice by passing only the
lower triangle, but Octave in this case seems to only factor the
diagonal:
octave:5> chol(tril(a),'lower')
ans =
3.46410 0.00000 0.00000 0.00000
0.00000 3.74166 0.00000 0.00000
0.00000 0.00000 4.00000 0.00000
0.00000 0.00000 0.00000 4.00000
Having been an LAPACK user for a long time this is a bit surprising.
Note that "that other package" (i.e. Matlab) does what I expect it to
do, i.e. with choll(tril(a),'lower') produces the transpose of the
output from chol(triu(a))
Bug or feature? If the latter, then maybe the documentation ought to
be more specific.
Thanks a lot
Salvatore
- chol(a),
Salvatore Filippone <=