Fix quotient digit estimate in division (fixes #123)

[pietroppeter]

fix #123

Fix was produced by Claude, tests are added, Claude made additional fuzzing (not reported here). My effort was to review, understand and validate.

Here is a quick explanation from first principles:

• division in bigints is implemented using Algorithm D of Knuth which is a variation of standard long division
• as in long division the quotient is obtained digit by digit, in bigints case digits are in base 2^32, in standard long division base is 10
• you might not aware (I was not) but standard long division involves guessing each digit of the quotient, current guess sometimes called q_hat. The guess can be higher than the correct quotient digit, so you need to lower it
• knuth algorithm is a way of making a good guess and adjust efficiently if the guess is wrong (he proves the guess he gets is at most two higher than the correct digit)
• what can happen is that the guess can be higher than the base and current algorithm does not reduce enough the guess in order to have it fit the base and the assertion raises the exception
• with this fix, we check the case where this can happen and make sure that qhat is reduced enough, also making sure it is correct

In the details below (claude generated) you might also want to know more about knuth algorithm (aka Algorithm D):

• algo "normalizes" the dividend and divisor scaling it up in case the divisor "digit" (aka limb) is big
• q hat guess is produced by using top two limbs of dividend (or one if one is enough) and first digit of divisor
• the guess has higher error if the divisor digit is small, but if you scale up the divisor digit you can reduce the error
• for an example take 80 and 19. The guess of quotient will use 8 and 1 and use 8 as guess (8/1) which is quite wrong since 4 is the correct quotient (19 is almost 20 which fits 4 times in 80). If you multiply both by 4 you get 320 / 76 and in this case the guess would be 32 / 7 or 4 (4*7=28) which is already correct

details

When the top limb of the running remainder equals the top limb of the normalized divisor, vv div wm1 >= 2^32 and the two-limb correction test cannot always bring the estimate back below 2^32 (it stops at exactly 2^32 when the next limbs are equal too), tripping assert q1 <= uint32.high.

Restore the clamp from Knuth, TAOCP vol. 2, Algorithm D, step D3: if the top limbs are equal, set q1 = 2^32 - 1 and r1 = next limb + wm1, and skip the correction test once r1 >= 2^32.

Adds regression tests with the numbers from the issue and minimal cases.

[pietroppeter]

posting also here the initial conversation with the fix I shared in the issue thread https://claude.ai/share/84e37c1a-fc9c-4117-85a6-e2a975247ea0

with respect to the original patch, this PR fix adds an additional test case.

Also of interest might be a recent article that actually finds a bug in Knuth's original algorithm (and collects the check!). The bug apply only to odd bases but the article is very interesting to have a better understanding and appreciation of the algo https://kolja.rs/algorithm-d/