diff --git a/.github/workflows/dotnet.yml b/.github/workflows/dotnet.yml index fd26737..e6c4e35 100644 --- a/.github/workflows/dotnet.yml +++ b/.github/workflows/dotnet.yml @@ -352,7 +352,7 @@ jobs: '/d:sonar.host.url=https://sonarcloud.io' '/d:sonar.projectBaseDir=${{ github.workspace }}' '/d:sonar.cs.vscoveragexml.reportsPaths=coverage/**/coverage.xml' - '/d:sonar.coverage.exclusions=**/*Test*.cs,**/*.Tests.cs,**/*.Tests/**/*,**/*.Benchmarks/**/*,**/obj/**/*,**/*.dll,**/NativeExports.cs' + '/d:sonar.coverage.exclusions=**/*Test*.cs,**/*.Tests.cs,**/*.Tests/**/*,**/obj/**/*,**/*.dll,**/NativeExports.cs' '/d:sonar.cs.vstest.reportsPaths=coverage/**/*.trx' '/d:sonar.exclusions=**/NativeExports.cs' ) diff --git a/CLAUDE.md b/CLAUDE.md index e4872c2..141b003 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -33,7 +33,8 @@ dotnet run -c Release --project PreciseNumber.Benchmarks -- --filter '*' --job s ### Key Design Patterns - Factory methods `CreateFromInteger()` and `CreateFromFloatingPoint()` handle type-specific conversion logic -- Arithmetic operations use `MakeCommonized()` to align exponents before calculation +- Addition, subtraction and modulus align exponents before calculating; multiplication and division work on the significands directly +- `Divide` is exact when the quotient terminates, and otherwise rounds to a precision that never falls below the wider operand or `MinimumDivisionPrecision`. `Exp` and non-integer `Pow` still route through `double` - The `sanitize` constructor parameter controls whether trailing zeros are removed (default: true) - Constants (`Zero`, `One`, `Pi`, `E`, `Tau`) are pre-computed static instances diff --git a/PreciseNumber.Benchmarks/PreciseNumber.Benchmarks.csproj b/PreciseNumber.Benchmarks/PreciseNumber.Benchmarks.csproj index 54f9cc8..c6aa6ec 100644 --- a/PreciseNumber.Benchmarks/PreciseNumber.Benchmarks.csproj +++ b/PreciseNumber.Benchmarks/PreciseNumber.Benchmarks.csproj @@ -16,6 +16,12 @@ internal copy here through InternalsVisibleTo - referencing it again would make every polyfilled member ambiguous. --> $(NoWarn);KTSU0001 + + true diff --git a/PreciseNumber.Benchmarks/README.md b/PreciseNumber.Benchmarks/README.md index cce8188..f9dd876 100644 --- a/PreciseNumber.Benchmarks/README.md +++ b/PreciseNumber.Benchmarks/README.md @@ -57,11 +57,11 @@ operation returns a new instance, so a change that avoids an intermediate value `Allocated` column before it shows up in `Mean`. A comparison that allocates at all is a regression — none of them should. -One entry is easy to misread. `Divide` barely moves between 8 and 200 digits, and at 200 digits -it comes out faster than `Add`. That is not division being efficient — it computes its fractional -part through a `double`, so it discards everything past roughly the 17th significant digit and -the result it constructs is small no matter how large the operands were. Read it as a measure of -how much precision the operation throws away, not how fast it is. +`Divide` is the one row whose cost is not driven by the operands alone. It produces a terminating +quotient exactly and a repeating one to a chosen precision, so the `Digits` column moves it twice +over: wider operands are more work to divide, and they also raise the precision the quotient is +taken to. A number that looks expensive at 200 digits is doing proportionally more work, not doing +the same work badly. Benchmark operands come from a fixed digit pattern rather than a random source, so two runs on the same machine measure the same work. Numbers are still only comparable within a single run on diff --git a/PreciseNumber.Test/PreciseNumberTests.cs b/PreciseNumber.Test/PreciseNumberTests.cs index e35d0c8..4ef448c 100644 --- a/PreciseNumber.Test/PreciseNumberTests.cs +++ b/PreciseNumber.Test/PreciseNumberTests.cs @@ -1031,8 +1031,10 @@ public void TestOperatorDivide() PreciseNumber number1 = PreciseNumber.CreateFromComponents(-2, 12345); PreciseNumber number2 = PreciseNumber.CreateFromComponents(-3, 678); PreciseNumber result = number1 / number2; - Assert.AreEqual(BigInteger.Parse("18207964601769911504"), result.Significand); - Assert.AreEqual(-17, result.Exponent); + // 123.45 / 0.678 does not terminate, so it comes out at the default precision. + Assert.AreEqual(BigInteger.Parse("18207964601769911504424778761061946902654867256637"), result.Significand); + Assert.AreEqual(-47, result.Exponent); + Assert.AreEqual(PreciseNumber.MinimumDivisionPrecision, result.SignificantDigits); } [TestMethod] @@ -1603,7 +1605,11 @@ public void TestExpWithPositivePower() public void TestExpWithNegativePower() { PreciseNumber result = PreciseNumber.Exp(-1.ToPreciseNumber()); - PreciseNumber expected = PreciseNumber.One / PreciseNumber.E; // e^-1 = 1/e + + // Exp routes through a double, so 1/e is taken to the precision Exp actually delivers + // rather than the full precision Divide is now capable of. + PreciseNumber expected = PreciseNumber.Divide(PreciseNumber.One, PreciseNumber.E, result.SignificantDigits); + Assert.AreEqual(expected, result); } @@ -2000,6 +2006,122 @@ public void TestCountDigitsMatchesDecimalText() } } + [TestMethod] + public void TestDivideIsExactWhenTheQuotientTerminates() + { + // A quotient terminates exactly when the reduced denominator is a product of twos and + // fives, and those come out exact however many digits that takes. + (int Numerator, int Denominator, string Expected)[] cases = + [ + (1, 2, "0.5"), + (1, 4, "0.25"), + (1, 5, "0.2"), + (1, 8, "0.125"), + (1, 10, "0.1"), + (1, 16, "0.0625"), + (1, 20, "0.05"), + (1, 25, "0.04"), + (3, 8, "0.375"), + (-3, 8, "-0.375"), + (7, 1, "7"), + ]; + + foreach ((int numerator, int denominator, string expected) in cases) + { + PreciseNumber quotient = numerator.ToPreciseNumber() / denominator.ToPreciseNumber(); + Assert.AreEqual(expected, quotient.ToString(CultureInfo.InvariantCulture), $"{numerator}/{denominator}"); + } + } + + [TestMethod] + public void TestDivideKeepsEveryDigitOfALongTerminatingQuotient() + { + // 2^-64 terminates, but only after 64 decimal places - far more than the default + // precision would allow, and far more than a double could carry. + PreciseNumber quotient = PreciseNumber.One / PreciseNumber.CreateFromComponents(0, BigInteger.Pow(2, 64)); + + Assert.AreEqual( + "0.0000000000000000000542101086242752217003726400434970855712890625", + quotient.ToString(CultureInfo.InvariantCulture)); + + // Exact means exact: multiplying back reproduces the dividend. + Assert.AreEqual(PreciseNumber.One, quotient * PreciseNumber.CreateFromComponents(0, BigInteger.Pow(2, 64))); + } + + [TestMethod] + public void TestDivideProducesTheRequestedPrecisionWhenTheQuotientRepeats() + { + PreciseNumber one = PreciseNumber.One; + PreciseNumber three = 3.ToPreciseNumber(); + + Assert.AreEqual("0.3", PreciseNumber.Divide(one, three, 1).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("0.333", PreciseNumber.Divide(one, three, 3).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual(20, PreciseNumber.Divide(one, three, 20).SignificantDigits); + Assert.AreEqual(200, PreciseNumber.Divide(one, three, 200).SignificantDigits); + } + + [TestMethod] + public void TestDivideRoundsHalfAwayFromZero() + { + PreciseNumber two = 2.ToPreciseNumber(); + PreciseNumber three = 3.ToPreciseNumber(); + PreciseNumber five = 5.ToPreciseNumber(); + PreciseNumber nine = 9.ToPreciseNumber(); + + // 0.666... rounds up, 0.333... rounds down, and the sign does not change which way. + Assert.AreEqual("0.667", PreciseNumber.Divide(two, three, 3).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("0.333", PreciseNumber.Divide(PreciseNumber.One, three, 3).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("-0.667", PreciseNumber.Divide(-two, three, 3).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("-0.333", PreciseNumber.Divide(-PreciseNumber.One, three, 3).ToString(CultureInfo.InvariantCulture)); + + // 0.555... cut after one digit sits just above the halfway mark, so it rounds away. + Assert.AreEqual("0.6", PreciseNumber.Divide(five, nine, 1).ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("-0.6", PreciseNumber.Divide(-five, nine, 1).ToString(CultureInfo.InvariantCulture)); + } + + [TestMethod] + public void TestDivideNeverReducesTheOperandsPrecision() + { + // An operand carrying more digits than the floor pulls the quotient up to match it. + PreciseNumber wide = PreciseNumber.Parse(new string('7', 120), CultureInfo.InvariantCulture); + PreciseNumber three = 3.ToPreciseNumber(); + + Assert.AreEqual(120, (wide / three).SignificantDigits); + Assert.AreEqual(120, (three / wide).SignificantDigits); + Assert.AreEqual(PreciseNumber.MinimumDivisionPrecision, (three / 7.ToPreciseNumber()).SignificantDigits); + } + + [TestMethod] + public void TestDivideRejectsNonPositivePrecision() + { + PreciseNumber one = PreciseNumber.One; + PreciseNumber three = 3.ToPreciseNumber(); + + Assert.ThrowsExactly(() => PreciseNumber.Divide(one, three, 0)); + Assert.ThrowsExactly(() => PreciseNumber.Divide(one, three, -1)); + } + + [TestMethod] + public void TestDivideZeroDividend() + { + Assert.AreEqual(PreciseNumber.Zero, PreciseNumber.Zero / 3.ToPreciseNumber()); + Assert.AreEqual(PreciseNumber.Zero, PreciseNumber.Zero / PreciseNumber.NegativeOne); + } + + [TestMethod] + public void TestDivideRoundTripsThroughMultiplication() + { + // For a terminating quotient the round trip is exact, which is the strongest statement + // that can be made about a division. + foreach (int denominator in new[] { 2, 4, 5, 8, 10, 16, 20, 25, 32, 50, 64, 100, 125, 128 }) + { + PreciseNumber divisor = denominator.ToPreciseNumber(); + PreciseNumber quotient = 123456789.ToPreciseNumber() / divisor; + + Assert.AreEqual(123456789.ToPreciseNumber(), quotient * divisor, $"123456789/{denominator}"); + } + } + [TestMethod] public void TestPow10IsCorrectAcrossCacheBoundaries() { @@ -2057,7 +2179,7 @@ public void TestDivideWithExponentsBeyondDoubleRange() PreciseNumber result = left / right; - Assert.AreEqual("2.3333333333333333", result.ToString(CultureInfo.InvariantCulture)); + Assert.AreEqual("2.3333333333333333333333333333333333333333333333333", result.ToString(CultureInfo.InvariantCulture)); } [TestMethod] diff --git a/PreciseNumber/PreciseNumber.cs b/PreciseNumber/PreciseNumber.cs index 912389c..f7512fb 100644 --- a/PreciseNumber/PreciseNumber.cs +++ b/PreciseNumber/PreciseNumber.cs @@ -44,6 +44,12 @@ public record PreciseNumber /// private const string InvalidFormatMessage = "Input string was not in a correct format."; + /// + /// The fewest significant digits produces + /// when a quotient does not terminate. + /// + public const int MinimumDivisionPrecision = 50; + /// /// Pre-computed powers of ten, grown on demand. Declared before any other static state so /// that the static constants below can rely on it while they are being initialized. @@ -1291,30 +1297,149 @@ public static PreciseNumber Multiply(PreciseNumber left, PreciseNumber right) /// The number to divide. /// The number to divide by. /// The result of the division. + /// Thrown when is zero. + /// + /// A quotient whose decimal expansion terminates is produced exactly, however many digits that + /// takes. One that repeats is produced to the precision of the wider operand, and never fewer + /// than significant digits, with the last digit rounded + /// half away from zero. Use to choose + /// that precision. + /// public static PreciseNumber Divide(PreciseNumber left, PreciseNumber right) { Ensure.NotNull(left); Ensure.NotNull(right); + // Dividing must not silently discard precision the operands already carry. + int significantDigits = Math.Max( + Math.Max(left.SignificantDigits, right.SignificantDigits), + MinimumDivisionPrecision); + + return Divide(left, right, significantDigits); + } + + /// + /// Divides one number by another, to a chosen number of significant digits. + /// + /// The number to divide. + /// The number to divide by. + /// + /// The number of significant digits to produce when the quotient does not terminate. A quotient + /// that does terminate is exact regardless of this value. + /// + /// The result of the division. + /// Thrown when is zero. + /// Thrown when is less than one. + public static PreciseNumber Divide(PreciseNumber left, PreciseNumber right, int significantDigits) + { + Ensure.NotNull(left); + Ensure.NotNull(right); + + if (significantDigits < 1) + { + throw new ArgumentOutOfRangeException(nameof(significantDigits), significantDigits, "At least one significant digit is required."); + } + if (right.Significand.IsZero) { throw new DivideByZeroException(); } - if (Compare(left, right) == 0) + if (left.Significand.IsZero) { - return One; + return Zero; } - (BigInteger commonLeft, BigInteger commonRight, _) = CommonizeSignificands(left, right); + BigInteger numerator = left.Significand; + BigInteger denominator = right.Significand; + int exponent = left.Exponent - right.Exponent; - BigInteger integerComponent = BigInteger.DivRem(commonLeft, commonRight, out BigInteger remainder); + // Carry the sign on the numerator so the denominator can be factorized as a positive value. + if (denominator.Sign < 0) + { + numerator = -numerator; + denominator = -denominator; + } - // The common power of ten cancels between the remainder and the divisor, so it is left out - // entirely; including it would overflow to infinity for large exponents. - double fractionalComponent = double.CreateTruncating(remainder) / double.CreateTruncating(commonRight); + return TryDivideExactly(numerator, denominator, exponent, out PreciseNumber? exact) + ? exact + : DivideToPrecision(numerator, denominator, exponent, significantDigits); + } + + /// + /// Divides exactly, when the quotient has a terminating decimal expansion. + /// + /// The numerator, carrying the sign of the quotient. + /// The denominator, which must be positive. + /// The exponent the quotient's significand sits at. + /// The exact quotient, when there is one. + /// true if the quotient terminates and is exact; otherwise false. + private static bool TryDivideExactly(BigInteger numerator, BigInteger denominator, int exponent, [NotNullWhen(true)] out PreciseNumber? result) + { + // A fraction terminates in base ten exactly when its denominator is 2^twos * 5^fives. Most + // denominators are rejected by the first remainder test, which is why this is worth trying + // before falling back to a rounded quotient. + int twos = (int)BigInteger.TrailingZeroCount(denominator); + BigInteger remaining = denominator >> twos; + + int fives = 0; + while ((remaining % 5).IsZero) + { + remaining /= 5; + fives++; + } + + if (!remaining.IsOne) + { + result = null; + return false; + } + + // 1 / (2^p * 5^q) == (2^(k-p) * 5^(k-q)) / 10^k, where k is the larger of p and q. + int scale = Math.Max(twos, fives); + BigInteger significand = numerator * BigInteger.Pow(2, scale - twos) * BigInteger.Pow(5, scale - fives); + + result = new PreciseNumber(exponent - scale, significand); + return true; + } + + /// + /// Divides to a fixed number of significant digits, rounding the last of them half away from zero. + /// + /// The numerator, carrying the sign of the quotient. + /// The denominator, which must be positive. + /// The exponent the quotient's significand sits at. + /// The number of significant digits to produce. + /// The rounded quotient. + private static PreciseNumber DivideToPrecision(BigInteger numerator, BigInteger denominator, int exponent, int significantDigits) + { + // A quotient has either digits(numerator) - digits(denominator) digits or one more, so + // scaling by this much leaves at least one digit past the ones being kept: the digit the + // rounding decision is made on. + int scale = significantDigits + 1 - CountDigits(numerator) + CountDigits(denominator); + BigInteger scaled = scale > 0 ? numerator * Pow10(scale) : numerator; + int scaledExponent = exponent - Math.Max(scale, 0); + + BigInteger quotient = scaled / denominator; + int excess = CountDigits(quotient) - significantDigits; + + if (excess <= 0) + { + return new PreciseNumber(scaledExponent, quotient); + } + + BigInteger divisor = Pow10(excess); + BigInteger kept = BigInteger.DivRem(quotient, divisor, out BigInteger dropped); + + // Round half away from zero. Whatever the division above discarded is worth less than one + // unit of the dropped digits, so it can never carry the comparison across the halfway mark; + // at most it turns an exact tie into something above it, which rounds the same way. + if (BigInteger.Abs(dropped) * 2 >= divisor) + { + kept += quotient.Sign; + } - return new PreciseNumber(0, integerComponent) + fractionalComponent.ToPreciseNumber(); + return new PreciseNumber(scaledExponent + excess, kept); } /// diff --git a/README.md b/README.md index 8b9ac84..f2d1abf 100644 --- a/README.md +++ b/README.md @@ -351,9 +351,17 @@ You can control precision using: - **ReduceSignificance()**: Reduces to a specific number of significant digits +- **Divide(left, right, significantDigits)**: Chooses the precision of a quotient + +Division produces a terminating quotient exactly, however many digits that takes — `1 / 8` is +`0.125`, and `1 / 2^64` keeps all 64 decimal places. A repeating quotient is produced to the +precision of the wider operand, never fewer than `MinimumDivisionPrecision` (50) significant +digits, with the last digit rounded half away from zero. Pass an explicit precision to the +three-argument overload when you want something other than that. + ## Limitations -- Operations that inherently require approximation (like certain roots or logarithms) fall back to `double` precision for calculation +- `Exp()`, and `Pow()` with a non-integer power, are computed through `double` and are therefore limited to its precision. Addition, subtraction, multiplication and division are not - Conversion to standard types may throw `OverflowException` if the value is too large