Project Euler 66

Consider a Diophantine equation of the form:
 * $$x^2 - Dy^2 = 1$$

For example, if D=13, the minimal solution in x is $$649^2 - 13 * 180^2 = 1$$.

There are no solutions in positive integers when D is square.

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


 * https://projecteuler.net/problem=66

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


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

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

Explanation
This is specifically Pell's equation. Solutions can be obtained with continued fractions converging on $$\sqrt{D}$$, but there is also a fast solution discovered by Bhāskara II in 1150.

Algorithm:
 * 1) Start with $$( x_0 = \lceil \sqrt D \rceil, y_0 = 1, k_0 = \lceil \sqrt D \rceil^2 - D ) $$
 * 2) Find first m, where $$\frac{a_n +m b_n}{k_n}$$ is an integer; m>0
 * 3) Minimize $$\frac{a_n + m b_n}{k_n}$$ - can add or subtract $$k_n$$; m>0
 * 4) Next step $$( x_{n+1} = \frac{m x_n + D*y_n}{k_n}, y_{n+1} = \frac{x_n + m *y_n}{k_n}, k_{n+1} = \frac{m^2 - D}{k_n} )$$
 * 5) If $$k_n=1$$, solution found.

Gotchas

 * Requires a big number library.

Implementations
Notes/Hints on actual implementation here.

Optimizations
Optimizations here.