sumaDivisores :: Integer -> Integer

   sumaDivisores 12                 ==  28
   sumaDivisores 25                 ==  31
   sumaDivisores (product [1..25])  ==  93383273455325195473152000
   length (show (sumaDivisores (product [1..30000])))  ==  121289
   maximum (map sumaDivisores [1..2*10^6])             ==  8851392

import Data.List (foldl', genericLength, group, inits)
import Data.Set (toList)
import Data.Numbers.Primes (primeFactors)
import Math.NumberTheory.ArithmeticFunctions (divisors, sigma)
import Test.QuickCheck

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

sumaDivisores1 :: Integer -> Integer
sumaDivisores1 n = sum (divisores1 n)

-- (divisores x) es la lista de los divisores de x. Por ejemplo,
--    divisores 60  ==  [1,5,3,15,2,10,6,30,4,20,12,60]
divisores1 :: Integer -> [Integer]
divisores1 n = [x | x <- [1..n], n `rem` x == 0]

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

-- Sustituyendo la definición de divisores de la solución anterior por
-- cada una de las del ejercicio [Divisores de un número](https://bit.ly/3S1HYwi)
-- Se obtiene una nueva definición de sumaDivisores. La usada en la
-- definición anterior es la menos eficiente y la que se usa en la
-- siguiente definición es la más eficiente.

sumaDivisores2 :: Integer -> Integer
sumaDivisores2 = sum . divisores2

divisores2 :: Integer -> [Integer]
divisores2 = toList . divisors

-- 3ª solución
-- ===========

-- La solución anterior se puede simplificar

sumaDivisores3 :: Integer -> Integer
sumaDivisores3 = sum . divisors

-- 4ª solución
-- ===========

sumaDivisores4 :: Integer -> Integer
sumaDivisores4 = foldl' (+) 0 . divisores2

-- 5ª solución
-- ===========

sumaDivisores5 :: Integer -> Integer
sumaDivisores5 n = aux [1..n]
  where aux [] = 0
        aux (x:xs) | n `rem` x == 0 = x + aux xs
                   | otherwise      = aux xs

-- 6ª solución
-- ===========

sumaDivisores6 :: Integer -> Integer
sumaDivisores6 = sum
               . map (product . concat)
               . mapM inits
               . group
               . primeFactors

-- 7ª 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

sumaDivisores7 :: Integer -> Integer
sumaDivisores7 x =
  product [(p^(e+1)-1) `div` (p-1) | (p,e) <- factorizacion 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)

-- 8ª solución
-- ===========

sumaDivisores8 :: Integer -> Integer
sumaDivisores8 = sigma 1

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

-- La propiedad es
prop_sumaDivisores :: Positive Integer -> Bool
prop_sumaDivisores (Positive x) =
  all (== sumaDivisores1 x)
      [ sumaDivisores2 x
      , sumaDivisores3 x
      , sumaDivisores4 x
      , sumaDivisores5 x
      , sumaDivisores6 x
      , sumaDivisores7 x
      , sumaDivisores8 x
      ]

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

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

-- La comparación es
--    λ> sumaDivisores1 5336100
--    21386001
--    (2.25 secs, 1,067,805,248 bytes)
--    λ> sumaDivisores2 5336100
--    21386001
--    (0.01 secs, 659,112 bytes)
--    λ> sumaDivisores3 5336100
--    21386001
--    (0.01 secs, 635,688 bytes)
--    λ> sumaDivisores4 5336100
--    21386001
--    (0.01 secs, 648,992 bytes)
--    λ> sumaDivisores5 5336100
--    21386001
--    (2.44 secs, 1,323,924,176 bytes)
--    λ> sumaDivisores6 5336100
--    21386001
--    (0.01 secs, 832,104 bytes)
--    λ> sumaDivisores7 5336100
--    21386001
--    (0.01 secs, 571,040 bytes)
--    λ> sumaDivisores8 5336100
--    21386001
--    (0.00 secs, 558,296 bytes)
--
--    λ> sumaDivisores2 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (2.30 secs, 1,130,862,080 bytes)
--    λ> sumaDivisores3 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (1.83 secs, 896,386,232 bytes)
--    λ> sumaDivisores4 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (1.52 secs, 997,992,328 bytes)
--    λ> sumaDivisores6 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (2.35 secs, 5,719,848,600 bytes)
--    λ> sumaDivisores7 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (0.00 secs, 628,136 bytes)
--    λ> sumaDivisores8 251888923423315469521109880000000
--    1471072204661054993275791673480320
--    (0.00 secs, 591,352 bytes)
--
--    λ> length (show (sumaDivisores7 (product [1..30000])))
--    121289
--    (2.76 secs, 4,864,576,304 bytes)
--    λ> length (show (sumaDivisores8 (product [1..30000])))
--    121289
--    (1.65 secs, 3,173,319,312 bytes)

from operator import mul
from functools import reduce
from timeit import Timer, default_timer
from sys import setrecursionlimit
from sympy import divisors, divisor_sigma, factorint
from hypothesis import given, strategies as st

setrecursionlimit(10**6)

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

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

def sumaDivisores1(n: int) -> int:
    return sum(divisores(n))

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

# Sustituyendo la definición de divisores de la solución anterior por
# cada una de las del ejercicio Divisores de un número https://bit.ly/3S1HYwi)
# Se obtiene una nueva definición de sumaDivisores. La usada en la
# definición anterior es la menos eficiente y la que se usa en la
# siguiente definición es la más eficiente.
def sumaDivisores2(n: int) -> int:
    return sum(divisors(n))

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

def sumaDivisores3(n: int) -> int:
    def aux(xs: list[int]) -> int:
        if xs:
            if n % xs[0] == 0:
                return xs[0] + aux(xs[1:])
            return aux(xs[1:])
        return 0

    return aux(list(range(1, n + 1)))

# 4ª 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

def sumaDivisores4(n: int) -> int:
    return reduce(mul, [(p ** (e + 1) - 1) // (p - 1)
                        for (p, e) in factorint(n).items()])

# 5ª solución
# ===========

def sumaDivisores5(n: int) -> int:
    x = 1
    r1 = 0
    r2 = 0
    while x * x < n:
        if n % x == 0:
            r1 += x
            r2 += n // x
        x += 1
    if x * x == n:
        r1 += x
    return r1 + r2

# 6ª solución
# ===========

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

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

# La propiedad es
@given(st.integers(min_value=2, max_value=1000))
def test_sumaDivisores(n):
    r = sumaDivisores1(n)
    assert sumaDivisores2(n) == r
    assert sumaDivisores3(n) == r
    assert sumaDivisores4(n) == r
    assert sumaDivisores5(n) == r
    assert sumaDivisores6(n) == r

# La comprobación es
#    src> poetry run pytest -q suma_de_divisores.py
#    1 passed in 0.90s

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

def tiempo(e):
    """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('sumaDivisores1(5336100)')
#    0.29 segundos
#    >>> tiempo('sumaDivisores2(5336100)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores3(5336100)')
#    Process Python terminado (killed)
#    >>> tiempo('sumaDivisores4(5336100)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores5(5336100)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores6(5336100)')
#    0.00 segundos
#
#    >>> tiempo('sumaDivisores1(2**9 * 3**8 * 5**2)')
#    4.52 segundos
#    >>> tiempo('sumaDivisores2(2**9 * 3**8 * 5**2)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores4(2**9 * 3**8 * 5**2)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores5(2**9 * 3**8 * 5**2)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores6(2**9 * 3**8 * 5**2)')
#    0.00 segundos
#
#    >>> tiempo('sumaDivisores2(2**9 * 3**8 * 5**7 * 7**4)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores4(2**9 * 3**8 * 5**7 * 7**4)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores5(2**9 * 3**8 * 5**7 * 7**4)')
#    3.24 segundos
#    >>> tiempo('sumaDivisores6(2**9 * 3**8 * 5**7 * 7**4)')
#    0.00 segundos
#
#    >>> tiempo('sumaDivisores2(251888923423315469521109880000000)')
#    1.13 segundos
#    >>> tiempo('sumaDivisores4(251888923423315469521109880000000)')
#    0.00 segundos
#    >>> tiempo('sumaDivisores6(251888923423315469521109880000000)')
#    0.00 segundos
#
#    >>> tiempo('sumaDivisores4(reduce(mul, list(range(1, 30000))))')
#    1.89 segundos
#    >>> tiempo('sumaDivisores6(reduce(mul, list(range(1, 30000))))')
#    1.88 segundos