Menu Close

Etiqueta: primeFactors

Teorema de Liouville sobre listas CuCu

Una lista CuCu es una lista de números enteros positivos tales que la suma de sus Cubos es igual al Cuadrado de su suma. Por ejemplo, [1, 2, 3, 2, 4, 6] es una lista CuCu ya que

   1³ + 2³ + 3³ + 2³ + 4³ + 6³ = (1 + 2 + 3 + 2 + 4 + 6)²

La lista de Liouville correspondiente al número entero positivo n es la lista formada por el número de divisores de cada divisor de n. Por ejemplo, para el número 20 se tiene que sus divisores son

   1, 2, 4, 5, 10, 20

puesto que el número de sus divisores es

  • El 1 tiene 1 divisor (el 1 solamente).
  • El 2 tiene 2 divisores (el 1 y el 2).
  • El 4 tiene 3 divisores (el 1, el 2 y el 4).
  • El 5 tiene 2 divisores (el 1 y el 5).
  • El 10 tiene 4 divisores (el 1, el 2, el 5 y el 10).
  • El 20 tiene 6 divisores (el 1, el 2, el 4, el 5, el 10 y el 20).

la lista de Liouville de 20 es [1, 2, 3, 2, 4, 6] que, como se comentó anteriormente, es una lista CuCu.

El teorema de Lioville afirma que todas las lista de Lioville son CuCu.

Definir las funciones

   esCuCu :: [Integer] -> Bool
   liouville :: Integer -> [Integer]

tales que

  • (esCuCu xs) se verifica si la lista xs es CuCu; es decir, la suma de los cubos de sus elementos es igual al cuadrado de su suma. Por ejemplo,
     esCuCu [1,2,3]        ==  True
     esCuCu [1,2,3,2]      ==  False
     esCuCu [1,2,3,2,4,6]  ==  True
  • (liouville n) es la lista de Lioville correspondiente al número n. Por ejemplo,
     liouville 20  ==  [1,2,3,2,4,6]
     liouville 60  ==  [1,2,2,3,2,4,4,6,4,6,8,12]
     length (liouville (product [1..25]))  ==  340032

Comprobar con QuickCheck

  • que para todo entero positivo n, (liouville (2^n)) es la lista [1,2,3,…,n+1] y
  • el teorema de Lioville; es decir, para todo entero positivo n, (liouville n) es una lista CuCu.

Nota: Este ejercicio está basado en Cómo generar conjuntos CuCu de Gaussianos.

Soluciones

import Data.List (genericLength, group, inits, sort)
import Data.Numbers.Primes (primeFactors)
import Test.QuickCheck
 
esCuCu :: [Integer] -> Bool
esCuCu xs = sum (map (^3) xs) == (sum xs)^2
 
-- 1ª definición de liouville
-- ==========================
 
liouville :: Integer -> [Integer]
liouville n = map numeroDivisores (divisores n)
 
-- (divisores x) es el conjunto de divisores de los x. Por ejemplo, 
--   divisores 30  ==  [1,2,3,5,6,10,15,30]
divisores :: Integer -> [Integer]
divisores n = [x | x <- [1..n], n `mod` x == 0]
 
-- (numeroDivisores x) es el número de divisores de x. Por ejemplo, 
--    numeroDivisores 12  ==  6
--    numeroDivisores 25  ==  3
numeroDivisores :: Integer -> Integer
numeroDivisores n = genericLength (divisores n) 
 
  -- 2ª definición de liouville
-- ============================
 
liouville2 :: Integer -> [Integer]
liouville2 n = map numeroDivisores2 (divisores2 n)
 
-- Se usan las funciones
-- + divisores de "Conjunto de divisores" http://bit.ly/2OtbFIj
-- + numeroDivisores de "Número de divisores" http://bit.ly/2DgVh74
 
-- (divisores2 x) es el conjunto de divisores de los x. Por ejemplo, 
--   divisores2 30  ==  [1,2,3,5,6,10,15,30]
divisores2 :: Integer -> [Integer]
divisores2 = sort
           . map (product . concat)
           . sequence
           . map inits
           . group
           . primeFactors
 
-- (numeroDivisores2 x) es el número de divisores de x. Por ejemplo, 
--    numeroDivisores2 12  ==  6
--    numeroDivisores2 25  ==  3
numeroDivisores2 :: Integer -> Integer
numeroDivisores2 =
  product . map ((+1) . genericLength) . group . primeFactors
 
-- Comparación de eficiencia
-- =========================
 
