sumasDeDosAbundantes :: [Integer]

   take 10 sumasDeDosAbundantes  ==  [24,30,32,36,38,40,42,44,48,50]

module Sumas_de_dos_abundantes where

import Data.List (genericLength, group)
import Data.Numbers.Primes (primeFactors)
import Test.Hspec (Spec, describe, hspec, it, shouldBe)
import Test.QuickCheck

-- 1ª solución
-- ===========

sumasDeDosAbundantes1 :: [Integer]
sumasDeDosAbundantes1 = [n | n <- [1..], esSumaDeDosAbundantes n]

-- (esSumaDeDosAbundantes n) se verifica si n es suma de dos números
-- abundantes. Por ejemplo,
--    esSumaDeDosAbundantes 24           ==  True
--    any esSumaDeDosAbundantes [1..22]  ==  False
esSumaDeDosAbundantes :: Integer -> Bool
esSumaDeDosAbundantes n = (not . null) [x | x <- xs, n-x `elem` xs]
  where xs = takeWhile (<n) abundantes

-- abundantes es la lista de los números abundantes. Por ejemplo,
--    take 10 abundantes  ==  [12,18,20,24,30,36,40,42,48,54]
abundantes :: [Integer]
abundantes = [n | n <- [2..], abundante n]

-- (abundante n) se verifica si n es abundante. Por ejemplo,
--    abundante 12  ==  True
--    abundante 11  ==  False
abundante :: Integer -> Bool
abundante n = sum (divisores n) > n

-- (divisores n) es la lista de los divisores propios de n. Por ejemplo,
--    divisores 12  ==  [1,2,3,4,6]
divisores :: Integer -> [Integer]
divisores n = [x | x <- [1..n `div` 2], n `mod` x == 0]

-- 2ª solución
-- ===========

sumasDeDosAbundantes2 :: [Integer]
sumasDeDosAbundantes2 = filter esSumaDeDosAbundantes2 [1..]

esSumaDeDosAbundantes2 :: Integer -> Bool
esSumaDeDosAbundantes2 n = (not . null) [x | x <- xs, n-x `elem` xs]
  where xs = takeWhile (<n) abundantes2

abundantes2 :: [Integer]
abundantes2 = filter abundante2 [2..]

abundante2 :: Integer -> Bool
abundante2 n = sumaDivisores n > n

sumaDivisores :: Integer -> Integer
sumaDivisores x =
  product [(p^(e+1)-1) `div` (p-1) | (p,e) <- factorizacion x] - x

-- (factorizacion x) es la lista de las bases y exponentes de la
-- descomposición prima de x. Por ejemplo,
--    factorizacion 600  ==  [(2,3),(3,1),(5,2)]
factorizacion :: Integer -> [(Integer,Integer)]
factorizacion = map primeroYlongitud . group . primeFactors

-- (primeroYlongitud xs) es el par formado por el primer elemento de xs
-- y la longitud de xs. Por ejemplo,
--    primeroYlongitud [3,2,5,7] == (3,4)
primeroYlongitud :: [a] -> (a,Integer)
primeroYlongitud (x:xs) =
  (x, 1 + genericLength xs)

-- Verificación
-- ============

verifica :: IO ()
verifica = hspec spec

specG :: [Integer] -> Spec
specG sumasDeDosAbundantes = do
  it "e1" $
    take 10 sumasDeDosAbundantes `shouldBe` [24,30,32,36,38,40,42,44,48,50]

spec :: Spec
spec = do
  describe "def. 1" $ specG sumasDeDosAbundantes1
  describe "def. 2" $ specG sumasDeDosAbundantes2

-- La verificación es
--    λ> verifica
--    2 examples, 0 failures

-- Comprobación de equivalencia
-- ============================

-- La propiedad es
prop_sumasDeDosAbundantes :: Positive Int -> Bool
prop_sumasDeDosAbundantes (Positive n) =
  sumasDeDosAbundantes1 !! n == sumasDeDosAbundantes2 !! n

-- La comprobación es
--    λ> quickCheck prop_sumasDeDosAbundantes
--    +++ OK, passed 100 tests.

-- Comparación de eficiencia
-- =========================

-- La comparación es
--    λ> sumasDeDosAbundantes1 !! (2*10^3)
--    2887
--    (2.54 secs, 516,685,168 bytes)
--    λ> sumasDeDosAbundantes2 !! (2*10^3)
--    2887
--    (1.43 secs, 141,606,136 bytes)

from functools import reduce
from itertools import count, islice, takewhile
from operator import mul
from timeit import Timer, default_timer
from typing import Iterator

from sympy import divisor_sigma, factorint

# 1ª solución
# ===========

