You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by bi...@apache.org on 2009/10/19 01:29:37 UTC
svn commit: r826550 - in /commons/proper/math/trunk/src:
main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
Author: billbarker
Date: Sun Oct 18 23:29:37 2009
New Revision: 826550
URL: http://svn.apache.org/viewvc?rev=826550&view=rev
Log:
Fix for possible zero divide on an indefinite matrix
Fix for: MATH-297
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
Modified: commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java?rev=826550&r1=826549&r2=826550&view=diff
==============================================================================
--- commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java (original)
+++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/linear/EigenDecompositionImpl.java Sun Oct 18 23:29:37 2009
@@ -1699,19 +1699,21 @@
// perform an initial non-shifted LDLt decomposition
final double[] d = new double[m];
final double[] l = new double[m - 1];
- double di = main[0];
+ // avoid zero divide on indefinite matrix
+ final double mu = realEigenvalues[m-1] <= 0 && realEigenvalues[0] > 0 ? 0.5-realEigenvalues[m-1] : 0;
+ double di = main[0]+mu;
d[0] = di;
for (int i = 1; i < m; ++i) {
final double eiM1 = secondary[i - 1];
final double ratio = eiM1 / di;
- di = main[i] - eiM1 * ratio;
+ di = main[i] - eiM1 * ratio + mu;
l[i - 1] = ratio;
d[i] = di;
}
// compute eigenvectors
for (int i = 0; i < m; ++i) {
- eigenvectors[i] = findEigenvector(realEigenvalues[i], d, l);
+ eigenvectors[i] = findEigenvector(realEigenvalues[i]+mu, d, l);
}
}
Modified: commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java
URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java?rev=826550&r1=826549&r2=826550&view=diff
==============================================================================
--- commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java (original)
+++ commons/proper/math/trunk/src/test/java/org/apache/commons/math/linear/EigenDecompositionImplTest.java Sun Oct 18 23:29:37 2009
@@ -244,6 +244,24 @@
}
/**
+ * Verifies operation on indefinite matrix
+ */
+ public void testZeroDivide() {
+ RealMatrix indefinite = MatrixUtils.createRealMatrix(new double [][] {
+ { 0.0, 1.0, -1.0 },
+ { 1.0, 1.0, 0.0 },
+ { -1.0,0.0, 1.0 }
+ });
+ EigenDecomposition ed = new EigenDecompositionImpl(indefinite, MathUtils.SAFE_MIN);
+ checkEigenValues((new double[] {2, 1, -1}), ed, 1E-12);
+ double isqrt3 = 1/Math.sqrt(3.0);
+ checkEigenVector((new double[] {isqrt3,isqrt3,-isqrt3}), ed, 1E-12);
+ double isqrt2 = 1/Math.sqrt(2.0);
+ checkEigenVector((new double[] {0.0,-isqrt2,-isqrt2}), ed, 1E-12);
+ double isqrt6 = 1/Math.sqrt(6.0);
+ checkEigenVector((new double[] {2*isqrt6,-isqrt6,isqrt6}), ed, 1E-12);
+ }
+ /**
* Verifies that the given EigenDecomposition has eigenvalues equivalent to
* the targetValues, ignoring the order of the values and allowing
* values to differ by tolerance.
@@ -257,6 +275,7 @@
}
}
+
/**
* Returns true iff there is an entry within tolerance of value in
* searchArray.