import Data.Array
import Data.Numbers.Primes
import Graphics.Gnuplot.Simple
import Test.QuickCheck
-- 1ª definición
-- =============
indiceGoldbach :: Int -> Int
indiceGoldbach n =
minimum (map length (particiones n))
particiones :: Int -> [[Int]]
particiones n = v ! n where
v = array (0,n) [(i,f i) | i <- [0..n]]
where f 0 = [[]]
f m = [x:y | x <- xs,
y <- v ! (m-x),
[x] >= take 1 y]
where xs = reverse (takeWhile (<= m) primes)
-- 2ª definición
-- =============
indiceGoldbach2 :: Int -> Int
indiceGoldbach2 x =
head [n | n <- [1..], esSumaDe x n]
-- (esSumaDe x n) se verifica si x se puede escribir como la suma de n
-- primos. Por ejemplo,
-- esSumaDe 2 1 == True
-- esSumaDe 4 1 == False
-- esSumaDe 4 2 == True
-- esSumaDe 27 2 == False
-- esSumaDe 27 3 == True
esSumaDe :: Int -> Int -> Bool
esSumaDe x 1 = isPrime x
esSumaDe x n = or [esSumaDe (x-p) (n-1) | p <- takeWhile (<= x) primes]
-- 3ª definición
-- =============
indiceGoldbach3 :: Int -> Int
indiceGoldbach3 x =
head [n | n <- [1..], esSumaDe3 x n]
esSumaDe3 :: Int -> Int -> Bool
esSumaDe3 x n = a ! (x,n) where
a = array ((2,1),(x,9)) [((i,j),f i j) | i <- [2..x], j <- [1..9]]
f i 1 = isPrime i
f i j = or [a!(i-k,j-1) | k <- takeWhile (<= i) primes]
-- 4ª definición
-- =============
indiceGoldbach4 :: Int -> Int
indiceGoldbach4 n = v ! n where
v = array (2,n) [(i,f i) | i <- [2..n]]
f i | isPrime i = 1
| otherwise = 1 + minimum [v!(i-p) | p <- takeWhile (< (i-1)) primes]
-- Comparación de eficiencia
-- =========================
-- λ> sum (map indiceGoldbach [2..80])
-- 142
-- (2.66 secs, 1,194,330,496 bytes)
-- λ> sum (map indiceGoldbach2 [2..80])
-- 142
-- (0.01 secs, 1,689,944 bytes)
-- λ> sum (map indiceGoldbach3 [2..80])
-- 142
-- (0.03 secs, 27,319,296 bytes)
-- λ> sum (map indiceGoldbach4 [2..80])
-- 142
-- (0.03 secs, 47,823,656 bytes)
--
-- λ> sum (map indiceGoldbach2 [2..1000])
-- 2030
-- (0.10 secs, 200,140,264 bytes)
-- λ> sum (map indiceGoldbach3 [2..1000])
-- 2030
-- (3.10 secs, 4,687,467,664 bytes)
-- Gráfica
-- =======
graficaGoldbach :: Int -> IO ()
graficaGoldbach n =
plotList [ Key Nothing
, XRange (2,fromIntegral n)
, PNG ("Conjetura_de_Goldbach_" ++ show n ++ ".png")
]
[indiceGoldbach2 k | k <- [2..n]]
-- Comprobación de la conjetura de Goldbach
-- ========================================
-- La propiedad es
prop_Goldbach :: Int -> Property
prop_Goldbach x =
x >= 2 ==> indiceGoldbach2 x < 4
-- La comprobación es
-- λ> quickCheck prop_Goldbach
-- +++ OK, passed 100 tests.