数字幂的整数除法可以优化吗?

Can integer division of numer power be optimized?

在Python中有一个内置函数pow优化了a**b%c的计算。为什么没有计算a**b//c的函数?

我会在这里尝试总结为什么 a ** b // c 的计算无法优化,不太可能 a ** b % c

先决条件:计算a ** b % c

让我们举个例子:说 a=2b=11。如果我们假设它很慢,那么我们可以推断出 b = 1 + 2 + 8 = 2**0 + 2**1 + 0*2**2 + 2**3。之后这个推导可以作为结果相乘的规则aa**2a**4a**8。每个结果在对前一个结果求平方后分配。最后,a**11 = a*(a**2)*(a**8),这个过程只需要 3 次平方。

如果我们概括这个过程,可以这样完成:

a, b, r = 2, 11 , []

while b>0:
    if b % 2: r.append(a)
    b = b//2
    a = a*a
else:
    if b % 2: r.append(a)

print(r)

一个输出是r=[2, 4, 256]。接下来,我们需要将这些乘数相乘。可以使用 from functools import reduce 和命令 reduce(lambda x,y: x*y, r) 来完成。

最后,如果乘数变得很大,乘法就会变得很慢,所以我们需要用它的模数 m%c 替换每个乘数 m,并在 reduce 函数中做同样的事情。最后,我们有:

from functools import reduce
def pow(a, b, c):
    # returns a ** b % c
    r = []
    while b > 0:
        if b % 2: r.append(a)
        b = b // 2
        a = (a * a) % c
    else:
        if b % 2: r.append(a)

    return reduce(lambda x, y: (x * y) % c, r)

输出是 4 因为 2 ** 11 % 74

我在我的电脑上也测试了结果 2046457 ** 1103207 % 71872。输出是 18249,计算需要 9 秒,而 pow(2046457, 1103207, 71872) 立即给出相同的结果。

更新:将 a ** b // c 插入计算

按照上述思路,我将尝试对a**b // c的计算进行类似的优化。我假设平方过程保持不变,这里的主要区别是我们在计算平方时需要同时考虑积分部分和残差部分(之前很容易,因为积分部分并不重要)。如果x是一个积分部分,y是残差部分,我们有一个关系:

我们还需要为两个不同的乘数引入类似的计算:

我的脚本现在看起来像这样:

from functools import reduce
def pow(a, b, c):
    #returns tuple (a ** b // c,  a ** b % c)
    print(f'calculating: {a}**{b} = ', end='')
    r = []
    ir = (a//c, a%c) # we keep integral and residual part of a instead of a
    while b > 0:
        if b % 2: r.append(ir)
        b = b // 2
        ir = (ir[0]*ir[0]*c + 2*ir[0]*ir[1]+ (ir[1]*ir[1])//c, (ir[1]*ir[1]) % c)
    else:
        if b % 2: r.append(ir)
    out = reduce(lambda x, y: (c*x[0]*y[0] + x[0]*y[1] + x[1]*y[0] + (x[1] * y[1])//c, (x[1] * y[1]) % c), [(2, 2)]+[r[-1]])

    print(' * '.join(str(n[0]*c+n[1]) for n in r), end=' = ')
    print(' * '.join(str(n) for n in r),'=', out)
    return out

pow(2,7,3)

输出

calculating: 2**7 = 2 * 4 * 16 = (0, 2) * (1, 1) * (5, 1) = (42, 2)

备注

为什么还没有优化?我们可以看到每个因子中的第二项始终很小,但这不是第一项的规则,例如 pow(26,31,35):

calculating: 26**31 = 26 * 676 * 456976 * 208827064576 * 43608742899428874059776 = 
(0, 26) * (19, 11) * (13056, 16) * (5966487559, 11) * (1245964082840824973136, 16) = 
(89709413964539398065824, 32)

在这种情况下,我们无法避免 a%b // c 的指数增长。这就是为什么 a%b // c 的现有内置函数对我来说似乎不合理。