We'll study here only the case D positive not a square and first of all
the case k = 1 :
x² - Dy² = 1 (proper Pell equation).
We then have the following results :
xn+1 = x0xn + Dy0yn, yn+1 = y0xn + x0yn
We can then deduce recurence relation for the sequences xn and yn :
xn+2 = 2x0xn+1-xn yn+2 = 2x0yn+1-yn
We can also write : xn + yn√D = (x0 + y0√ D)n+1 or the conjugate form, replacing + by -
This gives closed forms for xn and yn as a function of n.
These formulaes seems to be of little use because of irrational terms in √D, although giving
integral xn and yn !
(They are easy to write from the previous relation combined with its conjugate)
|Let U0 = 0, V0 = 1, a0 = E[√D], P0 = a0, Q0 = 1
and by convention P-1 = 1 and Q-1 = 0
Noting E[x] = integer part of x. Then repeat :
Un+1 = anVn - Un,
Vn+1 = (D - Un+1²)/Vn,
an+1 = E[ (a0 + Un+1)/Vn+1 ],
until ak > a0 (or Vk = 1).
If k even it is the required solution.
This algorithm is a specific case of PQa algorithm giving the expansion as a continued fraction of a quadratic number (u + √D)/v with here u = 0 and v = 1.
Apply this to x² - 13y² = 1 :
|U0 = 0||V0 = 1||a0 = 3||P0 = 3||Q0 = 1|
|U1 = 3||V1 = (13 - 9)/1 = 4||a1 = E[(3 + 3)/4] = 1||P1 = 4||Q1 = 1|
|U2 = 4 - 3 = 1||V2 = (13 - 1)/4 = 3||a2 = E[(3 + 1)/3] = 1||P2 = 7||Q2 = 2|
|U3 = 3 - 1 = 2||V3 = (13 - 4)/3 = 3||a3 = E[(3 + 2)/3] = 1||P3 = 11||Q3 = 3|
|U4 = 3 - 2 = 1||V4 = (13 - 1)/1 = 4||a4 = E[(3 + 1)/4] = 1||P4 = 18||Q4 = 5|
|U5 = 4 - 1 = 3||V5 = (13 - 9)/4 = 1||a5 = E[(3 + 3)/1] = 6 > 3 => stop|
5 being odd, (18, 5) is solution of x² - 13y² = -1.
The fundamental solution of x² - 13y² = 1 is (18² + 13×5², 2×18×5) = (649, 180) which was not obvious without this algorithm.
Let's go on with the general case k ≠ 1, and first with
or the conjugate form by replacing + by -
Note : solutions of equation x² - Dy² = 1 are obtained from the fundamental solution
(x0, y0) of x² - Dy² = -1 by : xn + yn√D = (x0 + y0√D)2n
The fundamental solution of x² - Dy²=1 is then (x0² + Dy0², 2x0y0).
We also have the recurrence relations :
xn+1 = (x0² + Dy0²)xn + 2Dx0y0yn, yn+1 = 2x0y0xn + (x0² + Dy0²)yn
As well as relations in xn+2 = axn+1 + bxn that are left as an exercice.
As for the searching of the fundamental solution, just use the previous algorithm, and note that when k is even,
then x² - Dy² = -1 has no solution.
There is no easy criterion to know it in advance. Note the case when D is a prime p :
Equation x² - py² = -1 has solutions if and anly if p = 4k + 1.
Equation x² - 2y² = -1
x = 1, y = 1 is an obvious solution, and is fundamental as any smaller solutions would have x and y < 1.
The fundamental solution of x² - 2y² = +1 is then x = 1² + 2*1² = 3 and y = 1² + 1² = 2 .
All other solutions are then xn+1 = 3xn + 4yn, yn+1 = 2xn + 3yn, that is :
(1,1) (7,5) (41,29) (239,169) (1393,985) ...
xn + yn√D = (x0 + y0√D)n(ξ + η√D) and also the corresponding recurrence relations.
The problem is now to find these independant primitive solutions. A first method is to scan systematically all values of x or y, the theory giving just an upper bound. (remind the equation may have no solutions). This is the most efficient method if the upper bound is not too big.
But at first we must explain what are independant solutions,
or more precisely when they are not, that is when they are said to be "equivallent".
Two solutions (x,y) and (x',y') are equivalent if there is a solution (X,Y) of x² - Dy² = 1 with :
x + y√D = (X + Y√D)(x' + y'√D)
This criterion being not very usefull, add the following :
|Two solutions (x,y) and (x',y') of x² - Dy² = K are equivallent
if and only if :
xx' - Dyy' and x'y - xy' are multiple of K.
Remains to find bounds for the primitive solutions :
|If K > 0, the primitive solutions are bounded by
0 < x ≤ √( (x0 + 1)K/2 ) or 0 ≤ y ≤ y0√ ( K/2(x0 + 1) )
If K < 0, the primitive solutions are bounded by
If there is no solutions in the bounds, there is no solutions at all.
We would like to find a simple criterion for x² - Dy² = K to have solutions by writing
x² ≡ K mod (D), that is K must be a quadratic residue of D.
The only thing we can say then is that when K is not a quadratic residue of D, then the Pell equation has no solutions. For instance equation x² - 97y² = 51 has no solutions Legendre symbol (51/97) = -1 hence 51 is not a quadratic residue of 97.
But this necessary condition is not sufficient...
Hence there are 4 primitive solutions, but not all are independant
(±3, 0) are equivalent, more generally (x,y) ~ (-x,-y)
(4,±1) are independant : 16 + 7 = 26 coprime with 9.
(3,0) and (4,1) are independant
(3,0) and (4,-1) are independant
Finally three independant primitive solutions : (3, 0), (4, 1) and (4,-1)
If we just care of x,y > 0, (4,-1) ~ (11,4), is then smallest >0 solution in this family
All solutions are built from the three families :
xn + yn√7 = (8 + 3√7)n(3 + 0√7)
xn + yn√7 = (8 + 3√7)n(4 + √7)
xn + yn√7 = (8 + 3√7)n(4 - √7)
Or the recurrence relations :
xn+1 = 8xn + 21yn, yn+1 = 3xn + 8yn, from (x, y) = (3,0), (4,1) and (4,-1)
that is the three families :
(3,0) (24,9) (381,144) ...
(11,4) (172,65) (2741,1036) ...
(4,1) (53,20) (844,319) ...
Scanning for solutions becomes tedious if the bound is too big. In that case we have two efficient methods.
Pi² - DQi² = (-1)i+1Vi+1V0
This allows to find, with the same algorithm, the primitive solutions of x² - Dy² = K
First of all note that if (x,y) is solution of x² - Dy² = K, (mx,my) is a solution of x² - Dy² = m²K.
If we get a (-1)i+1Vi+1 = K/m², then x = mPi, y = mQi is solution of x² - Dy² = K
If K² < D, this gives all primitive solutions.
Example x² - 19y² = -3
|4||3||2||9||2||-3||4||1||is solution of x² - 19y² = -3|
|2||3||2||170||39||-3||61||14||is solution of x² - 19y² = -3|
|4||1||8||1||170||39||is solution of x² - 19y² = 1|
1) Equation x² - 19y² = K hence has no solutions for K = -1, +2, +3, and -4.
2) x² - 19y² = 5 has primitive solutions (9,2) and (48,11) but we don't know if there are others because 5² > 19.
|List the f with f² dividing K, m=K/f²
For each such f, find all z : -|m|/2 < z≤|m|/2 and z² ≡ D modulo |m|
For each z, apply the PQA algorithm with U0 = z, V0 = |m| up to end of period or a Vi = ±1
- If Vi = 1, then (Pi-1,Qi-1) is solution of x² - Dy² = m, add this solution to the list.
- If Vi = -1, then (Pi-1,Qi-1) is solution of x² - Dy² = -m. If x² - Dy² = -1 has a solution (r,s) (primitive)
then (Pi-1,Qi-1)o(r,s) is solution of x² - Dy²=m. Add (Pi-1,Qi-1)o(r,s) = (Pi-1r + Qi-1sD, Pi-1s + Qi-1r)
All these solutions being independant, hence result into 6 solution families by combining them with the solution (649,180)n of x² - 13y² = 1