Merge pull request #69 from ktsu-dev/claude/library-optimization-9wvabg

Divide exactly, or to a chosen precision

Author: matt-edmondson

.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'
           )

CLAUDE.md

@@ -33,7 +33,8 @@ dotnet run -c Release --project PreciseNumber.Benchmarks -- --filter '*' --job s
 ### Key Design Patterns
 
 - Factory methods `CreateFromInteger<T>()` and `CreateFromFloatingPoint<T>()` 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
 

PreciseNumber.Benchmarks/PreciseNumber.Benchmarks.csproj

@@ -16,6 +16,12 @@
          internal copy here through InternalsVisibleTo - referencing it again would make
          every polyfilled member ambiguous. -->
     <NoWarn>$(NoWarn);KTSU0001</NoWarn>
+    <!-- Benchmarks are tooling: not shipped, and not something unit tests cover. Telling the
+         Sonar scanner so here keeps the coverage gate off this project's back without any
+         workflow needing to know benchmarks exist. The equivalent workflow-level exclusion
+         would have to be repeated in every repository that gains a benchmark project, and
+         re-applied whenever the shared .NET workflow is resynced. -->
+    <SonarQubeExclude>true</SonarQubeExclude>
   </PropertyGroup>
 
   <ItemGroup>

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

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<ArgumentOutOfRangeException>(() => PreciseNumber.Divide(one, three, 0));
+		Assert.ThrowsExactly<ArgumentOutOfRangeException>(() => 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]

PreciseNumber/PreciseNumber.cs

@@ -44,6 +44,12 @@ public record PreciseNumber
 	/// </summary>
 	private const string InvalidFormatMessage = "Input string was not in a correct format.";
 
+	/// <summary>
+	/// The fewest significant digits <see cref="Divide(PreciseNumber, PreciseNumber)"/> produces
+	/// when a quotient does not terminate.
+	/// </summary>
+	public const int MinimumDivisionPrecision = 50;
+
 	/// <summary>
 	/// 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)
 	/// <param name="left">The number to divide.</param>
 	/// <param name="right">The number to divide by.</param>
 	/// <returns>The result of the division.</returns>
+	/// <exception cref="DivideByZeroException">Thrown when <paramref name="right"/> is zero.</exception>
+	/// <remarks>
+	/// 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 <see cref="MinimumDivisionPrecision"/> significant digits, with the last digit rounded
+	/// half away from zero. Use <see cref="Divide(PreciseNumber, PreciseNumber, int)"/> to choose
+	/// that precision.
+	/// </remarks>
 	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);
+	}
+
+	/// <summary>
+	/// Divides one number by another, to a chosen number of significant digits.
+	/// </summary>
+	/// <param name="left">The number to divide.</param>
+	/// <param name="right">The number to divide by.</param>
+	/// <param name="significantDigits">
+	/// The number of significant digits to produce when the quotient does not terminate. A quotient
+	/// that does terminate is exact regardless of this value.
+	/// </param>
+	/// <returns>The result of the division.</returns>
+	/// <exception cref="DivideByZeroException">Thrown when <paramref name="right"/> is zero.</exception>
+	/// <exception cref="ArgumentOutOfRangeException">Thrown when <paramref name="significantDigits"/> is less than one.</exception>
+	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);
+	}
+
+	/// <summary>
+	/// Divides exactly, when the quotient has a terminating decimal expansion.
+	/// </summary>
+	/// <param name="numerator">The numerator, carrying the sign of the quotient.</param>
+	/// <param name="denominator">The denominator, which must be positive.</param>
+	/// <param name="exponent">The exponent the quotient's significand sits at.</param>
+	/// <param name="result">The exact quotient, when there is one.</param>
+	/// <returns><c>true</c> if the quotient terminates and <paramref name="result"/> is exact; otherwise <c>false</c>.</returns>
+	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;
+	}
+
+	/// <summary>
+	/// Divides to a fixed number of significant digits, rounding the last of them half away from zero.
+	/// </summary>
+	/// <param name="numerator">The numerator, carrying the sign of the quotient.</param>
+	/// <param name="denominator">The denominator, which must be positive.</param>
+	/// <param name="exponent">The exponent the quotient's significand sits at.</param>
+	/// <param name="significantDigits">The number of significant digits to produce.</param>
+	/// <returns>The rounded quotient.</returns>
+	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);
 	}
 
 	/// <summary>

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