PFH: La semana en Exercitium (del 16 al 20 de mayo de 2022)
Esta semana he publicado en Exercitium las soluciones de los siguientes problemas:
- 1. Suma de divisores
 - 2. Sucesión de sumas de dos números abundantes
 - 3. Sumas de 4 primos
 - 4. Parejas de números y divisores
 - 5. Sumas de divisores propios
 
A continuación se muestran las soluciones.
1. Suma de divisores
Definir la función
| 
					 1  | 
						   sumaDivisores :: Integer -> Integer  | 
					
tal que (sumaDivisores x) es la suma de los divisores de x. Por ejemplo,
| 
					 1 2 3 4 5  | 
						   sumaDivisores 12  ==  28    sumaDivisores 25  ==  31    sumaDivisores (product [1..25])  ==  93383273455325195473152000    length (show (sumaDivisores (product [1..30000])))  ==  121289    maximum (map sumaDivisores [1..10^5])  ==  403200  | 
					
Soluciones
| 
					 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95  | 
						import Data.List (genericLength, group, inits) import Data.Numbers.Primes (primeFactors) import Test.QuickCheck -- 1ª solución -- =========== sumaDivisores1 :: Integer -> Integer sumaDivisores1 = sum . divisores -- (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] divisores :: Integer -> [Integer] divisores = map (product . concat)           . productoCartesiano           . map inits           . group           . primeFactors -- (productoCartesiano xss) es el producto cartesiano de los conjuntos -- xss. Por ejemplo, --    λ> producto [[1,3],[2,5],[6,4]] --    [[1,2,6],[1,2,4],[1,5,6],[1,5,4],[3,2,6],[3,2,4],[3,5,6],[3,5,4]] productoCartesiano :: [[a]] -> [[a]] productoCartesiano []       = [[]] productoCartesiano (xs:xss) =   [x:ys | x <- xs, ys <- productoCartesiano xss] -- 2ª solución -- =========== sumaDivisores2 :: Integer -> Integer sumaDivisores2 = sum                . map (product . concat)                . sequence                . map inits                . group                . primeFactors -- 3ª 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 sumaDivisores3 :: Integer -> Integer sumaDivisores3 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) -- Comprobación de equivalencia -- ============================ -- La propiedad es prop_sumaDivisores :: Positive Integer -> Bool prop_sumaDivisores (Positive x) =   all (== sumaDivisores1 x)       [ sumaDivisores2 x       , sumaDivisores3 x       ] -- La comprobación es --    λ> quickCheck prop_sumaDivisores --    +++ OK, passed 100 tests. -- Comparación de eficiencia -- ========================= -- La comparación es --   λ> sumaDivisores1 251888923423315469521109880000000 --   1471072204661054993275791673480320 --   (10.63 secs, 10,614,618,080 bytes) --   λ> sumaDivisores2 251888923423315469521109880000000 --   1471072204661054993275791673480320 --   (2.51 secs, 5,719,399,056 bytes) --   λ> sumaDivisores3 251888923423315469521109880000000 --   1471072204661054993275791673480320 --   (0.01 secs, 177,480 bytes)  | 
					
El código se encuentra en GitHub.
2. Sucesión de sumas de dos números abundantes
Un número n es abundante si la suma de los divisores propios de n es mayor que n. El primer número abundante es el 12 (cuyos divisores propios son 1, 2, 3, 4 y 6 cuya suma es 16). Por tanto, el menor número que es la suma de dos números abundantes es el 24.
Definir la sucesión
| 
					 1  | 
						   sumasDeDosAbundantes :: [Integer]  | 
					
cuyos elementos son los números que se pueden escribir como suma de dos números abundantes. Por ejemplo,
| 
					 1  | 
						   take 10 sumasDeDosAbundantes  ==  [24,30,32,36,38,40,42,44,48,50]  | 
					
Soluciones
| 
					 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89  | 
						import Data.List (genericLength, group) import Data.Numbers.Primes (primeFactors) 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) -- 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)  | 
					
El código se encuentra en GitHub.
3. Sumas de 4 primos
La conjetura de Waring sobre los números primos establece que todo número impar es primo o la suma de tres primos. La conjetura de Goldbach afirma que todo par mayor que 2 es la suma de dos números primos. Ambos ha estado abiertos durante más de 200 años. En este problema no se propone su solución, sino una tarea más simple: buscar una manera de expresar los enteros mayores que 7 como suma de exactamente cuatro números primos; es decir, definir la función
| 
					 1  | 
						   suma4primos :: Integer -> [(Integer,Integer,Integer,Integer)]  | 
					