# divisores(n) es la lista de los divisores propios de n. Por ejemplo,
#    divisores(12)  ==  [1,2,3,4,6]
def divisores(n: int) -> list[int]:
    return [x for x in range(1, 1 + n // 2) if n % x == 0]

# abundante(n) se verifica si n es abundante. Por ejemplo,
#    abundante(12)  ==  True
#    abundante(11)  ==  False
def abundante(n: int) -> bool:
    return sum(divisores(n)) > n

# abundantes es la lista de los números abundantes. Por ejemplo,
#    >>> list(islice(abundantes(), 10))
#    [12, 18, 20, 24, 30, 36, 40, 42, 48, 54]
def abundantes() -> Iterator[int]:
    return (n for n in count(2) if abundante(n))

# esSumaDeDosAbundantes(n) se verifica si n es suma de dos números
# abundantes. Por ejemplo,
#    esSumaDeDosAbundantes(24)                           ==  True
#    any(esSumaDeDosAbundantes(n) for n in range(1, 23)) == False
def esSumaDeDosAbundantes(n: int) -> bool:
    xs = list(takewhile(lambda x: x <= n, abundantes()))
    return any(n - x in xs for x in xs)

def sumasDeDosAbundantes1() -> Iterator[int]:
    return (n for n in count(1) if esSumaDeDosAbundantes(n))

# 2ª solución
# ===========

# Si la descomposición de x en factores primos es
#    x = p(1)^e(1) . p(2)^e(2) . .... . p(n)^e(n)
# entonces la suma de los divisores de x es
#    p(1)^(e(1)+1) - 1     p(2)^(e(2)+1) - 1       p(n)^(e(2)+1) - 1
#   ------------------- . ------------------- ... -------------------
#        p(1)-1                p(2)-1                  p(n)-1
# Ver la demostración en http://bit.ly/2zUXZPc

# sumaDivisores(n) es la suma de los divisores propios de n. Por
# ejemplo,
#    sumaDivisores(12)  ==  16
def sumaDivisores(n: int) -> int:
    return reduce(mul, [(p ** (e + 1) - 1) // (p - 1)
                        for (p, e) in factorint(n).items()]) - n

def abundante2(n: int) -> bool:
    return sumaDivisores(n) > n

def abundantes2() -> Iterator[int]:
    return (n for n in count(2) if abundante2(n))

def esSumaDeDosAbundantes2(n: int) -> bool:
    xs = list(takewhile(lambda x: x <= n, abundantes2()))
    return any(n - x in xs for x in xs)

def sumasDeDosAbundantes2() -> Iterator[int]:
    return (n for n in count(1) if esSumaDeDosAbundantes2(n))

# 3ª solución
# ===========

def sumaDivisores2(n: int) -> int:
    return divisor_sigma(n, 1) - n

def abundante3(n: int) -> bool:
    return sumaDivisores2(n) > n

def abundantes3() -> Iterator[int]:
    return (n for n in count(2) if abundante3(n))

def esSumaDeDosAbundantes3(n: int) -> bool:
    xs = list(takewhile(lambda x: x <= n, abundantes3()))
    return any(n - x in xs for x in xs)

def sumasDeDosAbundantes3() -> Iterator[int]:
    return (n for n in count(1) if esSumaDeDosAbundantes3(n))

# Verificación
# ============

def test_sumasDeDosAbundantes() -> None:
    for sumasDeDosAbundantes in [sumasDeDosAbundantes1,
                                 sumasDeDosAbundantes2,
                                 sumasDeDosAbundantes3]:
        assert list(islice(sumasDeDosAbundantes(), 10)) ==\
            [24, 30, 32, 36, 38, 40, 42, 44, 48, 50]
    print("Verificado")

# La verificación es
#    >>> test_sumasDeDosAbundantes()
#    Verificado

# Comprobación de equivalencia
# ============================

# nth(i, n) es el n-ésimo elemento del iterador i. Por ejemplo,
#    nth(sumasDeDosAbundantes1(), 4) == 38
def nth(i: Iterator[int], n: int) -> int:
    return list(islice(i, n, n+1))[0]

# La propiedad es
def test_sumasDeDosAbundantes_equiv(n: int) -> bool:
    return list(islice(sumasDeDosAbundantes1(), n)) ==\
           list(islice(sumasDeDosAbundantes2(), n)) ==\
           list(islice(sumasDeDosAbundantes3(), n))

# La comprobación es
#    >>> test_sumasDeDosAbundantes_equiv(400)
#    True

# Comparación de eficiencia
# =========================

def tiempo(e: str) -> None:
    """Tiempo (en segundos) de evaluar la expresión e."""
    t = Timer(e, "", default_timer, globals()).timeit(1)
    print(f"{t:0.2f} segundos")

# La comparación es
#    >>> tiempo('nth(sumasDeDosAbundantes1(), 500)')
#    3.83 segundos
#    >>> tiempo('nth(sumasDeDosAbundantes2(), 500)')
#    2.30 segundos
#    >>> tiempo('nth(sumasDeDosAbundantes3(), 500)')
#    1.92 segundos