双子素数の予想、ブルン定数、ハーディ・リトルウッド定数

1. 双子素数

双子素数

　$1$ とその数自身以外に約数をもたない数のことを素数といいます：
　
\[2,\ 3,\ 5,\ 7,\ 11,\ 13,\ 19,\ 23,\ 29,\ 31,\ 37,\ 41,\ ......\]
　$2$ を除いたすべての素数は奇数です。隣り合う素数の差が $2$ となるような素数のペア $(p,\ p+2)$ を 双子素数 (twin prime) とよびます。$(2,\ 3)$ は、その差が $1$ なので双子素数には含まれません。双子素数を具体的に書き並べてみると
　
\[(3,\ 5),\ (5,\ 7),\ (11,\ 13),\ (17,\ 19),\ (29,\ 31),\ (41,\ 43),\ \cdots\]
のようになります。

双子素数生成関数

　Python で 100 未満の双子素数を生成してみましょう。
　sympy.primerange() で 3 以上 100 未満の素数リストを生成して NumPy の配列に格納します。

# PYTHON_TWIN_PRIME_01-1

import numpy as np
from sympy import primerange

# 素数イテレータで素数リストを作成
pn = list(primerange(3, 100))

# 素数リストをNumPy配列に格納
pn = np.array(pn)

print(pn)

[ 3  5  7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97]

　位置をスライドさせた素数リストを生成します。

# PYTHON_TWIN_PRIME_01-2

# 先頭～末尾の1つ手前までをスライス
x = pn[0:pn.shape[0]-1]

# 先頭の1つ後ろ～末尾までをスライス
y = pn[1:pn.shape[0]]

print("x:\n{}".format(x))
print("y:\n{}".format(y))

x:
[ 3  5  7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89]
y:
[ 5  7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97]

　マスク操作によって、x と y の要素同士の差が 2 となるような配列を生成します。これは双子素数のペアの片側となります。

# PYTHON_TWIN_PRIME_01-3

# xとyの要素同士の差が2の配列を生成
z=x[y-x==2]

print(z)

[ 3  5 11 17 29 41 59 71]

　この配列に 2 を加えると、ペアのもう片方が生成されます。

# PYTHON_TWIN_PRIME_01-4

# zの各要素に2を加える
w = z + 2

print(w)

[ 5  7 13 19 31 43 61 73]

　numpy.vstack() で z と w を結合して転置すると、双子素数を格納した配列が完成します。

# PYTHON_TWIN_PRIME_01-5

# zとwを結合して転置
tp = np.vstack((z, w)).T

print(tp)

[[ 3  5]
 [ 5  7]
 [11 13]
 [17 19]
 [29 31]
 [41 43]
 [59 61]
 [71 73]]

　イテレータとして双子素数を格納することもできます。

# PYTHON_TWIN_PRIME_01-6

# 双子素数のイテレータ
itp = iter(tp.tolist())

next(itp)

[3, 5]

　以上の処理を関数にまとめておきます。

# PYTHON_TWIN_PRIME_02-1

# a以上b未満の双子素数と個数
def twin_prime(a, b, inv=False):
    pn = list(primerange(a, b))
    pn = np.array(pn)
    x = pn[0 : pn.shape[0]-1]
    y = pn[1 : pn.shape[0]]
    z = x[y - x == 2]
    w = z + 2
    tp = np.vstack((z, w)).T
    if inv == True:
        tp = 1/tp
    ct = tp.shape[0]
    return tp, ct

# twin_prime() をユニバーサル関数に変換
twin_prime = np.vectorize(twin_prime)

　戻り値は双子素数と双子素数の個数です。オプション引数 inv を True に設定すると双子素数の逆数を返します。twin_prime() を使って 100 以上 200 未満の双子素数を探してみましょう。

# PYTHON_TWIN_PRIME_02-2

# 100以上200未満の双子素数
tp = twin_prime(100, 200)

print(tp)
print(ct)

[[101 103]
 [107 109]
 [137 139]
 [149 151]
 [179 181]
 [191 193]
 [197 199]]
7

双子素数の予想

　「双子素数は無数に存在するのか、あるいは有限個しかないのか」という問題 (双子素数の予想) は、いまだに数学上の未解決問題です。

双子素数の逆数和

　双子素数の逆数和
　
\[\left(\frac{1}{3}+\frac{1}{5}\right)+\left(\frac{1}{5}+\frac{1}{7}\right)+\left(\frac{1}{11}+\frac{1}{13}\right)+\cdots+\left(\frac{1}{p}+\frac{1}{p+2}\right)+\cdots\]
は収束することは証明されています。値はおよそ $1.9$ ぐらいで、ヴィーゴ・ブルンによって証明されたので ブルン定数 とよばれています。素数の逆数和は発散するので、双子素数は素数に比べて少ないことになりますが、無限級数が収束することもあるので、これだけで一概に双子素数が無限個あるとは断言できません。
　
　双子素数の逆数和の収束は非常に遅いので、数値計算でブルン定数 $1.9$ 強に達するのは大変です。特に私の激安ノート PC だと本当にきついですが、ぎりぎり頑張って 1 億以下の双子素数の逆数和を計算してみました。

# PYTHON_TWIN_PRIME_02-3

