## Pell Equations

Suppose we wish to solve

$x^2 - D y^2 = 1$

where $D$ is not a square. First we compute the continued fraction expansion of $\sqrt{D}$.

$\array { P_0 &=& 0 \\ Q_0 &=& 1 \\ a_0 &=& {\lfloor\sqrt{D}\rfloor} \\ P_1 &=& a_0 \\ Q_1 &=& D - a_0^2 \\ a_n &=& {\lfloor \frac{a_0 + P_n}{Q_{n}} \rfloor} \\ P_n &=& a_{n-1} Q_{n-1} - P_{n-1} \\ Q_n &=& \frac{D - P_n^2}{Q_{n-1}} }$

Then the continued fraction expansion of $\sqrt{D}$ is

$\sqrt{D} = [a_0; a_1, ..., a_{k+1}, a_1, ..., a_{k+1}, ... ]$

It turns out $a_{k+1} = 2 a_0$.

The convergents are computed by

$\array { p_0 &=& a_0 \\ p_1 &=& a_0 a_1 + 1 \\ p_n &=& a_n p_{n-1} + p_{n-2} \\ q_0 &=& 1 \\ q_1 &=& a_1 \\ q_n &=& a_n q_{n-1} + q_{n-2} \\ }$

The convergents satisfy

$p_n^2 - D q_n^2 = (-1)^{n+1}Q_{n+1}$

and $(x, y) = (p_k, q_k)$ is the smallest integer solution of the Pell equation for odd $k$, and $(x, y) = (p_{2k + 1}, q_{2k+1})$ if $k$ is even. From a minimal positive solution $(t,u)$, we may generate the other positive solutions $(x,y)$ via

$x+y\sqrt{D} = (t + u\sqrt{D})^n$

for all positive $n$.

## Generalized Pell Equations

Now consider the equation

$x^2 - D y^2 = N$

If $N^2 \lt D$ we may do the following. Compute the convergents $p_n, q_n$ until the smallest integer solution of the Pell equation $x^2 - D y^2 = 1$ is found. In the meantime, check if each $p_n, q_n$ satisfies $p_n^2 - D q_n^2 = N / f^2$ for some $f \gt 0$. If so, append $(f p_n , f q_n)$ to the list of solutions.

We now use this list of solutions to generate all other solutions. If $(r,s)$ is on the list, and $(t,u)$ is a minimal positive solution of the corresponding Pell equation, then we have a family of solutions $(x,y)$ given by

$x + y \sqrt{D} = (r + s \sqrt{D})(t + u \sqrt{D})^n$

for all positive $n$.

If $N^2 \ge D$ then one possibility is brute force. If $N \gt 0$ set $L_1 = 0, L_2 = \sqrt{N(t-1)/{2 D})}$, otherwise set $L_1 = \sqrt{-N/D}, L_2 = \sqrt{-N(t+1)/(2D)}$. For $L_1 \le y \le L_2$ check if $N + D y^2 = x^2$ for some integer $x$. If so, then append $(x,y)$ to the list of fundamental solutions. Also append $(-x,y)$ if not equivalent to $(x,y)$. (TODO: can avoid since we don’t want negative solutions?) Then proceed as above to generate all solutions.