Taking square roots modulo a prime.

The following algorithm, given an odd prime `p` and a quadratic residue `a` (mod `p`) -- which means a nonzero residue `a` (mod `p`) such that the equation

`x2 mod p = a`
has a solution -- returns an `x` such that `x2 mod p = a`. It comes from Leonard Adleman, Kenneth Manders, and Gary Miller ("On taking roots in finite fields", in 18th Annual Symposium on Foundations of Computer Science, IEEE (1977), pp. 175--178.)

This algorithm relies on Euler's criterion: if `p` is an odd prime and `a` is a nonzero residue mod `p`, then

• `a(p-1)/2 mod p = +1` if and only if `a` is a quadratic residue,
• `a(p-1)/2 mod p = -1` if and only if `a` is not a quadratic residue.
Furthermore, it relies on the fact that precisely one half of the nonzero residues mod `p` are quadratic residues.

```if   p mod 4 = 3
then return a(p+1)/4 mod p;
else if   p mod 8 = 5
then if   a(p-1)/4 mod p = 1
then return a(p+3)/8 mod p;
else /* now a(p-1)/4 mod p = -1 */
repeat choose a nonzero residue `b` (mod `p`) at random;
until  b(p-1)/2 mod p = -1
return a(p+3)/8·b(p-1)/4 mod p;
end
else /* now p mod 8 = 1 */
s = (p-1)/4;
while as mod p = 1
do    if   s is odd
then return a(s+1)/2 mod p;
else s = s/2;
end
end
/* now as mod p = -1 */
repeat choose a nonzero residue `b` (mod `p`) at random;
until  b(p-1)/2 mod p = -1
t=(p-1)/2;
while s is even do
do    s = s/2, t=t/2;
if as·bt mod p = -1 then t = t+(p-1)/2 end
/* now as·bt mod p = 1 */
end
return a(s+1)/2·bt/2 mod p;
end
end

```