日: 2024年10月27日

ペル方程式を解くコード

連分数の計算 を元にペル方程式を計算する。ペル方程式の連分数を用いた魔法の解法 – tsujimotterのノートブックを参考にした。

ペル方程式は以下のように表わされる。特にDが素数の時に興味がある。

\[ x^{2} – Dy^{2} = 1 \]

まず、D=61で計算してみる。計算手順は以下のようになる。

  • \(\sqrt{61}\) の連分数展開を求める
  • 連分数の周期を確認し、奇数なら2周期とする
  • 近似分数 \(\frac{P}{Q}\) を求める
  • \(P^{2} – 61Q^{2} = 1\) となる

cfrac = get_continued_fraction(61)
if 1 == (len(cfrac[1:]) % 2):
    cfrac += cfrac[1:]
get_approx_fraction(cfrac[:-1])

つぎに115以下の素数を ペル方程式 – Wikipedia にある一覧表で確認する。 

import sympy

lst_pq = [['D(prime)', 'P', 'Q'], None]
for i, p in enumerate(sympy.primerange(115)):
    cfrac = get_continued_fraction(p)
    if 1 == (len(cfrac[1:]) % 2):
        cfrac += cfrac[1:]
    pq = get_approx_fraction(cfrac[:-1])
    lst_pq.append([p, pq.numerator, pq.denominator])
lst_pq
D(prime) P Q
2 3 2
3 2 1
5 9 4
7 8 3
11 10 3
13 649 180
17 33 8
19 170 39
23 24 5
29 9801 1820
31 1520 273
37 73 12
41 2049 320
43 3482 531
47 48 7
53 66249 9100
59 530 69
61 1766319049 226153980
67 48842 5967
71 3480 413
73 2281249 267000
79 80 9
83 82 9
89 500001 53000
97 62809633 6377352
101 201 20
103 227528 22419
107 962 93
109 158070671986249 15140424455100
113 1204353 113296