# commons-dev mailing list archives

##### Site index · List index
Message view
Top
From Peter <ppl...@gmail.com>
Subject [PATCH] Cholesky decomposition
Date Mon, 13 Dec 2004 21:31:15 GMT
```This method implements the Cholesky decomposition of a RealMatrix. It
can be easily modified to be incorporated under the RealMatrix class and
to decompose BigMatrixes as well.
I hope you place it in the next release.

Peter P. Lupo
Computer Science barchelor student. UFRJ-Brazil
http://www.dcc.ufrj.br
http://www.ufrj.br

/**
* Computes the cholesky decomposition of matrix. The return is a lower triangular matrix
that may be empty and also coincide with the original matrix. The method only uses the upper
right part of the matrix, which must be square, symmetric and positive-definite by definition.
* @param matrix
*            The matrix to decompose.
* @return
*            The matrix with the result, which is lower triangular and has the same dimensions
of the original matrix.
* @exception InvalidMatrixException
*                if matrix is not square, symmetric or positive-definite.
*/
public static RealMatrix choleskyDecomposition(RealMatrix matrix) {
if (!matrix.isSquare()) throw new InvalidMatrixException("Not a
square/symmetric matrix."); //"matrix" is not square
if (!matrix.equals(matrix.transpose())) throw new
InvalidMatrixException("Not a symmetric matrix."); //"matrix" is not square
double[][] matrixData = matrix.copy().getData();
double[][] cholesky = new
double[matrix.getRowDimension()][matrix.getColumnDimension()];
for (int i = 0; i < matrix.getColumnDimension(); i++) {
double[] ithRowOfL = cholesky[i];
double[] ithRowOfA = matrixData[i];
for (int j = i; j < matrix.getRowDimension(); j++) {
double[] jthRowOfL = cholesky[j];
double sum = ithRowOfA[j];
for (int k = i - 1; k >= 0; k--){
sum -= ithRowOfL[k] * jthRowOfL[k];
}
if (i == j) {
if (sum <= 0)throw new InvalidMatrixException("Not a
positive definite matrix."); //"matrix", with rounding errors is not
positive definite
if (cholesky != null){
cholesky[j][i] = Math.sqrt(sum);
}
} else {
if (cholesky != null){
cholesky[j][i] = sum / cholesky[i][i];
}
if (i < j){
cholesky[i][j] = 0;
}
}
}
}
return new RealMatrixImpl(cholesky);
}

---------------------------------------------------------------------
To unsubscribe, e-mail: commons-dev-unsubscribe@jakarta.apache.org