# Solving Pell's equation: algorithm to converge $\sqrt n$

I'm trying to come up with an algorithm to solve the Diophantine equation

$$x^2 - ny^2 = 1$$

for minimum values of $x$ when $n$ is given. This equation is also known as Pell's Equation.

The only accepted solutions are positive integer values of $x$ and $y$, and it has been proven that no solutions exists for any $n$ that is a square.

I know that I can solve this by looking at all convergents for $\sqrt n$. For example, for $n = 2$, the solution is $x = 3$ and $y = 2$, which is the first convergent of $\sqrt 2$: $3 \over 2$.

However, I'm struggling with how to determine these convergents. I'm using this recursive formula:

$$\sqrt n = 1 + {n - 1 \over 1 + \sqrt n}$$

To calculate a $k ^{th}$ convergent, I apply the formula with depth $k$ and fill in $1$ for the base case.

This works for the first few cases, but for $n =6$, this gives me $7 \over 2$ for the first convergent, whereas I know this should be $5 \over 2$, because $x = 5$ and $y = 2$ is the solution for $x^2 - 6y^2 = 1$.

I arrive at $7 \over 2$ like this: $\sqrt 6 \approx 1 + {6 - 1 \over 1 + 1} = {7 \over 2 }$

I finally figured this out. There are remarkably few good resources on this topic. I finally stumbled upon this video, which (unexpectedly) does a great job explaining the convergence of square roots: https://www.youtube.com/watch?v=49CKAPq8g-w

I got the formula used in the question from wikipedia and although it kind of works, it's not the correct way of converging square roots and doesn't lead to correct solutions to the Pell's equation.

The correct way of converging is defining the term to be converged as an infinite formula of the following form:

$$T = b_0 + {1 \over b_1 + {1 \over b_2 + {1 \over ... }} }$$

In this formula, $b_0,b_1,...,b_n$ are the partial denominators that we can use to get the convergents. This answer describes an algorithms to calculate these values for the convergence of square roots of a known $n$.

Let's call the term to be converged $T$. $T$ can be expressed like this: $$T = {a_0 \sqrt x + a_1 \over a_2}$$ For the first iteration, $T = \sqrt x$,

so: $a_0, a_1, a_2 = (1,0,1)$

• Find the highest integer $b$ for which $b < T$, or: $b = \lfloor T \rfloor$. This is our first partial denominator.
• Separate $T$ and $b$: $$T = b + T - b$$
• Replacing $T$ gives: $$T = b + {a_0 \sqrt x + a_1 \over a_2} - b$$
• This can be unified: $$T = b + {a_0 \sqrt x + a_1 - ba_2 \over a_2}$$

• Flip this to $$T = b + {1 \over {a_2 \over a_0 \sqrt x - (ba_2-a_1)}}$$

• Rationalize the lower term by multiplying with $${a_0 \sqrt x + (ba_2-a_1) \over a_0 \sqrt x + (ba_2-a_1)}$$

• Like this: $$T = b + {1 \over {a_2 \over a_0 \sqrt x - (ba_2-a_1)}{a_0 \sqrt x + (ba_2-a_1) \over a_0 \sqrt x + (ba_2-a_1)} }$$

• This gives: $$T = b + {1 \over {a_2 a_0 \sqrt x + a_2(ba_2-a_1) \over a_0^2x - (ba_2-a_1)^2 } }$$

• Now this may seem like we didn't gain much, but we actually have a term here that is of the same form as what we started with. So now we can do the next iteration by getting new values for $a_0,a_1,a_2$:

$$a'_0 = a_2a_0$$ $$a'_1 = a_2(ba_2-a_1)$$ $$a'_2 = a_0^2x - (ba_2-a_1)^2$$

• And now simply get the next partial denominator with the floor function! (step 1):

$$b' = \left\lfloor {a'_0 \sqrt x + a'_1 \over a'_2} \right\rfloor$$

Final note: when implementing this as an algorithm, these terms $a_0,a_1,a_2$ should be simplified after each iteration by dividing them with their greatest common divisor. Otherwise the numbers might grow too big to handle.