Friday, April 30, 2021

[accudfui] Pell equation notes

regular (positive) Pell equation is x^2 - D*y^2 = +1.

negative Pell equation is x^2 - D*y^2 = -1.  if you have a solution (x,y) to the negative Pell equation, run Newton's method computing sqrt(D) for one iteration, starting with the guess old=x/y.  that is, new = (old + D/old)/2 (the Babylonian method is equivalent to Newton's method).  the components of the resulting fraction newx/newy will be a solution to the regular Pell equation: newx^2 - D*newy^2 = +1.

solving the Pell equation with continued fractions is very straightforward.  compute the continued fraction of sqrt(D), try all truncations in order of length.  for each truncation, convert to a regular fraction x/y and check if x and y solve the regular Pell equation.  there is no need to systematically modify the last term before the truncation point as you do when using continued fractions to find the sequence of best fractions that approximate a given number.

as an optimization, also check whether each (x,y) solves the negative Pell equation, then use the Newton's method trick above to convert the negative equation solution to a regular equation solution.  the smallest negative solution (if one exists) will convert to the smallest regular solution.

the following sequence of D require increasing numbers of continued fraction terms to solve.  the Newton's method optimization is included, so we are counting the number of continued fraction terms needed to reach the first regular or negative solution.  these D represent challenges of increasing difficulty.  the famous value 61 is not on the list.  61 requires only 11 continued fraction terms to find a negative equation solution (x=29718, y=3805).  it is eclipsed by the smaller 46 which requires 12 terms for its smallest solution (x=24335, y=3588), which happens to solve the regular equation.

2 3 7 13 19 31 43 46 94 139 151 166 211 331 421 526 571 604 631 751 886 919 1291 1324 1366 1516 1621 1726 2011 2311 2566 2671 3004 3019 3334 3691 3931 4174 4846 5119 6211 6451 6679 6694 7606 8254 8779 8941 9739 9949 10399 10651 10774 12541 12919 13126

the sequence below is a subsequence of the sequence above.  for these record-setting D values, the solution is to the negative equation, so you need to do one Newton's method iteration afterward, so the smallest regular solution unusually large.

2 13 421 1621 8941 9949 12541 17341 39901 40429 43261

some disorganized Haskell source code is here.  it includes numerical testing of the conjectures stated above.

we used the quadratic-irrational package to compute the continued fraction representation of square roots.  we used the continued-fraction package to convert continued fractions to Rational.  below is code to convert a quadratic irrational into an infinite list of continued fraction terms, using the cycle function.  Haskell handles infinite data structures elegantly.

qitocf :: Quadraticirrational.QI -> [Integer];
qitocf = Quadraticirrational.qiToContinuedFraction >>> \case { (integerpart, Quadraticirrational.CycList start repeatend) -> integerpart:(start ++ cycle repeatend);};

having constructed an infinite list, use the inits function to construct all truncations, then for each truncation, test whether it solves the Pell equation.  (but our code does not use inits because we were also testing whether the final term needs to be modified.)

No comments :