tal que (suma4primos n) es la lista de las cuádruplas crecientes (a,b,c,d) de números primos cuya suma es n (que se supone mayor que 7). Por ejemplo,
| 
					 1 2  | 
						   suma4primos 18             == [(2,2,3,11),(2,2,7,7),(3,3,5,7),(3,5,5,5)]    head (suma4primos (10^14)) == (2,2,23,99999999999973)  | 
					
Comprobar con QuickCheck que todo entero mayor que 7 se puede escribir como suma de exactamente cuatro números primos.
Soluciones
| 
					 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73  | 
						import Data.Numbers.Primes (isPrime, primes) import Test.QuickCheck -- 1ª solución -- =========== suma4primos1 :: Integer -> [(Integer, Integer, Integer, Integer)] suma4primos1 n =   [(a,b,c,d) | let as = takeWhile (< n) primes,                a <- as,                let bs = takeWhile (< n-a) as,                b <- bs, a <= b,                let cs = takeWhile (< n-a-b) bs,                c <- cs, b <= c,                let d = n-a-b-c, c <= d,                isPrime d] -- 2ª solución -- =========== suma4primos2 :: Integer -> [(Integer, Integer, Integer, Integer)] suma4primos2 n =   [(a,b,c,d) | let as = takeWhile (< n) primes,                a <- as,                let bs = takeWhile (< n-a) (dropWhile (< a) as),                b <- bs,                let cs = takeWhile (<n-a-b) (dropWhile (< b) bs),                c <- cs,                let d = n-a-b-c,                c <= d,                isPrime d] -- Comprobación de equivalencia -- ============================ -- La propiedad es prop_suma4primos :: Positive Integer -> Bool prop_suma4primos (Positive n) =   suma4primos1 n == suma4primos2 n -- La comprobación es --    λ> quickCheck prop_suma4primos --    +++ OK, passed 100 tests; 526 discarded. -- Comparación de eficiencia -- ========================= -- La comparación es --    λ> length (suma4primos1 2000) --    90219 --    (2.98 secs, 4,517,620,744 bytes) --    λ> length (suma4primos2 2000) --    90219 --    (2.22 secs, 4,223,251,928 bytes) -- --    λ> head (suma4primos1 (10^14)) --    (2,2,23,99999999999973) --    (1.67 secs, 5,963,327,168 bytes) --    λ> head (suma4primos2 (10^14)) --    (2,2,23,99999999999973) --    (1.70 secs, 5,963,326,848 bytes) -- Comprobación de la propiedad -- ============================ -- La propiedad es prop_suma4primos2 :: Integer -> Property prop_suma4primos2 n =   n > 7 ==> not (null (suma4primos1 n)) -- La comprobación es --    λ> quickCheck prop_suma4primos2 --    +++ OK, passed 100 tests; 582 discarded.  | 
					
El código se encuentra en GitHub.
4. Parejas de números y divisores
Definir la función
| 
					 1  | 
						   divisoresHasta :: Int -> [(Int,Int)]  | 
					
tal que (divisoresHasta n) es la lista de los pares (a,b) tales que a es un número entre 2 y n y b es un divisor propio de a. Por ejemplo,
| 
					 1 2 3 4 5 6  | 
						   λ> divisoresHasta 6    [(2,1),(3,1),(4,1),(5,1),(6,1),(4,2),(6,2),(6,3)]    λ> divisoresHasta 8    [(2,1),(3,1),(4,1),(5,1),(6,1),(7,1),(8,1),(4,2),(6,2),(8,2),(6,3),(8,4)]    λ> length (divisoresHasta 1234567)    16272448  | 
					
Soluciones
| 
					 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73  | 
						import Data.List (sort, sortBy) import Data.Ord (comparing) import Data.Tuple (swap) import Test.QuickCheck -- 1ª solución -- =========== divisoresHasta1 :: Int -> [(Int,Int)] divisoresHasta1 n =   ordena (concat [[(a,b) | b <- divisoresPropios a] | a <- [2..n]])   where ordena ps = [intercambia p | p <- sort (map intercambia ps)]         intercambia (x,y) = (y,x) divisoresPropios :: Int -> [Int] divisoresPropios n = [x | x <- [1..n `div` 2], n `mod` x == 0] -- 2ª solución -- =========== divisoresHasta2 :: Int -> [(Int,Int)] divisoresHasta2 n =   ordena (concat [[(a,b) | b <- divisoresPropios a] | a <- [2..n]])   where ordena           = sortBy comp         comp (x,y) (u,v) = compare (y,x) (v,u) -- 3ª solución -- =========== divisoresHasta3 :: Int -> [(Int,Int)] divisoresHasta3 n =   ordena (concat [[(a,b) | b <- divisoresPropios a] | a <- [2..n]])   where ordena = sortBy (comparing swap) -- 4ª solución -- =========== divisoresHasta4 :: Int -> [(Int,Int)] divisoresHasta4 n =   [(a,b) | b <- [1..n `div` 2], a <- [2*b, 3*b..n]] -- Comprobación de equivalencia -- ============================ -- La propiedad es prop_divisoresHasta :: Positive Int -> Bool prop_divisoresHasta (Positive n) =   all (== divisoresHasta1 n)       [ divisoresHasta2 n       , divisoresHasta3 n       , divisoresHasta4 n       ] -- La comprobación es --    λ> quickCheck prop_divisoresHasta --    +++ OK, passed 100 tests. -- Comparación de eficiencia -- ========================= -- La comparación es --    λ> length (divisoresHasta1 4000) --    29805 --    (1.78 secs, 843,584,472 bytes) --    λ> length (divisoresHasta2 4000) --    29805 --    (1.78 secs, 879,421,056 bytes) --    λ> length (divisoresHasta3 4000) --    29805 --    (1.70 secs, 876,475,584 bytes) --    λ> length (divisoresHasta4 4000) --    29805 --    (0.02 secs, 6,332,016 bytes)  | 
					
