import Data.List (genericLength, group, nub, sort, subsequences)
import Data.Numbers.Primes (primeFactors)
import Graphics.Gnuplot.Simple
-- 1ª definición
-- =============
esPractico1 :: Integer -> Bool
esPractico1 n =
takeWhile (<n) (sumas (divisores n)) == [0..n-1]
-- (divisores n) es la lista de los divisores de n. Por ejemplo,
-- divisores 12 == [1,2,3,4,6]
-- divisores 14 == [1,2,7]
divisores :: Integer -> [Integer]
divisores n = [k | k <- [1..n-1], n `mod` k == 0]
-- (sumas xs) es la lista ordenada de números que se pueden obtener como
-- sumas de elementos de xs sin usar ningún elemento más de una vez en
-- cada suma. Por ejemplo,
-- sumas [1,2,3] == [0,1,2,3,4,5,6]
-- sumas [1,2,7] == [0,1,2,3,7,8,9,10]
sumas :: [Integer] -> [Integer]
sumas xs = sort (nub (map sum (subsequences xs)))
-- 2ª definición
-- =============
esPractico2 :: Integer -> Bool
esPractico2 n = all (esSumable (divisores n)) [1..n-1]
-- (esSumable xs n) se verifica si n se puede escribir como una suma de
-- elementos distintos de la lista creciente xs. Por ejemplo,
-- esSumable [1,2,7] 8 == True
-- esSumable [1,2,7] 6 == False
-- esSumable [1,2,7] 4 == False
-- esSumable [1,2,7] 2 == True
-- esSumable [1,2,7] 0 == True
esSumable :: [Integer] -> Integer -> Bool
esSumable _ 0 = True
esSumable [] _ = False
esSumable (x:xs) n = x <= n && (esSumable xs (n-x) || esSumable xs n)
-- 3ª definición
-- =============
-- Usando la caracterización de Stewart y Sierpiński: un entero n >= 2
-- es práctico syss para su factorización prima
-- n = p(1)^e(1) * p(2)*e(2) *...* p(k)^e(k)
-- se cumple que p(1) = 2 y, para cada i de 2 a k se cumple que
-- 1+e(j)
-- i-1 p(j) - 1
-- p(i) <= 1 + ∏ ----------------
-- j=1 p(j) - 1
esPractico3 :: Integer -> Bool
esPractico3 1 = True
esPractico3 n =
x == 2 &&
and [p <= 1 + c | (p,c) <- zip bases cotas]
where xss = factorizacion n
(x:bases) = map fst xss
cotas = scanl1 (*) [(p^(1+e)-1) `div` (p-1) | (p,e) <- xss]
-- (factorizacion n) es la factorización de n. Por ejemplo,
-- factorizacion 600 == [(2,3),(3,1),(5,2)]
-- factorizacion 1400 == [(2,3),(5,2),(7,1)]
factorizacion :: Integer -> [(Integer,Integer)]
factorizacion n =
[(head xs,genericLength xs) | xs <- group (primeFactors n)]
-- Comparación de eficiencia
-- =========================
-- λ> length [n | n <- [1..400], esPractico1 n]
-- 92
-- (40.21 secs, 8,378,539,464 bytes)
-- λ> length [n | n <- [1..400], esPractico2 n]
-- 92
-- (8.29 secs, 1,109,669,760 bytes)
-- λ> length [n | n <- [1..400], esPractico3 n]
-- 92
-- (0.02 secs, 0 bytes)