-- La comparación es
--    λ> length (liouville (product [1..11]))
--    540
--    (13.66 secs, 7,983,550,640 bytes)
--    λ> length (liouville2 (product [1..11]))
--    540
--    (0.01 secs, 1,255,328 bytes)
 
-- Propiedad
-- =========
 
-- La propiedad es
prop_Liouville :: Integer -> Property
prop_Liouville n =
  n > 0 ==> liouville2 (2^n) == [1..n+1]
 
-- La comprobación es
--    λ> quickCheck prop_Liouville
--    +++ OK, passed 100 tests.
 
-- Teorema de Liouville
-- ====================
 
-- La propiedad es
teorema_Liouville :: Integer -> Property
teorema_Liouville n =
  n > 0 ==> esCuCu (liouville n)
 
-- La comprobación es
--    λ> quickCheck teorema_Liouville
--    +++ OK, passed 100 tests.

Pensamiento

¡Oh, tarde viva y quieta
que opuso al panta rhei su nada corre.

Antonio Machado

Conjetura de Grimm

La conjetura de Grimm establece que a cada elemento de un conjunto de números compuestos consecutivos se puede asignar un número primo que lo divide, de forma que cada uno de los números primos elegidos es distinto de todos los demás. Más formalmente, si n+1, n+2, …, n+k son números compuestos, entonces existen números primos p(i), distintos entre sí, tales que p(i) divide a n+i para 1 ≤ i ≤ k.

Diremos que la lista ps = [p(1),…,p(k)] es una sucesión de Grim para la lista xs = [x(1),…,x(k)] si p(i) son números primos distintos y p(i) divide a x(i), para 1 ≤ i ≤ k. Por ejemplo, 2, 5, 13, 3, 7 es una sucesión de Grim de 24, 25, 26, 27, 28.

Definir las funciones

   compuestos :: Integer -> [Integer]
   sucesionesDeGrim :: [Integer] -> [[Integer]]

tales que

  • (compuestos n) es la mayor lista de números enteros consecutivos empezando en n. Por ejemplo,
     compuestos 24  ==  [24,25,26,27,28]
     compuestos  8  ==  [8,9,10]
     compuestos 15  ==  [15,16]
     compuestos 16  ==  [16]
     compuestos 17  ==  []
  • (sucesionesDeGrim xs) es la lista de las sucesiones de Grim de xs. Por ejemplo,
     sucesionesDeGrim [15,16]          == [[3,2],[5,2]]
     sucesionesDeGrim [8,9,10]         == [[2,3,5]]
     sucesionesDeGrim [9,10]           == [[3,2],[3,5]]
     sucesionesDeGrim [24,25,26,27,28] == [[2,5,13,3,7]]
     sucesionesDeGrim [25,26,27,28]    == [[5,2,3,7],[5,13,3,2],[5,13,3,7]]

Comprobar con QuickCheck la conjetura de Grim; es decir, para todo número n > 1, (sucesionesDeGrim (compuestos n)) es una lista no vacía.

Soluciones

import Data.List (nub)
import Data.Numbers.Primes (isPrime, primeFactors)
import Test.QuickCheck
 
compuestos :: Integer -> [Integer]
compuestos n = takeWhile (not . isPrime) [n..]
 
sucesionesDeGrim :: [Integer] -> [[Integer]]
sucesionesDeGrim [] = [[]]
sucesionesDeGrim (x:xs) =
  [y:ys | y <- divisoresPrimos x
        , ys <- sucesionesDeGrim xs
        , y `notElem` ys]
 
-- (divisoresPrimos n) es la lista de los divisores primos de n. Por
-- ejemplo, 
--    divisoresPrimos 60  ==  [2,3,5]
divisoresPrimos :: Integer -> [Integer]
divisoresPrimos = nub . primeFactors
 
-- La propiedad es
conjeturaDeGrim :: Integer -> Property
conjeturaDeGrim n =
  n > 1 ==> not (null (sucesionesDeGrim (compuestos n))) 
 
-- La comprobación es
--    λ> quickCheck conjeturaDeGrim
--    +++ OK, passed 100 tests.

Pensamiento

De encinar en encinar
se va fatigando el día.

Antonio Machado

Teorema de Carmichael

La sucesión de Fibonacci, F(n), es la siguiente sucesión infinita de números naturales:

   0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, ...

La sucesión comieanza con los números 0 y 1. A partir de estos, cada término es la suma de los dos anteriores.

El teorema de Carmichael establece que para todo n mayor que 12, el n-ésimo número de Fibonacci F(n) tiene al menos un factor primo que no es factor de ninguno de los términos anteriores de la sucesión.

