Return the exact floor from 3D and 4D Hypot for 32-bit and narrower integers [patch]

Semantics.Quantities/StorageMath.cs

@@ -265,11 +265,12 @@ public static T RootN<T>(T value, int n)
 	/// therefore bounded by the type exactly as the caller's own <c>x * x + y * y</c> would be.
 	/// </para>
 	/// <para>
-	/// The primitives take the <see cref="double"/> round trip, so <c>Hypot&lt;double&gt;</c> is exactly
-	/// <see cref="double.Hypot(double, double)"/>. The 64-bit and 128-bit integers do not: they sum their
-	/// squares in <see cref="BigInteger"/> when they could overflow <typeparamref name="T"/>, and return
-	/// the exact floor, which throws <see cref="OverflowException"/> only when <typeparamref name="T"/>
-	/// cannot hold it (#520).
+	/// The floating point primitives take the <see cref="double"/> round trip, so <c>Hypot&lt;double&gt;</c>
+	/// is exactly <see cref="double.Hypot(double, double)"/>. The integer primitives do not: the 64-bit and
+	/// 128-bit ones sum their squares in <see cref="BigInteger"/> when they could overflow
+	/// <typeparamref name="T"/> (#520), the narrower ones sum theirs in <see cref="Int128"/> (#537), and
+	/// both return the exact floor, which throws <see cref="OverflowException"/> only when
+	/// <typeparamref name="T"/> cannot hold it.
 	/// </para>
 	/// </remarks>
 	public static T Hypot<T>(T x, T y)
@@ -280,6 +281,11 @@ public static T Hypot<T>(T x, T y)
 			return HypotOfWideInteger(x, y, T.Zero, T.Zero);
 		}
 
+		if (IsNarrowInteger<T>())
+		{
+			return HypotOfNarrowInteger(x, y, T.Zero, T.Zero);
+		}
+
 		if (IsRoundedThroughDouble<T>())
 		{
 			return T.CreateChecked(double.Hypot(double.CreateChecked(x), double.CreateChecked(y)));
@@ -319,6 +325,11 @@ internal static T Hypot<T>(T x, T y, T z)
 			return HypotOfWideInteger(x, y, z, T.Zero);
 		}
 
+		if (IsNarrowInteger<T>())
+		{
+			return HypotOfNarrowInteger(x, y, z, T.Zero);
+		}
+
 		if (IsRoundedThroughDouble<T>())
 		{
 			double dx = double.CreateChecked(x);
@@ -361,6 +372,11 @@ internal static T Hypot<T>(T x, T y, T z, T w)
 			return HypotOfWideInteger(x, y, z, w);
 		}
 
