You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by lu...@apache.org on 2009/02/15 18:53:33 UTC
svn commit: r744708 - in /commons/proper/math/trunk/src:
java/org/apache/commons/math/ java/org/apache/commons/math/linear/
site/xdoc/ site/xdoc/userguide/ test/org/apache/commons/math/linear/
Author: luc
Date: Sun Feb 15 17:53:32 2009
New Revision: 744708
URL: http://svn.apache.org/viewvc?rev=744708&view=rev
Log:
added Cholesky decomposition
Added:
commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java (with props)
commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java (with props)
commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java (with props)
commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java (with props)
commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java (with props)
Modified:
commons/proper/math/trunk/src/java/org/apache/commons/math/MessagesResources_fr.java
commons/proper/math/trunk/src/site/xdoc/changes.xml
commons/proper/math/trunk/src/site/xdoc/userguide/linear.xml
Modified: commons/proper/math/trunk/src/java/org/apache/commons/math/MessagesResources_fr.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/java/org/apache/commons/math/MessagesResources_fr.java?rev=744708&r1=744707&r2=744708&view=diff
==============================================================================
--- commons/proper/math/trunk/src/java/org/apache/commons/math/MessagesResources_fr.java (original)
+++ commons/proper/math/trunk/src/java/org/apache/commons/math/MessagesResources_fr.java Sun Feb 15 17:53:32 2009
@@ -68,10 +68,14 @@
{ "dimension mismatch {0} != {1}",
"dimensions incompatibles {0} != {1}" },
- // org.apache.commons.math.random.NotPositiveDefiniteMatrixException
+ // org.apache.commons.math.linear.NotPositiveDefiniteMatrixException
{ "not positive definite matrix",
"matrice non d\u00e9finie positive" },
+ // org.apache.commons.math.linear.NotSymmetricMatrixException
+ { "not symmetric matrix",
+ "matrice non symm\u00e9trique" },
+
// org.apache.commons.math.fraction.FractionConversionException
{ "Unable to convert {0} to fraction after {1} iterations",
"Impossible de convertir {0} en fraction apr\u00e8s {1} it\u00e9rations" },
Added: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java?rev=744708&view=auto
==============================================================================
--- commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java (added)
+++ commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java Sun Feb 15 17:53:32 2009
@@ -0,0 +1,72 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.linear;
+
+import java.io.Serializable;
+
+/**
+ * An interface to classes that implement an algorithm to calculate the
+ * Cholesky decomposition of a real symmetric positive-definite matrix.
+ * <p>This interface is based on the class with similar name from the now defunct
+ * <a href="http://math.nist.gov/javanumerics/jama/">JAMA</a> library, with the
+ * following changes:</p>
+ * <ul>
+ * <li>a {@link #getLT() getLT} method has been added,</li>
+ * <li>the <code>isspd</code> method has been removed, the constructors of
+ * implementation classes being expected to throw {@link
+ * NotPositiveDefiniteMatrixException} when a matrix cannot be decomposed,</li>
+ * <li>a {@link #getDeterminant() getDeterminant} method has been added,</li>
+ * <li>the <code>solve</code> method has been replaced by a {@link
+ * #getSolver() getSolver} method and the equivalent method provided by
+ * the returned {@link DecompositionSolver}.</li>
+ * </ul>
+ *
+ * @see <a href="http://mathworld.wolfram.com/CholeskyDecomposition.html">MathWorld</a>
+ * @see <a href="http://en.wikipedia.org/wiki/Cholesky_decomposition">Wikipedia</a>
+ * @version $Revision$ $Date$
+ * @since 2.0
+ */
+public interface CholeskyDecomposition extends Serializable {
+
+ /**
+ * Returns the matrix L of the decomposition.
+ * <p>L is an lower-triangular matrix</p>
+ * @return the L matrix
+ */
+ RealMatrix getL();
+
+ /**
+ * Returns the transpose of the matrix L of the decomposition.
+ * <p>L<sup>T</sup> is an upper-triangular matrix</p>
+ * @return the transpose of the matrix L of the decomposition
+ */
+ RealMatrix getLT();
+
+ /**
+ * Return the determinant of the matrix
+ * @return determinant of the matrix
+ */
+ double getDeterminant();
+
+ /**
+ * Get a solver for finding the A × X = B solution in least square sense.
+ * @return a solver
+ */
+ DecompositionSolver getSolver();
+
+}
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecomposition.java
------------------------------------------------------------------------------
svn:keywords = Author Date Id Revision
Added: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java?rev=744708&view=auto
==============================================================================
--- commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java (added)
+++ commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java Sun Feb 15 17:53:32 2009
@@ -0,0 +1,359 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.linear;
+
+import org.apache.commons.math.MathRuntimeException;
+
+
+/**
+ * Calculates the Cholesky decomposition of a matrix.
+ * <p>The Cholesky decomposition of a real symmetric positive-definite
+ * matrix A consists of a lower triangular matrix L with same size that
+ * satisfy: A = LL<sup>T</sup>Q = I). In a sense, this is the square root of A.</p>
+ *
+ * @see <a href="http://mathworld.wolfram.com/CholeskyDecomposition.html">MathWorld</a>
+ * @see <a href="http://en.wikipedia.org/wiki/Cholesky_decomposition">Wikipedia</a>
+ * @version $Revision$ $Date$
+ * @since 2.0
+ */
+public class CholeskyDecompositionImpl implements CholeskyDecomposition {
+
+ /** Serializable version identifier. */
+ private static final long serialVersionUID = -2036131698031167221L;
+
+ /** Default threshold above which off-diagonal elements are considered too different
+ * and matrix not symmetric. */
+ public static final double DEFAULT_RELATIVE_SYMMETRY_THRESHOLD = 1.0e-15;
+
+ /** Default threshold below which diagonal elements are considered null
+ * and matrix not positive definite. */
+ public static final double DEFAULT_ABSOLUTE_POSITIVITY_THRESHOLD = 1.0e-10;
+
+ /** Row-oriented storage for L<sup>T</sup> matrix data. */
+ private double[][] lTData;
+
+ /** Cached value of L. */
+ private RealMatrix cachedL;
+
+ /** Cached value of LT. */
+ private RealMatrix cachedLT;
+
+ /**
+ * Calculates the Cholesky decomposition of the given matrix.
+ * <p>
+ * Calling this constructor is equivalent to call {@link
+ * #CholeskyDecompositionImpl(RealMatrix, double, double)} with the
+ * thresholds set to the default values {@link
+ * #DEFAULT_RELATIVE_SYMMETRY_THRESHOLD} and {@link
+ * #DEFAULT_ABSOLUTE_POSITIVITY_THRESHOLD}
+ * </p>
+ * @param matrix the matrix to decompose
+ * @exception NonSquareMatrixException if matrix is not square
+ * @exception NotSymmetricMatrixException if matrix is not symmetric
+ * @exception NotPositiveDefiniteMatrixException if the matrix is not
+ * strictly positive definite
+ * @see #CholeskyDecompositionImpl(RealMatrix, double, double)
+ * @see #DEFAULT_RELATIVE_SYMMETRY_THRESHOLD
+ * @see #DEFAULT_ABSOLUTE_POSITIVITY_THRESHOLD
+ */
+ public CholeskyDecompositionImpl(final RealMatrix matrix)
+ throws NonSquareMatrixException,
+ NotSymmetricMatrixException, NotPositiveDefiniteMatrixException {
+ this(matrix, DEFAULT_RELATIVE_SYMMETRY_THRESHOLD,
+ DEFAULT_ABSOLUTE_POSITIVITY_THRESHOLD);
+ }
+
+ /**
+ * Calculates the Cholesky decomposition of the given matrix.
+ * @param matrix the matrix to decompose
+ * @param relativeSymmetryThreshold threshold above which off-diagonal
+ * elements are considered too different and matrix not symmetric
+ * @param absolutePositivityThreshold threshold below which diagonal
+ * elements are considered null and matrix not positive definite
+ * @exception NonSquareMatrixException if matrix is not square
+ * @exception NotSymmetricMatrixException if matrix is not symmetric
+ * @exception NotPositiveDefiniteMatrixException if the matrix is not
+ * strictly positive definite
+ * @see #CholeskyDecompositionImpl(RealMatrix)
+ * @see #DEFAULT_RELATIVE_SYMMETRY_THRESHOLD
+ * @see #DEFAULT_ABSOLUTE_POSITIVITY_THRESHOLD
+ */
+ public CholeskyDecompositionImpl(final RealMatrix matrix,
+ final double relativeSymmetryThreshold,
+ final double absolutePositivityThreshold)
+ throws NonSquareMatrixException,
+ NotSymmetricMatrixException, NotPositiveDefiniteMatrixException {
+
+ if (!matrix.isSquare()) {
+ throw new NonSquareMatrixException(matrix.getRowDimension(),
+ matrix.getColumnDimension());
+ }
+
+ final int order = matrix.getRowDimension();
+ lTData = matrix.getData();
+ cachedL = null;
+ cachedLT = null;
+
+ // check the matrix before transformation
+ for (int i = 0; i < order; ++i) {
+
+ final double[] lI = lTData[i];
+
+ // check diagonal element
+ if (lTData[i][i] < absolutePositivityThreshold) {
+ throw new NotPositiveDefiniteMatrixException();
+ }
+
+ // check off-diagonal elements (and reset them to 0)
+ for (int j = i + 1; j < order; ++j) {
+ final double[] lJ = lTData[j];
+ final double lIJ = lI[j];
+ final double lJI = lJ[i];
+ final double maxDelta =
+ relativeSymmetryThreshold * Math.max(Math.abs(lIJ), Math.abs(lJI));
+ if (Math.abs(lIJ - lJI) > maxDelta) {
+ throw new NotSymmetricMatrixException();
+ }
+ lJ[i] = 0;
+ }
+ }
+
+ // transform the matrix
+ for (int i = 0; i < order; ++i) {
+
+ final double[] ltI = lTData[i];
+ ltI[i] = Math.sqrt(ltI[i]);
+ final double inverse = 1.0 / ltI[i];
+
+ for (int q = order - 1; q > i; --q) {
+ ltI[q] *= inverse;
+ final double[] ltQ = lTData[q];
+ for (int p = q; p < order; ++p) {
+ ltQ[p] -= ltI[q] * ltI[p];
+ }
+ }
+
+ }
+
+ }
+
+ /** {@inheritDoc} */
+ public RealMatrix getL() {
+ if (cachedL == null) {
+ cachedL = getLT().transpose();
+ }
+ return cachedL;
+ }
+
+ /** {@inheritDoc} */
+ public RealMatrix getLT() {
+
+ if (cachedLT == null) {
+ cachedLT = MatrixUtils.createRealMatrix(lTData);
+ }
+
+ // return the cached matrix
+ return cachedLT;
+
+ }
+
+ /** {@inheritDoc} */
+ public double getDeterminant() {
+ double determinant = 1.0;
+ for (int i = 0; i < lTData.length; ++i) {
+ double lTii = lTData[i][i];
+ determinant *= lTii * lTii;
+ }
+ return determinant;
+ }
+
+ /** {@inheritDoc} */
+ public DecompositionSolver getSolver() {
+ return new Solver(lTData);
+ }
+
+ /** Specialized solver. */
+ private static class Solver implements DecompositionSolver {
+
+ /** Serializable version identifier. */
+ private static final long serialVersionUID = -7288829864732555901L;
+
+ /** Row-oriented storage for L<sup>T</sup> matrix data. */
+ private final double[][] lTData;
+
+ /**
+ * Build a solver from decomposed matrix.
+ * @param lData row-oriented storage for L<sup>T</sup> matrix data
+ */
+ private Solver(final double[][] lTData) {
+ this.lTData = lTData;
+ }
+
+ /** {@inheritDoc} */
+ public boolean isNonSingular() {
+ // if we get this far, the matrix was positive definite, hence non-singular
+ return true;
+ }
+
+ /** {@inheritDoc} */
+ public double[] solve(double[] b)
+ throws IllegalArgumentException, InvalidMatrixException {
+
+ final int m = lTData.length;
+ if (b.length != m) {
+ throw MathRuntimeException.createIllegalArgumentException(
+ "vector length mismatch: got {0} but expected {1}",
+ new Object[] { b.length, m });
+ }
+
+ final double[] x = b.clone();
+
+ // Solve LY = b
+ for (int j = 0; j < m; j++) {
+ final double[] lJ = lTData[j];
+ x[j] /= lJ[j];
+ final double xJ = x[j];
+ for (int i = j + 1; i < m; i++) {
+ x[i] -= xJ * lJ[i];
+ }
+ }
+
+ // Solve LTX = Y
+ for (int j = m - 1; j >= 0; j--) {
+ x[j] /= lTData[j][j];
+ final double xJ = x[j];
+ for (int i = 0; i < j; i++) {
+ x[i] -= xJ * lTData[i][j];
+ }
+ }
+
+ return x;
+
+ }
+
+ /** {@inheritDoc} */
+ public RealVector solve(RealVector b)
+ throws IllegalArgumentException, InvalidMatrixException {
+ try {
+ return solve((RealVectorImpl) b);
+ } catch (ClassCastException cce) {
+
+ final int m = lTData.length;
+ if (b.getDimension() != m) {
+ throw MathRuntimeException.createIllegalArgumentException(
+ "vector length mismatch: got {0} but expected {1}",
+ new Object[] { b.getDimension(), m });
+ }
+
+ final double[] x = b.getData();
+
+ // Solve LY = b
+ for (int j = 0; j < m; j++) {
+ final double[] lJ = lTData[j];
+ x[j] /= lJ[j];
+ final double xJ = x[j];
+ for (int i = j + 1; i < m; i++) {
+ x[i] -= xJ * lJ[i];
+ }
+ }
+
+ // Solve LTX = Y
+ for (int j = m - 1; j >= 0; j--) {
+ x[j] /= lTData[j][j];
+ final double xJ = x[j];
+ for (int i = 0; i < j; i++) {
+ x[i] -= xJ * lTData[i][j];
+ }
+ }
+
+ return new RealVectorImpl(x, false);
+
+ }
+ }
+
+ /** Solve the linear equation A × X = B.
+ * <p>The A matrix is implicit here. It is </p>
+ * @param b right-hand side of the equation A × X = B
+ * @return a vector X such that A × X = B
+ * @exception IllegalArgumentException if matrices dimensions don't match
+ * @exception InvalidMatrixException if decomposed matrix is singular
+ */
+ public RealVectorImpl solve(RealVectorImpl b)
+ throws IllegalArgumentException, InvalidMatrixException {
+ return new RealVectorImpl(solve(b.getDataRef()), false);
+ }
+
+ /** {@inheritDoc} */
+ public RealMatrix solve(RealMatrix b)
+ throws IllegalArgumentException, InvalidMatrixException {
+
+ final int m = lTData.length;
+ if (b.getRowDimension() != m) {
+ throw MathRuntimeException.createIllegalArgumentException(
+ "dimensions mismatch: got {0}x{1} but expected {2}x{3}",
+ new Object[] { b.getRowDimension(), b.getColumnDimension(), m, "n"});
+ }
+
+ final int nColB = b.getColumnDimension();
+ double[][] x = b.getData();
+
+ // Solve LY = b
+ for (int j = 0; j < m; j++) {
+ final double[] lJ = lTData[j];
+ final double lJJ = lJ[j];
+ final double[] xJ = x[j];
+ for (int k = 0; k < nColB; ++k) {
+ xJ[k] /= lJJ;
+ }
+ for (int i = j + 1; i < m; i++) {
+ final double[] xI = x[i];
+ final double lJI = lJ[i];
+ for (int k = 0; k < nColB; ++k) {
+ xI[k] -= xJ[k] * lJI;
+ }
+ }
+ }
+
+ // Solve LTX = Y
+ for (int j = m - 1; j >= 0; j--) {
+ final double lJJ = lTData[j][j];
+ final double[] xJ = x[j];
+ for (int k = 0; k < nColB; ++k) {
+ xJ[k] /= lJJ;
+ }
+ for (int i = 0; i < j; i++) {
+ final double[] xI = x[i];
+ final double lIJ = lTData[i][j];
+ for (int k = 0; k < nColB; ++k) {
+ xI[k] -= xJ[k] * lIJ;
+ }
+ }
+ }
+
+ return new RealMatrixImpl(x, false);
+
+ }
+
+ /** {@inheritDoc} */
+ public RealMatrix getInverse() throws InvalidMatrixException {
+ return solve(MatrixUtils.createRealIdentityMatrix(lTData.length));
+ }
+
+ }
+
+}
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/CholeskyDecompositionImpl.java
------------------------------------------------------------------------------
svn:keywords = Author Date Id Revision
Added: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java?rev=744708&view=auto
==============================================================================
--- commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java (added)
+++ commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java Sun Feb 15 17:53:32 2009
@@ -0,0 +1,42 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.linear;
+
+import org.apache.commons.math.MathException;
+
+/**
+ * This class represents exceptions thrown when a matrix expected to
+ * be symmetric is not
+ *
+ * @since 2.0
+ * @version $Revision$ $Date$
+ */
+
+public class NotSymmetricMatrixException extends MathException {
+
+ /** Serializable version identifier */
+ private static final long serialVersionUID = -7012803946709786097L;
+
+ /** Simple constructor.
+ * build an exception with a default message.
+ */
+ public NotSymmetricMatrixException() {
+ super("not symmetric matrix", null);
+ }
+
+}
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/java/org/apache/commons/math/linear/NotSymmetricMatrixException.java
------------------------------------------------------------------------------
svn:keywords = Author Date Id Revision
Modified: commons/proper/math/trunk/src/site/xdoc/changes.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/site/xdoc/changes.xml?rev=744708&r1=744707&r2=744708&view=diff
==============================================================================
--- commons/proper/math/trunk/src/site/xdoc/changes.xml (original)
+++ commons/proper/math/trunk/src/site/xdoc/changes.xml Sun Feb 15 17:53:32 2009
@@ -170,7 +170,8 @@
decompose a matrix as a product of several other matrices with predefined
properties and shapes depending on the algorithm. These algorithms allow to
solve the equation A * X = B, either for an exact linear solution
- (LU-decomposition) or an exact or least-squares solution (QR-decomposition).
+ (LU-decomposition, Cholesky decomposition) or an exact or least-squares
+ solution (QR-decomposition).
</action>
<action dev="luc" type="add" issue="MATH-219" due-to="Andrew Berry">
Added removeData methods for the SimpleRegression class. This allows
Modified: commons/proper/math/trunk/src/site/xdoc/userguide/linear.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/site/xdoc/userguide/linear.xml?rev=744708&r1=744707&r2=744708&view=diff
==============================================================================
--- commons/proper/math/trunk/src/site/xdoc/userguide/linear.xml (original)
+++ commons/proper/math/trunk/src/site/xdoc/userguide/linear.xml Sun Feb 15 17:53:32 2009
@@ -132,7 +132,8 @@
for X is such that the residual AX-B has minimal norm. If an exact solution
exist (i.e. if for some X the residual AX-B is exactly 0), then this exact
solution is also the solution in least square sense. Some solvers like
- the one obtained from <code>LUDecomposition</code> can only find the solution
+ the one obtained from <code>LUDecomposition</code> or
+ <code>CholeskyDecomposition</code>code> can only find the solution
for square matrices and when the solution is an exact linear solution. Other
solvers like the one obtained from <code>QRDecomposition</code>
are more versatile and can also find solutions with non-square matrix A or when
Added: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java?rev=744708&view=auto
==============================================================================
--- commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java (added)
+++ commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java Sun Feb 15 17:53:32 2009
@@ -0,0 +1,152 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.linear;
+
+import org.apache.commons.math.MathException;
+
+import junit.framework.Test;
+import junit.framework.TestCase;
+import junit.framework.TestSuite;
+
+public class CholeskyDecompositionImplTest extends TestCase {
+
+ private double[][] testData = new double[][] {
+ { 1, 2, 4, 7, 11 },
+ { 2, 13, 23, 38, 58 },
+ { 4, 23, 77, 122, 182 },
+ { 7, 38, 122, 294, 430 },
+ { 11, 58, 182, 430, 855 }
+ };
+
+ public CholeskyDecompositionImplTest(String name) {
+ super(name);
+ }
+
+ public static Test suite() {
+ TestSuite suite = new TestSuite(CholeskyDecompositionImplTest.class);
+ suite.setName("CholeskyDecompositionImpl Tests");
+ return suite;
+ }
+
+ /** test dimensions */
+ public void testDimensions() throws MathException {
+ CholeskyDecomposition llt =
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(testData));
+ assertEquals(testData.length, llt.getL().getRowDimension());
+ assertEquals(testData.length, llt.getL().getColumnDimension());
+ assertEquals(testData.length, llt.getLT().getRowDimension());
+ assertEquals(testData.length, llt.getLT().getColumnDimension());
+ }
+
+ /** test non-square matrix */
+ public void testNonSquare() {
+ try {
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(new double[3][2]));
+ } catch (NonSquareMatrixException ime) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ }
+
+ /** test non-symmetric matrix */
+ public void testNotSymmetricMatrixException() {
+ try {
+ double[][] changed = testData.clone();
+ changed[0][changed[0].length - 1] += 1.0e-5;
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(changed));
+ } catch (NotSymmetricMatrixException e) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ }
+
+ /** test non positive definite matrix */
+ public void testNotPositiveDefinite() {
+ try {
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(new double[][] {
+ { 14, 11, 13, 15, 24 },
+ { 11, 34, 13, 8, 25 },
+ { 13, 13, 14, 15, 21 },
+ { 15, 8, 15, 18, 23 },
+ { 24, 25, 21, 23, 45 }
+ }));
+ } catch (NotPositiveDefiniteMatrixException e) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ }
+
+ /** test A = LLT */
+ public void testAEqualLLT() throws MathException {
+ RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
+ CholeskyDecomposition llt = new CholeskyDecompositionImpl(matrix);
+ RealMatrix l = llt.getL();
+ RealMatrix lt = llt.getLT();
+ double norm = l.multiply(lt).subtract(matrix).getNorm();
+ assertEquals(0, norm, 1.0e-15);
+ }
+
+ /** test that L is lower triangular */
+ public void testLLowerTriangular() throws MathException {
+ RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
+ RealMatrix l = new CholeskyDecompositionImpl(matrix).getL();
+ for (int i = 0; i < l.getRowDimension(); i++) {
+ for (int j = i + 1; j < l.getColumnDimension(); j++) {
+ assertEquals(0.0, l.getEntry(i, j));
+ }
+ }
+ }
+
+ /** test that LT is transpose of L */
+ public void testLTTransposed() throws MathException {
+ RealMatrix matrix = MatrixUtils.createRealMatrix(testData);
+ CholeskyDecomposition llt = new CholeskyDecompositionImpl(matrix);
+ RealMatrix l = llt.getL();
+ RealMatrix lt = llt.getLT();
+ double norm = l.subtract(lt.transpose()).getNorm();
+ assertEquals(0, norm, 1.0e-15);
+ }
+
+ /** test matrices values */
+ public void testMatricesValues() throws MathException {
+ RealMatrix lRef = MatrixUtils.createRealMatrix(new double[][] {
+ { 1, 0, 0, 0, 0 },
+ { 2, 3, 0, 0, 0 },
+ { 4, 5, 6, 0, 0 },
+ { 7, 8, 9, 10, 0 },
+ { 11, 12, 13, 14, 15 }
+ });
+ CholeskyDecomposition llt =
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(testData));
+
+ // check values against known references
+ RealMatrix l = llt.getL();
+ assertEquals(0, l.subtract(lRef).getNorm(), 1.0e-13);
+ RealMatrix lt = llt.getLT();
+ assertEquals(0, lt.subtract(lRef.transpose()).getNorm(), 1.0e-13);
+
+ // check the same cached instance is returned the second time
+ assertTrue(l == llt.getL());
+ assertTrue(lt == llt.getLT());
+
+ }
+
+}
Propchange: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskyDecompositionImplTest.java
------------------------------------------------------------------------------
svn:keywords = Author Date Id Revision
Added: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java?rev=744708&view=auto
==============================================================================
--- commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java (added)
+++ commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java Sun Feb 15 17:53:32 2009
@@ -0,0 +1,133 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements. See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License. You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+
+package org.apache.commons.math.linear;
+
+import junit.framework.Test;
+import junit.framework.TestCase;
+import junit.framework.TestSuite;
+
+import org.apache.commons.math.MathException;
+
+public class CholeskySolverTest extends TestCase {
+
+ private double[][] testData = new double[][] {
+ { 1, 2, 4, 7, 11 },
+ { 2, 13, 23, 38, 58 },
+ { 4, 23, 77, 122, 182 },
+ { 7, 38, 122, 294, 430 },
+ { 11, 58, 182, 430, 855 }
+ };
+
+ public CholeskySolverTest(String name) {
+ super(name);
+ }
+
+ public static Test suite() {
+ TestSuite suite = new TestSuite(CholeskySolverTest.class);
+ suite.setName("LUSolver Tests");
+ return suite;
+ }
+
+ /** test solve dimension errors */
+ public void testSolveDimensionErrors() throws MathException {
+ DecompositionSolver solver =
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(testData)).getSolver();
+ RealMatrix b = MatrixUtils.createRealMatrix(new double[2][2]);
+ try {
+ solver.solve(b);
+ fail("an exception should have been thrown");
+ } catch (IllegalArgumentException iae) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ try {
+ solver.solve(b.getColumn(0));
+ fail("an exception should have been thrown");
+ } catch (IllegalArgumentException iae) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ try {
+ solver.solve(new RealVectorImplTest.RealVectorTestImpl(b.getColumn(0)));
+ fail("an exception should have been thrown");
+ } catch (IllegalArgumentException iae) {
+ // expected behavior
+ } catch (Exception e) {
+ fail("wrong exception caught");
+ }
+ }
+
+ /** test solve */
+ public void testSolve() throws MathException {
+ DecompositionSolver solver =
+ new CholeskyDecompositionImpl(MatrixUtils.createRealMatrix(testData)).getSolver();
+ RealMatrix b = MatrixUtils.createRealMatrix(new double[][] {
+ { 78, -13, 1 },
+ { 414, -62, -1 },
+ { 1312, -202, -37 },
+ { 2989, -542, 145 },
+ { 5510, -1465, 201 }
+ });
+ RealMatrix xRef = MatrixUtils.createRealMatrix(new double[][] {
+ { 1, 0, 1 },
+ { 0, 1, 1 },
+ { 2, 1, -4 },
+ { 2, 2, 2 },
+ { 5, -3, 0 }
+ });
+
+ // using RealMatrix
+ assertEquals(0, solver.solve(b).subtract(xRef).getNorm(), 1.0e-13);
+
+ // using double[]
+ for (int i = 0; i < b.getColumnDimension(); ++i) {
+ assertEquals(0,
+ new RealVectorImpl(solver.solve(b.getColumn(i))).subtract(xRef.getColumnVector(i)).getNorm(),
+ 1.0e-13);
+ }
+
+ // using RealVectorImpl
+ for (int i = 0; i < b.getColumnDimension(); ++i) {
+ assertEquals(0,
+ solver.solve(b.getColumnVector(i)).subtract(xRef.getColumnVector(i)).getNorm(),
+ 1.0e-13);
+ }
+
+ // using RealVector with an alternate implementation
+ for (int i = 0; i < b.getColumnDimension(); ++i) {
+ RealVectorImplTest.RealVectorTestImpl v =
+ new RealVectorImplTest.RealVectorTestImpl(b.getColumn(i));
+ assertEquals(0,
+ solver.solve(v).subtract(xRef.getColumnVector(i)).getNorm(),
+ 1.0e-13);
+ }
+
+ }
+
+ /** test determinant */
+ public void testDeterminant() throws MathException {
+ assertEquals(7290000.0, getDeterminant(MatrixUtils.createRealMatrix(testData)), 1.0e-15);
+ }
+
+ private double getDeterminant(RealMatrix m) throws MathException {
+ return new CholeskyDecompositionImpl(m).getDeterminant();
+ }
+
+}
Propchange: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java
------------------------------------------------------------------------------
svn:eol-style = native
Propchange: commons/proper/math/trunk/src/test/org/apache/commons/math/linear/CholeskySolverTest.java
------------------------------------------------------------------------------
svn:keywords = Author Date Id Revision