Pocklington's algorithm

Pocklington's algorithm is a technique for solving a congruence of the form

[math]\displaystyle{ x^2 \equiv a \pmod p, }[/math]

where x and a are integers and a is a quadratic residue.

The algorithm is one of the first efficient methods to solve such a congruence. It was described by H.C. Pocklington in 1917.

The algorithm

(Note: all [math]\displaystyle{ \equiv }[/math] are taken to mean [math]\displaystyle{ \pmod p }[/math], unless indicated otherwise.)

Inputs:

• p, an odd prime
• a, an integer which is a quadratic residue [math]\displaystyle{ \pmod p }[/math].

Outputs:

• x, an integer satisfying [math]\displaystyle{ x^2 \equiv a }[/math]. Note that if x is a solution, −x is a solution as well and since p is odd, [math]\displaystyle{ x \neq -x }[/math]. So there is always a second solution when one is found.

Solution method

Pocklington separates 3 different cases for p:

The first case, if [math]\displaystyle{ p=4m+3 }[/math], with [math]\displaystyle{ m \in \mathbb{N} }[/math], the solution is [math]\displaystyle{ x \equiv \pm a^{m+1} }[/math].

The second case, if [math]\displaystyle{ p=8m+5 }[/math], with [math]\displaystyle{ m \in \mathbb{N} }[/math] and

1. [math]\displaystyle{ a^{2m+1} \equiv 1 }[/math], the solution is [math]\displaystyle{ x \equiv \pm a^{m+1} }[/math].
2. [math]\displaystyle{ a^{2m+1} \equiv -1 }[/math], 2 is a (quadratic) non-residue so [math]\displaystyle{ 4^{2m+1} \equiv -1 }[/math]. This means that [math]\displaystyle{ (4a)^{2m+1} \equiv 1 }[/math] so [math]\displaystyle{ y \equiv \pm (4a)^{m+1} }[/math] is a solution of [math]\displaystyle{ y^2 \equiv 4a }[/math]. Hence [math]\displaystyle{ x \equiv \pm y/2 }[/math] or, if y is odd, [math]\displaystyle{ x \equiv \pm (p+y)/2 }[/math].

The third case, if [math]\displaystyle{ p=8m+1 }[/math], put [math]\displaystyle{ D \equiv -a }[/math], so the equation to solve becomes [math]\displaystyle{ x^2 + D \equiv 0 }[/math]. Now find by trial and error [math]\displaystyle{ t_1 }[/math] and [math]\displaystyle{ u_1 }[/math] so that [math]\displaystyle{ N = t_1^2 - D u_1^2 }[/math] is a quadratic non-residue. Furthermore, let

[math]\displaystyle{ t_n = \frac{(t_1 + u_1 \sqrt{D})^n + (t_1 - u_1 \sqrt{D})^n}{2}, \qquad u_n = \frac{(t_1 + u_1 \sqrt{D})^n - (t_1 - u_1 \sqrt{D})^n}{2 \sqrt{D}} }[/math].

The following equalities now hold:

[math]\displaystyle{ t_{m+n} = t_m t_n + D u_m u_n, \quad u_{m+n} = t_m u_n + t_n u_m \quad \mbox{and} \quad t_n^2 - D u_n^2 = N^n }[/math].

Supposing that p is of the form [math]\displaystyle{ 4m+1 }[/math] (which is true if p is of the form [math]\displaystyle{ 8m+1 }[/math]), D is a quadratic residue and [math]\displaystyle{ t_p \equiv t_1^p \equiv t_1, \quad u_p \equiv u_1^p D^{(p-1)/2} \equiv u_1 }[/math]. Now the equations

[math]\displaystyle{ t_1 \equiv t_{p-1} t_1 + D u_{p-1} u_1 \quad \mbox{and} \quad u_1 \equiv t_{p-1} u_1 + t_1 u_{p-1} }[/math]

give a solution [math]\displaystyle{ t_{p-1} \equiv 1, \quad u_{p-1} \equiv 0 }[/math].

