There are a couple of minor algorithmic improvements possible for the math.isqrt fast path (which is used for nonnegative integers smaller than 2**64). On my machine those improvements produce a little over a 10% speedup.

The current algorithm for values under 2**64 involves exactly four division instructions, corresponding to four Newton steps. The proposal is to:

• 1. Replace the first division with a table lookup. The necessary table is extremely small: 12 entries at one byte per entry.
• 2. Arrange for the return type of the _approximate_sqrt helper function to be uint32_t rather than uint64_t. That means that the correction step only involves a 32-bit-by-32-bit multiplication, not a 64-bit-by-64-bit multiplication.

The second part is a bit subtle: the input to _approximate_sqrt is a 64-bit integer n in the range [262, 264). Barring any overflow, the output u is guaranteed to satisfy (u-1)**2 < n < (u+1)**2. That implies that (u-1)**2 < 2**64, from which it follows that u <= 2**32. So the only possible case where u might overflow a uint32_t is when u == 2**32. But from the earlier inequality, that can only happen if n > (2**32 - 1)**2, and in that case the top 31 bits of n are completely determined and since the first steps of the algorithm only depend on the topmost bits of n, it's easy to follow through the algorithm and see that it's not possible for u to be 2**32 in that case. (We always get u = 128 from the lookup, followed by u = 255 after the first division, then u = 65536 after the second.)

New changeset d02c5e9 by Mark Dickinson in branch 'main':
bpo-46258: Streamline isqrt fast path (bpo-30333)
d02c5e9