( If V is diagonal, stop reading and say it !!! )
I initially wrote " is not diagonal " that was a typo. Corrected. I meant "is diagonal".
What I tried to say is that often V is diagonal SPD, but of course the problem is pretty much trivial.
If V is not diagonal, please read below.
1) Cholesky factorization of V. (Compute L such that V=LL^T ). Use DPOTRF.
2) Then scale the system. So overwrite A with ( L \ A ) (use DTRSM), overwrite B with ( L \ B ) (use DTRSM or DTRSV).
3) call DGELS on ( A, B )
Another way to do this would be to
1) call LAPACK DPOTRF on the matrix V to obtain its Cholesky factor, (upper version,)
2) zero out the lower part of the Cholesky factor
2) call LAPACK DGGGLM with B, the Cholesky factor of V.
I think this works. I never tried. I think this is more stable.
Last edited by Julien Langou
on Sat Apr 30, 2011 5:36 pm, edited 1 time in total.