External Issue: Cholesky decomposition returns upper/lower triangular part of a matrix instead of the decomposition

Notifications Issues
1

16487, “CaptainSharf”, “Cholesky decomposition returns upper/lower triangular part of a matrix instead of the decomposition”, “2020-09-27T17:55:37Z”

Summary of Problem

Cholesky decomposition returns the lower-triangular part of a given matrix instead of decomposing the matrix such that $$L^TL=A$$. One also needs to check if the matrix $$A$$ is symmetric positive definite before decomposing it into $$L,L^T$$. In the current implementation the return value of LAPACK.potrf is ignored.

Steps to Reproduce

Source Code:

use LinearAlgebra;
var dom = {1..3,1..3};
var arr:[dom] real;
for (i,j) in dom do
    arr[i,j] = i;
for i in 1..3 do
    arr[1,i] = 0;
var matr = Matrix(arr);
var x = cholesky(matr);
writeln(x);

Compile command:
chpl -L/usr/local/opt/lapack/lib -I/usr/local/opt/lapack/include -llapacke cholesky.chpl

Execution command:
./cholesky

Output

Matrix is 
0.0 0.0 0.0
2.0 2.0 2.0
3.0 3.0 3.0
L is 
0.0 0.0 0.0
2.0 2.0 0.0
3.0 3.0 3.0

The cholesky decomposition of the matrix does not exist as it is singular

Configuration Information

  • chpl --version: chpl version 1.23.0