commons-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From t.@apache.org
Subject svn commit: r1363105 - in /commons/proper/math/trunk/src: main/java/org/apache/commons/math3/linear/ test/java/org/apache/commons/math3/linear/
Date Wed, 18 Jul 2012 20:49:44 GMT
Author: tn
Date: Wed Jul 18 20:49:43 2012
New Revision: 1363105

URL: http://svn.apache.org/viewvc?rev=1363105&view=rev
Log:
[MATH-235] Added random data test for eigen decomposition, improved error handling.

Modified:
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/EigenDecomposition.java
    commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/SchurTransformer.java
    commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenDecompositionTest.java

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=1363105&r1=1363104&r2=1363105&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
Wed Jul 18 20:49:43 2012
@@ -18,6 +18,8 @@
 package org.apache.commons.math3.linear;
 
 import org.apache.commons.math3.complex.Complex;
+import org.apache.commons.math3.exception.MathArithmeticException;
+import org.apache.commons.math3.exception.MathUnsupportedOperationException;
 import org.apache.commons.math3.exception.MaxCountExceededException;
 import org.apache.commons.math3.exception.DimensionMismatchException;
 import org.apache.commons.math3.exception.util.LocalizedFormats;
@@ -102,9 +104,13 @@ public class EigenDecomposition {
 
     /**
      * Calculates the eigen decomposition of the given real matrix.
+     * <p>
+     * Supports decomposition of a general matrix since 3.1.
      *
      * @param matrix Matrix to decompose.
      * @throws MaxCountExceededException if the algorithm fails to converge.
+     * @throws MathArithmeticException if the decomposition of a general matrix
+     * results in a matrix with zero norm
      */
     public EigenDecomposition(final RealMatrix matrix)  {
         if (isSymmetric(matrix, false)) {
@@ -218,7 +224,6 @@ public class EigenDecomposition {
         }
         // return the cached matrix
         return cachedV;
-
     }
 
     /**
@@ -233,6 +238,7 @@ public class EigenDecomposition {
      * @see #getImagEigenvalues()
      */
     public RealMatrix getD() {
+
         if (cachedD == null) {
             // cache the matrix for subsequent calls
             cachedD = MatrixUtils.createRealDiagonalMatrix(realEigenvalues);
@@ -359,10 +365,20 @@ public class EigenDecomposition {
     /**
      * Gets a solver for finding the A &times; X = B solution in exact
      * linear sense.
-     *
-     * @return a solver.
+     * <p>
+     * Since 3.1, eigen decomposition of a general matrix is supported,
+     * but the {@link DecompositionSolver} only supports real eigenvalues.
+     *
+     * @return a solver
+     * @throws MathUnsupportedOperationException if the decomposition resulted in
+     * complex eigenvalues
      */
     public DecompositionSolver getSolver() {
+        for (int i = 0; i < imagEigenvalues.length; i++) {
+            if (imagEigenvalues[i] != 0.0) {
+                throw new MathUnsupportedOperationException();
+            }
+        }
         return new Solver(realEigenvalues, imagEigenvalues, eigenvectors);
     }
 
@@ -726,6 +742,7 @@ public class EigenDecomposition {
      * Find eigenvectors from a matrix transformed to Schur form.
      *
      * @param schur the schur transformation of the matrix
+     * @throws MathArithmeticException if the Schur form has a norm of zero
      */
     private void findEigenVectorsFromSchur(final SchurTransformer schur) {
         final double[][] matrixT = schur.getT().getData();
@@ -741,9 +758,9 @@ public class EigenDecomposition {
            }
         }
 
-        if (Precision.equals(norm, 0.0)) {
-            // TODO: we can not handle a zero matrix, what exception to throw?
-           return;
+        // we can not handle a matrix with zero norm
+        if (Precision.equals(norm, 0.0, epsilon)) {
+           throw new MathArithmeticException(LocalizedFormats.ZERO_NORM);
         }
 
         // Backsubstitute to find vectors of upper triangular form

Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/SchurTransformer.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/SchurTransformer.java?rev=1363105&r1=1363104&r2=1363105&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/SchurTransformer.java
(original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math3/linear/SchurTransformer.java
Wed Jul 18 20:49:43 2012
@@ -436,13 +436,17 @@ class SchurTransformer {
      * Contains variable names as present in the original JAMA code.
      */
     private static class ShiftInfo {
-        /** TODO: describe */
+        // CHECKSTYLE: stop all
+
+        /** x shift info */
         double x;
-        /** TODO: describe */
+        /** y shift info */
         double y;
-        /** TODO: describe */
+        /** w shift info */
         double w;
         /** Indicates an exceptional shift. */
         double exShift;
+
+        // CHECKSTYLE: resume all
     }
 }

Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenDecompositionTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenDecompositionTest.java?rev=1363105&r1=1363104&r2=1363105&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenDecompositionTest.java
(original)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math3/linear/EigenDecompositionTest.java
Wed Jul 18 20:49:43 2012
@@ -21,6 +21,7 @@ import java.util.Arrays;
 import java.util.Random;
 
 
+import org.apache.commons.math3.distribution.NormalDistribution;
 import org.apache.commons.math3.util.FastMath;
 import org.apache.commons.math3.util.Precision;
 import org.junit.After;
@@ -38,7 +39,7 @@ public class EigenDecompositionTest {
         RealMatrix matrix =
             MatrixUtils.createRealMatrix(new double[][] { { 1.5 } });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(1.5, ed.getRealEigenvalue(0), 1.0e-15);
     }
 
@@ -50,7 +51,7 @@ public class EigenDecompositionTest {
                     { 12.0, 66.0 }
             });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(75.0, ed.getRealEigenvalue(0), 1.0e-15);
         Assert.assertEquals(50.0, ed.getRealEigenvalue(1), 1.0e-15);
     }
@@ -64,7 +65,7 @@ public class EigenDecompositionTest {
                                    { -16560.0,  7920.0,  17300.0 }
                                });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(50000.0, ed.getRealEigenvalue(0), 3.0e-11);
         Assert.assertEquals(12500.0, ed.getRealEigenvalue(1), 3.0e-11);
         Assert.assertEquals( 3125.0, ed.getRealEigenvalue(2), 3.0e-11);
@@ -79,7 +80,7 @@ public class EigenDecompositionTest {
                     { 15,   30,   45 }
             });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(70.0, ed.getRealEigenvalue(0), 3.0e-11);
         Assert.assertEquals(0.0,  ed.getRealEigenvalue(1), 3.0e-11);
         Assert.assertEquals(0.0,  ed.getRealEigenvalue(2), 3.0e-11);
@@ -95,7 +96,7 @@ public class EigenDecompositionTest {
                                    {  0.000,  0.000, -0.048,  0.136 }
                                });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(1.0, ed.getRealEigenvalue(0), 1.0e-15);
         Assert.assertEquals(0.4, ed.getRealEigenvalue(1), 1.0e-15);
         Assert.assertEquals(0.2, ed.getRealEigenvalue(2), 1.0e-15);
@@ -112,7 +113,7 @@ public class EigenDecompositionTest {
                                    { -0.2976,  0.1152, -0.1344,  0.3872 }
                                });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(1.0, ed.getRealEigenvalue(0), 1.0e-15);
         Assert.assertEquals(0.4, ed.getRealEigenvalue(1), 1.0e-15);
         Assert.assertEquals(0.2, ed.getRealEigenvalue(2), 1.0e-15);
@@ -144,9 +145,7 @@ public class EigenDecompositionTest {
         };
 
         EigenDecomposition decomposition;
-        decomposition = new EigenDecomposition(mainTridiagonal,
-                                                   secondaryTridiagonal,
-                                                   Precision.SAFE_MIN);
+        decomposition = new EigenDecomposition(mainTridiagonal, secondaryTridiagonal);
 
         double[] eigenValues = decomposition.getRealEigenvalues();
         for (int i = 0; i < refEigenValues.length; ++i) {
@@ -189,9 +188,7 @@ public class EigenDecompositionTest {
 
         // the following line triggers the exception
         EigenDecomposition decomposition;
-        decomposition = new EigenDecomposition(mainTridiagonal,
-                                                   secondaryTridiagonal,
-                                                   Precision.SAFE_MIN);
+        decomposition = new EigenDecomposition(mainTridiagonal, secondaryTridiagonal);
 
         double[] eigenValues = decomposition.getRealEigenvalues();
         for (int i = 0; i < refEigenValues.length; ++i) {
@@ -236,9 +233,7 @@ public class EigenDecompositionTest {
 
         // the following line triggers the exception
         EigenDecomposition decomposition;
-        decomposition = new EigenDecomposition(mainTridiagonal,
-                                                   secondaryTridiagonal,
-                                                   Precision.SAFE_MIN);
+        decomposition = new EigenDecomposition(mainTridiagonal, secondaryTridiagonal);
 
         double[] eigenValues = decomposition.getRealEigenvalues();
         for (int i = 0; i < refEigenValues.length; ++i) {
@@ -268,9 +263,7 @@ public class EigenDecompositionTest {
         TriDiagonalTransformer t =
             new TriDiagonalTransformer(createTestMatrix(r, ref));
         EigenDecomposition ed;
-        ed = new EigenDecomposition(t.getMainDiagonalRef(),
-                                        t.getSecondaryDiagonalRef(),
-                                        Precision.SAFE_MIN);
+        ed = new EigenDecomposition(t.getMainDiagonalRef(), t.getSecondaryDiagonalRef());
         double[] eigenValues = ed.getRealEigenvalues();
         Assert.assertEquals(ref.length, eigenValues.length);
         for (int i = 0; i < ref.length; ++i) {
@@ -284,7 +277,7 @@ public class EigenDecompositionTest {
     public void testDimensions() {
         final int m = matrix.getRowDimension();
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         Assert.assertEquals(m, ed.getV().getRowDimension());
         Assert.assertEquals(m, ed.getV().getColumnDimension());
         Assert.assertEquals(m, ed.getD().getColumnDimension());
@@ -297,7 +290,7 @@ public class EigenDecompositionTest {
     @Test
     public void testEigenvalues() {
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         double[] eigenValues = ed.getRealEigenvalues();
         Assert.assertEquals(refValues.length, eigenValues.length);
         for (int i = 0; i < refValues.length; ++i) {
@@ -315,8 +308,7 @@ public class EigenDecompositionTest {
         }
         Arrays.sort(bigValues);
         EigenDecomposition ed;
-        ed = new EigenDecomposition(createTestMatrix(r, bigValues),
-                                        Precision.SAFE_MIN);
+        ed = new EigenDecomposition(createTestMatrix(r, bigValues));
         double[] eigenValues = ed.getRealEigenvalues();
         Assert.assertEquals(bigValues.length, eigenValues.length);
         for (int i = 0; i < bigValues.length; ++i) {
@@ -333,7 +325,7 @@ public class EigenDecompositionTest {
         });
 
         EigenDecomposition ed;
-        ed = new EigenDecomposition(symmetric, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(symmetric);
         
         RealMatrix d = ed.getD();
         RealMatrix v = ed.getV();
@@ -341,8 +333,6 @@ public class EigenDecompositionTest {
 
         double norm = v.multiply(d).multiply(vT).subtract(symmetric).getNorm();
         Assert.assertEquals(0, norm, 6.0e-13);
-
-//           check(A.times(V),V.times(D));
     }
 
     @Test
@@ -374,8 +364,53 @@ public class EigenDecompositionTest {
         checkUnsymmetricMatrix(MatrixUtils.createRealMatrix(randData2));
     }
     
+    @Test
+    public void testRandomUnsymmetricMatrix() {
+        for (int run = 0; run < 100; run++) {
+            Random r = new Random(System.currentTimeMillis());
+
+            // matrix size
+            int size = r.nextInt(20) + 4;
+
+            double[][] data = new double[size][size];
+            for (int i = 0; i < size; i++) {
+                for (int j = 0; j < size; j++) {
+                    data[i][j] = r.nextInt(100);
+                }
+            }
+
+            RealMatrix m = MatrixUtils.createRealMatrix(data);
+            checkUnsymmetricMatrix(m);
+        }        
+    }
+    
+    @Test
+    public void testNormalDistributionUnsymmetricMatrix() {
+        for (int run = 0; run < 100; run++) {
+            Random r = new Random(System.currentTimeMillis());
+            NormalDistribution dist = new NormalDistribution(0.0, r.nextDouble() * 5);
+
+            // matrix size
+            int size = r.nextInt(20) + 4;
+
+            double[][] data = new double[size][size];
+            for (int i = 0; i < size; i++) {
+                for (int j = 0; j < size; j++) {
+                    data[i][j] = dist.sample();
+                }
+            }
+
+            RealMatrix m = MatrixUtils.createRealMatrix(data);
+            checkUnsymmetricMatrix(m);
+        }
+    }
+
+    /**
+     * Checks that the eigen decomposition of a general (unsymmetric) matrix is valid by
+     * checking: A*V = V*D
+     */
     private void checkUnsymmetricMatrix(final RealMatrix m) {
-        EigenDecomposition ed = new EigenDecomposition(m, Precision.SAFE_MIN);
+        EigenDecomposition ed = new EigenDecomposition(m);
         
         RealMatrix d = ed.getD();
         RealMatrix v = ed.getV();
@@ -389,14 +424,14 @@ public class EigenDecompositionTest {
         
         RealMatrix invV = new LUDecomposition(v).getSolver().getInverse();
         double norm = v.multiply(d).multiply(invV).subtract(m).getNorm();
-        Assert.assertEquals(0.0, norm, 6.0e-13);
+        Assert.assertEquals(0.0, norm, 1.0e-10);
     }
 
     /** test eigenvectors */
     @Test
     public void testEigenvectors() {
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         for (int i = 0; i < matrix.getRowDimension(); ++i) {
             double lambda = ed.getRealEigenvalue(i);
             RealVector v  = ed.getEigenvector(i);
@@ -409,7 +444,7 @@ public class EigenDecompositionTest {
     @Test
     public void testAEqualVDVt() {
         EigenDecomposition ed;
-        ed = new EigenDecomposition(matrix, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(matrix);
         RealMatrix v  = ed.getV();
         RealMatrix d  = ed.getD();
         RealMatrix vT = ed.getVT();
@@ -420,7 +455,7 @@ public class EigenDecompositionTest {
     /** test that V is orthogonal */
     @Test
     public void testVOrthogonal() {
-        RealMatrix v = new EigenDecomposition(matrix, Precision.SAFE_MIN).getV();
+        RealMatrix v = new EigenDecomposition(matrix).getV();
         RealMatrix vTv = v.transpose().multiply(v);
         RealMatrix id  = MatrixUtils.createRealIdentityMatrix(vTv.getRowDimension());
         Assert.assertEquals(0, vTv.subtract(id).getNorm(), 2.0e-13);
@@ -432,7 +467,7 @@ public class EigenDecompositionTest {
         double[] diagonal = new double[] { -3.0, -2.0, 2.0, 5.0 };
         RealMatrix m = createDiagonalMatrix(diagonal, diagonal.length, diagonal.length);
         EigenDecomposition ed;
-        ed = new EigenDecomposition(m, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(m);
         Assert.assertEquals(diagonal[0], ed.getRealEigenvalue(3), 2.0e-15);
         Assert.assertEquals(diagonal[1], ed.getRealEigenvalue(2), 2.0e-15);
         Assert.assertEquals(diagonal[2], ed.getRealEigenvalue(1), 2.0e-15);
@@ -450,7 +485,7 @@ public class EigenDecompositionTest {
                 {4,  2,  3}
         });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(repeated, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(repeated);
         checkEigenValues((new double[] {8, -1, -1}), ed, 1E-12);
         checkEigenVector((new double[] {2, 1, 2}), ed, 1E-12);
     }
@@ -466,7 +501,7 @@ public class EigenDecompositionTest {
                 {-4, -4, 8}
         });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(distinct, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(distinct);
         checkEigenValues((new double[] {2, 0, 12}), ed, 1E-12);
         checkEigenVector((new double[] {1, -1, 0}), ed, 1E-12);
         checkEigenVector((new double[] {1, 1, 1}), ed, 1E-12);
@@ -484,7 +519,7 @@ public class EigenDecompositionTest {
                 { -1.0,0.0, 1.0 }
         });
         EigenDecomposition ed;
-        ed = new EigenDecomposition(indefinite, Precision.SAFE_MIN);
+        ed = new EigenDecomposition(indefinite);
         checkEigenValues((new double[] {2, 1, -1}), ed, 1E-12);
         double isqrt3 = 1/FastMath.sqrt(3.0);
         checkEigenVector((new double[] {isqrt3,isqrt3,-isqrt3}), ed, 1E-12);



Mime
View raw message