<![CDATA[
#include <numeric.h>
int main() {
int m = 3;
array double complex a[3][3] =
{ complex(2,0), complex(0,-1), complex(0,0),
complex(0,1), complex(2,0), complex(0,0),
complex(0,0), complex(0,0), complex(3,0)};
int status;
array double complex l[m][m];
status = choldecomp(a, l);
if (status == 0) {
printf("upper triangle l =\n%5.3f\n", l);
printf("transpose(l)*l - a =\n%5.3f\n", conj(transpose(l))*l-a);
}
else
printf("error: Matrix must be positive definite.\n"
"The %d order of matrix is not positive definite.\n",status);
status = choldecomp(a, l, 'l');
if (status == 0) {
printf("lower triangle l =\n%5.3f\n", l);
printf("l*transpose(l) - a =\n%5.3f\n", l*conj(transpose(l))-a);
}
else
printf("error: Matrix must be positive definite.\n"
"The %d order of matrix is not positive definite.\n",status);
}
]]>
<![CDATA[
upper triangle l =
complex(1.414,0.000) complex(0.000,-0.707) complex(0.000,0.000)
complex(0.000,0.000) complex(1.225,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(0.000,0.000) complex(1.732,0.000)
transpose(l)*l - a =
complex(0.000,0.000) complex(0.000,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(-0.000,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(0.000,0.000) complex(-0.000,0.000)
lower triangle l =
complex(1.414,0.000) complex(0.000,0.000) complex(0.000,0.000)
complex(0.000,0.707) complex(1.225,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(0.000,0.000) complex(1.732,0.000)
l*transpose(l) - a =
complex(0.000,0.000) complex(0.000,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(-0.000,0.000) complex(0.000,0.000)
complex(0.000,0.000) complex(0.000,0.000) complex(-0.000,0.000)
]]>