This is an automated email from the ASF dual-hosted git repository. asf-gitbox-commits pushed a commit to branch master in repository https://gitbox.apache.org/repos/asf/commons-math.git
commit 95cda75814c79e1af2fb414808486d45b612d508 Author: Kyon <[email protected]> AuthorDate: Mon Oct 19 12:14:24 2026 +0900 MATH-1424: Preserve converged real Schur boundaries Clear subdiagonals after accepted relative splits and real-pair rotations. Keep the existing convergence decisions and genuine complex blocks, with analytic, scaled and reconstruction regressions. Generated-by: OpenAI Codex --- .../math4/legacy/linear/SchurTransformer.java | 4 ++ .../legacy/linear/EigenDecompositionTest.java | 75 ++++++++++++++++++++++ src/changes/changes.xml | 3 + 3 files changed, 82 insertions(+) diff --git a/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/linear/SchurTransformer.java b/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/linear/SchurTransformer.java index 25373baff..f45059277 100644 --- a/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/linear/SchurTransformer.java +++ b/commons-math-legacy/src/main/java/org/apache/commons/math4/legacy/linear/SchurTransformer.java @@ -193,6 +193,8 @@ class SchurTransformer { matrixP[i][iu - 1] = q * z + p * matrixP[i][iu]; matrixP[i][iu] = q * matrixP[i][iu] - p * z; } + // The converged real pair is triangular. + matrixT[iu][iu - 1] = 0.0; } iu -= 2; iteration = 0; @@ -246,6 +248,8 @@ class SchurTransformer { s = norm; } if (JdkMath.abs(matrixT[l][l - 1]) < epsilon * s) { + // Make the accepted split explicit in the Schur form. + matrixT[l][l - 1] = 0.0; break; } l--; diff --git a/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/linear/EigenDecompositionTest.java b/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/linear/EigenDecompositionTest.java index 735dbc0fb..88aee9e95 100644 --- a/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/linear/EigenDecompositionTest.java +++ b/commons-math-legacy/src/test/java/org/apache/commons/math4/legacy/linear/EigenDecompositionTest.java @@ -60,6 +60,81 @@ public class EigenDecompositionTest { Assert.assertEquals(50.0, ed.getRealEigenvalue(1), 1.0e-15); } + @Test + public void testMath1424FiniteUnsymmetricScale() { + final double[][] data = {{1e7, -1e6}, {-1000000.1, 2e7}}; + // The two roots are (a + d) / 2 +/- sqrt(((a - d) / 2)^2 + b * c). + final double radius = JdkMath.sqrt(5e6 * 5e6 + 1e6 * 1000000.1); + final double[] expected = {15e6 - radius, 15e6 + radius}; + for (final double scale : new double[] {0x1.0p-20, 1.0}) { + final RealMatrix input = MatrixUtils.createRealMatrix(data).scalarMultiply(scale); + checkMath1424SchurResiduals(input); + final EigenDecomposition decomposition = new EigenDecomposition(input); + final double[] real = decomposition.getRealEigenvalues(); + final double[] imaginary = decomposition.getImagEigenvalues(); + Arrays.sort(real); + final String message = "scale=" + scale + ", real=" + Arrays.toString(real) + + ", imaginary=" + Arrays.toString(imaginary); + for (final double value : real) { + Assert.assertTrue(message, Double.isFinite(value)); + } + // Infinity norm of the input times a relative error allowance. + final double tolerance = 1e-12 * 21000000.1 * scale; + Assert.assertArrayEquals(message, + new double[] {expected[0] * scale, expected[1] * scale}, + real, tolerance); + Assert.assertArrayEquals(message, new double[] {0.0, 0.0}, imaginary, 0.0); + Assert.assertFalse(message, decomposition.hasComplexEigenvalues()); + } + } + + @Test + public void testMath1424ComplexBlockAndSmallSplit() { + for (final double boundary : new double[] {0.0, 1e-4}) { + final RealMatrix input = MatrixUtils.createRealMatrix(new double[][] { + {1e13, 0, 0}, {boundary, 2, 1}, {0, -1, 2} + }); + final RealMatrix t = checkMath1424SchurResiduals(input); + Assert.assertEquals(0.0, t.getEntry(1, 0), 0.0); + Assert.assertEquals(1.0, JdkMath.abs(t.getEntry(2, 1)), 1e-14); + final EigenDecomposition decomposition = new EigenDecomposition(input); + final double[] real = decomposition.getRealEigenvalues(); + final double[] imaginary = decomposition.getImagEigenvalues(); + Arrays.sort(real); + Arrays.sort(imaginary); + Assert.assertArrayEquals(new double[] {2, 2, 1e13}, real, 1e-12); + Assert.assertArrayEquals(new double[] {-1, 0, 1}, imaginary, 1e-14); + Assert.assertTrue(decomposition.hasComplexEigenvalues()); + } + } + + @Test + public void testMath1424SmallOneRootBoundary() { + final RealMatrix input = MatrixUtils.createRealMatrix(new double[][] { + {1e13, 1}, {1e-4, 2e13} + }); + Assert.assertEquals(0.0, checkMath1424SchurResiduals(input).getEntry(1, 0), 0.0); + final EigenDecomposition decomposition = new EigenDecomposition(input); + final double[] real = decomposition.getRealEigenvalues(); + Arrays.sort(real); + Assert.assertArrayEquals(new double[] {1e13, 2e13}, real, 1e-12 * input.getNorm()); + Assert.assertArrayEquals(new double[] {0, 0}, decomposition.getImagEigenvalues(), 0.0); + Assert.assertFalse(decomposition.hasComplexEigenvalues()); + } + + private RealMatrix checkMath1424SchurResiduals(final RealMatrix input) { + final SchurTransformer transformer = new SchurTransformer(input); + final RealMatrix p = transformer.getP(); + final RealMatrix t = transformer.getT(); + final double pNorm = p.getNorm(); + Assert.assertEquals(0.0, input.multiply(p).subtract(p.multiply(t)).getNorm() + / (input.getNorm() * pNorm), 1e-14); + final RealMatrix identity = MatrixUtils.createRealIdentityMatrix(input.getRowDimension()); + Assert.assertEquals(0.0, transformer.getPT().multiply(p).subtract(identity).getNorm() + / (pNorm * pNorm), 1e-14); + return t; + } + @Test public void testDimension3() { RealMatrix matrix = diff --git a/src/changes/changes.xml b/src/changes/changes.xml index dbf868654..84e021b55 100644 --- a/src/changes/changes.xml +++ b/src/changes/changes.xml @@ -96,6 +96,9 @@ Caveat: to support the whole codebase (it was one of the main reasons for creating more focused components). "> + <action type="fix" issue="MATH-1424"> + Preserve converged real eigenvalues by zeroing accepted Schur split boundaries. + </action> <action dev="aherbetrt" type="fix" issue="MATH-1692" due-to="Alex Herbert"> "AccurateMath": Avoid intermediate overflow in asinh and acosh. </action>
