更有效的python for循环(3个冗余循环)，在较大的迭代值下

Question

我正在尝试实现以下代码：

def foo(n, p):
    for i in range(1,n):
        for j in range(1,n):
            for k in range(1,n):
                if ((n-j)*i*k)==(j*(n-i)*(n-k)):
                    p=p-11

但是，n将接近10^10的值，这使得这是非常低效的。事实上，即使在n=1000中，这也是缓慢的。

有没有一种方法可以通过压缩for循环来加快速度，或者完全不需要for循环就可以做到这一点？

Answer 1

我要采取数学方法和计算机科学方法。减少这些循环显然有一些有趣的问题，但是数学方法可能会给你带来一个小小的错误。

我想知道这个序列是否有一个封闭的公式，因为它总是比任何循环都要快！在您提供的OEIS链接中，在公式中，有人提供了一个“经验”生成函数

x*(1+5*x+11*x^2+x^3+6*x^4)/(1-x)^3/(1+x)^2

我稍后会讲到“经验”部分。但是因为这是多项式的比率，所以如果你读到生成函数是如何工作的，得到一个封闭形式的解是相当容易的。如果这个方法最终是你喜欢的，我可以把代数加到我的答案中，但是现在，让我们直接切入公式：

def empirical(n):
    return ((-1)**n * (-1.5*n + 2.5)) + \
               (3.0*n**2 - 4.5*n + 3.5)

很干净很简单。这有多精确？我检查了前500个值。这两个函数通常是完美地排列在一起的，但有时empirical夸大了真正的顺序：

    correct  empirical  pct_diff
1         1        1.0  0.000000
2         6        6.0  0.000000
3        19       19.0  0.000000
4        30       30.0  0.000000
5        61       61.0  0.000000
6        78       78.0  0.000000
7       127      127.0  0.000000
8       150      150.0  0.000000
9       217      217.0  0.000000
10      246      246.0  0.000000
11      331      331.0  0.000000
12      366      366.0  0.000000
13      469      469.0  0.000000
14      510      510.0  0.000000
15      625      631.0  0.009600*
16      678      678.0  0.000000
17      817      817.0  0.000000
18      870      870.0  0.000000
19     1027     1027.0  0.000000
20     1080     1086.0  0.005556*
21     1261     1261.0  0.000000
22     1326     1326.0  0.000000

这种偶尔的差别几乎总是小于1%。现在，我不能保证这一模式将适用于n = 10**10 (即，经验几乎总是正确的，偶尔会出现轻微的多报)，但请查看OEIS页面上的另一条评论：

利用Ceva定理从朴素计数中扣除消失区域。第一个推论是在n=15 (n奇数)和n=20 (n偶数)。

15和20恰好是与empirical的第一次分歧！因此，经验生成函数在大多数情况下似乎都是正确的(“朴素计数”？)，但在某些地方，它是一个上限，需要进行推导。这就进入了特定领域，我对塞瓦定理还不太了解，也不知道什么时候以及如何做这些推论，所以我恐怕不能改进这个封闭的上界，因为我上面有它。

您最初的问题是要测试10**10。所以现在立即做int(empirical(10**10))：

299999999939999956992

这要么是完全正确的，要么是一个非常非常接近真实答案的上限。

我知道这是一个“替代”的解决方案，但希望这是一个信息转移。就像有人要求你找到(10**10)斐波那契数。您可以执行循环，但是如果存在封闭形式的公式，请使用它！

Answer 2

操纵(n-j)*i*k=j*(n-i)*(n-k)。我们有j=n/(((n-i)*(n-k))/(i*k) + 1)和j应该是介于1和n-1之间的整数，所以：

def foo(n, p):
    for i in range(1,n):
        for k in range(1,n):
            j=n/(((n-i)*(n-k))/(i*k) + 1)
            if n%(((n-i)*(n-k))/(i*k) + 1) == 0 and j > 0 and j < n:
                p=p-11

这将复杂度从O(n立方)降至O(n平方)。