Pell's Enormous Equation!

The continued fraction expansion of $\sqrt{6829}$ has period $161$ , which is odd, so the fundamental solution of the equation $x^2 - 6829y^2 = -1$ comes from the $160$ th convergent of the continued fraction expansion of $\sqrt{6829}$ . Thus, if this convergent is $\frac{p}{q}$ in lowest terms, the fundamental solution is $x=p$ , $y=q$ .

The fundamental solution of $x^2 - 6829y^2= 1$ is $x=P$ , $y=Q$ , where $P+Q\sqrt{6829} = (p + q\sqrt{6829})^2$ . This solution can also be derived from the $321$ st convergent $\frac{P}{Q}$ of $\sqrt{6829}$ . Since we only want to find $P+Q$ modulo $10000$ , we need to solve the recurrence relation $X_n \; \equiv \; a_n X_{n-1} + X_{n-2} \quad \pmod{10000} \qquad \qquad X_0 = 83\,,\, X_{-1} = 1$ (where $a_0,a_1,a_2,a_3,...$ are the coefficients of the continued fraction expansion of $\sqrt{6829}$ , so that $a_0,a_1,a_2,a_3,a_4,a_5,... = 82,1,1,1,3,5,...$ ) to find $X_{321} = 1789$ . These calculations can be performed readily in Excel, and do not require the ability to handle $162$ digit numbers like $P$ .

$\begin{aligned} x^2-6829y^2&=1 \\ \frac{x^2}{y^2} - 6829 &= \frac{1}{y^2} \\ \frac{x^2}{y^2} &= 6829 + \frac{1}{y^2} \end{aligned}$

Since $x$ and $y$ both differ by a constant factor, as one values grows, then so does the other for large values of $x$ , $y$ . With simple bruteforce, it can be checked that no solutions for $x$ and $y$ exist under 1000. So the solutions are expected to be large.

But as $x$ and $y$ grow, the value of $\frac{1}{y^2}$ gets insignificant and hence we can write it as

$\begin{aligned} \frac{x^2}{y^2} &\approx 6829 \\ \frac{x}{y} &\approx \sqrt{6829} \end{aligned}$

So $\frac{x}{y}$ are very close to $\sqrt{6829}$ . This fact can be exploited and successive rational approximations of $\sqrt{6829}$ can be found using continued fractions and the solution can be checked.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22

def Fundamental_Solution(a, P, k):          # x = floor(sqrt(n)), P = continued fraction period, x^2 - ky^2 = 1
  X = [1, a]                                # X_-1 = 1, X_0 = a, X_n = P_n+1 * X_n-1 + X_n-2
  Y = [0, 1]                                # Y_-1 = 0, Y_0 = 1, X_n = P_n+1 * Y_n-1 + Y_n-2

  i = 0                                     # ith partial continued fraction sum
  period = len(P)

  while (X[1] ** 2 - k * (Y[1] ** 2)) != 1:
    x_n = P[i] * X[1] + X[0]
    y_n = P[i] * Y[1] + Y[0]

    X[0] = X[1]
    X[1] = x_n

    Y[0] = Y[1]
    Y[1] = y_n

    i = (i + 1) % period

  return (X[1], Y[1])                       # (x,y) is the fundamental solution

print Fundamental_Solution(82, [long array containing the repeating continued fraction digits], 6829)

For record, the solutions are

1
2
3

329193470472193712949550776348318404515010693082347229388317940707942971956588011161803916649277857045856156953356297825072563425494288962207368895353723748671049

3983571862015289177268105541724820937760171200667222370948986501166335722772085747080516646954487089316624085365269517982453973452961546542626698572599406860740

Here's how I computed the continued fraction -

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22

D = 6829
R = int(D ** 0.5)

array = '[' + str(R) + '; '

a = R
P = 0
Q = 1

P = a * Q - P
Q = (D - P * P) / Q
a = (R + P) / Q

array += str(a) + ', '

while Q != 1:
    P = a * Q - P
    Q = (D - P * P) / Q
    a = (R + P) / Q
    array += str(a) + ', '

print array[:-2] + ']'

Output

1

[82; 1, 1, 1, 3, 5, 1, 1, 1, 2, 2, 1, 4, 54, 1, 7, 3, 1, 1, 4, 1, 3, 4, 1, 2, 1, 17, 1, 1, 1, 2, 10, 1, 1, 1, 3, 1, 14, 4, 5, 1, 7, 32, 1, 12, 1, 4, 12, 1, 1, 23, 10, 1, 40, 2, 2, 3, 1, 5, 7, 1, 2, 3, 3, 13, 2, 7, 1, 3, 1, 1, 2, 2, 4, 3, 3, 2, 4, 6, 2, 1, 1, 2, 6, 4, 2, 3, 3, 4, 2, 2, 1, 1, 3, 1, 7, 2, 13, 3, 3, 2, 1, 7, 5, 1, 3, 2, 2, 40, 1, 10, 23, 1, 1, 12, 4, 1, 12, 1, 32, 7, 1, 5, 4, 14, 1, 3, 1, 1, 1, 10, 2, 1, 1, 1, 17, 1, 2, 1, 4, 3, 1, 4, 1, 1, 3, 7, 1, 54, 4, 1, 2, 2, 1, 1, 1, 5, 3, 1, 1, 1, 164]