# 1億未満の双子素数と個数
tp = twin_prime(1, 1e8, inv=True)[0]

# 双子素数の逆数を全部足す
s = np.sum(tp)

print(s)

1.7588156210679744

　...... 全然足りていませんね。本気で計算してみたい人は、パワフルなデスクトップ PC を用意したほうがいいです。学生さんなら大学のコンピュータを使ってください。社会人のかたは会社の PC を ...... もちろん冗談です (会社の PC を私用に使うのはやめましょう)。あと、twin_prime() の内部で使用している SymPy の primerange() があまり速くないような気がするので、NumPy でもっと効率の良い素数生成関数を作成したほうがいいかもしれません (良いプログラムが書けたらコメントで教えてください)。

双子素数の漸近公式

　ハーディ (G.H.Hardy) とリトルウッド (J.E.Littlewood) は、$x$ 以下の双子素数の組の個数が漸近的に
　
\[N(x)=2C\int_{2}^{x}\frac{x}{(\log x)^2}\tag{1}\]
で与えられると予想しています。$C$ は $3$ 以上の素数についての $1-1/(p-1)$ の総積
　
\[C=\prod_{p\gt 2}\left\{1-\frac{1}{(p-1)^2}\right\}\tag{2}\]
で計算される定数であり、ハーディ・リトルウッド定数 とよばれています。$C$ は mpmath.twinprime から取得できます。

# PYTHON_TWIN_PRIME_02-4

from mpmath import mp, twinprime

# 精度を50桁に設定
mp.dps = 50; mp.pretty = True

# ハーディ・リトルウッド定数
c = twinprime

print(c)

0.66016181584686957392781211001455577843262336028473

　漸近式 (1) に含まれる積分
　
\[\int_{2}^{x}\frac{x}{(\log x)^2}\tag{3}\]
は普通の形で表すことのできない積分関数であり、一般に数値計算 (積分の近似計算) によって値を取得します。SciPy の integrate() を使えば、それほど難しくありません。

# PYTHON_TWIN_PRIME_02-5

import numpy as np
from scipy import integrate

# 被積分関数を定義
f = lambda x : 1 / np.log(x)**2

# fを2からxまで積分
g = lambda x : integrate.quad(f, 2, x)[0]

# gをユニバーサル化
g = np.frompyfunc(g, 1, 1)

　引数に配列を渡せるように、積分関数 g(x) をユニバーサル化しておきました。Matplotlib を使って、$100,000$ 未満の双子素数の個数の漸近曲線 (asymptotic curve) と真値 (true value) を確認しておきます。

# PYTHON_TWIN_PRIME_02-6

import matplotlib.pyplot as plt

# プロットするデータを用意
t = np.arange(2, 100000, 1000)
y = 2*c*g(t)
z = twin_prime(1, t)[1]

# FigureとAxesの設定
fig = plt.figure(figsize = (8, 6))
ax = fig.add_subplot(111)
ax.grid()
ax.set_title("Number of twin primes less than x", fontsize=16)
ax.set_xlabel("x", fontsize=14)
ax.set_ylabel("N(x)", fontsize=14)
ax.set_xlim(2, 100000)
ax.set_ylim(0, 1300)

# Axesに漸近曲線と真値をプロット
ax.plot(t, y, color = "blue", label="asymptotic curve")
ax.plot(t, z, color = "red", label=" true value")

# 凡例を表示
ax.legend()

双子素数 (twin prime) の個数とハーディ・リトルウッドの漸近線
　ハーディ・リトルウッドの漸近曲線を眺めると、双子素数は無限に存在するように思えますが、この漸近曲線自体があくまで予想なので、本当に正しい分布公式なのかどうかわかっていません。（参考：Wikipedia）

コメント（1件）

あとりえこばとより:

2023年5月12日 4:30 PM

ChatGPTに双子素数について講義してもらいました

　
双子素数は、素数のうち非常に特別な性質を持つ素数のペアです。具体的には、2つの素数が2つの差で隣り合っていることが特徴です。例えば、3と5、11と13、17と19などが双子素数のペアとなります。
　双子素数の存在と性質についての研究は、数学の分野で古くから行われています。有名な「ゴールドバッハ予想」という未解決問題と関連しており、この予想は、すべての偶数は2つの素数の和として表されるというものです。ゴールドバッハ予想が正しいとすると、双子素数は無限に存在することが示唆されます。双子素数の性質については、いくつかの重要な結果があります。例えば、双子素数のペアは無限に存在するという「ツインプライム予想」という予想がありますが、現在でも未解決のままです。また、双子素数は非常に希少であり、その出現パターンや性質についての詳細な理解はまだ不十分です。
　双子素数の研究は数学の興味深いトピックであり、素数理論や解析数論といった分野で重要な役割を果たしています。数学者たちは、双子素数の存在や性質についての洞察を深めるために、さまざまな手法と数学的なツールを使用しています。双子素数についての理解は、数学の研究のみならず、暗号学や情報セキュリティなどの応用分野でも重要です。素数の性質や分布を理解することで、暗号アルゴリズムの安全性や効率性を評価する上で役立ちます。双子素数は、数学の魅力的なテーマであり、その研究によって数学的な発展や新たな知見が生まれることが期待されています。

返信