如何提高/加快这个运行速度？

Question

我是一名Python初学者，正在努力提高自己的水平。我偶然发现了以下练习：

Let n be an integer greater than 1 and s(n) the sum of the dividors of n. For example,

s(12) 1 + 2 + 3 + 4 + 6 + 12 = 28

Also,

s(s(12)) = s(28) = 1 + 2 + 4 + 7 + 14 + 28 = 56

And

s(s(s(12))) = s(56) = 1 + 2 + 4 + 7 + 8 + 14 + 28 + 56 = 120

We use the notations:

s^1(n) = s(n)
s^2(n) = s(s(n))
s^3(n) = s(s(s(n)))
s^ m (n) = s(s(. . .s(n) . . .)), m times

For the integers n for which exists a positive integer k so that

 s^m(n) = k * n

are called (m, k)-perfect, for instance 12 is (3, 10)-perfect since s^3(12) = s(s(s(12))) = 120 = 10 * 12

Special categories:

For m =1 we have multiperfect numbers

A special case of the above exist for m = 1 and k = 2 which are called perfect numbers.

For m = 2 and k = 2 we have superperfect numbers.

Write a program which finds and prints all (m, k)-perfect numbers for m <= MAXM, which are less or equal to (<=) MAXNUM. If an integer belongs to one of the special categories mentioned above the program should print a relevant message. Also, the program has to print how many different (m, k)-perfect numbers were found, what percentage of the tested numbers they were, in how many occurrences for the different pairs of (m, k), and how many from each special category were found(perfect numbers are counted as multiperfect as well).

这是我的代码：

import time
start_time = time.time()
def s(n):
    tsum = 0
    i = 1
    con = n
    while i < con:
        if n % i == 0:
            temp = n / i
            tsum += i
            if temp != i:
                tsum += temp 
            con = temp
        i += 1                    
    return tsum
#MAXM
#MAXNUM

i = 2

perc = 0
perc1 = 0
perf = 0
multperf = 0
supperf = 0
while i <= MAXNUM:
    pert = perc1
    num = i
    for m in xrange(1, MAXM + 1):        
        tsum = s(num)                
        if tsum % i == 0:            
            perc1 += 1
            k = tsum / i            
            mes = "%d is a (%d-%d)-perfect number" % (i, m, k)
            if m == 1:
                multperf += 1
                if k == 2:
                    perf += 1 
                    print mes + ", that is a perfect number"
                else:
                    print mes + ", that is a multiperfect number"               
            elif m == 2 and k == 2:
                supperf += 1
                print mes + ", that is a superperfect number"
            else:
                print mes        
        num = tsum        
    i += 1
    if pert != perc1: perc += 1
print "Found %d distinct (m-k)-perfect numbers (%.5f per cent of %d ) in %d occurrences" % (
perc, float(perc) / MAXNUM * 100, MAXNUM, perc1)
print "Found %d perfect numbers" % perf
print "Found %d multiperfect numbers (including perfect numbers)" % multperf
print "Found %d superperfect numbers" % supperf  
print time.time() - start_time, "seconds"

这个程序可以正常运行，但我希望能够提高它的运行速度。比如说，使用什么方法会更快呢？

I = 1
while I <= MAXM:
    …..
    I += 1

取代

for I in xrange(1, MAXM + 1)

如果不把s(n)定义为函数，而是将代码放到主程序中，会更好吗？等等。 如果您有任何建议可以让我的程序运行得更快，我会很感激。 还有一件事，在最初的练习中，要求程序使用C语言（这是我不懂的），用Python编写后，将其转换成C语言会有多难？

- user2516856

---

这真的很慢吗？这只是一个建议，但作为一个初学者，你真的需要找到一种运行代码“更快”的方法吗？你应该首先关注编写更易读和可维护的代码。实际上，这是我曾经给（更有经验的）程序员的建议... - Sylvain Leroux

如果你会C语言，将其转换将非常容易...大部分都是C代码，你只需要加上大括号和分号 :) - Karoly Horvath

为了提高性能，有时候初学者会关注那些并不重要的东西。我使用的方法可以告诉你哪些代码片段值得关注。 - Mike Dunlavey

我认为提高性能的最大关键是尽量减少除法/取模运算。与您在此处采取的方法不同的方法是找到一个数字的质因数分解，然后重新组合这些因子来计算总和。记住因数分解，并注意到N=px（其中p是质数）的因子是x的因子列表和每个因子乘以p的乘积的组合列表将是好主意（例如6的因子是1、2、3；12=26的因子是1、2、3和2、4、6）。 - twalberg

也许你可以加入一些记忆化技术。 - Jan Matějka

Answer 1

最大的改进来自于使用更好的算法。像是否将s(n)定义为一个函数还是将代码放入主程序中，以及是否使用while循环代替for i in xrange(1, MAXM + 1):这些问题并没有太大的区别，因此在达到算法改进至少非常困难的状态之前，不应该考虑这些细小的事情。

