0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

2元2次不定方程式をfind_DNでペル型方程式に変換して解く

Last updated at Posted at 2021-12-19

前回の投稿でペル型方程式の解をSympyのdiop_DNで生成する方法について書きました。今回はこれを利用して以下のような2元2次不定方程式を解くプログラムを書いて行きます。

$\displaystyle \large \ \ \ ax^2+bxy+cy^2+dx+ey+f=0$

手順は以下のようになります。

  1. 2元2次不定方程式をペル型方程式に変換、および元の解への変換Matrixを取得
  2. 前回の投稿で作った関数pelEqを使ってペル型方程式の解を生成
  3. 1.で取得したA,Bで元の解に変換

例として以下の2元2次不定方程式を使います

$\displaystyle \large \ \ \ 2x^2+4x-y^2-y=0$

1. 2元2次不定方程式をペル型方程式に変換、および元の解への変換Matrixを取得

sympyのfind_DNを使うとあっけなく$X^2-DY^2=N$の$D,N$が求まります。

(D, N) = find_DN(2*x**2 + 4*x - y**2 -y)
print(D, N)
# (8, 56)

同時にtransformation_to_DNを使ってて$X^2-DY^2=N$の解から元の$x,y$に変換する$Matrix A, B$を取得します。

(A, B) = transformation_to_DN(2*x**2 + 4*x - y**2 -y)
print(A,B)
# Matrix([[1/8, 0], [0, 1/2]]) Matrix([[-1], [-1/2]])

2. 前回の投稿で作った関数pelEqを使ってペル型方程式の解を生成

ここは前回の投稿を参考にしてください。

3. 1.で取得したA,Bで元の解$x_n,y_n$に変換

1.で取得した$A,B$を使って$X_n,Yn$を元の2元2次不定方程式の解$x_n,y_n$に変換します。式は以下ようになります。

$\displaystyle \large \ \ \ \begin{pmatrix} x_n \\ y_n\end{pmatrix} = A \times \begin{pmatrix} X_n \\ Y_n\end{pmatrix}+ B$

(x, y) = A*sympy.Matrix([[abs(X)], [abs(Y)]])+B

これらを組み合わせて2元2次不定方程式を解くプログラムdiophQuadが以下となります。結果の最後の値が$0$になっているので正しい答えが生成されているのが確認できます。

import sympy
from sympy.solvers.diophantine.diophantine import find_DN, diop_DN, transformation_to_DN
from sympy.parsing.sympy_parser import parse_expr

def findDNAB(eqstr):
  sympy.var("x y")
  eq = parse_expr(eqstr)
  (D, N) = find_DN(eq)
  (A, B) = transformation_to_DN(eq)
  return (D, N, A, B)

# solve quadratic equation with interger solutions
def diophQuad(eq):
  (D, N, A, B) = findDNAB(eq)
  ansPelEq = pelEq(D, N)  
  while True:
    xyl = next(ansPelEq)
    for (X,Y) in xyl:
      (x, y) = A*sympy.Matrix([[abs(X)], [abs(Y)]])+B
      yield (x,y)

ansDQ = diophQuad(" 2*x**2 + 4*x - y**2 -y")
for i in range(10):
  (x,y) = next(ansDQ)
  print(x, y, 2*x**2 + 4*x - y**2 -y)
#0 0 0
#0 0 0
#1 2 0
#3 5 0
#10 15 0
#22 32 0
#63 90 0
#133 189 0
#372 527 0
#780 1104 0

(開発環境:Google Colab)

このアルゴリズムはProject Euler Problem 321: Swapping Countersを解くのに役に立ちます。

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?