+		if (IsNarrowInteger<T>())
+		{
+			return HypotOfNarrowInteger(x, y, z, w);
+		}
+
 		if (IsRoundedThroughDouble<T>())
 		{
 			double dx = double.CreateChecked(x);
@@ -433,6 +449,52 @@ private static bool IsWideInteger<T>()
 			|| typeof(T) == typeof(UInt128)
 			|| ((typeof(T) == typeof(nint) || typeof(T) == typeof(nuint)) && nint.Size == 8);
 
+	/// <summary>
+	/// Reports whether <typeparamref name="T"/> is an integer primitive of 32 bits or fewer, whose
+	/// squares and their sums fit in <see cref="Int128"/> however large the components are.
+	/// </summary>
+	/// <typeparam name="T">The numeric storage type.</typeparam>
+	/// <returns>
+	/// <see langword="true"/> for <see cref="byte"/>, <see cref="sbyte"/>, <see cref="short"/>,
+	/// <see cref="ushort"/>, <see cref="int"/>, <see cref="uint"/>, and <see cref="nint"/> and
+	/// <see cref="nuint"/> on a 32-bit process.
+	/// </returns>
+	private static bool IsNarrowInteger<T>()
+		where T : struct, INumber<T>
+		=> typeof(T) == typeof(byte)
+			|| typeof(T) == typeof(sbyte)
+			|| typeof(T) == typeof(short)
+			|| typeof(T) == typeof(ushort)
+			|| typeof(T) == typeof(int)
+			|| typeof(T) == typeof(uint)
+			|| ((typeof(T) == typeof(nint) || typeof(T) == typeof(nuint)) && nint.Size == 4);
+
+	/// <summary>
+	/// The exact floor of <c>sqrt(x² + y² + z² + w²)</c> for a narrow integer type.
+	/// </summary>
+	/// <typeparam name="T">A type <see cref="IsNarrowInteger{T}"/> accepts.</typeparam>
+	/// <param name="x">The first component.</param>
+	/// <param name="y">The second component.</param>
+	/// <param name="z">The third component, zero for fewer components.</param>
+	/// <param name="w">The fourth component, zero for fewer components.</param>
+	/// <returns>The floor of the length.</returns>
+	/// <exception cref="OverflowException">The length does not fit in <typeparamref name="T"/>.</exception>
+	/// <remarks>
+	/// Each component is at most 2^32 in magnitude, so four squares sum within 2^66, which
+	/// <see cref="Int128"/> holds exactly. Nesting <see cref="double.Hypot(double, double)"/> instead
+	/// rounds once per call, which can leave an exact length one ulp low, and truncating that back to
+	/// <typeparamref name="T"/> then answers one too low on a perfect square (#537).
+	/// </remarks>
+	private static T HypotOfNarrowInteger<T>(T x, T y, T z, T w)
+		where T : struct, INumber<T>
+	{
+		Int128 wx = Int128.CreateChecked(x);
+		Int128 wy = Int128.CreateChecked(y);
+		Int128 wz = Int128.CreateChecked(z);
+		Int128 ww = Int128.CreateChecked(w);
+		return T.CreateChecked(FloorSqrtOfWideInteger((wx * wx) + (wy * wy) + (wz * wz) + (ww * ww)));
+	}
+
 	/// <summary>
 	/// The exact floor of the square root of a non-negative wide integer.
 	/// </summary>

Semantics.Test/Quantities/StorageMathPreciseAndVectorTests.cs

@@ -59,6 +59,18 @@ public void IntegerVectorLengthsDoNotWrap()
 		Assert.AreEqual(50000, Position3D<int>.FromMeter(50000, 0, 0).Length());
 	}
 
