From: John Robertson Subject: Re: Using continued fractions to solve Pell's equations Date: Thu, 13 Apr 2000 23:00:08 GMT Newsgroups: sci.math To: rharper@intouchsurvey.com Summary: [missing] In article <38ECEF3B.40BD1E49@intouchsurvey.com>, Russell Harper wrote: > Does anyone know where I can find an online reference describing how to > use continued fractions to solve Pell's equations? > Thank you... > > Russell > rharper@intouchsurvey.com Let me give some references first, and then describe some algorithms. These are not online, but some good references are: Richard E. Mollin, Fundamental Number Theory with Applications, CRC Press, Boca Raton, 1998. Mollin gives the most complete set of methods for Pell’s equations of any modern reference. Chapter 5, pages 221 to 272, discusses the continued fractions generally. In particular, Section 5.3, pages 238 to 250 studies periodic continued fractions, and applies this study to find all solutions to the Pell equation x^2- dy^2=+/-1 for d>0, d not a square. The method below for computing the continued fraction expansion of a quadratic irrational is discussed in exercise 5.3.6, p. 251. While the method below for the +/-4 equation is not discussed explicitly, Mollin’s Theorem 5.3.4 on page 246 gives the main machinery needed to prove that that method is correct. Mollin also discusses the continued fraction expansion of (1+sqrt(d))/2 for d==1 (mod 4) in exercise 5.3.14 on page 252. For the general equation x^2-dy^2=m, for any positive nonsquare d and any m, limits to search on x or y are given correctly in Chapter 6 on pages 299 and 300. Theorem 5.2.5 on page 232 gives the main criterion for solving x^2-dy^2=m for d>0, d not a square, m^20, d not a square, m^2>d. Ivan Niven, Herbert S. Zuckerman, and Hugh L. Montgomery, An Introduction to the Theory of Numbers, Fifth Edition, John Wiley & Sons, Inc., New York, 1991. Sections 7.1 to 7.7, pages 325 to 351 cover continued fractions generally. Section 7.8, pages 351 to 356 covers the Pell +/-1 equation, and the case x^2-dy^2=n where n^2 d, but apart form the fact that the heading on page 482 is “Lagrange’s Chain of Reductions,” the section does not reference Lagrange. Section 19, on page 486, discusses the cases d<0 and d>0, d a square. Section 20, pages 486 to 488, reduces the general equation ax^2+2hxy+by^2+2gx+2fy+c=0 to a Pell equation covered previously. And now some methods: Solving Pell equations x^2 - dy^2 = m for m = +/-1, or +/-4, with d>0, d not a perfect square: Begin by setting a_0 = 0, b_0 = 1, p_-2 = 0, p_-1 = 2, q_-2 = 1, q_-1 = 0, but substitute as follows, depending on the case: If m=+/-1 then p_-1 = 1. If m=+/-4 and d==1 (mod 4) then a_0 = 1, b_0 = 2, p_-2 = -1. If m=+/-4 and d==0 (mod 4) then b_0 = 2. If m=+/-4 and d==2 or 3 (mod 4) then q_-2 = 2. Recursively define c_i = int((a_i + sqrt(d))/b_i) for i>=0, a_i = b_{i-1}*c_{i-1} - a_{i-1} for i>=1, b_i = (d - (a_i)^2)/b_{i-1} for i>=1, p_i = c_i * p_{i-1} + p_{i-2} for i>=0, and q_i = c_i * q_{i-1} + q_{i-2} for i>=0. [or b_i = b_{i-2} + c_{i-1}(a_{i-1} – a_i) for i>=2] (Where the formula for c_i has sqrt(d), one can use int(sqrt(d)).) If d is congruent to 1 modulo 4 and m = +/-4, let k be the minimal index such that c_k = 2c_0 - 1 (but if d=5, take k=1). If either d is not congruent to 1 modulo 4, or m = +/-1, then let k be the minimal index such that c_k = 2c_0. The minimal positive solution to the Pell equation is given by taking x_1 = p_{k-1} and y_1 = q_{k-1}. If k is odd, this is a solution to the "negative" equation, i.e. the equation with m = -1 or m = -4. If k is even, this is a solution to the "positive" equation, m = +1 or m = +4. If k is even, then the negative equation does not have solutions. Note that the sequence c_i is the continued fraction expansion of (a_0 + sqrt(d))/b_0, and at each stage (a_i + sqrt(d))/b_i is the inverse of the fractional part of (a_{i-1} + sqrt(d))/b_{i-1}. The remaining positive solutions are derived as follows: If m=+/-1, then the following three ways to generate all of the positive solutions are equivalent: x_n + y_n*sqrt(d)=(x_1 + y_1*sqrt(d))^n x_n = x_1*x_{n-1} + y_1*y_{n-1}*d; y_n = y_1*x_{n-1} + x_1*y_{n-1} x_n = 2*x_1*x_{n-1} +/- x_{n-2}, y_n = 2*x_1*y_{n-1} +/- y_{n-2} where the + signs are to be used if there is a solution to the -1 equation, and the - signs are to be used if there are only solutions to the +1 equation. If m=+/-4, then the following three ways to generate all of the positive solutions are equivalent: x_n + y_n*sqrt(d)=((x_1 + y_1*sqrt(d))^n)/(2^(n-1)) x_n = (x_1*x_{n-1} + y_1*y_{n-1}*d)/2; y_n = (y_1*x_{n-1} + x_1*y_{n- 1})/2 x_n = x_1*x_{n-1} +/- x_{n-2}, y_n = x_1*y_{n-1} +/- y_{n-2} where the + signs are to be used if there is a solution to the -4 equation, and the - signs are to be used if there are only solutions to the +4 equation. If the -1 or -4 equation has a solution, then the above formulas generate all of the +/-1 or +/-4 solutions, alternately -1 and +1, or - 4 and +4. If you only want the solutions to one of the +1, -1, +4, or - 4 equations, there are analogous formulas. Solving Pell equations x^2 - dy^2 = m for 0 < m^2 < d, m <> +/-1, or +/- 4, with d>0, d not a perfect square: First find the "fundamental" solutions, and then use these to find all solutions. To find fundamental solutions, first use the above method to solve x^2 - dy^2 = 1. Let k be the smallest index so that p_k^2 - dq_k^2 = 1. Compute p_i^2 - dq_i^2 for i from 1 to k. [Use p_i^2 - dq_i^2 = ((-1)^(i+1))*q_{i+1}.] For every positive f so that f^2 divides m, see whether m/f^2 is a value of one of these p_i^2 - dq_i^2. Each time this occurs, x=f*p_i, y=f*q_i is a fundamental solution. If this never occurs, the equation x^2 - dy^2 = m does not have any solutions. All solutions are generated from the fundamental solutions by applying the minimal positive solution to x^2 - dy^2 = 1. That is, if x, y is a solution then x*p_k +y*d*q_k, y*p_k+x*q_k is another solution. Repeatedly iterating gives all solutions. Solving Pell equations x^2 - dy^2 = m for m^2 > d, m <> +/-4, with d>0, d not a perfect square: The method below will find what I call the fundamental solutions. All solutions are generated from the fundamental solutions by repeatedly applying the minimal positive solution r, s to x^2 - dy^2 = 1. That is, if x, y is a solution to x^2 - dy^2 = m then x*r+y*d*s, y*r+x*s is another solution. Iterating this gives all solutions. This is a brute force search method. Let r, s be the minimal positive solution to x^2 - dy^2 = 1. Search on y as follows. If m > 0 then search on y from 0 to sqrt(m*(r - 1)/ (2*d)). If m < 0 then search on y from sqrt((-m)/d) to sqrt((-m)*(r + 1)/(2*d)). Search on y between these limits, and record any solutions found. If no solutions are found at this step, there are no solutions. Note that for m > 0 the above is equivalent to searching on x from sqrt (m) to sqrt((r + 1)*(m/2)). For m < 0, the above is equivalent to searching on x from 0 to sqrt((r - 1)*(-m/2)). If m > 0, then for each such solution x, y, append to the list of solutions r*x-d*s*y, s*x-r*y. If m < 0 then for each such solution x, y, append to the list of solutions d*s*y-r*x, r*y-s*x. In either case, eliminate any duplicates that might occur. The original list plus these added values constitute the fundamental solutions. For larger values of d, it might be better to use Lagrange's method of reduction, described in Mollin or Chrystal, above, or the cyclic method, described in Harold M. Edwards, Fermat's Last Theorem, Springer- Verlag, NY, 1977. In this, Chapters 7 (pp. 245 to 304) and 8 (pp. 305 to 341), especially sections 8.2 (pp. 313-318) and 8.7 (pp. 339-341) apply the cyclic method to solve equations of the form ax^2 + bxy + cy^2 +dx + ey + f = 0. That should get you started. John Robertson Sent via Deja.com http://www.deja.com/ Before you buy. ============================================================================== From: Jpr2718@aol.com (John Robertson) Subject: Re: Solving general pell's equation x^2-dy^2=n Date: 16 May 00 11:41:43 GMT Newsgroups: sci.math.numberthy In a message dated Mon, 15 May 2000 8:48:11 AM Eastern Daylight Time, Matt Herman writes: << Is there a better method of finding the fundamental solution than ... >> Yes, there are better methods to solve (changing notation on you) the Pell equation x^2 - dy^2 = M, for d>0. Note that there might be more than one fundamental solution. By a fundamental solution I mean a positive solution that cannot be generated from a smaller positive solution using the minimal positive solution to x^2 - dy^2 = 1 in the standard way. For M=+/-4, there are methods similar to the continued fraction methods for M=+/-1. If interested, let me know. For a brute-force search, you need only search on the variable y between the proscribed limits. For each y in that range see whether M+dy^2 is a square. Some care is needed to make sure you do get all the fundamental solutions, and not just half of them. Slightly better limits than the ones you gave are given below. If M^2 < D then all of the fundamental solutions can be found among the convergents to sqrt(D) (more details below). If M^2 > D then you can use Lagrange's method of reduction, as described in Mollin's Fundamental Number Theory with Applications, or Chrystal's Textbook of Algebra (sometimes titled just Algebra, recently republished by the AMS in their Chelsea book series). Or, you can use the cyclic method, as described in Edwards' Fermat's Last Theorem. A web page that solves these equations is at http://members.tripod.com/%7Ealpertron/ENGLISH.HTM Another link with some useful information is http://mathworld.wolfram.com/PellEquation.html Descriptions of the methods mentioned above, and more complete references follow. Solving Pell equations x^2 - dy^2 = M for 0 < M^2 < d, M <> +/-1, with d > 0, d not a perfect square: First find the fundamental solutions, and then use these to find all solutions. To find fundamental solutions, first use the continued fraction method to solve x^2 - dy^2 = 1 (let me know if you need more information on this). Let k be the smallest index so that A_k^2 - dB_k^2 = 1, where A_i/B_i is the i-th convergent to sqrt(d). Compute A_i^2 - dB_i^2 for i from 1 to k. [Use A_i^2 - dB_i^2 = ((-1)^(i+1))*Q_{i+1} - see Mollin for notation.] For every positive f so that f^2 divides M, see whether M/f^2 is a value of one or more of these A_i^2 - dB_i^2. Each time this occurs, x=f*A_i, y=f*B_i is a fundamental solution. If this never occurs, the equation x^2 - dy^2 = M does not have any solutions. All solutions are generated from the fundamental solutions by applying the minimal positive solution to x^2 - dy^2 = 1. Brute force search for x^2 - dy^2 = M for M^2 > d > 0, d not a perfect square: Let x1, y1 be the minimal positive solution to x^2 - dy^2 = 1. Search on y as follows. If M > 0 then search on y in the range 0 <= y <= sqrt(M*(x1 - 1)/(2*d)). If M < 0 then search on y in the range sqrt((-M)/d) <= y <= sqrt((-M)*(x1 + 1)/(2*d)). Search on y between these limits, and record any solutions found. If no solutions are found at this step, there are no solutions. Note that for M > 0 the above is equivalent to searching on x from sqrt(M) to sqrt((x1 + 1)*(M/2)). For M < 0, the above is equivalent to searching on x from 0 to sqrt((x1 - 1)*(-M/2)). If M > 0, then for each such solution x, y, append to the list of solutions x1*x-d*y1*y, y1*x-x1*y. If M < 0 then for each such solution x, y, append to the list of solutions d*y1*y-x1*x, x1*y-y1*x. In either case, eliminate any duplicates that might occur. The original list plus these added values constitute the fundamental solutions. Lagrange's Method of Reduction for x^2 - dy^2 = M for M^2 > d > 0, d not a perfect square: For each k 0 <= k <= abs(M)/2 so that (k^2-d)/M = h is an integer, solve x^2 - dy^2 = h. If (kx + dy)/h and (ky + d)/h are integers, they give a solution to the original equation. Also, the same applies if (kx - dy)/h and (ky - d)/h are integers. You might have to iterate if the reduced equation has h^2 > d. Note that there are efficient methods for solving k^2 == d mod M, although for small M, a brute force search is often sufficient. Some solutions found might not be the minimal positive solutions for the family, and so might need to be reduced. Also duplicates are possible. But, at least one member of each family of solutions will be found. Cyclic method for x^2 - dy^2 = M for M^2 > d > 0, d not a perfect square: For our purposes, a limited form of the general cyclic method is all that is needed. For a given a_0, r_0, and d, so that a_0 divides (r_0)^2-d, the limited form of the cyclic method that we need is as follows. We generate a_i, r_j by the following rules: a_i = ((r_{i-1})^2 - d)/a_{i-1} Choose r_i so that a_i divides r_i + r_{i-1} and r_i as large as possible so that r_i^2 < d. If r_i^2 > d for all possible r_i, then make abs(r_i) as small as possible. Continue until there is some pair i1, sr_{i+j+1-k}=tr_{k-1}. [We are marrying the arrays ma and mr up to the i-th index with the reversal of the arrays ta and tr, so that corresponding entries match up. Think of r_i as falling between a_{i-1} and a_i. The result is a cycle that starts 1, 0, and ends x, mm.] Let MM start as the 2x2 identity matrix. For k from 2 to i+j-1 replace MM with MM times [[0, -1], [1, (sr_{k-1} + sr_k)/sa_k]]. A solution to x^2 - dy^2 = m is given by x=f*MM_{1,2} and y=f*MM_{2,2}. The solutions found this way might not be the minimal positive members of their families, but at least one meber of each family will be found. References Richard E. Mollin, Fundamental Number Theory with Applications, CRC Press, Boca Raton, 1998. Mollin gives the most complete set of methods for Pell's equations of any modern reference. Chapter 5, pages 221 to 272, discusses continued fractions generally. In particular, Section 5.3, pages 238 to 250 studies periodic continued fractions, and applies this study to find all solutions to the Pell equation x^2-dy^2=+/-1 for d>0, d not a square. The "PQa" method for computing the continued fraction expansion of a quadratic irrational is discussed in exercise 5.3.6, p. 251. While the method alluded to above for the +/-4 equation is not discussed explicitly, Mollin's Theorem 5.3.4 on page 246 gives the main machinery needed to prove that such a method is possible. Mollin also discusses the continued fraction expansion of (1+sqrt(d))/2 for d==1 (mod 4) in exercise 5.3.14 on page 252. For the general equation x^2-dy^2=M, for any positive nonsquare d and any M, limits to search on x or y are given correctly in Chapter 6 on pages 299 and 300. Theorem 5.2.5 on page 232 gives the main criterion for solving x^2-dy^2=M for d>0, d not a square, M^20, d not a square, M^2>d. Ivan Niven, Herbert S. Zuckerman, and Hugh L. Montgomery, An Introduction to the Theory of Numbers, Fifth Edition, John Wiley & Sons, Inc., New York, 1991. Sections 7.1 to 7.7, pages 325 to 351 cover continued fractions generally. Section 7.8, pages 351 to 356 covers the Pell +/-1 equation, and the case x^2-dy^2=M where M^2 d, but apart form the fact that the heading on page 482 is "Lagrange's Chain of Reductions," the section does not reference Lagrange. Section 19, on page 486, discusses the cases d<0 and d>0, d a square. Section 20, pages 486 to 488, reduces the general equation ax^2+2hxy+by^2+2gx+2fy+c=0 to a Pell equation covered previously. Harold M. Edwards, Fermat's Last Theorem, Springer-Verlag, NY, 1977. Chapters 7 (pp. 245 to 304) and 8 (pp. 305 to 341), especially sections 8.2 (pp. 313-318) and 8.7 (pp. 339-341) apply the cyclic method to solve equations of the form ax^2 + bxy + cy^2 +dx + ey + f = 0. Andre Weil, Number Theory, An approach through history from Hammurapi to Legendre, Birkhauser, Boston, 1983, has some relevant discussion on approximatley pages 19-24, 89-99, and 229-239. Henri Cohen, A Course in Computational Algebraic Number Theory, Springer-Verlag, 1993, Section 1.5, Computing Square Roots Modulo p, pp. 31-36. Gives an efficient way to solve x^2==d mod m based on a method of Shanks. John Robertson ============================================================================== From: krm@maths.uq.edu.au (Keith Matthews) Subject: x^2-dy^2=n Date: 15 May 00 23:25:47 GMT Newsgroups: sci.math.numberthy 9am, Tuesday 16th May 2000. Dear Colleagues, I write in answer to the query of Matt Herman on 15th May 2000: Lagrange gave an algorithm for solving x^2-Dy^2=N in a forgotten paper Oeuvres de Lagrange II, (Ed. Serret 1868) pages 655-726, "Nouvelle me'thode pour re'soudre les probl`emes inde'termine's en nombres entiers", first published in Me'moires de l'Academie royale des Sciences et Belles-Lettres de Berlin, tome XXIV, 1770. Richard Mollin rediscovered the algorithm and one can find it his book "Fundamental Number Theory with Applications", pages 338-339. Unlike Lagrange's method, his uses ideals and semi-simple contnued fractions. I also rediscovered Lagrange's algorithm while studying a paper of W. Patz \"{U}ber die Gleichung $X^2-DY^2=\pm c\cdot (2^{31}-1)$, Bayer. Akad. Wiss. Math--Natur. Kl. S.--B (1948) 21--30. My version uses only simple continued fractions and can be considered as a generalisation of both the well--known method of solving Pell's equation $x^2-Dy^2=\pm 1$ using the simple continued fraction for $\sqrt{D}$ and that of solving $x^2-Dy^2=\pm 4$ using the continued fraction of $(\sqrt{D}-1)/2$. The preprint can be viewed at http://www.maths.uq.edu.au/~krm/patz.pdf (In my approach, there is an exceptional case D=2 or 3 and N<0 which requires separate treatment and it would be nice if this could be eliminated.) I have just finished coding the algorithm and it is now available in my CALC program at ftp://www.maths.uq.edu.au/pub/krm/calc/ One employs a simple necessary and sufficient condition for solubility and in the case of solubility, the fundamental solutions x+y\sqrt(D) with x positive are displayed. The cases |N| and -|N| are dealt with simultaneously. We only assume D>0 is not a perfect square. gcd(D,N) may be > 1. Lagrange also gave a recursive algorithm in an earlier paper Oeuvres II, pages 377--535, "Sur la solution des probl`emes inde'termine's du second degre' was published in the same journal, tome XXIII, 1769 and this algorithm is the one found in most textbooks. eg. Loo Keng Hua's book "Introduction to number theory", pages 283-285. Regards Keith Matthews Department of Mathematics University of Queensland Brisbane Email: krm@maths.uq.edu.au WWW: http://www.maths.uq.edu.au/~krm/ ============================================================================== From: john.cremona@maths.nottingham.ac.uk (John Cremona) Subject: Re: x^2-dy^2=n Date: 16 May 00 11:41:44 GMT Newsgroups: sci.math.numberthy For my 2 cents-worth on this and related equations, see the preprint at http://www.maths.nottingham.ac.uk/personal/jec/papers/conics2.ps. I recommend the LLL-based method as being very simple to implement (assuming that someone else has done LLL for you, as in Pari). John Cremona -- Prof. J. E. Cremona | University of Nottingham | Tel.: +44-115-9514920 School of Mathematical Sciences | Fax: +44-115-9514951 University Park | Email: John.Cremona@nottingham.ac.uk Nottingham NG7 2RD, UK |