First 30 solutions of Pell's equation.

Floating point precision is not sufficient for the larger solutions. You should probably look at using a recursive definition for the $k_n$ and $z_n$. They should satisfy $a_{n+2}=18a_{n+1}-a_{n}$ (derived from the minimum polynomial of $9\pm4\sqrt{5}$). Try to find a tail recursive solution, and you shouldn't run into problems until you overflow the ints.

One problem is that round-off error will gradually wreak havoc on any floating-point computation using $\sqrt{D}$. A surer way to generate solutions is to use Python's arbitrarily-long integers and do integer arithemetic, looking at the integer-coefficient recursion for the (integer!) coefficients of powers $(a+b\sqrt{D})^n$.

Any formula of the form $a_n=cx^{n}+dy^{n}$ can be computed recursively as:

$$a_{n+1}=(x+y)a_n - xya_{n-1}$$ if you know $a_0,a_1$.

In your case, you've got $x,y=9\pm 4\sqrt{5}$, and $x+y=18$, $xy=1$.

So you have $k_0=1,k_1=9,z_0=0,z_1=4$:

$$k_{n+1}=18k_n-k_{n-1}\\z_{n+1}=18z_n-z_{n-1}$$

If you only want to compute one pair, $(k_n,z_n)$ you can compute this in $O(\log n)$ time by computing:

$$\begin{pmatrix}0&1\\ -1&18 \end{pmatrix}^{n-1}\begin{pmatrix}1&0\\ 9&4\end{pmatrix}=\begin{pmatrix}k_{n-1}&z_{n-1}\\k_{n}&z_{n}\end{pmatrix}$$

The reason this is an $O(\log n)$ operation is that you can apply exponentiation by squaring. to raise the matrix to the $n$th power.

But if you need all values $k_1,\dots,k_n$ and $z_1,\dots,z_n$, you'll of course require $O(n)$ time (you can't compute a list of $2n$ values in less than $O(n)$ time, after all.)

---

Another matrix approach, probably more direct, is to note that $$k_{n+1}+z_{n+1}\sqrt{5}=(k_n+z_n\sqrt5)(9+4\sqrt{5})=(9k_n+20z_n)+(4k_n+9z_n)\sqrt{5}$$

So we can write:

$$\begin{pmatrix}k_n\\z_n\end{pmatrix}=\begin{pmatrix}9&20\\4&9\end{pmatrix}^n\begin{pmatrix}1\\0\end{pmatrix}$$

Again, exponentiation by squaring allows $O(\log n)$ time calculation. (It might be closer to $O(\log^2 n)$ in both cases because of the number of bits in the integer multiplications, but still much better than $O(n)$.)