You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by er...@apache.org on 2013/10/29 16:43:05 UTC

svn commit: r1536766 - in /commons/proper/math/trunk: pom.xml src/changes/changes.xml src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java

Author: erans
Date: Tue Oct 29 15:43:04 2013
New Revision: 1536766

URL: http://svn.apache.org/r1536766
Log:
MATH-1045
Singular matrices were considered non-singular due to strict comparison
with zero. Reported and fixed by Sean Owen.

Modified:
    commons/proper/math/trunk/pom.xml
    commons/proper/math/trunk/src/changes/changes.xml
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java

Modified: commons/proper/math/trunk/pom.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/pom.xml?rev=1536766&r1=1536765&r2=1536766&view=diff
==============================================================================
--- commons/proper/math/trunk/pom.xml (original)
+++ commons/proper/math/trunk/pom.xml Tue Oct 29 15:43:04 2013
@@ -253,6 +253,9 @@
       <name>Fredrik Norin</name>
     </contributor>
     <contributor>
+      <name>Sean Owen</name>
+    </contributor>
+    <contributor>
       <name>Sujit Pal</name>
     </contributor>
     <contributor>

Modified: commons/proper/math/trunk/src/changes/changes.xml
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/changes/changes.xml?rev=1536766&r1=1536765&r2=1536766&view=diff
==============================================================================
--- commons/proper/math/trunk/src/changes/changes.xml (original)
+++ commons/proper/math/trunk/src/changes/changes.xml Tue Oct 29 15:43:04 2013
@@ -51,6 +51,9 @@ If the output is not quite correct, chec
   </properties>
   <body>
     <release version="x.y" date="TBD" description="TBD">
+      <action dev="erans" type="fix" issue="MATH-1045" due-to="Sean Owen">
+        "EigenDecomposition": Using tolerance for detecting whether a matrix is singular.
+      </action>
       <action dev="luc" type="add" issue="MATH-1036" due-to="Ajo Fod">
         Added SparseGradient to deal efficiently with first derivatives when the number
         of variables is very large but most computations depend only on a few of the

Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java?rev=1536766&r1=1536765&r2=1536766&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java Tue Oct 29 15:43:04 2013
@@ -513,9 +513,16 @@ public class EigenDecomposition {
          * @return true if the decomposed matrix is non-singular.
          */
         public boolean isNonSingular() {
+            // The eigenvalues are sorted by size, descending
+            double largestEigenvalueNorm = eigenvalueNorm(0);
+            // Corner case: zero matrix, all exactly 0 eigenvalues
+            if (largestEigenvalueNorm == 0.0) {
+                return false;
+            }
             for (int i = 0; i < realEigenvalues.length; ++i) {
-                if (realEigenvalues[i] == 0 &&
-                    imagEigenvalues[i] == 0) {
+                // Looking for eigenvalues that are 0, where we consider anything much much smaller
+                // than the largest eigenvalue to be effectively 0.
+                if (Precision.equals(eigenvalueNorm(i) / largestEigenvalueNorm, 0, EPSILON)) {
                     return false;
                 }
             }
@@ -523,6 +530,16 @@ public class EigenDecomposition {
         }
 
         /**
+         * @param i which eigenvalue to find the norm of
+         * @return the norm of ith (complex) eigenvalue.
+         */
+        private double eigenvalueNorm(int i) {
+            final double re = realEigenvalues[i];
+            final double im = imagEigenvalues[i];
+            return FastMath.sqrt(re * re + im * im);
+        }
+
+        /**
          * Get the inverse of the decomposed matrix.
          *
          * @return the inverse matrix.

Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java?rev=1536766&r1=1536765&r2=1536766&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java (original)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenSolverTest.java Tue Oct 29 15:43:04 2013
@@ -27,6 +27,13 @@ import org.junit.Assert;
 
 public class EigenSolverTest {
 
+    private double[][] bigSingular = {
+        { 1.0, 2.0,   3.0,    4.0 },
+        { 2.0, 5.0,   3.0,    4.0 },
+        { 7.0, 3.0, 256.0, 1930.0 },
+        { 3.0, 7.0,   6.0,    8.0 }
+    }; // 4th row = 1st + 2nd
+
     /** test non invertible matrix */
     @Test
     public void testNonInvertible() {
@@ -86,6 +93,20 @@ public class EigenSolverTest {
         }
     }
 
+    @Test(expected=SingularMatrixException.class)
+    public void testNonInvertibleMath1045() {
+        EigenDecomposition eigen =
+            new EigenDecomposition(MatrixUtils.createRealMatrix(bigSingular));
+        eigen.getSolver().getInverse();
+    }
+
+    @Test(expected=SingularMatrixException.class)
+    public void testZeroMatrix() {
+        EigenDecomposition eigen =
+            new EigenDecomposition(MatrixUtils.createRealMatrix(new double[][] {{0}}));
+        eigen.getSolver().getInverse();
+    }
+
     /** test solve dimension errors */
     @Test
     public void testSolveDimensionErrors() {