多少个数字低于N的coprimes到N？ [英] How many numbers below N are coprimes to N?

问题描述

由于 是互质为 B 如果 GCD（A，B）= 1 （其中GCD代表的great公约数），以下N个正整数，有多少是互质为N？

Given that a is coprime to b if GCD(a,b) = 1 (where GCD stands for great common divisor), how many positive integers below N are coprime to N?

有一个聪明的办法？

---

这是最愚蠢的方式：

def count_coprime(N):
    counter = 0
    for n in xrange(1,N):
        if gcd(n,N) == 1:
            counter += 1
    return counter

它的工作原理，但它是缓慢的，和哑巴。我想用一个巧妙的和更快的算法。 我试图用质因数和N个约数，但我总是得到的东西，不与N越大的工作。

It works, but it is slow, and dumb. I'd like to use a clever and faster algorithm. I tried to use prime factors and divisors of N but I always get something that doesn't work with larger N.

我认为该算法应能指望他们不计算所有的人都喜欢的最愚蠢的算法呢：P

看来我已经找到了工作之一：

It seems I've found a working one:

def a_bit_more_clever_counter(N):
    result = N - 1
    factors = []
    for factor, multiplicity in factorGenerator(N):
        result -= N/factor - 1
        for pf in factors:
            if lcm(pf, factor) < N:
                result += N/lcm(pf, factor) - 1
        factors += [factor]
    return result

其中，LCM是最小公倍数。有没有人有一个更好的？

where lcm is least common multiple. Does anyone have a better one?

我使用的蟒蛇，我认为code应该是可读甚至谁不知道蟒蛇，如果你发现任何不明确的只是问的意见。我感兴趣的算法和数学的观念研究。

I'm using python, I think code should be readable even to who doesn't know python, if you find anything that is not clear just ask in the comments. I'm interested in the algorithm and the math, the idea.

推荐答案

最后一个想法，它（IMO）是足够重要，我会把它放在开头：如果你再搜罗了一堆totients的一次，就能避免很多重复的工作。不要打扰从大数开始找到自己的小因素 - 相反，迭代较小的因素，积累业绩较大的数字

One last thought, which (IMO) is important enough that I'll put it at the beginning: if you're collecting a bunch of totients at once, you can avoid a lot of redundant work. Don't bother starting from large numbers to find their smaller factors -- instead, iterate over the smaller factors and accumulate results for the larger numbers.

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)))

这只需8ms的我的桌面上。

This takes just 8ms on my desktop.

在 Euler函数的维基百科页面有一些很好的数学成绩。

The Wikipedia page on the Euler totient function has some nice mathematical results.

计数的数量互质并小于每除数：这有一个简单的*映射到从统计的整数到，所以总和是

counts the numbers coprime to and smaller than each divisor of : this has a trivial* mapping to counting the integers from to , so the sum total is .

*通过第二个定义的琐碎

这是完美的莫比乌斯反演公式的应用程序，一个聪明的把戏，反转款项这种具体形式。

This is perfect for an application of the Möbius inversion formula, a clever trick for inverting sums of this exact form.

这自然导致了code

This leads naturally to the code

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

不存在与莫比乌斯功能更好的实现，并且它可以memoized速度，但是这应该很容易理解。

There exist better implementations of the Möbius function, and it could be memoized for speed, but this should be easy enough to follow.

在Eul​​er函数的比较明显的计算是

The more obvious computation of the totient function is

在换句话说，完全因子的数量成独特的素数及指数，并做一个简单的乘法到了。

In other words, fully factor the number into unique primes and exponents, and do a simple multiplication from there.

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 的更好实现，但这是为了便于阅读。

Again, there exist better implementations of prime_factors, but this is meant for easy reading.

＃辅助函数

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的第一次尝试

def totient1(n):
    counter = 0
    for i in xrange(1, n):
        if gcd(i, n) == 1:
            counter += 1
    return counter

＃OP的第二次尝试

# I don't understand the algorithm, and just copying it yields inaccurate results

＃莫比乌斯反演

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

＃传统配方

def totient3(n):
    return int(reduce(mul, (1 - 1.0 / p for p, exp in factorize(n)), n))

＃传统配方，无师

def totient4(n):
    return reduce(mul, ((p-1) * p ** (exp-1) for p, exp in factorize(n)), 1)

使用这个code计算1所有号码的totients 9999我的桌面上，平均超过5奔跑，

Using this code to calculate the totients of all numbers from 1 to 9999 on my desktop, averaging over 5 runs,

• totient1 需要永远
• totient2 需要10秒
• totient3 需要1.3s
• totient4 需要1.3s

• totient1 takes forever
• totient2 takes 10s
• totient3 takes 1.3s
• totient4 takes 1.3s

这篇关于多少个数字低于N的coprimes到N？的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！