Let [math]\displaystyle{ p-1 = 2r }[/math]. Then [math]\displaystyle{ 0 \equiv u_{p-1} \equiv 2 t_r u_r }[/math]. This means that either [math]\displaystyle{ t_r }[/math] or [math]\displaystyle{ u_r }[/math] is divisible by p. If it is [math]\displaystyle{ u_r }[/math], put [math]\displaystyle{ r=2s }[/math] and proceed similarly with [math]\displaystyle{ 0 \equiv 2 t_s u_s }[/math]. Not every [math]\displaystyle{ u_i }[/math] is divisible by p, for [math]\displaystyle{ u_1 }[/math] is not. The case [math]\displaystyle{ u_m \equiv 0 }[/math] with m odd is impossible, because [math]\displaystyle{ t_m^2 - D u_m^2 \equiv N^m }[/math] holds and this would mean that [math]\displaystyle{ t_m^2 }[/math] is congruent to a quadratic non-residue, which is a contradiction. So this loop stops when [math]\displaystyle{ t_l \equiv 0 }[/math] for a particular l. This gives [math]\displaystyle{ -D u_l^2 \equiv N^l }[/math], and because [math]\displaystyle{ -D }[/math] is a quadratic residue, l must be even. Put [math]\displaystyle{ l = 2k }[/math]. Then [math]\displaystyle{ 0 \equiv t_l \equiv t_k^2 + D u_k^2 }[/math]. So the solution of [math]\displaystyle{ x^2 + D \equiv 0 }[/math] is got by solving the linear congruence [math]\displaystyle{ u_k x \equiv \pm t_k }[/math].

Examples

The following are 4 examples, corresponding to the 3 different cases in which Pocklington divided forms of p. All [math]\displaystyle{ \equiv }[/math] are taken with the modulus in the example.

Example 0

[math]\displaystyle{ x^2 \equiv 43 \pmod {47}. }[/math]

This is the first case, according to the algorithm, [math]\displaystyle{ x \equiv 43^{(47+1)/2} = 43^{12} \equiv 2 }[/math], but then [math]\displaystyle{ x^2=2^2=4 }[/math] not 43, so we should not apply the algorithm at all. The reason why the algorithm is not applicable is that a=43 is a quadratic non residue for p=47.

Example 1

Solve the congruence

[math]\displaystyle{ x^2 \equiv 18 \pmod {23}. }[/math]

The modulus is 23. This is [math]\displaystyle{ 23 = 4 \cdot 5 + 3 }[/math], so [math]\displaystyle{ m=5 }[/math]. The solution should be [math]\displaystyle{ x \equiv \pm 18^6 \equiv \pm 8\pmod {23} }[/math], which is indeed true: [math]\displaystyle{ (\pm 8)^2 \equiv 64 \equiv 18\pmod {23} }[/math].

Example 2

Solve the congruence

[math]\displaystyle{ x^2 \equiv 10 \pmod{13}. }[/math]

The modulus is 13. This is [math]\displaystyle{ 13 = 8 \cdot 1 + 5 }[/math], so [math]\displaystyle{ m=1 }[/math]. Now verifying [math]\displaystyle{ 10^{2m+1} \equiv 10^3 \equiv -1\pmod{13} }[/math]. So the solution is [math]\displaystyle{ x \equiv \pm y/2 \equiv \pm (4a)^{2}/2 \equiv \pm 800 \equiv \pm 7\pmod{13} }[/math]. This is indeed true: [math]\displaystyle{ (\pm 7)^2 \equiv 49 \equiv 10\pmod{13} }[/math].

Example 3

Solve the congruence [math]\displaystyle{ x^2 \equiv 13 \pmod {17} }[/math]. For this, write [math]\displaystyle{ x^2 -13 =0 }[/math]. First find a [math]\displaystyle{ t_1 }[/math] and [math]\displaystyle{ u_1 }[/math] such that [math]\displaystyle{ t_1^2 + 13u_1^2 }[/math] is a quadratic nonresidue. Take for example [math]\displaystyle{ t_1=3, u_1 = 1 }[/math]. Now find [math]\displaystyle{ t_8 }[/math], [math]\displaystyle{ u_8 }[/math] by computing

[math]\displaystyle{ t_2 = t_1 t_1 + 13u_1u_1 = 9-13 = -4 \equiv 13\pmod {17}, }[/math]

[math]\displaystyle{ u_2 = t_1u_1 + t_1u_1 = 3+3 \equiv 6\pmod {17}. }[/math]

And similarly [math]\displaystyle{ t_4 = -299 \equiv 7 \pmod {17}, u_4=156 \equiv 3\pmod {17} }[/math] such that [math]\displaystyle{ t_8=-68 \equiv 0\pmod {17}, u_8 = 42 \equiv 8\pmod {17}. }[/math]

Since [math]\displaystyle{ t_8 = 0 }[/math], the equation [math]\displaystyle{ 0 \equiv t_4^2 + 13u_4^2 \equiv 7^2 - 13 \cdot 3^2\pmod {17} }[/math] which leads to solving the equation [math]\displaystyle{ 3x \equiv \pm 7\pmod {17} }[/math]. This has solution [math]\displaystyle{ x \equiv \pm 8\pmod {17} }[/math]. Indeed, [math]\displaystyle{ (\pm 8)^2 = 64 \equiv 13\pmod {17} }[/math].

References

• Leonard Eugene Dickson, "History Of The Theory Of Numbers" vol 1 p 222, Chelsea Publishing 1952