On 09/13/2012 11:36 PM, Robert Dodier wrote:
> On 2012-09-13, Barton Willis <willisb at unk.edu> wrote:
>
>> Th Cholesky code doesn't work for symbolic cases, I think.
> Hmm ... there used to be a cholesky.mac in share which worked for
> symbolic cases, but it got nuked in favor of the one in the
> linearalgebra package. For the record here it is.
>
> load ("cholesky.mac");
> M : matrix ([1, p12, p13], [p12, 1, p23], [p13, p23, 1]);
> M1 : cholesky (M1);
> => matrix([1,0,0],[p12,sqrt(1-p12^2),0],
> [p13,(p23-p12*p13)/sqrt(1-p12^2),
> sqrt(-(p23-p12*p13)^2/(1-p12^2)-p13^2+1)])
> M1 . transpose (M1) - M;
> => matrix([0,0,0],[0,0,0],[0,0,0])
One could also use a L.D.transpose(L) decomposition
and avoid the usage of square roots.
Regards
Andre