$title Test cholesky utility (CHOLES02,SEQ=412)

$onText

Test factor returned by the cholesky utility

Contributor: Steve Dirkse

$offText

set I / i1 * i3 /;
alias (I,J,K);

parameters
  L0(I,J)  'lower triangular matrix'
  A(I,J)   'L * L-transpose'
  L1(I,J)  'computed cholesky factor'
  err(I,J) 'L0 - L1'
  ;

table L0
        i1      i2      i3
i1      4
i2      2       1
i3      1       4       2   ;

A(I,J) = sum {K, L0(I,K)*L0(J,K)};

execute_unload 'a.gdx', i, a;
executeTool.checkErrorLevel 'linalg.cholesky i a L1 -gdxin=a.gdx -gdxout=b.gdx';
execute_load 'b.gdx', L1;

err(I,J) = L0(I,J) - L1(I,J);
scalar sumerr; sumerr = sum{(I,J), abs(err(I,J))};
abort$[sumerr > 1e-8] 'bad Cholesky factor returned: expected L0, got L1', L0, L1, err;

option clear=L1;
executeTool.checkErrorLevel 'linalg.cholesky i a L1';

err(I,J) = L0(I,J) - L1(I,J);
sumerr = sum{(I,J), abs(err(I,J))};
abort$[sumerr > 1e-8] 'bad Cholesky factor returned: expected L0, got L1', L0, L1, err;