[编辑] 最后一个想法，哪个（IMO）足够重要，我将把它放在开头：如果你一次收集一堆 totients，你可以避免很多多余的工作。不必费心从大数开始寻找较小的因数，而是对较小的因数进行迭代并为较大的数累积结果。
class Totient:
def __init__(self, n):
self.totients = [1 for i in range(n)]
for i in range(2, n):
if self.totients[i] == 1:
for j in range(i, n, i):
self.totients[j] *= i - 1
k = j / i
while k % i == 0:
self.totients[j] *= i
k /= i
def __call__(self, i):
return self.totients[i]
if __name__ == '__main__':
from itertools import imap
totient = Totient(10000)
print sum(imap(totient, range(10000)))
这在我的桌面上只需要 8 毫秒。
Euler totient function 上的维基百科页面有一些很好的数学结果。
计算与 的每个除数互质且小于每个除数的数字：这有一个简单的* 映射来计算从 到 的整数，所以总和是。
* 由第二个定义 trivial
这对于 Möbius inversion formula 的应用来说是完美的，这是一个巧妙的技巧，用于反转这种精确形式的总和。
这自然会导致代码
def totient(n):
if n == 1: return 1
return sum(d * mobius(n / d) for d in range(1, n+1) if n % d == 0)
def mobius(n):
result, i = 1, 2
while n >= i:
if n % i == 0:
n = n / i
if n % i == 0:
return 0
result = -result
i = i + 1
return result
Möbius function 存在更好的实现，可以记忆它以提高速度，但这应该很容易理解。
totient函数的计算比较明显
换句话说，将数字完全分解为唯一的素数和指数，然后从那里做一个简单的乘法。
from operator import mul
def totient(n):
return int(reduce(mul, (1 - 1.0 / p for p in prime_factors(n)), n))
def primes_factors(n):
i = 2
while n >= i:
if n % i == 0:
yield i
n = n / i
while n % i == 0:
n = n / i
i = i + 1
同样，prime_factors 存在更好的实现，但这是为了便于阅读。
# helper functions
from collections import defaultdict
from itertools import count
from operator import mul
def gcd(a, b):
while a != 0: a, b = b % a, a
return b
def lcm(a, b): return a * b / gcd(a, b)
primes_cache, prime_jumps = [], defaultdict(list)
def primes():
prime = 1
for i in count():
if i < len(primes_cache): prime = primes_cache[i]
else:
prime += 1
while prime in prime_jumps:
for skip in prime_jumps[prime]:
prime_jumps[prime + skip] += [skip]
del prime_jumps[prime]
prime += 1
prime_jumps[prime + prime] += [prime]
primes_cache.append(prime)
yield prime
def factorize(n):
for prime in primes():
if prime > n: return
exponent = 0
while n % prime == 0:
exponent, n = exponent + 1, n / prime
if exponent != 0:
yield prime, exponent
# OP's first attempt
def totient1(n):
counter = 0
for i in xrange(1, n):
if gcd(i, n) == 1:
counter += 1
return counter
# OP's second attempt
# I don't understand the algorithm, and just copying it yields inaccurate results
# Möbius inversion
def totient2(n):
if n == 1: return 1
return sum(d * mobius(n / d) for d in xrange(1, n+1) if n % d == 0)
mobius_cache = {}
def mobius(n):
result, stack = 1, [n]
for prime in primes():
if n in mobius_cache:
result = mobius_cache[n]
break
if n % prime == 0:
n /= prime
if n % prime == 0:
result = 0
break
stack.append(n)
if prime > n: break
for n in stack[::-1]:
mobius_cache[n] = result
result = -result
return -result
# traditional formula
def totient3(n):
return int(reduce(mul, (1 - 1.0 / p for p, exp in factorize(n)), n))
# traditional formula, no division
def totient4(n):
return reduce(mul, ((p-1) * p ** (exp-1) for p, exp in factorize(n)), 1)
使用此代码计算我桌面上从 1 到 9999 的所有数字的总和，平均运行 5 次以上，
-
totient1 需要永远
-
totient2 需要 10 秒
-
totient3 需要 1.3 秒
-
totient4 需要 1.3 秒