因此，让我们来看一下您的算法，以及如何在不关心诸如while循环或for迭代哪个更快这样微小的事情的情况下，大幅度改进它。

def s(n):
    tsum = 0
    i = 1
    con = n
    while i < con:
        if n % i == 0:
            temp = n / i
            tsum += i
            if temp != i:
                tsum += temp 
            con = temp
        i += 1                    
    return tsum

这已经包含了一个好的思路，你知道 n 的因数成对出现，一旦找到较小的那个因数，就将两个因数加起来。你甚至正确地处理了平方。

对于像120这样的数字，它运作得非常好：当你找到因子2时，你将停止条件设置为60，当你找到3时，设置为40，...，当你找到8时，将其设置为15，当你找到10时，将其设置为12，然后只需除以11，当i增加到12时停止。不错。

但是当n是质数时，它的效果不太好，con永远不会被设置为小于n的值，你需要迭代到n才能找到因子和。对于形如n = 2*p的数字，其中p是质数，那么你需要循环到n/2，或者n = 3*p（n/3，除非p = 2）等。

根据素数定理，不超过x的质数数量渐近于x/log x（其中log是自然对数），并且您有一个下界。

Ω(MAXNUM² / log MAXNUM)

仅用于计算质数的除数和。这样不好。

既然您已经考虑了n的除数对d和n/d，请注意两者中较小的一个（暂时忽略当n是平方数时的情况d = n/d）小于n的平方根，因此一旦测试除数达到平方根，您就知道已经找到并添加了所有的除数，完成了任务。任何进一步的循环都是徒劳的浪费工作。

所以让我们考虑

def s(n):
    tsum = 0
    root = int(n**0.5)  # floor of the square root of n, at least for small enough n
    i = 1
    while i < root + 1:
        if n % i == 0:
            tsum += i + n/i
        i += 1
    # check whether n is a square, if it is, we have added root twice
    if root*root == n:
        tsum -= root
    return tsum

作为首次改进，然后您始终循环到平方根，并计算1 <= n <= MAXNUM的s(n)为Θ(MAXNUM^1.5)。这已经是相当大的改进了。（当然，您必须计算迭代除数总和，并且对于一些n <= MAXNUM，s(n)可能大于MAXNUM，因此您不能从中推断出总算法的复杂度界限。O(MAXM * MAXNUM^1.5) 但是，s(n)不会 much larger，因此复杂度也不会更糟糕。）

但是我们仍然可以在twalberg suggested提到的方法上进行改进，使用n的质因式分解来计算除数总和。

首先，如果 n = p^k 是一个质数的幂，那么 n 的因数是 1, p, p², ..., p^k，并且可以轻松计算出因数和（几何级数的闭式公式为

(p^(k+1) - 1) / (p - 1)

但是，无论是使用这个还是添加p的k+1次幂来除以n都不重要。

接下来，如果n = p^k * m，其中p是质数，m是一个不被p整除的数，则

s(n) = s(p^k) * s(m)

一个简单的方法是将n的每个因子d写成d = p^a * g的形式，其中p不被g整除。然后p^a必须整除p^k，即a <= k，而g必须整除m。反之，对于每个0 <= a <= k和每个整除m的g，p^a * g是n的一个因子。因此，我们可以列出n的因子（其中1 = g_1 < g_2 < ... < g_r = m是m的因子）。

 1*g_1   1*g_2  ...  1*g_r
 p*g_1   p*g_2  ...  p*g_r
  :        :           :
p^k*g_1 p^k*g_2 ... p^k*g_r

每一行的和为 p^a * s(m)。

如果我们有一份质数列表，那么我们可以写成：

def s(n):
    tsum = 1
    for p in primes:
        d = 1
        # divide out all factors p of n
        while n % p == 0:
            n = n//p
            d = p*d + 1
        tsum *= d
        if p*p > n: # n = 1, or n is prime
            break
    if n > 1:       # one last prime factor to account for
        tsum *= 1 + n
    return tsum

试除法会去到n的第二大质因数[如果n是合数]或者n最大质因数的平方根中较大的一个。对于最大试除因子，它有一个最坏情况的界限为n**0.5，这个界限只有在质数时才会达到，但对于大多数合数，试除的过程要早得多。

如果我们没有现成的质数列表，我们可以将for p in primes:这一行替换为for p in xrange(2, n):[上限不重要，因为如果它比n**0.5更大，则永远不会到达]，并获得一个不太慢的分解。（但是，通过避免使用大于2的偶数试除因子，即使用列表[2] + [3,5,7...] - 最好作为生成器 - 作为试除因子，可以轻松地加快速度，通过跳过3的倍数（除了3），[2,3] + [5,7, 11,13, 17,19, ...]以及如果需要的话再跳过一些其他小质数，速度可以进一步提高。）

