From 3ec7c559d2db8e66524435910ade24118a36ad42 Mon Sep 17 00:00:00 2001 From: Sam <40273116+Aweptimum@users.noreply.github.com> Date: Wed, 8 Nov 2023 13:21:08 -0600 Subject: [PATCH 1/2] Add successful rref case and its failing counterpart The matrix with the bigger magnitude passes, the same matrix scaled by its largest eigenvalue fails. Yet they should have the same rref. --- .../Reduction/ReducedRowEchelonFormTest.php | 38 +++++++++++++++++++ 1 file changed, 38 insertions(+) diff --git a/tests/LinearAlgebra/Reduction/ReducedRowEchelonFormTest.php b/tests/LinearAlgebra/Reduction/ReducedRowEchelonFormTest.php index bdb426758..240dba643 100644 --- a/tests/LinearAlgebra/Reduction/ReducedRowEchelonFormTest.php +++ b/tests/LinearAlgebra/Reduction/ReducedRowEchelonFormTest.php @@ -1045,6 +1045,44 @@ public function dataProviderForRref(): array [0, 0, 0, 0, 0], ], ], + // Large magnitude case + [ + [ + [-2.82878905765*10**6, 0, 0, 0, 0, 0], + [0, -2.82878905765*10**6, 0, 0, 0, 1729.7], + [0, 0, -2.82878905765*10**6, 0, -1729.7, 0], + [0, 0, 0, -2.82879105765*10**6, 0, 0], + [0, 0, -1729.7, 0, -1.05765, 0], + [0, 1729.7, 0, 0, 0, -1.05765] + ], + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [0, 0, 0, 0, 0, 1] + ] + ], + // Scaled version of large magnitude case + [ + [ + [-0.99999966687087, 0, 0, 0, 0, 0], + [0, -0.99999966687087, 0, 0, 0, 0.00061146285160793], + [0, 0, -0.99999966687087, 0, -0.00061146285160793, 0], + [0, 0, 0, -1.0000003738869, 0, 0], + [0, 0, -0.00061146285160793, 0, -3.7388694340557*10**-7, 0], + [0, 0.00061146285160793, 0, 0, 0, -3.7388694340557*10**-7] + ], + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [0, 0, 0, 0, 0, 1] + ] + ] ]; } From f3867713261d8ed7c011b13c25fd74cf1f1f95f5 Mon Sep 17 00:00:00 2001 From: Sam <40273116+Aweptimum@users.noreply.github.com> Date: Thu, 9 Nov 2023 13:51:05 -0600 Subject: [PATCH 2/2] REF: scale by smallest non-zero element Since we're technically manipulating a linear system, it seems that the property of the final result remaining the same despite scaling applies. If scaling by the inverse of the smallest element, the matrix error becomes relative to the magnitudes of the numbers rather than an absolute quantity. This intuition seems true because all of the tests pass --- src/LinearAlgebra/Reduction/RowEchelonForm.php | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/src/LinearAlgebra/Reduction/RowEchelonForm.php b/src/LinearAlgebra/Reduction/RowEchelonForm.php index 510c6d7d7..f2be68c8d 100644 --- a/src/LinearAlgebra/Reduction/RowEchelonForm.php +++ b/src/LinearAlgebra/Reduction/RowEchelonForm.php @@ -127,6 +127,19 @@ public static function gaussianElimination(NumericMatrix $A): array $swaps = 0; $ε = $A->getError(); + // Scale the matrix by its smallest non-zero element + $min = PHP_INT_MAX; + for ($i = 0; $i < $m; $i++) { + for ($j = 0; $j < $n; $j++) { + $elem = \abs($A[$i][$j]); + if (Support::isNotZero($elem, $ε) and $elem < $min) { + $min = $elem; + } + } + } + + $R = $A->scalarMultiply(1/$min)->getMatrix(); + for ($k = 0; $k < $size; $k++) { // Find column max $i_max = $k;