Skip to content

Test and fix Singular Value Decomposition (double) at different scales#1146

Open
lellid wants to merge 1 commit intomathnet:masterfrom
lellid:TestAndFixSvdAtDifferentScales
Open

Test and fix Singular Value Decomposition (double) at different scales#1146
lellid wants to merge 1 commit intomathnet:masterfrom
lellid:TestAndFixSvdAtDifferentScales

Conversation

@lellid
Copy link

@lellid lellid commented Jan 18, 2026

Currently, Singular value decomposition gives wrong results if the numbers in the matrix are very small (1E-16 or smaller).

This PR adds a test for SVD at different scales (1, 1E40, 1E-40) and compares both the reconstructed matrix with the original matrix as well as the results with high precision results generated by Mathematica.
Svd (only the double version) was fixed by using an own version of the function AlmostEqualRelative, which provides a really relative comparison method.
In contrast, the method AlmostEqualRelative in Precision.Equality calls AlmostEqualNormRelative, and in this method, there is a problematic section which decided on absolute differences instead of relative differences:

            // If one is almost zero, fall back to absolute equality
            if (Math.Abs(a) < DoublePrecision || Math.Abs(b) < DoublePrecision)
            {
                // The values are equal if the difference between the two numbers is smaller than
                // 10^(-numberOfDecimalPlaces). We divide by two so that we have half the range
                // on each side of the numbers, e.g. if decimalPlaces == 2,
                // then 0.01 will equal between 0.005 and 0.015, but not 0.02 and not 0.00
                return Math.Abs(diff) < Pow10(-decimalPlaces) * 0.5;
            }

The absolute comparsion is OK if the numbers are in the order of 1, which is the case for U and VT in SVD, but not for W.

Since AreAlmostRelative is abundantly used both in code and test, I suspect that more functions (and tests) are affected. If this PR is accepted, I can also add a test and fix the float version of SVD.

Dirk

@lellid
Copy link
Author

lellid commented Jan 18, 2026

One thing to add: I set the relative accuracy in the calls to AlmostEqualRelative to 4E-15 instead of 1E-15. This was because with 1E-15, one of the tests failed (but this test is already marked as questionable). But with 4E-15, we will sacrifice some of the accuracy of the SVD. Thus consider to set the relative accuracy to 1E-15 in the two calls of AlmostEqualRelative in the SVD, and have a look at the test that failed.

Dirk

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant