python - 求解大数的线性方程

Question

我有 a+b*n1=c+d*n2 形式的方程，其中 a、b、c、d 是已知数字，大约有 1000 位数字，我需要求解 n1。

我试过:

i=1
while True:
    a=(a+b)%d
    if(a==c):
        break
    i+=1
print(i)

，但是这种方法对于这么大的数字来说太慢了。在这种情况下有没有更好的方法可以使用？

Answer 1

您想找到 x 使得 x = a (mod b) 和 x = c (mod d)。然后，n1 = (x - a)/b 和 n2 = (x - c)/d。

如果 b 和 d 互质，则 x 的存在由 Chinese Remainder Theorem 保证——并且可以使用以下方法找到解决方案扩展欧几里德算法。

如果 b 和 d 不是互质的(也就是说，如果 gcd(b, d) != 1)，则 (注意 a = c (mod gcd(b, d))，我们可以从两边减去 a % gcd(b, d)，然后除以 gcd(b, d) 化简为上述问题。

将其写入代码

下面是使用此方法查找 n1 和 n2 的代码:

def egcd(a, b):
    if a == 0:
        return (b, 0, 1)
    else:
        g, y, x = egcd(b % a, a)
        return (g, x - (b // a) * y, y)

def modinv(a, m):
    return egcd(a, m)[1] % m

def solve(a, b, c, d):
    gcd = egcd(b, d)[0]
    if gcd != 1:
        if a % gcd != c % gcd:
            raise ValueError('no solution')
        a, c = a - a % gcd, c - c % gcd
        a //= gcd
        b //= gcd
        c //= gcd
        d //= gcd
    x = a * d * modinv(d, b) + c * b * modinv(b, d)
    return (x - a) // b, (x - c) // d

下面是一些测试代码，可以对 1000 位数字的输入进行 1000 次随机试验:

import sys
sys.setrecursionlimit(10000)

import random

digit = '0123456789'

def rn(k):
    return int(''.join(random.choice(digit) for _ in xrange(k)), 10)

k = 1000

for _ in xrange(1000):
    a, b, c, d, = rn(k), rn(k), rn(k), rn(k)
    print a, b, c, d
    try:
        n1, n2 = solve(a, b, c, d)
    except ValueError, exn:
        print 'no solution'
        print
        continue
    if a + b * n1 != c + d * n2:
        raise AssertionError('failed!')
    print 'found solution:', n1, n2
    print

(请注意，必须增加递归限制，因为实现扩展欧几里德算法的 egcd 函数是递归的，并且在 1000 位数字上运行它可能需要相当深的堆栈)。

另请注意，这会在返回解决方案时检查结果。但是当 a != c (mod gcd(b, d)) 并且引发异常时没有结果，则不进行检查。因此，当解决方案确实存在时，您需要仔细考虑是否无法找到结果。

这在我的机器上运行(1000 次试验)大约需要 7-8 秒，因此它的表现相当不错。