Si un número primo p es un factor de F(n) y no es factor de ningún F(m) con m < n, entonces se dice que p es un factor característico o un divisor primitivo de F(n).

Definir la función

   factoresCaracteristicos :: Int -> [Integer]

tal que (factoresCaracteristicos n) es la lista de los factores característicos de F(n). Por ejemplo,

   factoresCaracteristicos  4  ==  [3]
   factoresCaracteristicos  6  ==  []
   factoresCaracteristicos 19  ==  [37,113]
   factoresCaracteristicos 20  ==  [41]
   factoresCaracteristicos 37  ==  [73,149,2221]

Comprobar con QuickCheck el teorema de Carmichael; es decir, para todo número entero (factoresCaracteristicos (13 + abs n)) es una lista no vacía.

Soluciones

import Data.List (nub)
import Data.Numbers.Primes
import Test.QuickCheck
 
factoresCaracteristicos :: Int -> [Integer]
factoresCaracteristicos n =
  [x | x <- factoresPrimos (fib n)
     , and [fib m `mod` x /= 0 | m <- [1..n-1]]]
 
-- (fib n) es el n-ésimo término de la sucesión de Fibonacci. Por
-- ejemplo,
--    fib 6  ==  8
fib :: Int -> Integer
fib n = fibs !! n
 
-- fibs es la lista de términos de la sucesión de Fibonacci. Por ejemplo,
--    λ> take 20 fibs
--    [0,1,1,2,3,5,8,13,21,34,55,89,144,233,377,610,987,1597,2584,4181]
fibs :: [Integer]
fibs = 0 : 1 : zipWith (+) fibs (tail fibs)
 
-- (factoresPrimos n) es la lista de los factores primos de n. Por
-- ejemplo, 
--    factoresPrimos 600  ==  [2,3,5]
factoresPrimos :: Integer -> [Integer]
factoresPrimos 0 = []
factoresPrimos n = nub (primeFactors n)
 
-- Teorema
-- =======
 
