MATH-1045

Singular matrices were considered non-singular due to strict comparison
with zero. Reported and fixed by Sean Owen.


git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1536766 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Gilles Sadowski 2013-10-29 15:43:04 +00:00
parent 53f55d2ca0
commit a4ffd39333
4 changed files with 46 additions and 2 deletions

View File

@ -252,6 +252,9 @@
<contributor>
<name>Fredrik Norin</name>
</contributor>
<contributor>
<name>Sean Owen</name>
</contributor>
<contributor>
<name>Sujit Pal</name>
</contributor>

View File

@ -51,6 +51,9 @@ If the output is not quite correct, check for invisible trailing spaces!
</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

View File

@ -513,15 +513,32 @@ 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;
}
}
return true;
}
/**
* @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.
*

View File

@ -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() {