[
https://issues.apache.org/jira/browse/MATH-1424?page=com.atlassian.jira.plugin.system.issuetabpanels:comment-tabpanel&focusedCommentId=16076370#comment-16076370
]
Gilles commented on MATH-1424:
------------------------------
Fine, since it's unlikely that we'll be in a position to consider a release of
"Commons Math" before we finalize the initial release of ["Commons
Numbers"|http://commons.apache.org/proper/commons-numbers/modules.html].
> Wrong eigen values computed by EigenDecomposition when the input matrix has
> large values
> ----------------------------------------------------------------------------------------
>
> Key: MATH-1424
> URL: https://issues.apache.org/jira/browse/MATH-1424
> Project: Commons Math
> Issue Type: Bug
> Affects Versions: 3.6.1
> Environment: JDK 7.51 64 bits on Windows 7.
> Reporter: Jerome
> Labels: easyfix
>
> The following code gives a wrong result:
> RealMatrix m = [[10_000_000.0, -1_000_000.0],[-1_000_000.1, 20_000_000.0]];
> // pseudo code
> EigenDecomposition ed = new EigenDecomposition(m);
> double[] eigenValues = ed.getRealEigenvalues();
> Computed values: [1.57E13, 1.57E13].
> Expected values: [1.0E7, 2.0E7]
> The problem lies in method EigenDecomposition.transformToSchur(RealMatrix).
> At line 758, the value matT[i+1][i] is checked against 0.0 within an EPSILON
> margin.
> If the precision of the computation were perfect, matT[i+1][i] == 0.0 means
> that matT[i][i] is a solution of the characteristic polynomial of m. In the
> other case there are 2 complex solutions.
> But due to imprecisions, this value can be different from 0.0 while m has
> only real solutions.
> The else part assume that the solutions are complex, which is wrong in the
> provided example.
> To correct it, you should resolve the 2 degree polynomial without assuming
> the solutions are complex (that is: test whether p*p + matT[i+1][i] *
> matT[i][i+1] is negative for 2 complex solutions, or positive or null for 2
> real solutions).
> You should also avoid testing values against something within epsilon margin,
> because this method is almost always wrong in some cases. At least, check
> with a margin that depends on the amplitude of the value (ex: margin =
> highest absolute value of the matrix * EPSILON); this is still wrong but
> problems will occur less often.
> The problem occurs when the input matrix has large values because matT has
> values of magnitude E7. MatT[1, 0] is really low (E-10) and you can not
> expect a better precision due to the large values on the diagonal.
> The test within EPSILON margin fails, which does not occurs when the input
> matrix has lowest values.
> Testing the code with m2 = m / pow(2, 20) will work, because matT[1, 0] is
> now low enough.
--
This message was sent by Atlassian JIRA
(v6.4.14#64029)