现在，这很有帮助，但计算所有 n <= MAXNUM 的除数和仍需要 Ω(MAXNUM^1.5 / log MAXNUM) 的时间（我还没有分析，可能也是一个上界，或者 MAXNUM^1.5 仍然是一个下界，无论如何，对数因子很少有太大的区别 [超出常数因子])。

而且你会多次计算许多除数和（在你的示例中，当研究12时，你计算了 s(56)，再次研究28时，再次研究56）。为了缓解这种影响，记忆化 s(n) 是一个好主意。然后你只需要计算每个 s(n) 一次。

现在我们已经用空间换时间了，所以我们可以使用更好的算法一次性计算所有 1 <= n <= MAXNUM 的除数和，具有更好的时间复杂度（也更小的常数因子）。我们可以直接标记只有多个，而不是尝试每个足够小的（质数）数是否能够整除 n，从而避免所有留下余数的除法 - 这是绝大多数情况。

做到这一点的简单方法是

def divisorSums(n):
    dsums = [0] + [1]*n
    for k in xrange(2, n+1):
        for m in xrange(k, n+1, k):
            dsums[m] += k
    return dsums

使用时间复杂度为O(n * log n)的方法进行操作。如果使用质因数分解，可以稍微提高效率（时间复杂度为O(n * log log n)），但该方法较为复杂，现在不加入，或许以后会加入。

然后，您可以使用所有约数和的列表查找n小于等于MAXNUM的s(n)，并使用上述实现来计算大于MAXNUM的值的约数和[或者您可能想将值备忘录化到更大的限制]。

dsums = divisorSums(MAXNUM)
def memo_s(n):
    if n <= MAXNUM:
        return dsums[n]
    return s(n)

---

还不错嘛，

Found 414 distinct (m-k)-perfect numbers (0.10350 per cent of 400000 ) in 496 occurrences
Found 4 perfect numbers
Found 8 multiperfect numbers (including perfect numbers)
Found 7 superperfect numbers
12.709428072 seconds

为

import time
start_time = time.time()


def s(n):
    tsum = 1
    for p in xrange(2,n):
        d = 1
        # divide out all factors p of n
        while n % p == 0:
            n = n//p
            d = p*d + 1
        tsum *= d
        if p*p > n: # n = 1, or n is prime
            break
    if n > 1:       # one last prime factor to account for
        tsum *= 1 + n
    return tsum
def divisorSums(n):
    dsums = [0] + [1]*n
    for k in xrange(2, n+1):
        for m in xrange(k, n+1, k):
            dsums[m] += k
    return dsums

MAXM = 6
MAXNUM = 400000

dsums = divisorSums(MAXNUM)
def memo_s(n):
    if n <= MAXNUM:
        return dsums[n]
    return s(n)

i = 2

perc = 0
perc1 = 0
perf = 0
multperf = 0
supperf = 0
while i <= MAXNUM:
    pert = perc1
    num = i
    for m in xrange(1, MAXM + 1):
        tsum = memo_s(num)
        if tsum % i == 0:
            perc1 += 1
            k = tsum / i
            mes = "%d is a (%d-%d)-perfect number" % (i, m, k)
            if m == 1:
                multperf += 1
                if k == 2:
                    perf += 1
                    print mes + ", that is a perfect number"
                else:
                    print mes + ", that is a multiperfect number"
            elif m == 2 and k == 2:
                supperf += 1
                print mes + ", that is a superperfect number"
            else:
                print mes
        num = tsum
    i += 1
    if pert != perc1: perc += 1
print "Found %d distinct (m-k)-perfect numbers (%.5f per cent of %d ) in %d occurrences" % (
perc, float(perc) / MAXNUM * 100, MAXNUM, perc1)
print "Found %d perfect numbers" % perf
print "Found %d multiperfect numbers (including perfect numbers)" % multperf
print "Found %d superperfect numbers" % supperf
print time.time() - start_time, "seconds"

通过记忆化也需要为 n > MAXNUM 记录所需的除数和：

d = {}
for i in xrange(1, MAXNUM+1):
    d[i] = dsums[i]
def memo_s(n):
    if n in d:
        return d[n]
    else:
        t = s(n)
        d[n] = t
        return t

时间下降到

3.33684396744 seconds

- Daniel Fischer

---

大脑=炸裂，非常感谢。 - user2516856

Answer 2

from functools import lru_cache

...

@lru_cache
def s(n):
    ...

这应该会让它显着地更快。

[更新]哦，不好意思，根据文档，在 3.2 中添加了缓存。但是任何缓存都可以。请参见是否有装饰器来简单缓存函数返回值？

- andrew cooke

---

很遗憾，lru_cache无法导入，我猜测是因为我正在使用Python 2.7。不管怎样，感谢您的建议。 - user2516856