如何在代码中实现除数功能？

Question

总体问题：Project Euler 12 - 第一个具有超过 500 个除数的三角形数的值是多少？

问题焦点：除数函数

语言：Python

描述：我使用的函数是粗暴的，程序找到一个除数比 x 多的数字所需的时间几乎随着数字每增加 10 或 20 个而呈指数增长。我需要达到 500 个或更多除数。我已经确定除数函数是占用程序的原因。我所做的研究将我引向了除数函数，特别是 除数函数，它应该是一个可以计算任何整数的所有除数的函数。我看过的每一页似乎都是针对数学专业的，而我只有高中数学。虽然我确实遇到过一些提到关于素数和阿特金斯筛的分配的页面，但我无法在素数之间建立联系并找到任何整数的所有除数，也无法在网上找到任何关于它的信息。

主要问题：有人可以解释如何编写除数函数甚至提供示例吗？当我用代码查看数学概念时，它们对我来说更有意义。非常感谢。

蛮力除数函数：

def countdiv(a):
    count = 0
    for i in range(1,(a/2)+1):
        if a % i == 0:
            count += 1
    return count + 1    # +1 to account for number itself as a divisor

Answer 1

如果你需要一个 bruteforce 函数来计算除数的数量（也称为 tau(n)）

这就是它的样子

def tau(n):
        sqroot,t = int(n**0.5),0
        for factor in range(1,sqroot+1):
                if n % factor == 0:
                        t += 2 # both factor and N/factor
        if sqroot*sqroot == n: t = t - 1 # if sqroot is a factor then we counted it twice, so subtract 1
        return t

第二种方法涉及将n 分解为其主要因子（及其指数）。

tau(n) = (e1+1)(e2+1)....(em+1) 其中n = p1^e1 * p2^e2 .... pm^em 和p1,p2..pm are primes

更多信息here

第三种方法更容易理解，就是简单地使用筛子计算tau。

def sieve(N):
        t = [0]*(N+1)
        for factor in range(1,N+1):
                for multiple in range(factor,N+1,factor):
                        t[multiple]+=1
        return t[1:]

这是ideone的实际操作

Answer 2

我同意这里提交的其他两个答案，因为您只需要搜索到数字的平方根。不过，我还有一件事要补充。所提供的解决方案将在合理的时间内为您提供正确的答案。但是当问题开始变得更加棘手时，您将需要一个更强大的功能。

看看Euler's Totient function。虽然它只是间接地适用于这里，但它在以后的问题中非常有用。另一个相关的概念是Prime Factorization。

改进算法的一种快速方法是找到数字的素因数分解。在维基百科的文章中，他们以 36 为例，其素数分解为 2^2 * 3^2。因此，知道了这一点，您可以使用组合学来找到 36 的因数。这样，​​您实际上不会计算每个因数，而且您只需要在完成之前检查除数 2 和 3。

Answer 3

在搜索n 的除数时，您不必搜索超出数字n 的平方根。每当您找到一个小于 sqrt(n) 的除数时，恰好有一个匹配的除数大于根，因此您可以将您的 count 增加 2（如果您找到除数 d 的 n 然后 n/d将是对应的）。

不过，请注意平方数。 :) 当然，根将是一个不计两次的除数。

Answer 4

如果您要解决 Project Euler 问题，您需要一些处理素数和整数分解的函数。这是我的普通库，提供primes(n)、is_prime(n) 和factors(n)；重点在于以牺牲速度为代价的简单、清晰和简洁，尽管这些功能对于 Project Euler 来说应该足够了：

def primes(n):
    """
        list of primes not exceeding n in ascending
        order; assumes n is an integer greater than
        1; uses Sieve of Eratosthenes
    """
    m = (n-1) // 2
    b = [True] * m
    i, p, ps = 0, 3, [2]
    while p*p < n:
        if b[i]:
            ps.append(p)
            j = 2*i*i + 6*i + 3
            while j < m:
                b[j] = False
                j = j + 2*i + 3
        i += 1; p += 2
    while i < m:
        if b[i]:
            ps.append(p)
        i += 1; p += 2
    return ps

def is_prime(n):
    """
        False if n is provably composite, else
        True if n is probably prime; assumes n
        is an integer greater than 1; uses
        Miller-Rabin test on prime bases < 100
    """
    ps = [2,3,5,7,11,13,17,19,23,29,31,37,41,
         43,47,53,59,61,67,71,73,79,83,89,97]
    def is_spsp(n, a):
        d, s = n-1, 0
        while d%2 == 0:
            d /= 2; s += 1
        if pow(a,d,n) == 1:
            return True
        for r in xrange(s):
            if pow(a, d*pow(2,r), n) == n-1:
                return True
        return False
    if n in ps: return True
    for p in ps:
        if not is_spsp(n,p):
            return False
    return True

def factors(n):
    """
        list of prime factors of n in ascending
        order; assumes n is an integer, may be
        positive, zero or negative; uses Pollard's
        rho algorithm with Floyd's cycle finder
    """
    def gcd(a,b):
        while b: a, b = b, a%b
        return abs(a)
    def facts(n,c,fs):
        f = lambda(x): (x*x+c) % n
        if is_prime(n): return fs+[n]
        t, h, d = 2, 2, 1
        while d == 1:
            t = f(t); h = f(f(h))
            d = gcd(t-h, n)
        if d == n:
            return facts(n, c+1, fs)
        if is_prime(d):
            return facts(n//d, c+1, fs+[d])
        return facts(n, c+1, fs)
    if -1 <= n <= 1: return [n]
    if n < -1: return [-1] + factors(-n)
    fs = []
    while n%2 == 0:
        n = n//2; fs = fs+[2]
    if n == 1: return fs
    return sorted(facts(n,1,fs))

一旦您知道如何分解一个数字，就很容易计算除数的数量。考虑 76576500 = 2^2 * 3^2 * 5^3 * 7^1 * 11^1 * 13^1 * 17^1。忽略底数并查看指数，即 2、2、3、1、1、1 和 1。将每个指数加 1，得到 3、3、4、2、2、2 和 2。现在相乘该列表得到原始数字 76576500 的除数：3 * 3 * 4 * 2 * 2 * 2 * 2 = 576。函数如下：

def numdiv(n):
    fs = factors(n)
    f = fs.pop(0); d = 1; x = 2
    while fs:
        if f == fs[0]:
            x += 1
        else:
            d *= x; x = 2
        f = fs.pop(0)
    return d * x

您可以在http://codepad.org/4j8qp60u 看到这些功能的工作原理，并在my blog 了解更多关于它们如何工作的信息。问题 12 的解决方案就交给你了。