+	[TestMethod]
+	public void IntegerVectorLengthsOfPerfectSquaresAreExact()
+	{
+		// 36 + 1681 + 1764 = 3481 = 59², and 9 + 9 + 9 + 9 = 36 = 6² (#537).
+		Position3D<int> position = Position3D<int>.FromMeter(6, 41, 42);
+		Assert.AreEqual(59, position.Length());
+		Assert.AreEqual(59, position.Magnitude().Value);
+		Assert.AreEqual(59, position.DistanceTo(Position3D<int>.Zero).Value);
+		Assert.AreEqual((short)59, Position3D<short>.FromMeter(6, 41, 42).Length());
+		Assert.AreEqual(6, new Acceleration4D<int> { X = 3, Y = 3, Z = 3, W = 3 }.Length());
+	}
+
 	[TestMethod]
 	public void DecimalVectorLengthPastTheSquareOfItsRange()
 	{

Semantics.Test/Quantities/StorageMathTests.cs

@@ -180,7 +180,7 @@ public void TheRootsArePublicAndTheirWorkingsAreNot()
 		Assert.IsTrue(typeof(StorageMath).IsPublic, "StorageMath itself has to be reachable from outside the assembly.");
 
 		string[] contract = [nameof(StorageMath.Sqrt), nameof(StorageMath.Cbrt), nameof(StorageMath.RootN), nameof(StorageMath.Hypot)];
-		string[] workings = ["IsRoundedThroughDouble", "IsBinaryFloatingPoint", "HasFloorDivision", "Seed", "SeedForRoot", "TryRootThroughDouble", "RootThroughDouble", "RootByNewton", "RootByDescent", "RootBySettling", "TryPower", "IsWideInteger", "FloorSqrtOfWideInteger", "HypotOfWideInteger"];
+		string[] workings = ["IsRoundedThroughDouble", "IsBinaryFloatingPoint", "HasFloorDivision", "Seed", "SeedForRoot", "TryRootThroughDouble", "RootThroughDouble", "RootByNewton", "RootByDescent", "RootBySettling", "TryPower", "IsWideInteger", "FloorSqrtOfWideInteger", "HypotOfWideInteger", "IsNarrowInteger", "HypotOfNarrowInteger"];
 
 		foreach (string name in contract)
 		{
@@ -410,6 +410,77 @@ public void WideIntegerHypotenuseIsTheExactFloor()
 		Assert.AreEqual(2L, StorageMath.Hypot(1L, 1L, 1L, 1L));
 	}
 
+	/// <summary>
+	/// A 32-bit-or-narrower integer vector whose length is a perfect square gets that length, not one
+	/// less: nesting <see cref="double.Hypot(double, double)"/> rounds twice and can land one ulp low,
+	/// which truncation then turned into an off-by-one (#537).
+	/// </summary>
+	[TestMethod]
+	public void NarrowIntegerHypotenuseIsTheExactFloor()
+	{
+		Assert.AreEqual(59, StorageMath.Hypot(6, 41, 42));
+		Assert.AreEqual(59, StorageMath.Hypot(-6, 41, -42));
+		Assert.AreEqual(59, StorageMath.Hypot(14, 39, 42));
+		Assert.AreEqual((short)59, StorageMath.Hypot((short)6, (short)41, (short)42));
+		Assert.AreEqual(59u, StorageMath.Hypot(6u, 41u, 42u));
+		Assert.AreEqual(6, StorageMath.Hypot(3, 3, 3, 3));
+		Assert.AreEqual((short)6, StorageMath.Hypot((short)3, (short)3, (short)3, (short)3));
+		Assert.AreEqual((byte)15, StorageMath.Hypot((byte)2, (byte)3, (byte)4, (byte)14));
+		Assert.AreEqual(15, StorageMath.Hypot(1, 4, 8, 12));
+		Assert.AreEqual(5, StorageMath.Hypot(3, -4));
+
+		// Squares of the extremes sum far past 64 bits, and the length still fits the type.
+		Assert.AreEqual(uint.MaxValue, StorageMath.Hypot(uint.MaxValue, 0u, 0u, 0u));
+		Assert.AreEqual(int.MaxValue, StorageMath.Hypot(int.MinValue + 1, 0, 0));
+		Assert.ThrowsExactly<OverflowException>(static () => StorageMath.Hypot(int.MinValue, 0, 0));
+		Assert.ThrowsExactly<OverflowException>(static () => StorageMath.Hypot(uint.MaxValue, uint.MaxValue, uint.MaxValue, uint.MaxValue));
+	}
+
+	/// <summary>
+	/// A narrow integer gets the exact floor wherever the components lie: every 3D vector with
+	/// components up to 60, and every 4D vector up to 20, agrees with an integer square root of the
+	/// integer sum of squares (#537).
+	/// </summary>
+	[TestMethod]
+	public void NarrowIntegerHypotenuseMatchesTheIntegerRootEverywhere()
+	{
+		for (int x = 0; x <= 60; x++)
+		{
+			for (int y = 0; y <= 60; y++)
+			{
+				for (int z = 0; z <= 60; z++)
+				{
+					Assert.AreEqual(FloorSqrt((x * x) + (y * y) + (z * z)), StorageMath.Hypot(x, y, z), $"({x}, {y}, {z})");
+				}
+			}
+		}
+
+		for (int x = 0; x <= 20; x++)
+		{
+			for (int y = 0; y <= 20; y++)
+			{
+				for (int z = 0; z <= 20; z++)
+				{
+					for (int w = 0; w <= 20; w++)
+					{
+						Assert.AreEqual(FloorSqrt((x * x) + (y * y) + (z * z) + (w * w)), StorageMath.Hypot(x, y, z, w), $"({x}, {y}, {z}, {w})");
+					}
+				}
+			}
+		}
+
+		static int FloorSqrt(int value)
+		{
+			int root = 0;
+			while ((root + 1) * (root + 1) <= value)
+			{
+				root++;
+			}
+
+			return root;
+		}
+	}
+
 	/// <summary>
 	/// A wide-integer hypotenuse the type cannot hold is an overflow, as it was before (#520).
 	/// </summary>