Fast remainder

[Albert Chan]

Conclusion part of Algorithm 1 Wrote:On this basis, Algorithm 2 gives the minimal number of fractional bits F. It is sufficient to pick F ≥ N + log2(d)

We can intuitively show above is true, by thinking floating point.

Let M = 2^F - 1, c = ⌊M/d⌋ + 1, s = "significant bits" of c = F - log2(d)

x % d == (x*c & M) * d >> 32 == ⌊ (x*c & M) * d / 2**32 ⌋

(x*c & M) removed quotient of x/d, leaving "fractional part" of s significant bits.

If x has less than s significnat bits, result is accurate.
If x is too big, returning s significant bits is not enough.

Example:

>>> M = (1<<32) - 1 # F = 32
>>> fastmod = lambda x: int((x*c & M) * d >> 32)
>>> test = lambda x: (fastmod(x), x%d, len(bin(x))-2)

>>> d = 17
>>> c = M//d + 1 # = 0xf0f0f10, 4*7 = 28 bits
>>> test((1<<28) - 1) # 28 bits x, ok
(15, 15, 28)
>>> test((1<<28)) # 29 bits x, failed
(0, 16, 29)

>>> d = 29
>>> c = M//d + 1 # = 0x8d3dcb1, 4*7 = 28 bits
>>> test(0x13b13b13) # 29 bits x, ok
(5, 5, 29)
>>> test(0x13b13b14) # 29 bits x, failed
(7, 6, 29)