Fermat's theorem on sums of two squares

1. Introduction

• wiki: https://en.wikipedia.org/wiki/Fermat%27s_theorem_on_sums_of_two_squares
• A prime p is a sum of two integer squares if and only if p = 2 or p ≡ 1 (mod 4).

2. Brahmagupta–Fibonacci identity

• wiki: https://en.wikipedia.org/wiki/Brahmagupta%E2%80%93Fibonacci_identity
• A positive integer n is a sum of two squares if and only if ord_p(n) is even for all primes p ≡ 3 (mod 4).

3. Algorithm

• Reference: 
• http://www.ams.org/journals/mcom/1972-26-120/S0025-5718-1972-0314745-6/S0025-5718-1972-0314745-6.pdf
• http://math.stackexchange.com/questions/5877/efficiently-finding-two-squares-which-sum-to-a-prime/5883#5883
• Hermite-Serret Algorithm:
• Find x^2 ≡ -1 (mod p) where 0 < x < p/2.
• If a is a quadratic nonresidue modulo p then a^((p-1)/2) ≡ -1(mod p), so we can take x ≡ a^((p-1)/4) (mod p). 
• How does one find quadratic nonresidues? Well just over half of all numbers a such that 1< a < p-1 are quadratic nonresidues. 
• Taking a at random will need at most two goes on average.
• Carry out the Euclidean algorithm on p/x, producing the sequence of remainders R_1, R_2, ..., to the point where R_k is first less than Sqrt[p]. Then 
• p = (R_k)^2 + (R_{k+1})^2 if R_1 > 1, 
• p = x^2 + 1^2 if R_1 = 1.

4. Counting:

• Reference: 
• http://matwbn.icm.edu.pl/ksiazki/aa/aa40/aa4029.pdf
• http://mathworld.wolfram.com/SumofSquaresFunction.html