Re: CholeskyDecomposition

["Roldan Pozo" <roldanpozo@hotmail.com>]

Brian,

Thanks for your note.  The Cholesky solver is a known problem and we are in 
the process of releasing an update for this.

Best regards,
Roldan

----- Original Message ----- 
From: Suchomel, Brian (Mission Systems)
To: Multiple recipients of list
Sent: Thursday, March 03, 2005 11:09 AM
Subject: CholeskyDecomposition


Hi:

I recently downloaded version 1.0.1 of the Jama linear algebra package.
I ran into some problems with the CholeskyDecomposition class.  In
particular, the solve() method was giving the wrong solution on a symmetric
positive-definite matrix example.  Meanwhile, Matrix::solve() was giving the
correct solution.  The solve() routine does not contain standard forward
and backward solve routines.  I rewrote the code and was able to get
consistent solutions when the right hand side was a vector.  I've included
the code below.

// Solve L*Y = B;
for (int k=0; k<n; k++) {
   for (int i=0; i<k; i++) {
      for (int j=0; j<nx; j++) {
         X[k][j] -= X[i][j]*L[k][i];
      }
   }
   for (int j=0; j<nx; j++) {
      X[k][j] /= L[k][k];
   }
}

// Solve L'*X = Y
for (int k=n-1; k>=0; k--) {
   for (int i=k+1; i<n; i++) {
      for (int j=0; j<nx; j++) {
         X[k][j] -= X[i][j]*L[i][k];
      }
   }
   for (int j=0; j<nx; j++) {
      X[k][j] /= L[k][k];
   }
}

Otherwise, I found the package easy to use.  It's nice to have a simplified
version of LaPack available for Java.

I thought the documentation was OK.  Colleagues of mine who do not have
much experience with linear algebra thought the Documentation could be
clearer.  For example, in Cholesky, the constructor lists a matrix as "Arg".
Then, the solve routine solves A*X=B, given B.  Of course, it's obvious that
A = Arg and that you find X s.t. X = A^{-1}*B, but maybe to reach a larger 
audience . .

Thanks for making the package available,
Brian Suchomel