Integer and Rational Solutions of a Binary Quadratic Equation (4): Pell equation and generalized Pell equation
In this post, we are interested in finding all integer solutions of Pell equation \[ \begin{equation} x^2 - d y^2 = 1, \label{eq:pell} \tag{Pell Equation} \end{equation} \] and generalized Pell equation \[ \begin{equation} x^2 - d y^2 = n, \label{eq:generalized-pell} \tag{Generalized Pell Equation} \end{equation} \] for a non-square integer \(d > 0\) and \(n \neq 0\).
Connection to \(\ZZ[\sqrt{ d }]\)
For any element \(s = x + y \sqrt{ d }\) in \(\ZZ[\sqrt{ d }]\), we define \(\bar{s} := x - y \sqrt{ d }\). It is easy to show that \(\bar{s} \bar{t} = \overline{st}\) for any \(s, t \in \ZZ[\sqrt{d}]\).
Note that for \(s = x + y \sqrt{ d }\), \(x^2 - dy^2 = s \bar{s}.\) Now we define \[ S_{n} := \{ s \in \ZZ[\sqrt{ d }]: s \bar{s} = n \}. \] so solving Pell equations is equivalent to calculating \(S_1\) and solving generalized Pell equations is equivalent to calculating \(S_n\).
Most references I found online use the \((x, y)\) notation, but in certain cases, it might be cleaner to present the results using operations in \(\ZZ[\sqrt{ d }]\), especially when presenting the structure of solutions.
Primitive solutions
(Primitive solutions) A solution \((x, y)\) is primitive if \(\gcd(x, y) = 1\). If \((x, y)\) is primitive, \(\gcd(x, n) = \gcd(y, n) = 1\).
Similar to solving \(ax^2 + by^2 = n\), once we have an algorithm to find all primitive solutions, we can derive all solutions. Let \(S'_n\) be the set of all primitive solutions of \(x^2 - d y^2 = n\): \[ S'_{n}:= \left\{ s\in \ZZ[\sqrt{ d }]: s \bar{s} = n, \quad s = x + y \sqrt{ d }, \quad \gcd(x, y) = 1 \right\}, \] then it is straightforward to see that \[ S_{n} = \left\{ gs: g \in \NN, \quad g^2 | n, \quad s \in S'_{n / g^2} \right\} . \] Clearly, \(S_1 = S'_1\).
Pell equation \(x^2 - d y^2 = 1\)
The Pell equation has been studied extensively over time. See [Lenstra02] for a comprehensive overview. There are so many tutorials (e.g. [Conrad22A] and [Conrad22B]) on solving Pell equations, so here I only list a few key conclusions and the algorithm.
(Structure of solutions) It follows from Dirichlet's unit theorem that \(S_{1}\) has the structure \[ S_{1} = \{ \pm s^{*k}: k \in \ZZ \}, \] for a unique fundamental solution \(s^* := x^* + y^* \sqrt{d}\) where \(x^* > 0, y^* > 0\). Note \(n\) can be negative, as \(\bar{s} = s^{*-1}\). See [Conrad22A, Theorem 5.3] for proof.
(Fundamental solution). It is also shown that \(x^*/y^*\) must be a convergent of \(\sqrt{d}\). We can use the PQa algorithm below to find all convergents. See [Conrad22B, Theorem 5.1] for proof.
(PQa algorithm) Given \(p_0, q_0\) and \(d\) such that \(p_0^2 \equiv d \pmod {q_{0}}\), the PQa algorithm computes the continued fraction \([ a_0, a_1, a_2, \dots]\) of \((p_0 + \sqrt{d}) / q_0\). It works as follows: \[ \begin{aligned} a_i & = \floor{(p_i + \sqrt{d}) / q_i}, \\ p_{i+1} & = a_i q_i - p_i, \\ q_{i+1} & = (d -p_{i+1}^2)/q_i. \\ \end{aligned} \] Note that the division for calculating \(q_{i+1}\) is always exact.
Below is an implementation of the PQa algorithm optimized for solving the Pell equation. See [Robertson04] for more details.
def pell(d): |
Generalized Pell equation \(x^2 -d y^2 = n\)
(The idea in this section is taken from [Matthews00].)
(Nagell equivalence) For any \(s \in S'_{n}\) and \(t \in S'_1\), we observe \[ (s t) \overline{(st)} = (s \overline{s})(t \overline{t}) = n, \quad (st) \bar{t} = s \in S'_{n}, \] implying that \(s t \in S'_{n}\). This motivates us to define equivalence relation, the Nagell equivalence, for elements in \(S'_{n}\): Two solutions \(s \sim s'\) if and only if there exists an element \(t \in S'_{1}\), such that \(st = s'\).
(Partition of solution space) Let \(s = x + y \sqrt{ d }\) and \(s' = x' + y' \sqrt{ d }\) be two primitive solutions, and define \(\lambda(s) := -x y^{-1} \bmod |n|\). Since \(\lambda(s)^2 \equiv d \pmod {|n|}\), it follows that \[ \begin{aligned} s \sim s' \quad & \Longleftrightarrow \quad s s'^{-1} \in \ZZ[\sqrt{ d }] \\ & \Longleftrightarrow \quad x x' \equiv yy' d \pmod{|n|}, \quad xy' \equiv x' y \pmod {|n|} \\ & \Longleftrightarrow \quad \tau(s) = \tau(s'). \end{aligned} \] So an alternative definition of Nagall equivalence for primitive solutions is that \(s \sim s'\) if and only if \(\lambda(s) = \lambda(s')\).
(Structure of solutions) The alternative definition of Nagall equivalence directly implies the number of equivalent classes in \(S'_n\) is finite, that is, there exists a finite set \(F_{n} \subset \ZZ[\sqrt{ d }]\) such that \[ S'_{n} = \left\{f t: f \in F_{n}, t \in S'_{1} \right\} = \left\{ \pm f s^{*n}: F \in F_{n, }n \in \ZZ \right\} , \] where \(s^*\) is the fundamental solution of \(x^2 - dy^2 = 1\).
(Finding a minimum positive solution) Now we focus on one equivalence class, represented by \(\lambda \equiv -x y^{-1} \pmod{|n|}\) (with abuse of notation) where \(-|n| / 2 < \lambda \leq |n| / 2\). If we write \(x = x' |n| - \lambda y\), it is proved in [Matthews00, Theorem 1] that \(x' / y\) is a convergent of \((\lambda + \sqrt{ d }) / |n|\). Thus, we can find a minimum positive solution \(s = x + y \sqrt{ d }\) with smallest positive \(x\) and \(y\) using the PQa algorithm.
(Finding a fundamental solution) It might be more convenient to find a fundamental solution for an equivalence class, with smallest positive \(y\). If multiple solutions have the same smallest positive \(y\), choose the one with positive \(x\). A fundamental solution can be found as follows:
- If \(\lambda \in \{ 0, n / 2\}\), the minimum positive solution is the fundamental solution;
- Otherwise, let \((x, y)\) and \((x', y')\) be the minimum positive solutions of \(\lambda\) and \(-\lambda\). If \(y < y'\), the fundamental solutions are \((x, y)\) and \((-x, y)\); otherwise they're \((-x', y')\) and \((x', y')\).
The algorithm is almost the same as [Mollin01, Algorithm 4.1]. The Lagrange-Matthews-Mollin (LMM) algorithm, proposed in [Matthews00, Section 5], improves the process by considering properties of the continued fraction of \((\lambda + \sqrt{ d }) / |n|\).
Here is the pseudocode:
def solve_generalized_pell_primitive(d, n): |
A special generalized Pell equation \(x^2 - dy^2 = 4\)
We also discuss special generalized Pell equation \(x^2 - dy^2 = 4\), as it also shares many properties as the Pell equation \(x^2 - dy^2 = 1\).
As we have mentioned above, the set \(S_1\) has a special structure \(S_1 = \{s_{1}^{*n}: n \in \ZZ\}\) for a fundamental solution \(s_{1}^*\). Similarly, the solutions for \(x^2 - dy^2 = 4\) have a similar structure \[ T_{4} = \left\{ \pm t_{4}^{*k}: k \in \ZZ \right\}, \quad T_{4} := \left\{ \frac{x+y\sqrt{ d }}{2} : x^2 - dy^2 = 4 \right\}. \] for a fundamental solution \(t_4^* = \frac{x^* + y^* \sqrt{ d }}{2}\) where \((x^*, y^*)\) is minimum positive solution to \(x^2 - dy^2 = 4\) [Stolt57, p.4]. There is actually a connection between \(s_{1}^{*}\) and \(t_4^{*}\): \(s_1^{*} = t_{4}^{*k}\) for \(k \in \left\{ 1,2,3 \right\}\) ([Stolt57, p.4], [Matthews23, Theorem 1.1]).
(Equivalence relation) Observing that Negall equivalence is the associated equivalence relation of the group action of group \(\{\pm s_1^{*k}: k \in \ZZ\}\) on the set \(S'_n\). As \(T_4\) is also a group, we can define a similar equivalence relation, the Stolt equivalence. In some sense, the Stolt equivalence is more fundamental than the Negall equivalence, in the sense that \(s_1^*\) can be generated by \(t_4^*\) but \(t_4^*\) might not be generated by \(s_1^*\). As a consequence, the Stolt equivalence might have fewer equivalence classes than the Negall equivalence.
To find the fundamental solution \(t_4^*\), check out [Johnson04] for an efficient algorithm implemented below:
def pell4(d): |
References
Rathnayake13: Rathnayake T., 2013. Quadratic Diophantine equation – I [online]
Lenstra02: Lenstra Jr, H.W., 2002. Solving the Pell equation. Notices of the AMS, 49(2), pp.182-192.
Conrad22A: Conrad, K., 2022. Pell’s equation, I [online]
Conrad22B: Conrad, K., 2022. Pell’s equation, II [online]
Matthews00: Matthews, K., 2000. The Diophantine Equation \(x^ 2-Dy^ 2= N\), \(D> 0\). Expositiones Mathematicae, 18(4), pp.323-332.
- slides.
Tamang21: Tamang, B.B., 2021. On the study of quadratic Diophantine equations (Doctoral dissertation, Tribhuvan University Kathmandu).
- Discussed this topic in more details and shows more examples.
Mollin01: Mollin, R.A., 2001. Simple continued fraction solutions for Diophantine equations. Expositiones Mathematicae, 19(1), pp.55-73.
Robertson04: Robertson, J.P., 2004. Solving the generalized Pell equation \(x^2− Dy^2= N\). Unpublished manuscript.
- It gives a complete description of the LMM algorithm.
Stolt57: Stolt, B., 1957. On a Diophantine equation of the second degree. Arkiv för Matematik, 3(4), pp.381-390.
Update @2025-03-20: Added the discussion of \(x^2 - dy^2 = 4\) and the optimized algorithm for the Pell equation.