El código se encuentra en GitHub.
5. Sumas de divisores propios
Definir la función
| 
					 1  | 
						   sumaDivisoresHasta :: Integer -> [(Integer,Integer)]  | 
					
tal que (sumaDivisoresHasta n) es la lista de los pares (a,b) tales que a es un número entre 1 y n y b es la suma de los divisores propios de a. Por ejemplo,
| 
					 1 2 3 4  | 
						   λ> sumaDivisoresHasta 12    [(1,0),(2,1),(3,1),(4,3),(5,1),(6,6),(7,1),(8,7),(9,4),(10,8),(11,1),(12,16)]    λ> last (sumaDivisoresHasta2 (10^7))    (10000000,14902280)  | 
					
Soluciones
| 
					 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87  | 
						import Data.Array (accumArray, assocs) import Data.List (genericLength, group) import Data.Numbers.Primes (primeFactors) import Test.QuickCheck -- 1ª solución -- =========== sumaDivisoresHasta1 :: Integer -> [(Integer,Integer)] sumaDivisoresHasta1 n = [(x, sum (divisores x)) | x <- [1..n]] divisores :: Integer -> [Integer] divisores n = [x | x <- [1..n `div` 2], n `mod` x == 0] -- 2ª solución -- =========== sumaDivisoresHasta2 :: Integer -> [(Integer,Integer)] sumaDivisoresHasta2 n = [(x, sumaDivisores x) | x <- [1..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) -- 3ª solución -- =========== sumaDivisoresHasta3 :: Integer -> [(Integer,Integer)] sumaDivisoresHasta3 n = assocs (accumArray (+) 0 (1,n) (divisoresHasta n)) -- (divisoresHasta n) es la lista de los pares (a,b) tales que a es -- un número entre 2 y n y b es un divisor propio e x. Por ejemplo, --    λ> divisoresHasta 6 --    [(2,1),(3,1),(4,1),(5,1),(6,1),(4,2),(6,2),(6,3)] --    λ> divisoresHasta 8 --    [(2,1),(3,1),(4,1),(5,1),(6,1),(7,1),(8,1),(4,2),(6,2),(8,2),(6,3),(8,4)] divisoresHasta :: Integer -> [(Integer,Integer)] divisoresHasta n = [(a,b) | b <- [1..n `div` 2], a <- [b*2, b*3..n]] -- Comprobación de equivalencia -- ============================ -- La propiedad es prop_sumaDivisoresHasta :: Positive Integer -> Bool prop_sumaDivisoresHasta (Positive n) =   all (== sumaDivisoresHasta1 n)       [ sumaDivisoresHasta2 n       , sumaDivisoresHasta3 n       ] -- La comprobación es --    λ> quickCheck prop_sumaDivisoresHasta --    +++ OK, passed 100 tests. -- Comparación de eficiencia -- ========================= -- La comparación es --    λ> last (sumaDivisoresHasta1 (10^6)) --    (1000000,1480437) --    (0.47 secs, 308,542,392 bytes) --    λ> last (sumaDivisoresHasta2 (10^6)) --    (1000000,2480437) --    (0.26 secs, 208,548,944 bytes) --    λ> last (sumaDivisoresHasta3 (10^6)) --    (1000000,1480437) --    (6.65 secs, 3,249,831,856 bytes) -- --    λ> last (sumaDivisoresHasta1 (5*10^6)) --    (5000000,7402312) --    (2.27 secs, 1,540,543,352 bytes) --    λ> last (sumaDivisoresHasta2 (5*10^6)) --    (5000000,12402312) --    (1.19 secs, 1,040,549,800 bytes)  | 
					
El código se encuentra en GitHub.