前回の記事「PythonとRustで！友愛数の組を見つけてみよう！」では、友愛数の定義からPythonとRustを使ってコードを実装し、友愛数の組を探しました。
今回は「友愛数入門! 数の面白シリーズ!」の記事の中で紹介したサービト・イブン=クッラの方法を使用して友愛数の組を探します。
友愛数とサービト・イブン=クッラの方法
まずは友愛数の定義を思い出します。
２つの自然数aとbが友愛数であるとは、aの真の約数の和がbに等しく、かつbの真の約数の和がaに等しいことである。
以下は友愛数の組の例です。
- (220,284)
- (1184, 1210)
- (2620, 2924)
- (5020, 5564)
- (6232, 6368)
友愛数を具体的に探そうとするとかなり大変です。
友愛数を探す手段として、次のサービト・イブン＝クッラの方法があります。
ただし、この方法で探すことができる友愛数の組はかなり限定的です。
サービト・イブン=クッラの方法で探すことができない友愛数の組の例として、(6232, 6368) があります。
整数n \geq 2に対して、次のようにp,q,rを定める:
p = 3 \times 2^{n-1} – 1
q = 3 \times 2^n -1
r = 9 \times 2^{2n-1} – 1
この時、p,q,rが素数であるならば、(2^npq,2^nr)は友愛数の組となる。
Pythonでサービト・イブン=クッラの方法を実装しよう
Pythonでサービト・イブン=クッラの方法を実装します。
ソースコードは以下になります。
# 素数判定を行う関数
def is_prime(n):
# 2未満の数は素数ではない
if n < 2:
return False
# (*1) 2から√nまでの数で割り切れるか調べる
for i in range(2,int(n**0.5) + 1):
if n % i == 0:
# 割り切れたら素数ではない
return False
# 割り切れる数がなければ素数
return True
# 友愛数を探す関数
# limit : サービト・イブン=クッラの方法に出てくる整数nの範囲の最大値
def thabit_ibn_qurra_method(limit):
"""友愛数のペアをサービト・イブン=クッラの方法で探す"""
result = [] # 友愛数の組のリストを初期化
for n in range(2, limit):
# サービト・イブン=クッラの方法に従って計算
p = 3 * (2 ** (n - 1)) - 1
q = 3 * (2 ** n) - 1
r = 9 * (2 ** (2 * n - 1)) - 1
# p,q,rがすべて素数である必要がある
if is_prime(p) and is_prime(q) and is_prime(r):
# 友愛数の候補を計算
a = 2 ** n * p * q
b = 2 ** n * r
# 友愛数の組をリストに追加
result.append((a, b))
return result
if __name__ == '__main__':
amicable_pairs = thabit_ibn_qurra_method(20)
for a, b in amicable_pairs:
print(f"{a} と {b}")
実行すると次の結果が表示されます。
220 と 284 17296 と 18416 9363584 と 9437056
ソースコードについて、解説します。
is_prime関数は完全数のところでも出てきましたが、与えられた引数nが素数かどうかを判定するものです。
thabit_ibn_qurra_method関数がサービト・イブン=クッラの方法を使う関数です。
limitがサービト・イブン=クッラの方法に出てくる整数nの動く範囲の最大値になります。
この関数は次のような処理を行っています。
- 2からlimit-1までの各整数nをサービト・イブン=クッラの方法に出てくる3つの式に代入
- その3つの式により算出される値全てが素数であるかをis_prime関数でチェック
- 3つの値がすべて素数ならば、(2^npq,2^nr)は友愛数の組になるので、その２つの値を計算して友愛数の組のリストに追加
- limit-1までの全てのnを調べたら、友愛数の組のリストを返却
Rustでもサービト・イブン=クッラの方法を実装しよう
Rust版のソースコードは以下になります。
// 素数判定を行う関数
fn is_prime(n: u128) -> bool {
// 2未満の数は素数ではない
if n < 2 {
return false;
}
// 2から√nまでの数で割り切れるか調べる
let sqrt_n = (n as f64).sqrt() as u128;
for i in 2..sqrt_n {
if n % i == 0 {
// 割り切れたら素数ではない
return false;
}
}
// 割り切れる数がなければ素数
true
}
fn thabit_ibn_qurra_method(limit: u32) -> Vec<(u128, u128)> {
let mut result = Vec::new(); //完全数を格納するベクタ
for n in 2..limit {
// サービト・イブン=クッラの方法に従って計算
let p = 3 * 2u128.pow(n - 1) - 1;
let q = 3 * 2u128.pow(n) - 1;
let r = 9 * 2u128.pow(2 * n - 1) - 1;
// p,q,rがすべて素数である必要がある
if is_prime(p) && is_prime(q) && is_prime(r) {
// 友愛数の候補を計算
let a: u128 = 2u128.pow(n) * p * q;
let b: u128 = 2u128.pow(n) * r;
// 友愛数の組をリストに追加
result.push((a, b));
}
}
result
}
fn main() {
let limit: u32 = 60;
let amicable_pairs = thabit_ibn_qurra_method(limit);
for (a, b) in amicable_pairs {
println!("{} と {}", a, b);
}
}
やっていることはPythonと同じです。
limitの値が60を超えるとサービト・イブン=クッラの方法に出てくる計算式がu128の値を超えるのでオーバーフローして実行途中でエラーが吐かれます。
こちらも結果は以下になります。
220 と 284 17296 と 18416 9363584 と 9437056
まとめ
サービト・イブン=クッラの方法を使って、PythonとRustで友愛数の組を探しました。
この方法でn=2からn=60までで探すことができた組は
- (220, 284)
- (17296, 18416)
- (9363584, 9437056)
でした。
友愛数の組の例を上で紹介していますが、(1184, 1210)や(2620, 2924)など、この方法で探すことができていません。
そういう意味でサービト・イブン=クッラの方法は友愛数の組を探すという意味でかなり限定的な方法と言えます。
これまで
- 「友愛数入門! 数の面白シリーズ!」
- 「PythonとRustで！友愛数の組を見つけてみよう！」
- そして今回の「Python×Rustで数学体験！サービト・イブン=クッラの方法で友愛数を見つけよう」
という３本の記事で、友愛数の定義や、友愛数の組を探す方法を解説して、プログラムで実際に定義とサービト・イブン=クッラの方法を使用して友愛数の組を探してみました。
次回は、素数の「並び方」に注目してみます。
素数の出現には、実は不思議な“リズム”があるのです。
Pythonで素数の間隔を可視化し、さらには音に変換して聴いてみる──
数学の世界を目と耳で体験する「素数ビート」の記事をお楽しみに！
