Pythonで素数判定(1)

どの素数判定法を使うべきか迷うことがあるのでまとめてみます。

試し割り法

ある整数が素数であるか2から順に割っていきます。nが素数かどうかの判定には√nまでの素数で割り、全て割り切れなかったら素数です。途中で割り切れたらその時点で合成数と分かります。
素数を列挙するコードを書いてみましょう。まず、素数を用意するのはめんどうなので2, 3, 4, …と割っていきます。

def is_prime1(n):
    return all(n % p != 0 for p in takewhile(lambda p: p * p <= n, count(2)))

def make_primes1(N):
    return [ n for n in xrange(2, N + 1) if is_prime1(n) ]

t0 = time.clock()
N = 10 ** 6
primes = make_primes1(N)
print len(primes)
print time.clock() - t0
78498
25.4130005749

takewhileを使って自乗して終了判定しています。√nまでの整数で割っていけばよさそうですが、そうではありません。たいていの場合すぐに割り切れてしまうためsqrtの計算が相対的に重くなってしまうからです。試してみましょう。

def is_prime1a(n):
    max_p = int(sqrt(n))
    return all(n % p != 0 for p in xrange(2, max_p + 1))
78498
9.92075910043

あ、あれ?

def is_prime1b(n):
    max_p = int(sqrt(n))
    return all(n % p != 0 for p in takewhile(lambda p: p <= max_p, count(2)))
78498
24.1319481871

どうやらtakewhileを使うと遅いようですね。こう書くと、

def is_prime1c(n):
    for p in count(2):
        if p * p > n:
            return True
        elif n % p == 0:
            return False
78498
12.3978489139

こうすると、まあまあ速いですね。しかし、C++だと自乗で判定したほうが速いんですが。
2, 3, 4, 5と割っていくより2, 3, 5, 7と割っていく方がもちろん速いです。

def is_prime2(n):
    if n == 2:
        return True
    elif n % 2 == 0:
        return False
    max_p = int(sqrt(n))
    return all(n % p != 0 for p in xrange(3, max_p + 1, 2))
78498
5.44453929431

あらかじめ素数を用意する方が速いはずです。

def is_prime3(n):
    max_p = int(sqrt(n))
    for p in ps:    # psは必要な素数の列
        if p > max_p:
            return True
        elif n % p == 0:
            return False
    return True
78498
2.54170416017

takewhileを使うと遅い以上こう書くしかないですね。関数型の書き方を覚えると許せない書き方ですが。
さて、この試し割りは計算量がわからないので調べてみます。10000〜10000000まで調べます。

O(n1.18)くらいということになっていますが、ちょっと下に凸っぽいですね。適当にいじってみたら、O(n1.5/log n)のほうが近そうな気がします。