TOPICS
Search

Pell Equation


A special case of the quadratic Diophantine equation having the form

 x^2-Dy^2=1,
(1)

where D>0 is a nonsquare natural number (Dickson 2005). The equation

 x^2-Dy^2=+/-4
(2)

arising in the computation of fundamental units is sometimes also called the Pell equation (Dörrie 1965, Itô 1987), and Dörrie calls the positive form of (2) the Fermat difference equation. While Fermat deserves credit for being the first to extensively study the equation, the erroneous attribution to Pell was perpetrated by none other than Euler himself (Nagell 1951, p. 197; Burton 1989; Dickson 2005, p. 341). The Pell equation was also solved by the Indian mathematician Bhaskara. Pell equations are extremely important in number theory, and arise in the investigation of numbers which are figurate in more than one way, for example, simultaneously square and triangular.

The equation has an obvious generalization to the Pell-like equation

 ax^2+/-by^2=c,
(3)

as well as the general second-order bivariate Diophantine equation

 ax^2+bxy+cy^2+dx+ey+f=0.
(4)

However, several different techniques are required to solve this equation for arbitrary values of a, b, and c. The Wolfram Language command Reduce[f[x, y] && Element[x|y, Integers]] finds solutions to the general equation (4) when they exist.

Pell equations of the form (1), as well as certain cases of the analogous equation with a minus sign on the right,

 x^2-Dy^2=-1,
(5)

can be solved by finding the continued fraction [a_0,a_1,...] of sqrt(D). Note that although the equation (5) is solvable for only certain values of D, the continued fraction technique provides solutions when they exist, and always in the case of (1), for which a solution always exists. A necessary condition that (5) be solvable is that all odd prime factors of D be of the form 4n+1, and that D cannot be doubly even (i.e., divisible by 4). However, these conditions are not sufficient for a solution to exist, as demonstrated by the equation x^2-34y^2=-1, which has no solutions in integers (Nagell 1951, pp. 201 and 204).

In all subsequent discussion, ignore the trivial solution x=1, y=0. Let p_n/q_n denote the nth convergent [a_0,a_1,...,a_n], then we will have solved (1) or (5) if we can find a convergent which obeys the identity

 p_n^2-Dq_n^2=(-1)^(n+1).
(6)

Amazingly, this turns out to always be possible as a result of the fact that the continued fraction of a quadratic surd always becomes periodic at some term a_(r+1), where a_(r+1)=2a_0, i.e.,

 sqrt(D)=[a_0,a_1,...,a_r,2a_0^_].
(7)

To compute the continued fraction convergents to sqrt(D), use the usual recurrence relations

a_0=|_sqrt(D)_|
(8)
p_0=a_0
(9)
p_1=a_0a_1+1
(10)
p_n=a_np_(n-1)+p_(n-2)
(11)
q_0=1
(12)
q_1=a_1
(13)
q_n=a_nq_(n-1)+q_(n-2),
(14)

where |_x_| is the floor function. For reasons to be explained shortly, also compute the two additional quantities P_n and Q_n defined by

P_0=0
(15)
P_1=a_0
(16)
P_n=a_(n-1)Q_(n-1)-P_(n-1)
(17)
Q_0=1
(18)
Q_1=D-a_0^2
(19)
Q_n=(D-P_n^2)/(Q_(n-1))
(20)
a_n=|_(a_0+P_n)/(Q_n)_|.
(21)

Now, two important identities satisfied by continued fraction convergents are

 p_nq_(n-1)-p_(n-1)q_n=(-1)^(n+1)
(22)
 p_n^2-Dq_n^2=(-1)^(n+1)Q_(n+1)
(23)

(Beiler 1966, p. 262), so both linear

 ax-by=+/-1
(24)

and quadratic

 x^2-Dy^2=+/-c
(25)

equations are solved simply by finding an appropriate continued fraction.

Let a_(r+1)=2a_0 be the term at which the continued fraction becomes periodic (which will always happen for a quadratic surd). For the Pell equation

 x^2-Dy^2=1
(26)

