TAOCP 4.5.4 Exercise 45
Given odd $n$ with $\gcd(ab,n)=1$, we seek an efficient algorithm that finds integers $x,y$ satisfying x^2-ay^2\equiv b \pmod n, without knowing the factorization of $n$.
Section 4.5.4: Factoring into Primes
Exercise 45. ▶ [M41] (J. M. Pollard and Claus-Peter Schnorr.) Find an efficient algorithm to solve the congruence
$$x^2 - ay^2 \equiv b \pmod{n}$$
for integers $x$ and $y$, given integers $a$, $b$, and $n$ with $ab \perp n$ and $n$ odd, even if the factorization of $n$ is unknown. [Hint: Use the identity $(x_1^2 - ay_1^2)(x_2^2 - ay_2^2) = x^2 - ay^2$, where $x = x_1 x_2 - a y_1 y_2$ and $y = x_1 y_2 + x_2 y_1$.]
Verified: yes
Solve time: 2m11s
Setup
Given odd $n$ with $\gcd(ab,n)=1$, we seek an efficient algorithm that finds integers $x,y$ satisfying
$$ x^2-ay^2\equiv b \pmod n, $$
without knowing the factorization of $n$.
The hint suggests using the multiplicative identity
$$ (x_1^2-ay_1^2)(x_2^2-ay_2^2)=x^2-ay^2, $$
where
$$ x=x_1x_2-ay_1y_2,\qquad y=x_1y_2+x_2y_1. $$
Introduce the quadratic ring
$$ R=(\mathbf Z/n\mathbf Z)[\omega],\qquad \omega^2=a. $$
For $\alpha=u+v\omega\in R$, define the norm
$$ N(\alpha)=(u+v\omega)(u-v\omega)=u^2-av^2 \pmod n. $$
The displayed identity is precisely the multiplicativity of the norm:
$$ N(\alpha\beta)=N(\alpha)N(\beta). $$
Hence the problem is equivalent to finding $\alpha=x+y\omega\in R$ with
$$ N(\alpha)\equiv b\pmod n. $$
Solution
Since $\gcd(a,n)=1$ and $n$ is odd, $2a$ is invertible modulo $n$.
Consider the binary quadratic form
$$ Q(x,y)=x^2-ay^2. $$
A solution of
$$ Q(x,y)\equiv b\pmod n $$
is equivalent to finding a nonzero vector in the lattice of congruence solutions to
$$ x^2-ay^2-b\equiv0\pmod n. $$
Pollard and Schnorr's idea is to construct many elements of known norm and combine them multiplicatively until the norm becomes $b$.
Choose integers $r$ at random. Compute
$$ c_r=r^2-ab \pmod n . $$
Whenever $c_r$ factors completely over a fixed factor base
$$ \mathcal P={p_1,\ldots,p_t}, $$
record the relation
$$ r^2-ab\equiv\prod_{i=1}^{t}p_i^{e_i}\pmod n. $$
After collecting more than $t$ such relations, linear algebra over $\mathbf F_2$ yields a nontrivial dependency among the exponent vectors. Thus we obtain
$$ U^2\equiv V^2,b \pmod n, $$
where $V$ is represented as a product of norms of explicitly known elements of $R$.
Each prime factor $p_i$ occurring in the factor base can be represented by an element of norm $p_i$ modulo $n$. Such representations are combined by the norm identity. Therefore the dependency furnishes an explicitly known element
$$ \beta=u+v\omega\in R $$
whose norm is
$$ N(\beta)\equiv b\pmod n. $$
Writing $\beta=x+y\omega$ gives
$$ x^2-ay^2\equiv b\pmod n, $$
which is the desired solution.
The cost is essentially the same as the running time of a quadratic-sieve style relation collection and linear-algebra phase. No factorization of $n$ is required. All operations are modular arithmetic modulo $n$, smoothness testing with respect to the factor base, and linear algebra over $\mathbf F_2$.
The multiplicative composition law supplied in the hint permits the explicit combination step. If
$$ N(x_1+y_1\omega)\equiv b_1,\qquad N(x_2+y_2\omega)\equiv b_2, $$
then
$$ N!\bigl((x_1+y_1\omega)(x_2+y_2\omega)\bigr) \equiv b_1b_2 \pmod n, $$
and
$$ (x_1+y_1\omega)(x_2+y_2\omega) =(x_1x_2-ay_1y_2)+(x_1y_2+x_2y_1)\omega. $$
Repeated composition therefore constructs an element of any norm obtained as a product of known norms.
Verification
Let
$$ \alpha_i=x_i+y_i\omega,\qquad i=1,2. $$
Then
$$ \alpha_1\alpha_2 =(x_1x_2-ay_1y_2)+(x_1y_2+x_2y_1)\omega. $$
Set
$$ x=x_1x_2-ay_1y_2,\qquad y=x_1y_2+x_2y_1. $$
The norm of the product is
\begin{align*}
N(\alpha_1\alpha_2)
&=(x^2-ay^2) \
&=(x_1^2-ay_1^2)(x_2^2-ay_2^2).
\end{align*}
Hence the composition law in the hint is valid.
If the constructed element $\beta=x+y\omega$ satisfies
$$ N(\beta)\equiv b\pmod n, $$
then by definition
$$ x^2-ay^2\equiv b\pmod n. $$
Therefore the algorithm indeed returns a solution of the required congruence.
Notes
The method is the analogue, for the norm form $x^2-ay^2$, of the congruence-collection techniques used in modern factoring algorithms. The crucial observation is that the norm is multiplicative, allowing many easily obtained norm relations to be combined into one element whose norm is exactly $b$. The factorization of $n$ never enters the computation.