http://blog.practical-scheme.net/shiro/20110929-primes に触発されて、Squeak Smalltalk で無限素数列ジェネレーター版と配列の特異メソッド版で書き、それぞれ Python（ジェネレーター版）と Ruby（特異メソッド版）とで対決させてみました。Cog VM、盤石です。
きっかけとなったジェネレーター版で爆速の Gauche が同じ環境でどのくらいかも調べたかったのですが、なぜかここのところのバージョンが手元の環境でうまくビルドできないのでやむなく断念。
使用マシンは MacBook Air (Mid 2011; 1.8GHz Core i7, Win7) 。Squeak Smalltalk は結果の第一要素が計測時間（ミリ秒）です。
追記： id:coze さんのご指摘で当方のケアレスミスでコードが元の Gauche のものと違って不必要なループを回していることが判明しましたので、念のためすべて書き直して再度計測しました。結果、全体的にずいぶんと速くなりました。同時に、Gauche 爆速の謎も解けました。ありがとうございます。なお差し替える前のコードは最後に移動してありますのであしからず。
ジェネレーター版
| primesGen result | primesGen := Generator on: [:g | | primes next | primes := OrderedCollection withAll: #(2 3 5 7). primes do: [:prime | g value: (next := prime)]. [ [:exit | [ next := next + 2. primes anySatisfy: [:prime | prime * prime > next ifTrue: [exit value]. next isDivisibleBy: prime] ] whileTrue] valueWithExit. g value: (primes add: next) ] repeat]. {[5000 timesRepeat: [result := primesGen next]] timeToRun. result} "=> #(23 48611) "
Python 3
def primesgen(): primes = [2, 3, 5, 7] for prime in primes: yield(prime) nxt = primes[-1] while True: nxt += 2 done = False for prime in primes: if prime * prime > nxt: done = True break if nxt % prime == 0: break if done: primes.append(nxt) yield(nxt) g = primesgen() n = 0 result = None while n < 5000: n += 1 result = next(g) print(result)
$ python3 -V Python 3.1.2 $ time python3 primes-gen3.py 48611 real 0m0.241s user 0m0.155s sys 0m0.077s
特異メソッド版
| primes result | primes := OrderedCollection withAll: #(2 3 5 7). primes assureUniClass. primes class compile: 'at: index | next | index <= self size ifTrue: [^super at: index]. next := self last. [self size < index] whileTrue: [ [:exit | [ next := next + 2. self anySatisfy: [:pn | pn * pn > next ifTrue: [exit value]. next isDivisibleBy: pn] ] whileTrue ] valueWithExit. self add: next]. ^next'. {[result := primes at: 5000] timeToRun. result} "=> #(23 48611) "
Ruby 1.9.3dev
primes = [2, 3, 5, 7] def primes.[](idx) return super if idx < size nxt = last + 2 while size <= idx while any?{ |pn| break if pn * pn > nxt; nxt % pn == 0 } nxt += 2 end self << nxt end nxt end result = nil 5000.times{ |i| result = primes[i] } p result
$ time ruby -v primes-gen2.rb ruby 1.9.3dev (2010-11-19 trunk 29837) [i386-cygwin] 48611 real 0m0.216s user 0m0.078s sys 0m0.077s
参考：Squeak 組み込みの #primesUpTo:
無限素数列の文脈からは離れてしまいますが、これぞまさに爆速です。
| result | {[result := (Integer primesUpTo: 48611+1) last] timeToRun. result} "=> #(4 48611) "
エラトステネスのふるいで実装されています。
Integer class >> primesUpTo: max "Return a list of prime integers up to the given integer." "Integer primesUpTo: 100" ^Array streamContents:[:s| self primesUpTo: max do:[:prime| s nextPut: prime]]
Integer class >> primesUpTo: max do: aBlock "Compute aBlock with all prime integers up to the given integer." "Integer primesUpTo: 100" | limit flags prime k | limit := max asInteger - 1. "Fall back into #largePrimesUpTo:do: if we'd require more than 100k of memory; the alternative will only requre 1/154th of the amount we need here and is almost as fast." limit > 25000 ifTrue:[^self largePrimesUpTo: max do: aBlock]. flags := (Array new: limit) atAllPut: true. 1 to: limit - 1 do: [:i | (flags at: i) ifTrue: [ prime := i + 1. k := i + prime. [k <= limit] whileTrue: [ flags at: k put: false. k := k + prime]. aBlock value: prime]].
Integer class >> largePrimesUpTo: max do: aBlock "Evaluate aBlock with all primes up to maxValue. The Algorithm is adapted from http://www.rsok.com/~jrm/printprimes.html It encodes prime numbers much more compactly than #primesUpTo: 38.5 integer per byte (2310 numbers per 60 byte) allow for some fun large primes. (all primes up to SmallInteger maxVal can be computed within ~27MB of memory; the regular #primesUpTo: would require 4 *GIGA*bytes). Note: The algorithm could be re-written to produce the first primes (which require the longest time to sieve) faster but only at the cost of clarity." | limit flags maskBitIndex bitIndex maskBit byteIndex index primesUpTo2310 indexLimit | limit := max asInteger - 1. indexLimit := max sqrt truncated + 1. "Create the array of flags." flags := ByteArray new: (limit + 2309) // 2310 * 60 + 60. flags atAllPut: 16rFF. "set all to true" "Compute the primes up to 2310" primesUpTo2310 := self primesUpTo: 2310. "Create a mapping from 2310 integers to 480 bits (60 byte)" maskBitIndex := Array new: 2310. bitIndex := -1. "for pre-increment" maskBitIndex at: 1 put: (bitIndex := bitIndex + 1). maskBitIndex at: 2 put: (bitIndex := bitIndex + 1). 1 to: 5 do:[:i| aBlock value: (primesUpTo2310 at: i)]. index := 6. 2 to: 2309 do:[:n| [(primesUpTo2310 at: index) < n] whileTrue:[index := index + 1]. n = (primesUpTo2310 at: index) ifTrue:[ maskBitIndex at: n+1 put: (bitIndex := bitIndex + 1). ] ifFalse:[ "if modulo any of the prime factors of 2310, then could not be prime" (n \\ 2 = 0 or:[n \\ 3 = 0 or:[n \\ 5 = 0 or:[n \\ 7 = 0 or:[n \\ 11 = 0]]]]) ifTrue:[maskBitIndex at: n+1 put: 0] ifFalse:[maskBitIndex at: n+1 put: (bitIndex := bitIndex + 1)]. ]. ]. "Now the real work begins... Start with 13 since multiples of 2,3,5,7,11 are handled by the storage method; increment by 2 for odd numbers only." 13 to: limit by: 2 do:[:n| (maskBit := maskBitIndex at: (n \\ 2310 + 1)) = 0 ifFalse:["not a multiple of 2,3,5,7,11" byteIndex := n // 2310 * 60 + (maskBit-1 bitShift: -3) + 1. bitIndex := 1 bitShift: (maskBit bitAnd: 7). ((flags at: byteIndex) bitAnd: bitIndex) = 0 ifFalse:["not marked -- n is prime" aBlock value: n. "Start with n*n since any integer < n has already been sieved (e.g., any multiple of n with a number k < n has been cleared when k was sieved); add 2 * i to avoid even numbers and mark all multiples of this prime. Note: n < indexLimit below limits running into LargeInts -- nothing more." n < indexLimit ifTrue:[ index := n * n. (index bitAnd: 1) = 0 ifTrue:[index := index + n]. [index <= limit] whileTrue:[ (maskBit := maskBitIndex at: (index \\ 2310 + 1)) = 0 ifFalse:[ byteIndex := (index // 2310 * 60) + (maskBit-1 bitShift: -3) + 1. maskBit := 255 - (1 bitShift: (maskBit bitAnd: 7)). flags at: byteIndex put: ((flags at: byteIndex) bitAnd: maskBit). ]. index := index + (2 * n)]. ]. ]. ]. ].
追記：コメント欄で あんちもん2 さんに、require "prime" で使える Prime がけっこう速いとの情報をいただきました。ありがとうございます。手元の環境での結果はこんな感じでした。
$ time ruby -v primes-gen-req.rb ruby 1.9.3dev (2010-11-19 trunk 29837) [i386-cygwin] 48611 real 0m0.270s user 0m0.109s sys 0m0.109s
そういえば以前、前述 Squeak 組み込みの #largePrimesUpTo:do: を Ruby に移植したことがあったなぁと思い出したので、こちらも試したところやはりそれなりに少しだけ高速でした。参考まで。
$ time ruby -v primes-gen-sq.rb ruby 1.9.3dev (2010-11-19 trunk 29837) [i386-cygwin] 4999 real 0m0.167s user 0m0.046s sys 0m0.078s
改訂前掲載のコード
| primesGen result | primesGen := Generator on: [:g | | primes next | primes := OrderedCollection withAll: #(2 3 5 7). next := primes last. 1 to: Float infinity do: [:idx | idx <= primes size ifTrue: [g value: (primes at: idx)] ifFalse: [ [ next := next + 2. primes anySatisfy: [:prime | prime * prime <= next and: [next isDivisibleBy: prime]] ] whileTrue. g value: (primes add: next)]]]. {[5000 timesRepeat: [result := primesGen next]] timeToRun. result} "=> #(879 48611) "
def primesgen(): primes = [2, 3, 5, 7] nxt = primes[-1] idx = 0 while True: if idx < len(primes): yield(primes[idx]) else: nxt += 2 while any(nxt % prime == 0 for prime in primes if prime * prime <= nxt): nxt += 2 primes.append(nxt) yield(nxt) idx += 1 g = primesgen() n = 0 result = None while n < 5000: n += 1 result = next(g) print(result)
$ python3 -V Python 3.1.2 $ time python3 primes-gen.py 48611 real 0m2.182s user 0m2.059s sys 0m0.109s
| primes result | primes := OrderedCollection withAll: #(2 3 5 7). primes assureUniClass. primes class compile: 'at: index | next | index <= self size ifTrue: [^super at: index]. next := self last. [self size < index] whileTrue: [ [ next := next + 2. self anySatisfy: [:pn | pn * pn <= next and: [next isDivisibleBy: pn]] ] whileTrue. self add: next]. ^next'. {[result := primes at: 5000] timeToRun. result} "=> #(875 48611) "
primes = [2, 3, 5, 7] def primes.[](idx) return super if idx < size nxt = last + 2 while size <= idx while any?{ |pn| pn * pn <= nxt and nxt % pn == 0 } nxt += 2 end self << nxt end nxt end result = nil 5000.times{ |i| result = primes[i] } p result
$ time ruby -v primes-gen.rb ruby 1.9.3dev (2010-11-19 trunk 29837) [i386-cygwin] 48611 real 0m2.458s user 0m2.340s sys 0m0.061s
