1.プログラム


題記の不定方程式について、-100 ≤ n ≤ 100 の範囲で解を探してみる。
まず、1 ≤ |x| ≤ |y| ≤ |z| としても、一般性を失わない。

(x,y,z) が 1 より大きい最大公約数を持つ場合、i.e. d=gcd(x,y,z)>1 のとき、x=dx', y=dy', z=dz' とおくと、

  (x+y+z)(1/x+1/y+1/z)
= (dx'+dy'+dz')(1/dx'+1/dy'+1/dz')
= (x'+y'+z')(1/x'+1/y'+1/z')

となるので、最初から gcd(x,y,z)=1 としてよい。

また、l=lcm(x,y,z), x'=l/x, y'=l/y, z'=l/z とおくと、

  (x+y+z)(1/x+1/y+1/z)
= (lx'+ly'+lz')(1/lx'+1/ly'+1/lz')
= (x'+y'+z')(1/x'+1/y'+1/z')

となり、(x',y',z') も解となる。これを双対な解と呼ぶことにする。

この式が、d= -1 についても成り立つことを考慮すると、調べるべきパターンは、

( x, y, z)
(-x, y, z)
( x,-y, z)
( x, y,-z)

の 4通り。プログラムは以下のようになる。

 10   ' xyz.ub
 20   M=1000
 30   for X=1 to M
 40     for Y=X to M:D=gcd(X,Y)
 50       for Z=Y to M:D=gcd(D,Z):if D>1 then 110
 60         if and{X=Y,Y=Z} then 110
 70         R=fnSub(X,Y,Z)
 80         R=fnSub(-X,Y,Z)
 90         R=fnSub(X,-Y,Z)
100         R=fnSub(X,Y,-Z)
110       next Z
120     next Y
130   next X
140   end
150   fnSub(X,Y,Z)
160   local N
170   N=(X+Y+Z)*(1//X+1//Y+1//Z)
180   if or{N=0,N=1} then 220
190   if den(N)>1 then 220
200   if abs(N)>100 then 220
210   print N;":";X;",";Y;",";Z
220   return(0)

UBASICでは、何も修飾子の付かない変数は多倍長変数を表す
多倍長変数は整数型に比べてメモリも喰うし、計算時間もかかる。
にもかかわらず、これまで、変数に多倍長型のみを使ってきたのは、
この方がタイプ数が少なく、タイプミスを防げるからである。
例えば、nという変数を使っていて、途中で nは整数型で十分と気付いて、n%に打ち直したとする。
このとき修正漏れがあれば、UBASICでは n と n% という変数を別物として扱うので、
予期しない結果が得られることとなり、かつ、この間違いは発見しにくい。

したがって、まず全て修飾子のない多倍長型で打ち込んでしまう。
次にvchgという、変数置き換えのコマンドを使い、整数型で充分なものは、整数型に置き換える。
このようにすると、タイプミスを防ぐことができる。

例えば上のプログラムでは、n以外は全て整数型でよい。
よって実行するときは、n以外の変数を vchg で整数型に置き換えて欲しい。
何をやっているのかわからない、という部分は、特にないと思う。
190行目は、整数かどうかの判定(分母 den(n) が 1 より大きい場合は分数。符号は分子に付く)。
また、

x+y+z=0 または 1/x+1/y+1/z=0 のとき、n=0
x+y=0, y+z=0, z+x=0 のどれかが成り立つとき、n=1

これらの解は大量に出てくるので、180行目の条件ではじいている。

とりあえず、1 ≤ |x| ≤ |y| ≤ |z| ≤ 1000 の範囲で探してみる。
このアルゴリズムは、O(n3) のオーダだから、
探査範囲を広げるときは、何かやり方を考えなければならない。

結果は次頁のとおり。


この章の目次

E-mail : kc2h-msm@asahi-net.or.jp
三島 久典