素数の計算（３）

 素数の計算（２）からの続き

JuliaのPrime.jlのprime()がおかしいのであれば，n番目までの素数が登場する範囲の整数の範囲がわかれば，その整数をすべてisprime()でチェックしたほうが早いのではないかという考えで，n番目の素数から整数範囲を求める方法を調べた。

素数定理によれば，素数の数$n=\pi(x)$は${\rm Li}(x)=\int_2^x \frac{1}{log t}dt$で与えられるということだ。対数積分は，Pythonのライブラリにはあるが，JuliaのSpecialFunctions.jlには入っていない。いずれにせよ近似値でよいのだから，素数定理の$\pi(x)$の表をながめて，えいやっと逆関数の近似式を求めた。

$x=\pi^{-1}(n)=10^{\dfrac{9.11*\log_{10}n+6.5}{8.5}}$

残念ながらこれで計算したprime()の代替関数はオリジナルの15倍時間がかかるという情けない結果に終わってしまった。猿の浅知恵の典型的パターンである。

using Primes

function logint(n)
  x=10.0^((9.11*log10(n)+6.5)/8.5)
  m = BigInt(floor(x))
  println(m)
  return m
end

function test(n)
  m = logint(n)
  i = 0
  for p in 1:m
    if(isprime(p))
      i = i+1
      if(i == n)
        println(p)
        break
      end
    end
  end
end

@time test(10^6)
@time println(prime(10^6))

15678124
15485863
32.740030 seconds (124.81 M allocations: 4.688 GiB, 12.73% gc time)
15485863
2.493369 seconds (2.53 M allocations: 38.664 MiB, 0.99% gc time)