import Data.Numbers.Primes (primes)
import Test.QuickCheck
-- 1ª definición de representaciones
-- =================================
representaciones :: Integer -> [(Integer,Integer)]
representaciones n =
[(x,y) | x <- [0..n], y <- [x..n], n == x*x + y*y]
-- 2ª definición de representaciones
-- =================================
representaciones2 :: Integer -> [(Integer,Integer)]
representaciones2 n =
[(x,raiz z) | x <- [0..raiz (n `div` 2)]
, let z = n - x*x
, esCuadrado z]
-- (esCuadrado x) se verifica si x es un número al cuadrado. Por
-- ejemplo,
-- esCuadrado 25 == True
-- esCuadrado 26 == False
esCuadrado :: Integer -> Bool
esCuadrado x = x == y * y
where y = raiz x
-- (raiz x) es la raíz cuadrada entera de x. Por ejemplo,
-- raiz 25 == 5
-- raiz 24 == 4
-- raiz 26 == 5
raiz :: Integer -> Integer
raiz 0 = 0
raiz 1 = 1
raiz x = aux (0,x)
where aux (a,b) | d == x = c
| c == a = a
| d < x = aux (c,b)
| otherwise = aux (a,c)
where c = (a+b) `div` 2
d = c^2
-- 3ª definición de representaciones
-- =================================
representaciones3 :: Integer -> [(Integer,Integer)]
representaciones3 n =
[(x,raiz3 z) | x <- [0..raiz3 (n `div` 2)]
, let z = n - x*x
, esCuadrado3 z]
-- (esCuadrado x) se verifica si x es un número al cuadrado. Por
-- ejemplo,
-- esCuadrado3 25 == True
-- esCuadrado3 26 == False
esCuadrado3 :: Integer -> Bool
esCuadrado3 x = x == y * y
where y = raiz3 x
-- (raiz3 x) es la raíz cuadrada entera de x. Por ejemplo,
-- raiz3 25 == 5
-- raiz3 24 == 4
-- raiz3 26 == 5
raiz3 :: Integer -> Integer
raiz3 x = floor (sqrt (fromIntegral x))
-- 4ª definición de representaciones
-- =================================
representaciones4 :: Integer -> [(Integer, Integer)]
representaciones4 n = aux 0 (floor (sqrt (fromIntegral n)))
where aux x y
| x > y = []
| otherwise = case compare (x*x + y*y) n of
LT -> aux (x + 1) y
EQ -> (x, y) : aux (x + 1) (y - 1)
GT -> aux x (y - 1)
-- Equivalencia de las definiciones de representaciones
-- ====================================================
-- La propiedad es
prop_representaciones_equiv :: (Positive Integer) -> Bool
prop_representaciones_equiv (Positive n) =
representaciones n == representaciones2 n &&
representaciones2 n == representaciones3 n &&
representaciones3 n == representaciones4 n
-- La comprobación es
-- λ> quickCheck prop_representaciones_equiv
-- +++ OK, passed 100 tests.
-- Comparación de eficiencia de las definiciones de representaciones
-- =================================================================
-- λ> representaciones 3025
-- [(0,55),(33,44)]
-- (2.86 secs, 1,393,133,528 bytes)
-- λ> representaciones2 3025
-- [(0,55),(33,44)]
-- (0.00 secs, 867,944 bytes)
-- λ> representaciones3 3025
-- [(0,55),(33,44)]
-- (0.00 secs, 173,512 bytes)
-- λ> representaciones4 3025
-- [(0,55),(33,44)]
-- (0.00 secs, 423,424 bytes)
--
-- λ> length (representaciones2 (10^10))
-- 6
-- (3.38 secs, 2,188,903,544 bytes)
-- λ> length (representaciones3 (10^10))
-- 6
-- (0.10 secs, 62,349,048 bytes)
-- λ> length (representaciones4 (10^10))
-- 6
-- (0.11 secs, 48,052,360 bytes)
--
-- λ> length (representaciones3 (4*10^12))
-- 7
-- (1.85 secs, 1,222,007,176 bytes)
-- λ> length (representaciones4 (4*10^12))
-- 7
-- (1.79 secs, 953,497,480 bytes)
-- Definición de primosImparesConRepresentacionUnica
-- =================================================
primosImparesConRepresentacionUnica :: [Integer]
primosImparesConRepresentacionUnica =
[x | x <- tail primes
, length (representaciones4 x) == 1]
-- Definición de primos4nM1
-- ========================
primos4nM1 :: [Integer]
primos4nM1 = [x | x <- primes
, x `mod` 4 == 1]
-- Teorema de Navidad de Fermat
-- ============================
-- La propiedad es
prop_teoremaDeNavidadDeFermat :: Positive Int -> Bool
prop_teoremaDeNavidadDeFermat (Positive n) =
primosImparesConRepresentacionUnica !! n == primos4nM1 !! n
-- La comprobación es
-- λ> quickCheck prop_teoremaDeNavidadDeFermat
-- +++ OK, passed 100 tests.