From 51f446bf0c09a0d051341ccec6ae50a0daa1728b Mon Sep 17 00:00:00 2001 From: Thomas Neidhart Date: Fri, 14 Sep 2012 21:55:46 +0000 Subject: [PATCH] [MATH-789] Fixed rank calculation in case of dependant columns, added additional constructor that repaces small parameter. git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1384945 13f79535-47bb-0310-9956-ffa450edef68 --- .../RectangularCholeskyDecomposition.java | 24 ++++++++++++++++--- .../RectangularCholeskyDecompositionTest.java | 20 +++++++++------- .../CorrelatedRandomVectorGeneratorTest.java | 2 +- 3 files changed, 33 insertions(+), 13 deletions(-) diff --git a/src/main/java/org/apache/commons/math3/linear/RectangularCholeskyDecomposition.java b/src/main/java/org/apache/commons/math3/linear/RectangularCholeskyDecomposition.java index aba7b9806..b870800f0 100644 --- a/src/main/java/org/apache/commons/math3/linear/RectangularCholeskyDecomposition.java +++ b/src/main/java/org/apache/commons/math3/linear/RectangularCholeskyDecomposition.java @@ -52,9 +52,28 @@ public class RectangularCholeskyDecomposition { /** * Decompose a symmetric positive semidefinite matrix. + *

+ * Note: this constructor follows the linpack method to detect dependent + * columns by proceeding with the Cholesky algorithm until a nonpositive diagonal + * element is encountered. + * + * @see + * Analysis of the Cholesky Decomposition of a Semi-definite Matrix * * @param matrix Symmetric positive semidefinite matrix. - * @param small Diagonal elements threshold under which column are + * @exception NonPositiveDefiniteMatrixException if the matrix is not + * positive semidefinite. + */ + public RectangularCholeskyDecomposition(RealMatrix matrix) + throws NonPositiveDefiniteMatrixException { + this(matrix, 0); + } + + /** + * Decompose a symmetric positive semidefinite matrix. + * + * @param matrix Symmetric positive semidefinite matrix. + * @param small Diagonal elements threshold under which columns are * considered to be dependent on previous ones and are discarded. * @exception NonPositiveDefiniteMatrixException if the matrix is not * positive semidefinite. @@ -97,7 +116,7 @@ public class RectangularCholeskyDecomposition { // check diagonal element int ir = index[r]; - if (c[ir][ir] < small) { + if (c[ir][ir] <= small) { if (r == 0) { throw new NonPositiveDefiniteMatrixException(c[ir][ir], ir, small); @@ -114,7 +133,6 @@ public class RectangularCholeskyDecomposition { // all remaining diagonal elements are close to zero, we consider we have // found the rank of the symmetric positive semidefinite matrix - ++r; loop = false; } else { diff --git a/src/test/java/org/apache/commons/math3/linear/RectangularCholeskyDecompositionTest.java b/src/test/java/org/apache/commons/math3/linear/RectangularCholeskyDecompositionTest.java index 6fbfdfd3c..fa43eac33 100644 --- a/src/test/java/org/apache/commons/math3/linear/RectangularCholeskyDecompositionTest.java +++ b/src/test/java/org/apache/commons/math3/linear/RectangularCholeskyDecompositionTest.java @@ -81,9 +81,7 @@ public class RectangularCholeskyDecompositionTest { {0.009881156, 0.008196856, 0.019023866, 0.009210099}, {0.010499559, 0.010732709, 0.009210099, 0.019107243} }); - RealMatrix root1 = new RectangularCholeskyDecomposition(m1, 1.0e-10).getRootMatrix(); - RealMatrix rebuiltM1 = root1.multiply(root1.transpose()); - Assert.assertEquals(0.0, m1.subtract(rebuiltM1).getNorm(), 1.0e-16); + composeAndTest(m1, 4); final RealMatrix m2 = MatrixUtils.createRealMatrix(new double[][]{ {0.0, 0.0, 0.0, 0.0, 0.0}, @@ -92,9 +90,7 @@ public class RectangularCholeskyDecompositionTest { {0.0, 0.009881156, 0.008196856, 0.019023866, 0.009210099}, {0.0, 0.010499559, 0.010732709, 0.009210099, 0.019107243} }); - RealMatrix root2 = new RectangularCholeskyDecomposition(m2, 1.0e-10).getRootMatrix(); - RealMatrix rebuiltM2 = root2.multiply(root2.transpose()); - Assert.assertEquals(0.0, m2.subtract(rebuiltM2).getNorm(), 1.0e-16); + composeAndTest(m2, 4); final RealMatrix m3 = MatrixUtils.createRealMatrix(new double[][]{ {0.013445532, 0.010394690, 0.0, 0.009881156, 0.010499559}, @@ -103,10 +99,16 @@ public class RectangularCholeskyDecompositionTest { {0.009881156, 0.008196856, 0.0, 0.019023866, 0.009210099}, {0.010499559, 0.010732709, 0.0, 0.009210099, 0.019107243} }); - RealMatrix root3 = new RectangularCholeskyDecomposition(m3, 1.0e-10).getRootMatrix(); - RealMatrix rebuiltM3 = root3.multiply(root3.transpose()); - Assert.assertEquals(0.0, m3.subtract(rebuiltM3).getNorm(), 1.0e-16); + composeAndTest(m3, 4); } + + private void composeAndTest(RealMatrix m, int expectedRank) { + RectangularCholeskyDecomposition r = new RectangularCholeskyDecomposition(m); + Assert.assertEquals(expectedRank, r.getRank()); + RealMatrix root = r.getRootMatrix(); + RealMatrix rebuiltMatrix = root.multiply(root.transpose()); + Assert.assertEquals(0.0, m.subtract(rebuiltMatrix).getNorm(), 1.0e-16); + } } diff --git a/src/test/java/org/apache/commons/math3/random/CorrelatedRandomVectorGeneratorTest.java b/src/test/java/org/apache/commons/math3/random/CorrelatedRandomVectorGeneratorTest.java index eba52a387..0496e9832 100644 --- a/src/test/java/org/apache/commons/math3/random/CorrelatedRandomVectorGeneratorTest.java +++ b/src/test/java/org/apache/commons/math3/random/CorrelatedRandomVectorGeneratorTest.java @@ -68,7 +68,7 @@ public class CorrelatedRandomVectorGeneratorTest { @Test public void testRank() { - Assert.assertEquals(3, generator.getRank()); + Assert.assertEquals(2, generator.getRank()); } @Test