Sum of squares forms

In the last post we saw that, for example, -3 is a square mod any prime $p = 3k+1$. This means we have $p | u^2 + 3$ or we can say we have a multiple of $p$ expressed as a sum of two squares if we write $p | u^2 + 3\cdot 1^2$. Similarly for primes $p = 4k+1$ we have $p | u^2 + 1^2$.

In some lecture notes ( http://math.bu.edu/people/kost/teaching/MA341/ - Lecture 6 - Sums of Squares ) I found a very nice proof that shows if we can write $p | x^2 + y^2$ then we can "descend" to a smaller $(x,y)$ pair with the same property. By iterating this - since you can only get smaller so many times before you reach 1 - we eventually reach $p = x^2 + y^2$. In fact this proof constitues an algorithm and I have it implemented.

Another key fact about these forms is the uniqueness of representation: Euler showed how to factor a number $n$ given two distinct representations $n = x^2 + y^2 = u^2 + v^2$. This can be found in ( https://archive.org/stream/NumberTheoryItsHistory/Ore-NumberTheoryItsHistory page 61 ). The consequence of this factoring algorithm is that any representation of a prime must be unique.

I was able to adapt both of these algorithms to the form $x^2 + 3y^2$. This proves that all primes $p = 3k + 1$ are represented, and in a unique way by this form.

Factoring

Let $x^2 + 3y^2 = u^2 + 3v^2$ be two distinct primitive representations of a number, by primitive I mean that $(x,y) = 1$ and $(u,v) = 1$ otherwise we would already have a factorization. Hence

$(x-u)(x+u)=3(v-y)(v+y)$

Now $x \equiv \pm u$, so by inverting the sign of $u$ assume $x \equiv u$ this means $3 | x - u$. Using the gcd we can resolve $x-u$ and $x+u$ into parts:

$3 \alpha \beta \gamma \delta = 3(\frac{x-u}{3},v-y)(\frac{x-u}{3},v+y)(x+u,v-y)(x+u,v+y)$

with $\alpha, \beta, \gamma, \delta$ all coprime to each other. Now

• $x-u = 3 \alpha \beta$
• $x+u = \gamma \delta$
• $v-y = \alpha \gamma$
• $v+y = \beta \delta$

from which we can recover

• $x = \frac{3 \alpha \beta + \gamma \delta}{2}$
• $u = \frac{\gamma \delta - 3 \alpha \beta}{2}$
• $v = \frac{\alpha \gamma + \beta \delta}{2}$
• $y = \frac{\beta \delta - \alpha \gamma}{2}$

multiplying out $x^2 + 3y^2$ in terms of these quantities, then considering the identity

$(x^2+3y^2)(u^2+3v^2) = (xu + 3yv)^2 + 3(xv - yu)^2$

leads to the factorization

$x^2 + 3y^2 = (\frac{\gamma^2 + 3\alpha^2}{2})(\frac{\delta^2 + 3\beta^2}{2})$

Descent

Let $dp = x^2 + 3y^2$, we wish to find $d'p = u^2 + 3v^2$ for some smaller $d' < d$ until $p = u^2 + 3v^2$

d=1

If $d=1$ we are done and have a representation using this form.

3|d

If $3|d$ then we can simply "divide by 3": $3d'p = x^2 + 3y^2$ implies $3 | x$ so put $3u = x$ and $3d'p = 3^2\cdot u^2 + 3y^2$ divide out the 3 to get: $d'p = y^2 + 3\cdot u^2$.

3 doesn't divide d

This is the really interesting case, by congruence considerations $d \equiv 1$. Pick $u \equiv x \,(\text{mod}\,d)$ and $v \equiv y \,(\text{mod}\,d)$ with $|u|,|v| < d/2$.

With the identity $(x^2+3y^2)(u^2+3v^2) = (xu + 3yv)^2 + 3(xv - yu)^2$ in mind note

• $a = xy + 3yv \equiv x^2 + 3y^2 \equiv 0 \,(\text{mod}\,d)$
• $b = xv - yu \equiv xy - yx \equiv 0 \,(\text{mod}\,d)$

so $pd^2$ divides the product $a^2 + 3b^2$, thus the pair $(a/d, b/d)$ leads to a smaller representation and lets us recurse.

Code

Here is a bit of code that implements this

{-# LANGUAGE BangPatterns, NoMonomorphismRestriction  #-}
import Prelude hiding ((/))

import System.Environment
--import System.Random
import Debug.Trace

x / y = if not test then error ("doesnt divide " ++ show (x,y)) else x `div` y
 where test = x `mod` y == 0

factorialm _ 0 = 1
factorialm !m !n = (n * factorialm m (n - 1)) `mod` m

--f :: Integer -> 

{-
   forall prime p, p = 4k + 1, exists n. n^2 = -1 mod p
-}

f !p = if not test then error "fail" else u where
    test = (u^2 + 1) `mod` p == 0
    u = factorialm p ((p - 1) / 2)

{-
f' p = do
    n <- randomRIO (2, p - 1)
    let k = (p - 1) / 4
    let g = (n^k) `mod` p
    if (g^2 + 1) `mod` p == 0 then return g else f' p
-}

pow p n 0 !r = r
pow p n 1 !r = n * r
pow p n k !r = pow p n (k - 1) ((n * r)`mod`p)

f3' p = do
     n <- [ 2 .. p - 1 ]
     let k = (p - 1) / 3
     let g = (2*pow p n k 1 + 1) `mod` p
     if (g^2 + 3) `mod` p == 0 then return g else [] -- f3' p

g !p = (f p, 1)

h !p !(a, b) = case () of
    _ | d == 1 -> if not test then error "really fail" else (a, b)
    _ | even d -> if even a then h p (a / 2, b / 2)
                            else h p ((a + b) / 2, (a - b) / 2)
    _ | odd d -> let (u, v) = (q a, q b) in h p ((a * u + b * v) / d, (a * v - b * u) / d)
    where
        d = (a^2 + b^2) / p
        test = a^2 + b^2 == p
        q !a = let a' = a `mod` d in if 2 * a' > d then a' - d else a'

h3 !p !(a, b) = {- trace (show [a,b,d, q a, q b]) $ -} case () of
     _ | d == 1 -> if not test then error "really fail" else (a, b)
     _ | (d`mod`3) == 0 -> h3 p (b, a/3)
     _ | otherwise -> let (u, v) = (q a, q b) in
           h3 p ((a * u + 3 * b * v) / d, (a * v - b * u) / d)
    where
        d = (a^2 + 3*b^2) / p
        test = a^2 + 3*b^2 == p
        q !a = let a' = a `mod` d in if 2 * a' > d then a' - d else a'
{-

*Main> :l n
[1 of 1] Compiling Main             ( n.hs, interpreted )
Ok, modules loaded: Main.
*Main> h3 10687 (7307, 1)
(98,19)
*Main> 98^2 + 3*19^2
10687
*Main> 3*3562210+1
10686631
*Main> head $ f3' (3*3562210+1)
10239214
*Main> h3 (3*3562210+1) (10239214, 1)
(2078,1457)
*Main> 2078^2 + 3*1457^2
10686631
*Main> 


-}

Consequences

The reason for all of this (besides being able to resolve primes into sum of squares forms) is that I want to prove something about $n^3 = a^2 + 3b^2$

supposing $(a,b) = 1$. This theorem comes from Weil - Number theory through history.

In this case any prime divisor $p$ of $n$ divides $a^2 + 3b^2$ and hence is of this form - therefore $n$ is a product of these forms and so we must have $n = h^2 + 3k^2$ and therefore $n^3 = (h^3 + 3k^2h)^2 + 3(kh^2 + 3k^3)^2$

by the following calculation

? f(x,y,u,v) = [x*u+3*y*v, x*v - y*u]
%1 = (x,y,u,v)->[x*u+3*y*v,x*v-y*u]
? f(h,k,h,k)
%2 = [h^2 + 3*k^2, 0]
? f(h^2 + 3*k^2,0,h,k)
%3 = [h^3 + 3*k^2*h, k*h^2 + 3*k^3]

we can infact choose $h,k$ to make

• $a = h^3 + 3k^2h$
• $b = kh^2 + 3k^3$
• $(h,k) = 1$

This will be useful for a proof of fermats last theorem in the case $n=3$.

All of this was completely elementary and did not make use of any irrational or imaginary numbers. The same results can be proved much quicker using algebraic number fields - they are of course what lurks behind all the mathematics here.