$J(x)$ と第2項について
文献[1, P22]より
$$ J(x)=\frac{1}{2}\Biggl[\sum_{p^n\lt x}\frac{1}{n}+\sum_{p^n\leq x}\frac{1}{n}\Biggr] $$と定義すると $J(x)$ は素数のべき乗 $p^n$ が現れるたびに $1/n$ だけ増加する関数となる. ただし $x$ がちょうど素数のべき乗になるときは $1/2n$ だけ増加する.
$J(x)$ はリーマンゼータ関数の零点を用いて表わせる[1, P48].
$$ J(x)=\text{Li}(x)-\sum_{\rho}\text{Li}(x^\rho)-\log 2+\int_{x}^{\infty}\frac{dt}{t(t^2-1)\log t}\ \ \ (x\gt 1) $$$J(x)$ の第2項は, 零点のペア $\rho,\overline\rho$ を足していくことによって収束が保証される. より厳密に書こうとすれば
$$ -\sum_{\rho}\text{Li}(x^\rho) = -\lim_{T\to\infty}\sum_{\text{Im}(\rho)\leq T}\{\text{Li}(x^\rho)+\text{Li}(x^{\overline\rho})\} $$となる. ここで, クリティカルライン上の零点を虚部の小さい順に並べてみよう.
$$ \rho_1=1/2+i14.1347...\\ \rho_2=1/2+i21.0220...\\ \rho_3=1/2+i25.0108...\\ \rho_4=1/2+i30.4248...\\ \rho_5=1/2+i32.9350...\\ \cdots $$クリティカルライン上には無限に零点が存在するから, 第2項を厳密に計算することはできない. 有限個の項で計算するために, 関数 $f(x)$ を以下で定義する.
$$ f(x)=-\sum_{i=1}^{500}\{\text{Li}(x^{\rho_i})+\text{Li}(x^{\overline\rho_i})\} $$次に関数 $f(x)$ の実装について述べる.
プログラム
- プログラムをPythonで作成する
- $x$ を $10$ から $50$ までとり, $f(x)$ をプロットする
- リーマンゼータ関数の零点は Odlyzko による数値計算の結果1を用いる
- 対数積分 $\text{Li}(x)$ の計算には指数積分との関係式を利用し, scipy の expi を活用する
- $x$ が素数または素数のべき乗のときの $f(x)$ の値を調べる
# -*- encoding:utf-8 -*-
import matplotlib.pyplot as plt
import numpy as np
from scipy.special import expi
# 計算に使用する零点の個数
nzeros = 500
# 零点の読み込み
# 零点のデータをダウンロードして適宜ファイルパスを修正してください
with open('/path/to/zeros.txt','r') as f:
lines = f.readlines()
zeros = [float(line.strip()) for line in lines[:nzeros]]
# 対数積分
def li(z):
return expi(np.log(z))
# J(x) の第2項を有限個の項で打ち切ったもの
def f(x):
s = 0.0
for t in zeros:
rho = 0.5+1j*t
s += (li(x**rho) + li(x**(1-rho))).real
return -s
# ラベル・目盛りの設定
plt.tick_params(labelbottom=True, labelleft=False, labelright=False, labeltop=False)
plt.title(r'$f(x)=-\sum_{i=1}^{500}\{Li(x^{\rho_i})+Li(x^{1-\rho_i})\}$')
plt.xlabel("x")
plt.ylabel(r'$f(x)$')
# 関数 f(x) のプロット
x = np.linspace(10,50,800)
y = [f(x_i) for x_i in x]
plt.plot(x, y)
# 素数と素数べきのプロット
primes = [11,13,17,19,23,29,31,37,41,43,47]
prime_powers = [16,25,27,32,49]
for p in primes:
plt.plot(p,f(p),'ro')
for p in prime_powers:
plt.plot(p,f(p),'go')
plt.show()
プロット
作成したプログラムを Python 3.11.9 で実行し, 下記の結果を得た.
赤丸は素数 (11,13,17,19,23,29,31,37,41,43,47) , 緑丸は素数のべき乗 (16,25,27,32,49) を表す. 零点は素数の位置を"知っている"ということが, このプロットから実感できる.
参考文献
[1] H.M.Edwards, Riemann’s Zeta Function, Dover, 2001.