-- El teorema es
teorema_de_Carmichael :: Int -> Bool
teorema_de_Carmichael n =
  not (null (factoresCaracteristicos n'))
  where n' = 13 + abs n
 
-- La comprobación es
--    λ> quickCheckWith (stdArgs {maxSize=50}) teorema_de_Carmichael
--    +++ OK, passed 100 tests.

Pensamiento

No puede ser
amor de tanta fortuna:
dos soledades en una.

Antonio Machado

Derivada aritmética

La derivada aritmética es una función definida sobre los números naturales por analogía con la regla del producto para el cálculo de las derivadas usada en análisis.

Para un número natural n su derivada D(n) se define por

   D(0)  = 0
   D(1)  = 0
   D(p)  = 1, si p es primo
   D(ab) = D(a)b + aD(b) (regla de Leibniz para derivar productos)

Por ejemplo,

   D(6)  = D(2*3) = D(2)*3 + 2*D(3) = 1*3 + 2*1 =  5
   D(12) = D(2*6) = D(2)*6 + 2*D(6) = 1*6 + 2*5 = 16

Definir la función

   derivada :: Integer -> Integer

tal que (derivada n) es la derivada aritmética de n. Por ejemplo,

   derivada  6  ==  5
   derivada 12  ==  16
   maximum [derivada n | n <- [1..60000]]  ==  380928

Comprobar con QuickCheck que si x es un número entero positivo y su descomposición en factores primos es

   x = p(1)^e(1) + p(2)^e(2) +...+ p(n)^e(n)

entonces la derivada de x es

   x * [e(1)/p(1) + e(2)/p(2) +...+ e(n)/p(n)]

Nota: No usar en la definición la propiedad que hay que comprobar.

Soluciones

import Data.List (genericLength, group)
import Data.Numbers.Primes (isPrime, primeFactors)
import Test.QuickCheck
 
-- 1ª solución
-- ===========
 
derivada :: Integer -> Integer
derivada 0 = 0
derivada 1 = 0
derivada n | esPrimo n = 1
           | otherwise = (derivada a) * b + a * (derivada b)
  where a = menorFactor n
        b = n `div` a
 
-- (esPrimo n) se verifica si n es primo. Por ejemplo,
--    esPrimo 5  ==  True
--    esPrimo 6  ==  False
esPrimo :: Integer -> Bool
esPrimo 0 = False
esPrimo 1 = False
esPrimo n = n == menorFactor n
 
-- (menorFactor n) es el menor divisor primo de n (con n >= 2). Por
-- ejemplo, 
--    menorFactor 6   ==  2
--    menorFactor 7   ==  7
--    menorFactor 15  ==  3
menorFactor :: Integer -> Integer
menorFactor n
  | even n = 2
  | otherwise = head [x | x <- [3,5..]
                        , n `mod` x == 0]
 
-- 2ª solución
-- ===========
 
derivada2 :: Integer -> Integer
derivada2 0 = 0
derivada2 1 = 0
derivada2 n | isPrime n = 1
            | otherwise = (derivada2 a) * b + a * (derivada2 b)
  where (a:_) = primeFactors n
        b     = n `div` a
 
-- Comparación de eficiencia
-- =========================
 
--    λ> maximum [derivada n | n <- [1..10000]]
--    53248
--    (1.59 secs, 1,091,452,552 bytes)
--    λ> maximum [derivada2 n | n <- [1..10000]]
--    53248
--    (0.17 secs, 457,819,120 bytes)
 
-- Propiedad
-- =========
 
-- La propiedad es
prop_derivada :: Integer -> Property
prop_derivada x =
  x > 0 ==>
  derivada x == sum [(x * e) `div` p | (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 n =
  [(head xs,genericLength xs) | xs <- group (primeFactors n)]
 
-- Su comprobación es
--    λ> quickCheck prop_derivada
--    +++ OK, passed 100 tests.

Referencias

Pensamiento

En ese jardín, Guiomar,
el mutuo jardín que inventan
dos corazones al par,
se funden y complementan
nuestras horas.

Antonio Machado

Árbol binario de divisores

El árbol binario de los divisores de 24 es

    90
    /\
   2  45
      /\
     3  15
        /\
       3  5

Se puede representar por

   N 90 (H 2) (N 45 (H 3) (N 15 (H 3) (H 5)))

usando el tipo de dato definido por

   data Arbol = H Int
              | N Int Arbol Arbol
     deriving (Eq, Show)

Análogamente se obtiene el árbol binario de cualquier número x: se comienza en x y en cada paso se tiene dos hijos (su menor divisor y su cociente) hasta obtener números primos en las hojas.

Definir las funciones

   arbolDivisores      :: Int -> Arbol
   hojasArbolDivisores :: Int -> [Int]

tales que

  • (arbolDivisores x) es el árbol binario de los divisores de x. Por ejemplo,
     λ> arbolDivisores 90
     N 90 (H 2) (N 45 (H 3) (N 15 (H 3) (H 5)))
     λ> arbolDivisores 24
     N 24 (H 2) (N 12 (H 2) (N 6 (H 2) (H 3)))
     λ> arbolDivisores 300
     N 300 (H 2) (N 150 (H 2) (N 75 (H 3) (N 25 (H 5) (H 5))))
  • (hojasArbolDivisores x) es la lista de las hohas del árbol binario de los divisores de x. Por ejemplo
     hojasArbolDivisores 90   ==  [2,3,3,5]
     hojasArbolDivisores 24   ==  [2,2,2,3]
     hojasArbolDivisores 300  ==  [2,2,3,5,5]

Soluciones

import Data.Numbers.Primes (isPrime, primeFactors)
 
data Arbol = H Int
           | N Int Arbol Arbol
  deriving (Eq, Show)
 
-- 1ª solución
-- ===========
 
arbolDivisores :: Int -> Arbol
arbolDivisores x
  | y == x    = H x
  | otherwise = N x (H y) (arbolDivisores (x `div` y))
  where y = menorDivisor x
 
-- (menorDivisor x) es el menor divisor primo de x. Por ejemplo,
--    menorDivisor 45  ==  3
--    menorDivisor 5   ==  5
menorDivisor :: Int -> Int
menorDivisor x =
  head [y | y <- [2..x], x `mod` y == 0]
 
hojasArbolDivisores :: Int -> [Int]
hojasArbolDivisores = hojas . arbolDivisores
 
-- (hojas a) es la lista de las hojas del árbol a. Por ejemplo,
--    hojas (N 3 (H 4) (N 5 (H 7) (H 2)))  ==  [4,7,2]
hojas :: Arbol -> [Int]
hojas (H x)     = [x]
hojas (N _ i d) = hojas i ++ hojas d
 
-- 2ª solución
-- ===========
 
arbolDivisores2 :: Int -> Arbol
arbolDivisores2 x
  | y == x    = H x
  | otherwise = N x (H y) (arbolDivisores (x `div` y))
  where (y:_) = primeFactors x
 
hojasArbolDivisores2 :: Int -> [Int]
hojasArbolDivisores2 = primeFactors

Pensamiento

Cuando el Ser que se es hizo la nada
y reposó que bien lo merecía,
ya tuvo el día noche, y compañía
tuvo el amante en la ausencia de la amada.

Antonio Machado