diff --git a/0-source/Just.PreciseMath/DDMath.cs b/0-source/Just.PreciseMath/DDMath.cs index 8b0abc5..810c69d 100644 --- a/0-source/Just.PreciseMath/DDMath.cs +++ b/0-source/Just.PreciseMath/DDMath.cs @@ -36,16 +36,7 @@ public static partial class DDMath return PreciseMathHelper.DivideBoundary(one, value); } - double quotient = one / value._high; - double remainder = Math.FusedMultiplyAdd(-quotient, value._high, one); - remainder = Math.FusedMultiplyAdd(-quotient, value._low, remainder); - double correction = remainder / value._high; - // For normalized input, using the high denominator in the correction - // adds only O(u^2) error, with u=2^-53. The denominator exponent guard - // keeps the quotient normal and finite; it dominates the correction, - // so QuickTwoSum is ordered and its sum remains finite. - (double high, double low) = PreciseMathHelper.TwoQuickAdd(quotient, correction); - return new DoubleDouble(high, low == 0.0 ? 0.0 : low); + return PreciseMathHelper.DivideScalarFinite(one, value); } /// diff --git a/0-source/Just.PreciseMath/DoubleDouble.Arithmetic.cs b/0-source/Just.PreciseMath/DoubleDouble.Arithmetic.cs index ca86a5f..37085ed 100644 --- a/0-source/Just.PreciseMath/DoubleDouble.Arithmetic.cs +++ b/0-source/Just.PreciseMath/DoubleDouble.Arithmetic.cs @@ -103,24 +103,28 @@ public readonly partial struct DoubleDouble : } /// Applies the expansion operation without discarding the low component. + /// Preserves the component bits of addition with the equivalent zero-low scalar expansion. public static DoubleDouble operator +(DoubleDouble left, double right) { return PreciseMathHelper.AddScalar(left._high, left._low, right); } /// Applies the expansion operation without discarding the low component. + /// Preserves the component bits of addition with the equivalent zero-low scalar expansion. public static DoubleDouble operator +(double left, DoubleDouble right) { return right + left; } /// Applies the expansion operation without discarding the low component. + /// Preserves the component bits of subtraction with the equivalent zero-low scalar expansion. public static DoubleDouble operator -(DoubleDouble left, double right) { return PreciseMathHelper.AddScalar(left._high, left._low, -right); } /// Applies the expansion operation without discarding the low component. + /// Preserves the component bits of subtraction with the equivalent zero-low scalar expansion. public static DoubleDouble operator -(double left, DoubleDouble right) { // Negate the components, not the result: exact cancellation must yield +0. @@ -153,7 +157,8 @@ public readonly partial struct DoubleDouble : return right * left; } - /// Applies the expansion operation without discarding the low component. + /// Divides by a scalar using a quotient estimate and two residual corrections. + /// Preserves the component bits of division by the equivalent zero-low expansion. public static DoubleDouble operator /(DoubleDouble left, double right) { if (!IsFinite(left) || !double.IsFinite(right) || left._high == 0.0 || right == 0.0) @@ -166,15 +171,24 @@ public readonly partial struct DoubleDouble : } double quotient = left._high / right; - double remainder = Math.FusedMultiplyAdd(-quotient, right, left._high); - double correction = (remainder + left._low) / right; - // The exponent guard keeps the quotient normal. The correction is - // O(u * quotient), so QuickTwoSum is ordered; one correction gives O(u^2) error. - (double high, double low) = PreciseMathHelper.TwoQuickAdd(quotient, correction); - return new DoubleDouble(high, low == 0.0 ? 0.0 : low); + // The initial product is near left.High, so the exponent guard makes + // TwoMultiply exact. Keep the complete remainder, including its low. + (double product, double productError) = PreciseMathHelper.TwoMultiply(right, quotient); + DoubleDouble remainder = PreciseMathHelper.AddFinite(left._high, left._low, -product, -productError); + double correction = remainder._high / right; + (double correctionProduct, double correctionError) = PreciseMathHelper.TwoMultiply(right, correction); + // This product can be subnormal. FMA rounds its residual just as the + // exact boundary product does; normalize and canonicalize its zero low. + DoubleDouble correctedProduct = PreciseMathHelper.NormalizeFinite(correctionProduct, correctionError); + double finalRemainder = PreciseMathHelper.SubtractDivisionCorrectionHigh(remainder, correctedProduct); + double finalCorrection = finalRemainder / right; + // The same finite bounds as DD/DD apply. Normalize the first correction + // before adding the last, preserving the DD/DD evaluation order. + return PreciseMathHelper.NormalizeFinite(quotient, correction) + finalCorrection; } - /// Applies the expansion operation without discarding the low component. + /// Divides a scalar using a quotient estimate and two residual corrections. + /// Preserves the component bits of division with the equivalent zero-low numerator. public static DoubleDouble operator /(double left, DoubleDouble right) { if (!double.IsFinite(left) || !IsFinite(right) || left == 0.0 || right._high == 0.0) @@ -186,13 +200,6 @@ public readonly partial struct DoubleDouble : return PreciseMathHelper.DivideBoundary(left, right); } - double quotient = left / right._high; - double remainder = Math.FusedMultiplyAdd(-quotient, right._high, left); - remainder = Math.FusedMultiplyAdd(-quotient, right._low, remainder); - double correction = remainder / right._high; - // Using the high denominator in the correction adds only O(u^2) error. - // As above, the guarded quotient dominates its correction in magnitude. - (double high, double low) = PreciseMathHelper.TwoQuickAdd(quotient, correction); - return new DoubleDouble(high, low == 0.0 ? 0.0 : low); + return PreciseMathHelper.DivideScalarFinite(left, right); } } diff --git a/0-source/Just.PreciseMath/PreciseMathHelper.cs b/0-source/Just.PreciseMath/PreciseMathHelper.cs index 28b4092..b31b79e 100644 --- a/0-source/Just.PreciseMath/PreciseMathHelper.cs +++ b/0-source/Just.PreciseMath/PreciseMathHelper.cs @@ -104,7 +104,25 @@ internal static class PreciseMathHelper return NormalizeFinite(sum, sumError + (middleError + lowError)); } - // Only for DD division's second remainder: normalized finite inputs with highs + // Finite, nonzero numerator and normalized denominator, both with high + // exponents in [-450, 450]. Shared with Reciprocal after its constant-numerator + // checks. Do not promote left or dispatch to the DD/DD division operator. + internal static DoubleDouble DivideScalarFinite(double left, DoubleDouble right) + { + double quotient = left / right._high; + DoubleDouble product = right * quotient; + // Preserve DD/DD's complete first remainder. A chained scalar FMA drops + // residual bits needed by the second correction; retain the expansion + // subtraction before computing the next quotient correction. + DoubleDouble remainder = AddFinite(left, 0.0, -product._high, -product._low); + double correction = remainder._high / right._high; + double finalRemainder = SubtractDivisionCorrectionHigh(remainder, right * correction); + double finalCorrection = finalRemainder / right._high; + // DD/DD's guard-derived finite bounds and normalization order apply. + return NormalizeFinite(quotient, correction) + finalCorrection; + } + + // Only for division's second remainder: normalized finite inputs with highs // bounded by 2^462, including the correction product's exact fallback results. // Same sign and normal binade imply Sterbenz-exact high subtraction. Its error // is +0; with canonical input lows, TwoAdd(+0, lowSum) is (lowSum, +0). @@ -164,11 +182,14 @@ internal static class PreciseMathHelper } (double sum, double error) = TwoAdd(high, value); - // Near high-component cancellation, Sterbenz makes the first sum exact, - // so error is zero and this retains low exactly. Otherwise its rounding - // contributes only O(u^2) relative error. The final TwoSum normalizes. - (double result, double residual) = TwoAdd(sum, error + low); - return new DoubleDouble(result, residual == 0.0 ? 0.0 : residual); + // Specialize AddFinite for a zero-low scalar: adding that zero to low + // has no residual. Retain the rounding error when combining the high + // sum's error with low, then fold it into the final normalization. + // TwoSum also absorbs a negated zero low from scalar-left subtraction. + (double middle, double middleError) = TwoAdd(error, low); + (double result, double residual) = TwoAdd(sum, middle); + // The unchanged addition guard gives AddFinite's finite-sum bounds. + return NormalizeFinite(result, residual + middleError); } // Keep the complete BigInteger expressions out of ordinary arithmetic bodies, diff --git a/1-tests/Directory.Build.props b/1-tests/Directory.Build.props index 41c5032..a9288b0 100644 --- a/1-tests/Directory.Build.props +++ b/1-tests/Directory.Build.props @@ -18,4 +18,9 @@ + + + + + diff --git a/1-tests/Just.PreciseMath.Tests/ArithmeticRangeTests.cs b/1-tests/Just.PreciseMath.Tests/ArithmeticRangeTests.cs index 59fe0a2..c500019 100644 --- a/1-tests/Just.PreciseMath.Tests/ArithmeticRangeTests.cs +++ b/1-tests/Just.PreciseMath.Tests/ArithmeticRangeTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class ArithmeticRangeTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleArithmeticTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleArithmeticTests.cs index d93e405..ce09c43 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleArithmeticTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleArithmeticTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleArithmeticTests @@ -332,6 +328,53 @@ public class DoubleDoubleArithmeticTests Check(6.0 / new DoubleDouble(2.0), 3.0, 0.0); } + [Theory] + [InlineData(2.0)] + [InlineData(3.0)] + public void ScalarAdditionAndSubtractionRetainTheMiddleSumResidual(double magnitude) + { + // Exact Fraction arithmetic on the stored E components confirms that + // all signed E +/- 2 and E +/- 3 results fit exactly in two components. + // Compare exact dyadic sums, not another overload or a rounded tolerance. + foreach (double valueSign in new[] { -1.0, 1.0 }) + { + foreach (double scalarSign in new[] { -1.0, 1.0 }) + { + DoubleDouble value = valueSign * DoubleDouble.E; + double scalar = scalarSign * magnitude; + BigInteger expectedSum = Units(value) + Units(scalar); + BigInteger expectedDifference = Units(value) - Units(scalar); + Units(value + scalar).ShouldBe(expectedSum); + Units(scalar + value).ShouldBe(expectedSum); + Units(value - scalar).ShouldBe(expectedDifference); + Units(scalar - value).ShouldBe(-expectedDifference); + } + } + } + + [Theory] + [InlineData("+")] + [InlineData("-")] + public void ScalarSumsAndDifferencesMatchExpansionBitsAcrossRanges(string operation) + { + // Bitwise compatibility; independent exact-rational accuracy checks are separate. + (double[] scalars, DoubleDouble[] values) = ScalarArithmeticCompatibilityCases(); + foreach (DoubleDouble value in values) + { + foreach (double scalar in scalars) + { + DoubleDouble expanded = new(scalar); + DoubleDouble expected = operation == "+" ? value + expanded : value - expanded; + DoubleDouble actual = operation == "+" ? value + scalar : value - scalar; + CheckBoundary(actual, expected.High, expected.Low); + + DoubleDouble expectedReverse = operation == "+" ? expanded + value : expanded - value; + DoubleDouble actualReverse = operation == "+" ? scalar + value : scalar - value; + CheckBoundary(actualReverse, expectedReverse.High, expectedReverse.Low); + } + } + } + [Fact] public void ExpansionCancellationRetainsExactComponentsAcrossTheAdditionGuard() { @@ -579,6 +622,25 @@ public class DoubleDoubleArithmeticTests } } + [Fact] + public void DivisionOverloadsPreserveComponentBitsAcrossRangesAndSpecialValues() + { + // Compatibility, not an accuracy oracle. GeneralArithmeticMeetsExactRationalBound + // separately checks complete results against exact rational bounds. + (double[] scalars, DoubleDouble[] values) = ScalarArithmeticCompatibilityCases(); + foreach (DoubleDouble value in values) + { + foreach (double scalar in scalars) + { + DoubleDouble expectedQuotient = value / new DoubleDouble(scalar); + CheckBoundary(value / scalar, expectedQuotient.High, expectedQuotient.Low); + + DoubleDouble expectedReverseQuotient = new DoubleDouble(scalar) / value; + CheckBoundary(scalar / value, expectedReverseQuotient.High, expectedReverseQuotient.Low); + } + } + } + [Fact] public void DivisionRetainsSubnormalCorrectionsWithOrdinaryHighComponents() { @@ -597,6 +659,38 @@ public class DoubleDoubleArithmeticTests } } + [Theory] + [InlineData(-1000)] + [InlineData(-1021)] + [InlineData(-1022)] + [InlineData(-1023)] + public void ScalarDivisionPreservesSparseCorrectionProductsNearUnderflow(int lowExponent) + { + // At -1000, the second product is normal but its exact residual is + // (962132647665515 / 1073741824) * epsilon, which rounds to 896056 * epsilon. + // Derived from exact binary64 rational products; the other rows bracket + // the correction product's normal/subnormal transition. + foreach (double numeratorSign in new[] { -1.0, 1.0 }) + { + foreach (double denominatorSign in new[] { -1.0, 1.0 }) + { + foreach (double lowSign in new[] { -1.0, 1.0 }) + { + double low = Math.ScaleB(lowSign, lowExponent); + double denominator = denominatorSign * 1.1; + DoubleDouble numerator = DoubleDouble.FromComponents(numeratorSign * 1.1, low); + DoubleDouble actual = numerator / denominator; + DoubleDouble expanded = numerator / new DoubleDouble(denominator); + CheckBoundary(actual, expanded.High, expanded.Low); + // Exact quotient = +/-1 + low/denominator. Its high cannot + // change here; binary64 division independently rounds the low. + CheckBoundary(actual, numeratorSign / denominatorSign, low / denominator); + AssertRelative(actual, Units(numerator) << 1074, Units(denominator)); + } + } + } + } + [Fact] public void SubnormalProductsWithLargeNormalsRetainExactResultsInBothOrders() { @@ -884,6 +978,39 @@ public class DoubleDoubleArithmeticTests } } + private static (double[] Scalars, DoubleDouble[] Values) ScalarArithmeticCompatibilityCases() + { + List scalars = [0.0, -0.0, 1.0, -1.0, 2.0, -2.0, 3.0, -3.0, 100.0, -100.0, + double.MaxValue, double.MinValue, double.PositiveInfinity, + double.NegativeInfinity, double.NaN]; + List values = [DoubleDouble.E, -DoubleDouble.E]; + // Include addition/division dispatch edges, exponent extremes, dense/sparse + // lows of either sign, and significands on either side of a binade. + int[] exponents = [-1074, -1022, -451, -450, -1, 0, 1, 450, 451, 1020, 1021, 1023]; + double[] significands = [1.0, Math.BitIncrement(1.0), 1.5, Math.BitDecrement(2.0)]; + foreach (int exponent in exponents) + { + foreach (double significand in significands) + { + foreach (double sign in new[] { -1.0, 1.0 }) + { + double high = Math.ScaleB(sign * significand, exponent); + scalars.Add(high); + double denseLow = Math.ScaleB(1.0, exponent - 54); + foreach (double low in new[] { 0.0, denseLow, -denseLow, double.Epsilon, -double.Epsilon }) + { + values.Add(DoubleDouble.FromComponents(high, low)); + } + } + } + } + foreach (double scalar in scalars) + { + values.Add(new DoubleDouble(scalar)); + } + return (scalars.ToArray(), values.ToArray()); + } + private static void AssertDivisionMatchesPrevious(DoubleDouble left, DoubleDouble right) { // Freeze the pre-specialization expression. The unchanged public operators diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleBoundaryTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleBoundaryTests.cs index 0f7e90d..fadfc95 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleBoundaryTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleBoundaryTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleBoundaryTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleCbrtTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleCbrtTests.cs index 7e3cd8e..f8d999f 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleCbrtTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleCbrtTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleCbrtTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleComparisonTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleComparisonTests.cs index 022696a..11a6b5a 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleComparisonTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleComparisonTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleComparisonTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleConversionTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleConversionTests.cs index 3867bcc..e7d3d9e 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleConversionTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleConversionTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleConversionTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleExponentialFunctionsTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleExponentialFunctionsTests.cs index fe45528..7e90f0e 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleExponentialFunctionsTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleExponentialFunctionsTests.cs @@ -1,8 +1,5 @@ using System.Globalization; -using System.Numerics; using Just.PreciseMath.Tests.ReferenceData; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleFormattingTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleFormattingTests.cs index 0efb243..fab29ed 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleFormattingTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleFormattingTests.cs @@ -1,6 +1,4 @@ using System.Globalization; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleHypotTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleHypotTests.cs index e455895..d228c6b 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleHypotTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleHypotTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleHypotTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleNumberStylesTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleNumberStylesTests.cs index c8a10bd..59e8a9e 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleNumberStylesTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleNumberStylesTests.cs @@ -1,7 +1,4 @@ using System.Globalization; -using System.Numerics; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleParsingTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleParsingTests.cs index 0140ec2..20e2dee 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleParsingTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleParsingTests.cs @@ -1,7 +1,4 @@ using System.Globalization; -using System.Numerics; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRepresentationTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRepresentationTests.cs index 695dc78..d4385f2 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRepresentationTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRepresentationTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleRepresentationTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootFunctionsTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootFunctionsTests.cs index d2929f0..54bbaee 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootFunctionsTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootFunctionsTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleRootFunctionsTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootNTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootNTests.cs index c7f996e..1a4d8ed 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootNTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleRootNTests.cs @@ -1,8 +1,5 @@ using System.Globalization; -using System.Numerics; using Just.PreciseMath.Tests.ReferenceData; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSignedNumberTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSignedNumberTests.cs index f3895f8..a0fef1d 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSignedNumberTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSignedNumberTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleSignedNumberTests diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpanFormattingTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpanFormattingTests.cs index 4d6ecc0..27aa5ae 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpanFormattingTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpanFormattingTests.cs @@ -1,6 +1,4 @@ using System.Globalization; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpecialValueTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpecialValueTests.cs index 7134170..8a1ba53 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpecialValueTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleSpecialValueTests.cs @@ -1,7 +1,4 @@ using System.Globalization; -using System.Numerics; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/DoubleDoubleTests.cs b/1-tests/Just.PreciseMath.Tests/DoubleDoubleTests.cs index 790f0e7..7be04c4 100644 --- a/1-tests/Just.PreciseMath.Tests/DoubleDoubleTests.cs +++ b/1-tests/Just.PreciseMath.Tests/DoubleDoubleTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class DoubleDoubleTests diff --git a/1-tests/Just.PreciseMath.Tests/GenericConversionTests.cs b/1-tests/Just.PreciseMath.Tests/GenericConversionTests.cs index 5c9d7b6..189306a 100644 --- a/1-tests/Just.PreciseMath.Tests/GenericConversionTests.cs +++ b/1-tests/Just.PreciseMath.Tests/GenericConversionTests.cs @@ -1,7 +1,4 @@ using System.Diagnostics.CodeAnalysis; -using System.Numerics; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathExpTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathExpTests.cs index 4e7d4fb..9b9a3dd 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathExpTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathExpTests.cs @@ -1,8 +1,5 @@ using System.Globalization; -using System.Numerics; using Just.PreciseMath.Tests.ReferenceData; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathInvSqrtTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathInvSqrtTests.cs index c029044..0f26469 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathInvSqrtTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathInvSqrtTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathInvSqrtTests diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathLogTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathLogTests.cs index 788b113..94c3308 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathLogTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathLogTests.cs @@ -1,8 +1,5 @@ using System.Globalization; -using System.Numerics; using Just.PreciseMath.Tests.ReferenceData; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathPowTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathPowTests.cs index 3670a2f..9f76dc6 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathPowTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathPowTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathPowTests diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowSpecialTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowSpecialTests.cs index 5a2afac..08af8a8 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowSpecialTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowSpecialTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathRealPowSpecialTests diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowTests.cs index d67dd17..607c5a4 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathRealPowTests.cs @@ -1,8 +1,5 @@ using System.Globalization; -using System.Numerics; using Just.PreciseMath.Tests.ReferenceData; -using Shouldly; -using Xunit; namespace Just.PreciseMath.Tests; diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathReciprocalTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathReciprocalTests.cs index ab63d4a..ded089b 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathReciprocalTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathReciprocalTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathReciprocalTests diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathSqrtTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathSqrtTests.cs index f48c1f2..94fbcb8 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathSqrtTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathSqrtTests.cs @@ -1,7 +1,3 @@ -using System.Numerics; -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathSqrtTests diff --git a/1-tests/Just.PreciseMath.Tests/PreciseMathTests.cs b/1-tests/Just.PreciseMath.Tests/PreciseMathTests.cs index 66ca109..ba030f8 100644 --- a/1-tests/Just.PreciseMath.Tests/PreciseMathTests.cs +++ b/1-tests/Just.PreciseMath.Tests/PreciseMathTests.cs @@ -1,6 +1,3 @@ -using Shouldly; -using Xunit; - namespace Just.PreciseMath.Tests; public class PreciseMathTests diff --git a/1-tests/Just.PreciseMath.Tests/ReferenceData/ExpReferenceData.cs b/1-tests/Just.PreciseMath.Tests/ReferenceData/ExpReferenceData.cs index 14e5238..132e5cc 100644 --- a/1-tests/Just.PreciseMath.Tests/ReferenceData/ExpReferenceData.cs +++ b/1-tests/Just.PreciseMath.Tests/ReferenceData/ExpReferenceData.cs @@ -1,5 +1,3 @@ -using Xunit; - namespace Just.PreciseMath.Tests.ReferenceData; // Generated by generate_exp.py; do not derive expected values from DD arithmetic. diff --git a/1-tests/Just.PreciseMath.Tests/ReferenceData/ExponentialFunctionsReferenceData.cs b/1-tests/Just.PreciseMath.Tests/ReferenceData/ExponentialFunctionsReferenceData.cs index 1de9ca4..ed22dd5 100644 --- a/1-tests/Just.PreciseMath.Tests/ReferenceData/ExponentialFunctionsReferenceData.cs +++ b/1-tests/Just.PreciseMath.Tests/ReferenceData/ExponentialFunctionsReferenceData.cs @@ -9,8 +9,6 @@ // absolute; reference rounding uncertainty is bounded by 2^-350 relative. // Tiny corrections to 1/-1 below 120 digits are not component-retention oracles. // Signed zeros, infinities, NaNs, and enormous arguments are tested separately. -using Xunit; - namespace Just.PreciseMath.Tests.ReferenceData; internal static class ExponentialFunctionsReferenceData diff --git a/1-tests/Just.PreciseMath.Tests/ReferenceData/LogReferenceData.cs b/1-tests/Just.PreciseMath.Tests/ReferenceData/LogReferenceData.cs index e120849..c4002e3 100644 --- a/1-tests/Just.PreciseMath.Tests/ReferenceData/LogReferenceData.cs +++ b/1-tests/Just.PreciseMath.Tests/ReferenceData/LogReferenceData.cs @@ -1,7 +1,5 @@ // Generated by generate_log.py; do not hand-edit reference literals. // Exact binary64 sums; Decimal.ln at 450/650 digits, rounded to 120 digits. -using Xunit; - namespace Just.PreciseMath.Tests.ReferenceData; internal static class LogReferenceData diff --git a/1-tests/Just.PreciseMath.Tests/ReferenceData/RealPowReferenceData.cs b/1-tests/Just.PreciseMath.Tests/ReferenceData/RealPowReferenceData.cs index 5719924..ce1c580 100644 --- a/1-tests/Just.PreciseMath.Tests/ReferenceData/RealPowReferenceData.cs +++ b/1-tests/Just.PreciseMath.Tests/ReferenceData/RealPowReferenceData.cs @@ -1,8 +1,6 @@ // Generated by generate_real_pow.py; do not hand-edit reference literals. // Exact binary64 sums; Decimal ln/exp at 450/650 digits, rounded to 120 digits. // Reference uncertainty is explicitly allowed as 2^-350 relative in the tests. -using Xunit; - namespace Just.PreciseMath.Tests.ReferenceData; internal static class RealPowReferenceData diff --git a/1-tests/Just.PreciseMath.Tests/ReferenceData/RootNReferenceData.cs b/1-tests/Just.PreciseMath.Tests/ReferenceData/RootNReferenceData.cs index 913bc5a..c868645 100644 --- a/1-tests/Just.PreciseMath.Tests/ReferenceData/RootNReferenceData.cs +++ b/1-tests/Just.PreciseMath.Tests/ReferenceData/RootNReferenceData.cs @@ -1,7 +1,5 @@ // Generated by generate_rootn.py; do not hand-edit reference literals. // Exact component sums; Decimal ln/exp at 450/650 digits, rounded to 120 digits. -using Xunit; - namespace Just.PreciseMath.Tests.ReferenceData; internal static class RootNReferenceData diff --git a/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleAdditionChecks.cs b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleAdditionChecks.cs new file mode 100644 index 0000000..9b993ea --- /dev/null +++ b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleAdditionChecks.cs @@ -0,0 +1,27 @@ +namespace Just.PreciseMath.Tests.SanityChecks; + +public class DoubleDoubleAdditionChecks +{ + public static IEnumerable> ReferenceCases => + [ + (new DoubleDouble(1.0), new DoubleDouble(1.0), new DoubleDouble(2.0)), + (new DoubleDouble(1.0), new DoubleDouble(2.0), new DoubleDouble(3.0)), + (new DoubleDouble(1.25), new DoubleDouble(2.5), new DoubleDouble(3.75)), + (new DoubleDouble(65535.0), new DoubleDouble(1.0), new DoubleDouble(65536.0)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), DoubleDouble.Parse("-812563124576179134504", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(42.0)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(127.0), DoubleDouble.Parse("812563124576179134673", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(27.0), DoubleDouble.Parse("812563124576179134573", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(2.70), DoubleDouble.Parse("812563124576179134548.7", System.Globalization.CultureInfo.InvariantCulture)), + ]; + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void AdditionShouldProduceExpectedResults(DoubleDouble left, DoubleDouble right, DoubleDouble expected) + { + DoubleDouble result = left + right; + + result.High.ShouldBe(expected.High); + result.Low.ShouldBe(expected.Low); + result.ShouldBe(expected); + } +} diff --git a/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleArithmeticEqualityChecks.cs b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleArithmeticEqualityChecks.cs new file mode 100644 index 0000000..c71e6f2 --- /dev/null +++ b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleArithmeticEqualityChecks.cs @@ -0,0 +1,108 @@ +namespace Just.PreciseMath.Tests.SanityChecks; + +public class DoubleDoubleArithmeticEqualityChecks +{ + public static IEnumerable> ReferenceCases => [1.0, 2.0, 3.0, 100.0, 0.1, 0.5, 0.333333333333333333, -1.0, 1553.0, 1.0e+10]; + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void AdditionShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = testValue + (new DoubleDouble(reference)); + DoubleDouble ddxd_Result = testValue + reference; + DoubleDouble dxdd_Result = reference + testValue; + + ddxdd_Result.High.ShouldBe(ddxd_Result.High); + ddxdd_Result.Low.ShouldBe(ddxd_Result.Low); + ddxdd_Result.ShouldBe(ddxd_Result); + + ddxdd_Result.High.ShouldBe(dxdd_Result.High); + ddxdd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxdd_Result.ShouldBe(dxdd_Result); + + ddxd_Result.High.ShouldBe(dxdd_Result.High); + ddxd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxd_Result.ShouldBe(dxdd_Result); + } + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void SubtractionFromShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = testValue - (new DoubleDouble(reference)); + DoubleDouble ddxd_Result = testValue - reference; + + ddxdd_Result.High.ShouldBe(ddxd_Result.High); + ddxdd_Result.Low.ShouldBe(ddxd_Result.Low); + ddxdd_Result.ShouldBe(ddxd_Result); + } + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void SubtractionShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = (new DoubleDouble(reference)) - testValue; + DoubleDouble dxdd_Result = reference - testValue; + + ddxdd_Result.High.ShouldBe(dxdd_Result.High); + ddxdd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxdd_Result.ShouldBe(dxdd_Result); + } + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void MultiplicationShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = testValue * (new DoubleDouble(reference)); + DoubleDouble ddxd_Result = testValue * reference; + DoubleDouble dxdd_Result = reference * testValue; + + ddxdd_Result.High.ShouldBe(ddxd_Result.High); + ddxdd_Result.Low.ShouldBe(ddxd_Result.Low); + ddxdd_Result.ShouldBe(ddxd_Result); + + ddxdd_Result.High.ShouldBe(dxdd_Result.High); + ddxdd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxdd_Result.ShouldBe(dxdd_Result); + + ddxd_Result.High.ShouldBe(dxdd_Result.High); + ddxd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxd_Result.ShouldBe(dxdd_Result); + } + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void DivisionShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = testValue / (new DoubleDouble(reference)); + DoubleDouble ddxd_Result = testValue / reference; + + ddxdd_Result.High.ShouldBe(ddxd_Result.High); + ddxdd_Result.Low.ShouldBe(ddxd_Result.Low); + ddxdd_Result.ShouldBe(ddxd_Result); + } + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void DivisionByShouldProduceIdenticalResults(double reference) + { + DoubleDouble testValue = DoubleDouble.E; // any value with full precision range used + + DoubleDouble ddxdd_Result = (new DoubleDouble(reference)) / testValue; + DoubleDouble dxdd_Result = reference / testValue; + + ddxdd_Result.High.ShouldBe(dxdd_Result.High); + ddxdd_Result.Low.ShouldBe(dxdd_Result.Low); + ddxdd_Result.ShouldBe(dxdd_Result); + } +} diff --git a/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleDivisionChecks.cs b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleDivisionChecks.cs new file mode 100644 index 0000000..7891033 --- /dev/null +++ b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleDivisionChecks.cs @@ -0,0 +1,29 @@ +namespace Just.PreciseMath.Tests.SanityChecks; + +public class DoubleDoubleDivisionChecks +{ + public static IEnumerable> ReferenceCases => + [ + (new DoubleDouble(6.0), new DoubleDouble(2.0), new DoubleDouble(3.0)), + (new DoubleDouble(6.0), new DoubleDouble(3.0), new DoubleDouble(2.0)), + (new DoubleDouble(5.0), new DoubleDouble(2.0), new DoubleDouble(2.5)), + (new DoubleDouble(1.0), new DoubleDouble(4.0), new DoubleDouble(0.25)), + (new DoubleDouble(3.75), new DoubleDouble(1.25), new DoubleDouble(3.0)), + (new DoubleDouble(65536.0), new DoubleDouble(1024.0), new DoubleDouble(64.0)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(2.0), DoubleDouble.Parse("406281562288089567273", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(4.0), DoubleDouble.Parse("203140781144044783636.5", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(8.0), DoubleDouble.Parse("101570390572022391818.25", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(1.0)), + ]; + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void DivisionShouldProduceExpectedResults(DoubleDouble left, DoubleDouble right, DoubleDouble expected) + { + DoubleDouble result = left / right; + + result.High.ShouldBe(expected.High); + result.Low.ShouldBe(expected.Low); + result.ShouldBe(expected); + } +} diff --git a/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleMultiplicationChecks.cs b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleMultiplicationChecks.cs new file mode 100644 index 0000000..2f8892b --- /dev/null +++ b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleMultiplicationChecks.cs @@ -0,0 +1,29 @@ +namespace Just.PreciseMath.Tests.SanityChecks; + +public class DoubleDoubleMultiplicationChecks +{ + public static IEnumerable> ReferenceCases => + [ + (new DoubleDouble(1.0), new DoubleDouble(1.0), new DoubleDouble(1.0)), + (new DoubleDouble(2.0), new DoubleDouble(3.0), new DoubleDouble(6.0)), + (new DoubleDouble(1.25), new DoubleDouble(2.5), new DoubleDouble(3.125)), + (new DoubleDouble(-1.25), new DoubleDouble(2.5), new DoubleDouble(-3.125)), + (new DoubleDouble(65535.0), new DoubleDouble(65536.0), new DoubleDouble(4294901760.0)), + (new DoubleDouble(65536.0), new DoubleDouble(65536.0), new DoubleDouble(4294967296.0)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(2.0), DoubleDouble.Parse("1625126249152358269092", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(0.5), DoubleDouble.Parse("406281562288089567273", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(27.0), DoubleDouble.Parse("21939204363556836632742", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(127.0), DoubleDouble.Parse("103195516821174750087342", System.Globalization.CultureInfo.InvariantCulture)), + ]; + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void MultiplicationShouldProduceExpectedResults(DoubleDouble left, DoubleDouble right, DoubleDouble expected) + { + DoubleDouble result = left * right; + + result.High.ShouldBe(expected.High); + result.Low.ShouldBe(expected.Low); + result.ShouldBe(expected); + } +} diff --git a/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleSubtractionChecks.cs b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleSubtractionChecks.cs new file mode 100644 index 0000000..ca8943b --- /dev/null +++ b/1-tests/Just.PreciseMath.Tests/SanityChecks/DoubleDoubleSubtractionChecks.cs @@ -0,0 +1,30 @@ +namespace Just.PreciseMath.Tests.SanityChecks; + +public class DoubleDoubleSubtractionChecks +{ + public static IEnumerable> ReferenceCases => + [ + (new DoubleDouble(2.0), new DoubleDouble(1.0), new DoubleDouble(1.0)), + (new DoubleDouble(1.0), new DoubleDouble(2.0), new DoubleDouble(-1.0)), + (new DoubleDouble(3.75), new DoubleDouble(2.5), new DoubleDouble(1.25)), + (new DoubleDouble(2.5), new DoubleDouble(3.75), new DoubleDouble(-1.25)), + (new DoubleDouble(65536.0), new DoubleDouble(1.0), new DoubleDouble(65535.0)), + (new DoubleDouble(65536.0), new DoubleDouble(65535.0), new DoubleDouble(1.0)), + (DoubleDouble.Parse("812563124576179134673", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(127.0), DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(127.0), DoubleDouble.Parse("812563124576179134419", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), DoubleDouble.Parse("812563124576179134504", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(42.0)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(0.5), DoubleDouble.Parse("812563124576179134545.5", System.Globalization.CultureInfo.InvariantCulture)), + (DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), DoubleDouble.Parse("812563124576179134546", System.Globalization.CultureInfo.InvariantCulture), new DoubleDouble(0.0)), + ]; + + [Theory] + [MemberData(nameof(ReferenceCases))] + public void SubtractionShouldProduceExpectedResults(DoubleDouble left, DoubleDouble right, DoubleDouble expected) + { + DoubleDouble result = left - right; + + result.High.ShouldBe(expected.High); + result.Low.ShouldBe(expected.Low); + result.ShouldBe(expected); + } +} diff --git a/AGENTS.md b/AGENTS.md index 87a3362..d0fef2c 100644 --- a/AGENTS.md +++ b/AGENTS.md @@ -16,12 +16,6 @@ Do not claim API completeness or accuracy beyond tested contracts. comparisons; use the Dry job for smoke validation, not performance conclusions. - Root `Directory.Build.props` holds shared settings. Each numbered directory's props explicitly imports it; preserve this import chain. -- `review-legacy/`: optional local reference material, excluded by `.gitignore`. - Start with `review-legacy/Review-Revisited-DoubleDouble.md` when present. The code - and proposed fixes contain known defects: reproduce findings against current - code rather than treating them as correctness oracles. Do not compile, copy - wholesale, or force-add these files. If absent, proceed without them; do not - invent their contents or make validation depend on them. - Use the ignored, repository-local `.hermes/` directory for plans, checklists, investigation notes, and handoff context. Create or update notes as useful; keep them concise and revalidate them against current files. Optional diff --git a/README.md b/README.md index 2e34be6..036e13e 100644 --- a/README.md +++ b/README.md @@ -24,11 +24,15 @@ not perform double-double arithmetic or allocate on the heap. a `double`. Addition retains residuals under cancellation; multiplication uses fused multiply-add; division uses residual corrections. Mixed `double` operators use specialized scalar paths rather than promoting the scalar to `DoubleDouble`. - Their finite fast paths normalize once with a final sum transform; scalar - division uses one compensated quotient correction within the error contract below. + Division retains the complete first remainder and applies two quotient corrections, + normalizing before the final correction. Mixed addition and subtraction also retain + intermediate sum residuals through final normalization. Mixed addition, subtraction, + and division produce the same component bits as their `DoubleDouble` operations + with the scalar represented as a zero-low pair. - Exponent boundaries: bounded `BigInteger` calculations avoid intermediate overflow and underflow on the exceptional finite path. Ordinary arithmetic uses - floating-point transforms without allocations. The stored value remains two + floating-point transforms without allocations, although sparse division correction + products can also reach a boundary path. The stored value remains two doubles; this is not an arbitrary-precision API. - Comparisons use both components. `Equals` treats NaNs as equal and signed zeros as equal for collections. `CompareTo` orders NaN before other values. Numerical @@ -107,8 +111,8 @@ The `DDMath` static class provides: It shares the existing `DoubleDouble.Abs` implementation. - `Reciprocal(DoubleDouble)`: returns exactly the same high and low component bits as `1.0 / value`. It specializes scalar/DD division for a numerator of one, - omitting only redundant numerator checks while preserving both divisions, - the FMA sequence, and normalization. Signed zeros map to signed infinities, + omitting redundant numerator checks and sharing the finite scalar-numerator kernel, + including both residual corrections and normalization. Signed zeros map to signed infinities, signed infinities to signed zeros, and NaN to canonical NaN. The allocating exact boundary path handles extreme exponents; finite overflow produces signed infinity. No speedup over scalar/DD division has been measured.