Project Euler 66: Investigate the Diophantine equation x² − Dy² = 1.

Problem Description

Consider quadratic Diophantine equations of the form:

x² – Dy² = 1

For example, when D=13, the minimal solution in x is 649² – 13×180² = 1.

It can be assumed that there are no solutions in positive integers when D is square.

By finding minimal solutions in x for D = {2, 3, 5, 6, 7}, we obtain the following:

3² – 2×2² = 1
2² – 3×1² = 1
9² – 5×4² = 1
5² – 6×2² = 1
8² – 7×3² = 1

Hence, by considering minimal solutions in x for D ≤ 7, the largest x is obtained when D=5.

Find the value of D ≤ 1000 in minimal solutions of x for which the largest value of x is obtained.

Analysis

Using the chakravala method for solving minimal solutions to Pell’s Equation. A detailed explanation is posted in the comments section.

Project Euler 66 Solution

Afterthoughts

• D=92821 for D ≤ 100000 in 15 seconds.
  

Solution of “Pell’s Equation” by Chakravala Method by Gopal Menon

It was discovered by Brahmagupta that a triple (a, b, k), representing the equation a² – Nb² = k, can be composed with the trivial triple (m, 1, m² – N) to get a new triple (am + Nb, a + bm, k(m² – N)). Scaling down by k, using Bhaskara’s lemma, we get,

a² – Nb² = k ==> {(am + Nb)/k}² – N{(a + bm)/k}² = (m² – N)/k

The following is the algorithm for solving “Pell’s equation”, x² – Ny² = 1, as given by the Indian mathematicians who were the first to solve this equation:

1. For the given value of N (which is not a perfect square) replace x by a, y by b and 1 by k and solve the equation a² – Nb² = k. Let the initial value of b be 1 and select the initial value of a such that the square of the number is as close to N as possible so as to get the least absolute value of k. Note that k may be either positive or negative but not zero.
2. Now select an integer m such that k divides a + bm and the absolute value of m² – N is minimum.
3. Now replace a by the absolute value of (am + Nb)/k and b by the absolute value of (a + bm)/k.
4. Now replace k by a² – Nb². Note that k may be either positive or negative but not zero.
5. If the new value of k is 1 we have got the solution. The value of x is a and that of y is b.
6. If the value of k is not equal to 1 we go for the first iteration and repeat steps 1 to 5. We keep on doing this until the value of k becomes 1 and we thus get the values of x and y which satisfies “Pell’s equation” for the given value of N.
7. Using Brahmagupta’s lemma, if the value of k assumes a value of -1 or ± 2 for any value of a and b or if the value of k assumes a value of ± 4 and either a or b is an even number, we can stop further iterations and compose the triple with itself to get the final solution.
   a² – Nb² = k ==> {(a² + Nb²)/k}² – N{(ab)/k}² = 1

An example with N = 53 is given below.

1. Since N = 53, we initially select a = 7 and b = 1 to get k = -4.
2. We have to now select an integer m such that 4 divides 7 + m and the absolute value of m² – 53 is minimum. Therefore m should be of the type 4t + 1 for integer values of t. To minimise the absolute value of m² – 53, m should be 5.
3. Hence, in the first iteration we get the new values of a and b as the absolute value of (am + Nb)/k and (a + bm)/k respectively. Thus we get a = (7×5 + 53×1)/4 = 22, b = (7 + 1×5)/4 = 3 and k = (5² – 53)/(-4) = 7.
4. In the second iteration we have to select an integer m such that 7 divides 22 + 3m and the absolute value of m² – 53 is minimum. Therefore m should be of the type 7t + 2 for integer values of t. To minimise the absolute value of m² – 53, m should be 9.
5. The new values of a and b are now the absolute value of (am + Nb)/k and (a + bm)/k respectively. Thus we get a = (22×9 + 53×3)/7 = 51, b = (22 + 3×9)/7 = 7 and k = (9² – 53)/7 = 4.
6. In the third iteration we select m = 7 to get the revised values of a = 182, b = 25 and k = -1.
7. Since the value of k is now -1, we can use Bhaskara’s lemma an compose the triple with itself. Since a = 182, b = 25 and k = -1, we get x = (182² + 53×25²)/1 = 66249, y = 2x182x25/1 = 9100 and k = 1.
8. Thus, the solution to the equation is 66249² – 53 x 9100² = 1.
9. From this we can also get the value of √53 as 7.28010989010989. This value differs from the value obtained from the square root function of the spreadsheet by just 0.000000011 %.
10. This method always terminates with a solution (when the value of k becomes 1) for all values of N (other than perfect squares).