with r odd, (-1)^(r+1) is positive and the solution in terms of smallest integers is x=p_r and y=q_r, where p_r/q_r is the rth convergent. If r is even, then (-1)^(r+1) is negative, but

 p_(2r+1)^2-Dq_(2r+1)^2=1,
(27)

so the solution in smallest integers is x=p_(2r+1), y=q_(2r+1). Summarizing,

 (x,y)={(p_r,q_r)   for r odd; (p_(2r+1),q_(2r+1))   for r even.
(28)

The equation

 x^2-Dy^2=-1
(29)

can be solved analogously to the equation with +1 on the right side iff r is even, but has no solution if r is odd,

 (x,y)={(p_r,q_r)   for r even; no solution   for r odd.
(30)

Given one solution (x,y)=(p,q) (which can be found as above), a whole family of solutions can be found by taking each side to the nth power,

 x^2-Dy^2=(p^2-Dq^2)^n=1.
(31)

Factoring gives

 (x+sqrt(D)y)(x-sqrt(D)y)=(p+sqrt(D)q)^n(p-sqrt(D)q)^n
(32)

and

x+sqrt(D)y=(p+sqrt(D)q)^n
(33)
x-sqrt(D)y=(p-sqrt(D)q)^n,
(34)

which gives the family of solutions

x=((p+qsqrt(D))^n+(p-qsqrt(D))^n)/2
(35)
y=((p+qsqrt(D))^n-(p-qsqrt(D))^n)/(2sqrt(D)).
(36)

These solutions also hold for

 x^2-Dy^2=-1,
(37)

except that n can take on only odd values.

The following table gives the smallest integer solutions (x,y) to the Pell equation with constant D<=102 (Beiler 1966, p. 254). Square D=d^2 are not included, since they would result in an equation of the form

 x^2-d^2y^2=x^2-(dy)^2=x^2-y^('2)=1,
(38)

which has no solutions (since the difference of two squares cannot be 1).

DxyDxy
2325448566
321558912
59456152
6525715120
78358196032574
8315953069
1019660314
11103611766319049226153980
127262638
136491806381
141546512916
154166658
1733867488425967
1817468334
1917039697775936
20927025130
215512713480413
221974272172
23245732281249267000
2451743699430
26511075263
2726576577996630
28127247735140
299801182078536
3011279809
3115202738091
321738216318
3323483829
3435684556
35618528576930996
37731286104051122
3837687283
392548819721
401938950000153000
41204932090192
42132911574165
433482531921151120
441993093121511260
4516124942143295221064
4624335358895394
4748796495
487197628096336377352
509914989910
5150799101
526499010120120
5366249910010210110

The first few minimal values of x and y for nonsquare D are 3, 2, 9, 5, 8, 3, 19, 10, 7, 649, ... (OEIS A033313) and 2, 1, 4, 2, 3, 1, 6, 3, 2, 180, ... (OEIS A033317), respectively. The values of D having x=2, 3, ... are 3, 2, 15, 6, 35, 12, 7, 5, 11, 30, ... (OEIS A033314) and the values of D having y=1, 2, ... are 3, 2, 7, 5, 23, 10, 47, 17, 79, 26, ... (OEIS A033318). Values of the incrementally largest minimal x are 3, 9, 19, 649, 9801, 24335, 66249, ... (OEIS A033315) which occur at D=2, 5, 10, 13, 29, 46, 53, 61, 109, 181, ... (OEIS A033316). Values of the incrementally largest minimal y are 2, 4, 6, 180, 1820, 3588, 9100, 226153980, ... (OEIS A033319), which occur at D=2, 5, 10, 13, 29, 46, 53, 61, ... (OEIS A033316).

The more complicated Pell-like equation

 x^2-Dy^2=c
(39)

with |c|<sqrt(D) has solution iff c is one of the values (-1)^kQ_k for k=1, 2, ..., r computed in the process of finding the convergents to sqrt(D) (where, as above, a_(r+1)=2a_0 is the term at which the continued fraction becomes periodic). If |c|>sqrt(D), the procedure is significantly more complicated (Beiler 1966, p. 265; Dickson 2005, pp. 387-388) and is discussed by Gérardin (1910) and Chrystal (1961).

Regardless of how it is found, if a single solution x=p, y=q to (39) is known, other solutions can be found. Let p and q be solutions to (39), and r and s solutions to the "unit" form

 x^2-Dy^2=1.
(40)

Then the identity

(p^2-Dq^2)(r^2-Ds^2)=(pr+/-Dqs)^2-D(ps+/-qr)^2
(41)
=c
(42)

allows larger solutions (x,y)=(pr+/-Dqs,ps+/-qr) to the c equation to be found by using incrementally larger values of the (r,s), which can be easily computed using the standard technique for the Pell equation. Such a family of solutions does not necessarily generate all solutions, however. For example, the equation

 x^2-10y^2=9
(43)

has three distinct sets of fundamental solutions, (x,y)=(7,2), (13, 4), and (57, 18). Using (42), these generate the solutions shown in the following table, from which the set of all solutions (7, 2), (13, 4), (57, 18), (253, 80), (487, 154), (2163, 684), (9607, 3038), ... can be generated.

fundamentalgenerated solutions
(7, 2)(253, 80), (9607, 3038), (364813, 115364), (13853287, 4380794), ...
(13, 4)(487, 154), (18493, 5848), (702247, 222070), (26666893, 8432812), ...
(57, 18)(2163, 684), (82137, 25974), (3119043, 986328), (118441497, 37454490), ...

The case

 ax^2-by^2=c
(44)

can be reduced to the one above by multiplying through by a,

 (ax)^2-(ab)y^2=ac,
(45)

finding solutions in (x^'=ax,y), and then selecting those for which x^'/a is an integer.

According to Dickson (2005, pp. 408 and 411), the equation

 ax^2+by^2=c
(46)

with a,b,c>0, which has either no solutions or a finite number of solutions, was solved by Gauss in 1863 using the method of exclusions and considered by Euler (1773) and Nasimoff (1885), although Euler's methods were incomplete (Smith 1965; Dickson 2005, p. 378). According to Itô (1987), this equation can be solved completely using solutions to Pell's equation. Nasimoff (1885) applied Jacobi elliptic functions to express the number of solutions of this equation for a,c odd (Dickson 2005, p. 411). Additional discussion including the connection with elliptic functions is given in Dickson (2005, pp. 387-391).

The special case of a=1 and c prime was solved by Cornacchia (Cornacchia 1908, Cox 1989, Wagon 1990). A deterministic algorithm for finding all primitive solutions to (46) for a,b,c>0 fixed relatively prime integers, c>=a+b+1, and (c,ab)=1 was given by Hardy et al. (1990). This algorithm generalizes those of Hermite (1848), Serret (1848), Brillhart (1972), Cornacchia (1908), and Wilker (1980). It requires factorization of c, and has worst case running time of O(c^(1/4)(lnc)^3(lnlnc)(lnlnlnc)), independent of a and b. An algorithmic method for finding all solutions is implemented in the Wolfram Language as Reduce[x^2 + d y^2 == n, {x, y}, Integers].


See also

Binary Quadratic Form, Diophantine Equation, Diophantine Equation--2nd Powers, Fundamental Unit, Hilbert Symbol, Lagrange Number, Monomorph, Polymorph

Explore with Wolfram|Alpha

References

Beiler, A. H. "The Pellian." Ch. 22 in Recreations in the Theory of Numbers: The Queen of Mathematics Entertains. New York: Dover, pp. 248-268, 1966.Brillhart, J. "Note on Representing a Prime as a Sum of Two Squares." Math. Comput. 26, 1011-1013, 1972.Burton, D. M. Elementary Number Theory, 4th ed. Boston, MA: Allyn and Bacon, pp. 379-382 and 392, 1989.Chrystal, G. Textbook of Algebra, 2nd ed., Vol. 2. New York: Chelsea, pp. 478-486, 1961.Cipolla, M. "Un metodo per la risoluzione della congruenza di secondo grado." Rend. Accad. Sci. Fis. Mat. Napoli 9, 154-163, 1903.Cohn, H. "Pell's Equation." §6.9 in Advanced Number Theory. New York: Dover, pp. 110-111, 1980.Cornacchia, G. "Su di un metodo per la risoluzione in numeri unteri dell' equazione sum_(h=0)^(n)c_hx^(n-h)y^h=P." Giornale di Matematiche di Battaglini 46, 33-90, 1908.Cox, D. A. Primes of the form x^2+ny^2. New York: Wiley, 1989.Degan, C. F. Canon Pellianus. Copenhagen, Denmark, 1817.Dickson, L. E. "Pell Equation: ax^2+bx+c Made Square." Ch. 12 in History of the Theory of Numbers, Vol. 2: Diophantine Analysis. New York: Dover, pp. 341-400, 2005.Dörrie, H. 100 Great Problems of Elementary Mathematics: Their History and Solutions. New York: Dover, 1965.Euler, L. Novi Comm. Acad. Petrop. 18, 218, 1773. Reprinted in Opera Omnia, Series Prima, Vol. 3, p. 310.Euler, L. Comm. Arith. 570. Reprinted in Opera Omnia, Series Prima, Vol. 3, p. 310.Gérardin, A. "Formules de récurrence." Sphinx-Oedipe 5, 17-29, 1910.Hardy, K.; Muskat, J. B.; and Williams, K. S. "A Deterministic Algorithm for Solving n=fu^2+gv^2 in Coprime Integers u and v." Math. Comput. 55, 327-343, 1990.Havil, J. Gamma: Exploring Euler's Constant. Princeton, NJ: Princeton University Press, pp. 91-92, 2003.Hermite, C. "Note au sujet de l'article précédent." J. Math. Pures Appl. 13, 15, 1848.Itô, K. (Ed.). Encyclopedic Dictionary of Mathematics, 2nd ed., Vol. 1. Cambridge, MA: MIT Press, p. 450, 1987.Lagarias, J. C. "On the Computational Complexity of Determining the Solvability or Unsolvability of the Equation X^2-DY^2=-1." Trans. Amer. Math. Soc. 260, 485-508, 1980.Lenstra, H. W. Jr. "Solving the Pell Equation." Not. Amer. Math. Soc. 49, 182-192, 2002.Nagell, T. "The Diophantine Equation x^2-Dy^2=1," "The Diophantine Equation x^2-Dy^2=-1," and "The Diophantine Equation u^2-Dv^2=C." §56-58 in Introduction to Number Theory. New York: Wiley, pp. 195-212, 1951.Nasimoff, P. S. Ch. 1 in Application of Elliptic Functions to the Theory of Numbers. Moscow, 1885. French summary in Ann. sci. de l'École normale supér. 5, 23-31, 1888.Robertson, J. "Home Page for John Robertson." http://www.jpr2718.org/.Serret, J. A. "Sur un théorème rélatif aux nombres enti'eres." J. Math. Pures Appl. 13, 12-14, 1848.Sloane, N. J. A. Sequences A033313, A033314, A033315, A033316, A033317, A033318, and A033319 in "The On-Line Encyclopedia of Integer Sequences."Smith, H. J. S. Collected Mathematical Papers I. New York: Chelsea, pp. 195-202, 1965.Smarandache, F. "Un metodo de resolucion de la ecuacion diofantica." Gaz. Math. 1, 151-157, 1988.Smarandache, F. "Method to Solve the Diophantine Equation ax^2-by^2+c=0." In Collected Papers, Vol. 1. Lupton, AZ: Erhus University Press, 1996.Stillwell, J. C. Mathematics and Its History. New York: Springer-Verlag, 1989.Wagon, S. "The Euclidean Algorithm Strikes Again." Amer. Math. Monthly 97, 124-125, 1990.Whitford, E. E. Pell Equation. New York: Columbia University Press, 1912.Wilker, P. "An efficient Algorithmic Solution of the Diophantine Equation u^2+5v^2=m." Math. Comput. 35, 1347-1352, 1980.

Referenced on Wolfram|Alpha

Pell Equation

Cite this as:

Weisstein, Eric W. "Pell Equation." From MathWorld--A Wolfram Web Resource. https://mathworld.wolfram.com/PellEquation.html

Subject classifications