gooブログはじめました!

写真付きで日記や趣味を書くならgooブログ

ディオファントス方程式を解く

2023-02-19 08:07:58 | ブログ
 参考文献などには、ユークリッドの互除法を用いて二つの自然数の最大公約数を求める計算が例示されている。計算法の一つは、二つの自然数のうち大きい方の自然数から小さい方の自然数を、大小関係が逆転するまで何回か引く。その結果の数について、同様に大小関係が逆転するか一方が0になるまで何回か引く。このようにして、一方の数が0になったとき他方に残った数が答えとなる最大公約数である。

 もう一つの計算法は、二つの自然数のうち大きい自然数を小さい自然数で割って商と余りを求める。次に小さい数を得られた余りで割って商と余りを求める。この演算を余りが0になるまで繰り返す。このようにして余りが0になったときの被除数が答えとなる最大公約数である。

 次に、参考文献に習って、ax+by=1の形式のディオファントス方程式を、ユークリッドの互除法を用いて解いてみよう。ここで、aとbは互いに素である負でない整数とする。

 もしaまたはbが1なら、すなわち(a,b)=(1,0)または(0,1)であれば、(x,y)が(1,0)または(0,1)であることはすぐ分かる。そこで、計算式にユークリッドの互除法を適用して、方程式がこの形の計算式になるように変形する。

 aがbより大とする。aをbで割った商をq1、余りをr1とする。式で書くと、a=q1b+r1である。この式を最初の式に代入すると、
   r1x+b(q1x+y)=1
であるから、x=x1,q1x+y=y1とおくと、r1x1+by1=1を得る。

 次にbをr1で割った商をq2、余りをr2とすると、b=q2r1+r2である。この式を上の式に代入すると、
   r1(x1+q2y1)+r2y1=1
であるから、x1+q2y1=x2,y1=y2とおくと、r1x2+r2y2=1を得る。

 この手順を繰り返し、xkの係数かykの係数が1になるまで続ける。そうすると、(x,y)=(1,0)か(0,1)が解になるから、その解を、手順を逆にたどって、元の式の解を得る。

 例として、方程式16x+7y=1の解を求めてみよう。
   16=7×2+2
2x+7(2x+y)=1
2x1+7y1=1; ここでx1=x,y1=2x+y
7=2×3+1
2(x1+3y1)+1×y1=1
2x2+1y2=1; ここでx2=x1+3y1,y2=y1
x2=0, y2=1
y1=1, x1+3=0
x=x1=-3, 1=2x+yからy=7
x=-3, y=7

 次に、x^2-Dy^2=1の形式のディオファントス方程式を解いてみよう。Dは自然数である。この方程式は、ペル方程式と呼ばれ、(x,y)が自然数となるような解を無限個もつことが知られている。

 D=2の場合の(x,y)が最小となる解は、(3,2)であることは式の形からすぐ分かる。(3,2)が求まれば、
   3^2-2×2^2=(3+2SQRT(2))(3-2SQRT(2))=1
であるから、両辺を2乗すると、17^2-2×12^2=1となり、次の解は(17,12)であることが分かる。

 D=7の場合には、x^2=7y^2+1から、yの平方数を7倍して1を加えたものは、xの平方数であるから、7y^2=63とおくと、y=3、したがってx^2=64からx=8、つまり(8,3)が求まる。

 D=10の場合には、x^2-1=10nとおいて、y=4,5,6の場合の平方数を求めると、n=16,25,36である。そうすると、10n+1は、各々161,251,361である。このうちxの平方数となっているのは361=19^2であるから、(19,6)が求まる。

 それでは、D=13の場合の解はどうなるのだろうか。解のxは大きな数になることが予想されるので、xの値はおよそySQRT(13)になるとみてよい。SQRT(13)の値は、電卓を用いて3.6055くらいまでは容易に計算できる。そうすると、y=10くらいから始めてySQRT(13)を計算し、その値が整数に近いものを選んで、x^2-13y^2=1になるような(x,y)を探索することが考えられる。しかし、y=100まで計算しても解が見つからないので、この方法の限界を悟ることになる。

 そこで、参考文献に習って、連分数を用いて解の探索計算を試みることにした。連分数を用いないで解を求める技法があるのかどうかは分からない。

 与えられた数を連分数で表現するには、次の操作を行えばよい。与えられた数を0または自然数とxとの和で表す。ここでxは0と1との間にある実数とする。まずxの逆数をとると、分子は1、分母は1より小さな実数で表せる。分母の実数をその最も近い整数a1と正の小数部分x1との和で表す。x1は分数で表す必要がある。多くの場合、x1の逆数をとるのが簡単であるが、一般的にはx1の正当な分数表現であればよく、分子が1でなくてもよい。次にx1に対して同様の操作を行い、以下この操作を繰り返すと、与えられた数の連分数展開ができる。

 与えられた数がSQRT(13)のような無理数であると、難しく感じられるかも知れないが、操作パターンが一本道に決まっているので、ワンパターンの機械的計算を繰り返すだけである。

 SQRT(13)-3=1/xとおくと、x=(SQRT(13)+3)/4が求まるから、次は(SQRT(13)+3)/4=1+1/xとおくと、x=(SQRT(13)+1)/3が求まる。以下同じ操作の繰り返しとなる。

 連分数は、下位の分数が何段も続いて無限回並ぶが、上位の段から始めて多くの下位段をとるほど近似の精度が良くなる。SQRT(13)の場合、最初の3を0近似とよび、次の3+1/1=4までとると1次近似、その次の段の3+1/(1+1)=7/2=3.5までとると2次近似、・・・という具合である。8次近似までとると、393/109=3.6055となり、電卓による計算の精度と一致してくる。

 SQRT(13)の場合、1次近似から4次近似まで1+が続くが、5次近似で6+が現れ、その後6次近似から9次近似まで1+が続き、10次近似で再び6+が現れる。また、4次近似は1+1/(SQRT(13)+3)という分母であるが、9次近似は同じ分母となっている。つまり、9次近似より下の段は4次近似から8次近似までのパターンを繰り返す循環連分数になっていることが分かる。

 SQRT(13)の5次近似ですでに119/33=3.606となっているので、5次近似以後の分子^2-13×分母^2を計算すると、
   5次近似: 119/33=3.606; 119^2-13×33^2=+4
   6次近似: 137/38=3.6052; 137^2-13×38^2=-3
   7次近似: 256/71=3.6056; 256^2-13×71^2=+3
   8次近似: 393/109=3.6055; 393^2-13×109^2=-4
   9次近似: 649/180=3.60556; 649^2-13×180^2=+1
となり、9次近似で方程式を満足する(x,y)=(649,180)の組が現れている。

 連分数の計算をしたことのない人には、連分数を使う計算を体験することをお勧めしたい。参考文献にあるように、円周率(パイ)や自然対数の底(e)の連分数展開は、自然数のみで構築される美しいパターン構成になっていることが分かる。また、連分数には、循環連分数にみられるように無理数の小数点以下の表示では現れない規則正しいパターンが潜んでおり、数学がもつ深遠な構造を知ることができる。

 参考文献
 木村俊一著「連分数のふしぎ」(ブルーバックス)
 マーク・カチ他著「数学の展開」(エンサイクロペディア ブリタニカ)