Test.QuickCheck

Posiciones de conjuntos finitos de naturales

PorJosé A. Alonso 27-diciembre-20193-enero-2020

En un ejercicio anterior se mostró que los conjuntos finitos de números naturales se pueden enumerar como sigue

    0: []
    1: [0]
    2: [1]
    3: [1,0]
    4: [2]
    5: [2,0]
    6: [2,1]
    7: [2,1,0]
    8: [3]
    9: [3,0]
   10: [3,1]
   11: [3,1,0]
   12: [3,2]
   13: [3,2,0]
   14: [3,2,1]
   15: [3,2,1,0]
   16: [4]
   17: [4,0]
   18: [4,1]
   19: [4,1,0]

0: []

1: [0]

2: [1]

3: [1,0]

4: [2]

5: [2,0]

6: [2,1]

7: [2,1,0]

8: [3]

9: [3,0]

10: [3,1]

11: [3,1,0]

12: [3,2]

13: [3,2,0]

14: [3,2,1]

15: [3,2,1,0]

16: [4]

17: [4,0]

18: [4,1]

19: [4,1,0]

en la que los elementos están ordenados de manera decreciente.

Además, se definió la constante

   enumeracionCFN :: [[Integer]]

1	enumeracionCFN :: [[Integer]]

tal que sus elementos son los conjuntos de los números naturales con la ordenación descrita anteriormente. Por ejemplo,

   λ> take 20 enumeracionCFN
   [[],
    [0],
    [1],[1,0],
    [2],[2,0],[2,1],[2,1,0],
    [3],[3,0],[3,1],[3,1,0],[3,2],[3,2,0],[3,2,1],[3,2,1,0],
    [4],[4,0],[4,1],[4,1,0]]

λ> take 20 enumeracionCFN

[[],

[0],

[1],[1,0],

[2],[2,0],[2,1],[2,1,0],

[3],[3,0],[3,1],[3,1,0],[3,2],[3,2,0],[3,2,1],[3,2,1,0],

[4],[4,0],[4,1],[4,1,0]]

Definir la función

   posicion :: [Integer] -> Integer

1	posicion :: [Integer] -> Integer

tal que (posicion xs) es la posición del conjunto finito de números naturales xs, representado por una lista decreciente, en enumeracionCFN. Por ejemplo,

   posicion [2,0]          ==  5
   posicion [2,1]          ==  6
   posicion [2,1,0]        ==  7
   posicion [0]            ==  1
   posicion [1,0]          ==  3
   posicion [2,1,0]        ==  7
   posicion [3,2,1,0]      ==  15
   posicion [4,3,2,1,0]    ==  31
   posicion [5,4,3,2,1,0]  ==  63
   length (show (posicion [3*10^7])) == 9030900

posicion [2,0] == 5

posicion [2,1] == 6

posicion [2,1,0] == 7

posicion [0] == 1

posicion [1,0] == 3

posicion [2,1,0] == 7

posicion [3,2,1,0] == 15

posicion [4,3,2,1,0] == 31

posicion [5,4,3,2,1,0] == 63

length (show (posicion [3*10^7])) == 9030900

Comprobar con QuickCheck que para todo número natural n,

   posicion [n,n-1..0] == 2^(n+1) - 1.

1	posicion [n,n-1..0] == 2^(n+1) - 1.

Soluciones

import Data.List (genericLength, nub, sort)
import Test.QuickCheck

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

posicion :: [Integer] -> Integer
posicion xs =
  genericLength (takeWhile (< xs) enumeracionCFN)

enumeracionCFN :: [[Integer]]
enumeracionCFN = concatMap enumeracionCFNHasta [0..]

-- (enumeracionCFNHasta n) es la lista de conjuntos con la enumeración
-- anterior cuyo primer elemento es n. Por ejemplo,
--    λ> enumeracionCFNHasta 1
--    [[1],[1,0]]
--    λ> enumeracionCFNHasta 2
--    [[2],[2,0],[2,1],[2,1,0]]
--    λ> enumeracionCFNHasta 3
--    [[3],[3,0],[3,1],[3,1,0],[3,2],[3,2,0],[3,2,1],[3,2,1,0]]
enumeracionCFNHasta :: Integer -> [[Integer]]
enumeracionCFNHasta 0 = [[],[0]]
enumeracionCFNHasta n =
  [n:xs | k <- [0..n-1], xs <- enumeracionCFNHasta k]

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

posicion2 :: [Integer] -> Integer
posicion2 []     = 0
posicion2 (x:xs) = 2^x + posicion2 xs

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

posicion3 :: [Integer] -> Integer
posicion3 = foldr (\x -> (+) (2^x)) 0 

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

posicion4 :: [Integer] -> Integer
posicion4 xs = sum [2^x | x <- xs ]

-- Equivalencia de las definiciones
-- ================================

-- La propiedad es
prop_equiv :: [Integer] -> Bool
prop_equiv xs =
  all (== posicion xs') [f xs' | f <- [ posicion2
                                      , posicion3
                                      , posicion4]]
  where xs' = reverse (sort (nub (map abs xs)))

-- La comprobación es
--    λ> quickCheckWith (stdArgs {maxSize=15}) prop_equiv
--    +++ OK, passed 100 tests.

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

-- La comparación es
--    λ> posicion [19,18,17,16,14,9,6]
--    1000000
--    (2.61 secs, 1,754,265,776 bytes)
--    λ> posicion2 [19,18,17,16,14,9,6]
--    1000000
--    (0.01 secs, 111,808 bytes)
--    λ> posicion3 [19,18,17,16,14,9,6]
--    1000000
--    λ> posicion4 [19,18,17,16,14,9,6]
--    1000000
--    (0.01 secs, 111,704 bytes)

--    λ> length (show (posicion2 [3*10^7]))
--    9030900
--    (2.06 secs, 571,911,304 bytes)
--    λ> length (show (posicion3 [3*10^7]))
--    9030900
--    (2.06 secs, 571,911,544 bytes)
--    λ> length (show (posicion4 [3*10^7]))
--    9030900
--    (2.05 secs, 571,911,464 bytes)

-- Propiedad
-- =========

-- La propiedad es
prop_posicion :: Integer -> Property
prop_posicion n =
  n >= 0 ==> posicion3 [n,n-1..0] == 2^(n+1) - 1

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

import Data.List (genericLength, nub, sort)

import Test.QuickCheck

-- 1ª solución

-- ===========

posicion :: [Integer] -> Integer

posicion xs =

genericLength (takeWhile (< xs) enumeracionCFN)

enumeracionCFN :: [[Integer]]

enumeracionCFN = concatMap enumeracionCFNHasta [0..]

-- (enumeracionCFNHasta n) es la lista de conjuntos con la enumeración

-- anterior cuyo primer elemento es n. Por ejemplo,

-- λ> enumeracionCFNHasta 1

-- [[1],[1,0]]

-- λ> enumeracionCFNHasta 2

-- [[2],[2,0],[2,1],[2,1,0]]

-- λ> enumeracionCFNHasta 3

-- [[3],[3,0],[3,1],[3,1,0],[3,2],[3,2,0],[3,2,1],[3,2,1,0]]

enumeracionCFNHasta :: Integer -> [[Integer]]

enumeracionCFNHasta 0 = [[],[0]]

enumeracionCFNHasta n =

[n:xs | k <- [0..n-1], xs <- enumeracionCFNHasta k]

-- 2ª solución

-- ===========

posicion2 :: [Integer] -> Integer

posicion2 [] = 0

posicion2 (x:xs) = 2^x + posicion2 xs

-- 3ª solución

-- ===========

posicion3 :: [Integer] -> Integer

posicion3 = foldr (\x -> (+) (2^x)) 0

-- 4ª solución

-- ===========

posicion4 :: [Integer] -> Integer

posicion4 xs = sum [2^x | x <- xs ]

-- Equivalencia de las definiciones

-- ================================

-- La propiedad es

prop_equiv :: [Integer] -> Bool

prop_equiv xs =

all (== posicion xs') [f xs' | f <- [ posicion2

, posicion3

, posicion4]]

where xs' = reverse (sort (nub (map abs xs)))

-- La comprobación es

-- λ> quickCheckWith (stdArgs {maxSize=15}) prop_equiv

-- +++ OK, passed 100 tests.

-- Comparación de eficiencia

-- =========================

-- La comparación es

-- λ> posicion [19,18,17,16,14,9,6]

-- 1000000

-- (2.61 secs, 1,754,265,776 bytes)

-- λ> posicion2 [19,18,17,16,14,9,6]

-- 1000000

-- (0.01 secs, 111,808 bytes)

-- λ> posicion3 [19,18,17,16,14,9,6]

-- 1000000

-- λ> posicion4 [19,18,17,16,14,9,6]

-- 1000000

-- (0.01 secs, 111,704 bytes)

-- λ> length (show (posicion2 [3*10^7]))

-- 9030900

-- (2.06 secs, 571,911,304 bytes)

-- λ> length (show (posicion3 [3*10^7]))

-- 9030900

-- (2.06 secs, 571,911,544 bytes)

-- λ> length (show (posicion4 [3*10^7]))

-- 9030900

-- (2.05 secs, 571,911,464 bytes)

-- Propiedad

-- =========

-- La propiedad es

prop_posicion :: Integer -> Property

prop_posicion n =

n >= 0 ==> posicion3 [n,n-1..0] == 2^(n+1) - 1

-- La comprobación es

-- λ> quickCheck prop_posicion

-- +++ OK, passed 100 tests.

Pensamiento

¡Volar sin alas donde todo es cielo!

Antonio Machado

Infinitud de primos gemelos

PorJosé A. Alonso 24-diciembre-201931-diciembre-2019

Un par de números primos (p,q) es un par de números primos gemelos si su distancia de 2; es decir, si q = p+2. Por ejemplo, (17,19) es una par de números primos gemelos.

La conjetura de los primos gemelos postula la existencia de infinitos pares de primos gemelos.

Definir la constante

   primosGemelos :: [(Integer,Integer)]

1	primosGemelos :: [(Integer,Integer)]

tal que sus elementos son los pares de primos gemelos. Por ejemplo,

   λ> take 7 primosGemelos
   [(3,5),(5,7),(11,13),(17,19),(29,31),(41,43),(59,61)]
   λ> primosGemelos !! (4*10^4)
   (6381911,6381913)

λ> take 7 primosGemelos

[(3,5),(5,7),(11,13),(17,19),(29,31),(41,43),(59,61)]

λ> primosGemelos !! (4*10^4)

(6381911,6381913)

Comprobar con QuickCheck la conjetura de los primos gemelos.

Soluciones

import Data.Numbers.Primes (primes, isPrime)
import Test.QuickCheck

-- 1ª solución
primosGemelos :: [(Integer,Integer)]
primosGemelos = [(x,x+2) | x <- primes, isPrime (x+2)]

-- 2ª solución
primosGemelos2 :: [(Integer,Integer)]
primosGemelos2 = [(x,y) | (x,y) <- zip primes (tail primes)
                        , y - x == 2]

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

--    λ> primosGemelos !! (2*10^4)
--    (2840447,2840449)
--    (3.93 secs, 12,230,474,952 bytes)
--    λ> primosGemelos2 !! (2*10^4)
--    (2840447,2840449)
--    (0.77 secs, 2,202,822,456 bytes)

-- Propiedad
-- =========

-- La propiedad es
prop_primosGemelos :: Integer -> Property
prop_primosGemelos n =
  n >= 0 ==> not (null [(x,y) | (x,y) <- primosGemelos2, x > n])

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

import Data.Numbers.Primes (primes, isPrime)

import Test.QuickCheck

-- 1ª solución

primosGemelos :: [(Integer,Integer)]

primosGemelos = [(x,x+2) | x <- primes, isPrime (x+2)]

-- 2ª solución

primosGemelos2 :: [(Integer,Integer)]

primosGemelos2 = [(x,y) | (x,y) <- zip primes (tail primes)

, y - x == 2]

-- Comparación de eficiencia

-- =========================

-- λ> primosGemelos !! (2*10^4)

-- (2840447,2840449)

-- (3.93 secs, 12,230,474,952 bytes)

-- λ> primosGemelos2 !! (2*10^4)

-- (2840447,2840449)

-- (0.77 secs, 2,202,822,456 bytes)

-- Propiedad

-- =========

-- La propiedad es

prop_primosGemelos :: Integer -> Property

prop_primosGemelos n =

n >= 0 ==> not (null [(x,y) | (x,y) <- primosGemelos2, x > n])

-- La comprobación es

-- λ> quickCheck prop_primosGemelos

-- +++ OK, passed 100 tests.

Pensamiento

El sentimiento ha de tener tanto de individual como de genérico; debe orientarse hacia valores universales, o que pretenden serlo.

Antonio Machado

Suma de números de Fibonacci con índice impar

PorJosé A. Alonso 23-diciembre-201930-diciembre-2019

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, ...

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

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

Definir la función

   sumaFibsIndiceImpar :: Int -> Integer

1	sumaFibsIndiceImpar :: Int -> Integer

tal que (sumaFibsIndiceImpar n) es la suma de los n primeros términos de la sucesión de Fibonacci no índice impar; es decir,

   sumaFibsIndiceImpar n = F(1) + F(3) + ... + F(2*n-1)

1	sumaFibsIndiceImpar n = F(1) + F(3) + ... + F(2*n-1)

Por ejemplo,

   sumaFibsIndiceImpar 1  ==  1
   sumaFibsIndiceImpar 2  ==  3
   sumaFibsIndiceImpar 3  ==  8
   sumaFibsIndiceImpar 4  ==  21
   sumaFibsIndiceImpar 5  ==  55
   sumaFibsIndiceImpar (10^4) `rem` (10^9)  ==  213093125

sumaFibsIndiceImpar 1 == 1

sumaFibsIndiceImpar 2 == 3

sumaFibsIndiceImpar 3 == 8

sumaFibsIndiceImpar 4 == 21

sumaFibsIndiceImpar 5 == 55

sumaFibsIndiceImpar (10^4) `rem` (10^9) == 213093125

En los ejemplos anteriores se observa que

   sumaFibsIndiceImpar 1  ==  F(2)
   sumaFibsIndiceImpar 2  ==  F(4)
   sumaFibsIndiceImpar 3  ==  F(6)
   sumaFibsIndiceImpar 4  ==  F(8)
   sumaFibsIndiceImpar 5  ==  F(10)

sumaFibsIndiceImpar 1 == F(2)

sumaFibsIndiceImpar 2 == F(4)

sumaFibsIndiceImpar 3 == F(6)

sumaFibsIndiceImpar 4 == F(8)

sumaFibsIndiceImpar 5 == F(10)

Comprobar con QuickCheck que (sumaFibsIndiceImpar n) es F(2n); es decir, el 2n-ésimo número de Fibonacci

Soluciones

import Test.QuickCheck

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

sumaFibsIndiceImpar :: Int -> Integer
sumaFibsIndiceImpar n = sum [fib (2*k-1) | k <- [1..n]]

-- (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)

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

sumaFibsIndiceImpar2 :: Int -> Integer
sumaFibsIndiceImpar2 n =
  sum [a | (a,b) <- zip fibs [0..2*n], odd b]

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

--    λ> sumaFibsIndiceImpar (10^4) `rem` (10^9)
--    213093125
--    (0.98 secs, 13,889,312 bytes)
--    λ> sumaFibsIndiceImpar2 (10^4) `rem` (10^9)
--    213093125
--    (0.05 secs, 18,047,720 bytes)

-- Comprobación
-- ============

-- La propiedad es
prop_sumaFibsIndiceImpar :: Int -> Property
prop_sumaFibsIndiceImpar n =
  n >= 0 ==> sumaFibsIndiceImpar n == fib (2*n)

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

import Test.QuickCheck

-- 1ª solución

-- ===========

sumaFibsIndiceImpar :: Int -> Integer

sumaFibsIndiceImpar n = sum [fib (2*k-1) | k <- [1..n]]

-- (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)

-- 2ª solución

-- ============

sumaFibsIndiceImpar2 :: Int -> Integer

sumaFibsIndiceImpar2 n =

sum [a | (a,b) <- zip fibs [0..2*n], odd b]

-- Comparación de eficiencia

-- =========================

-- λ> sumaFibsIndiceImpar (10^4) `rem` (10^9)

-- 213093125

-- (0.98 secs, 13,889,312 bytes)

-- λ> sumaFibsIndiceImpar2 (10^4) `rem` (10^9)

-- 213093125

-- (0.05 secs, 18,047,720 bytes)

-- Comprobación

-- ============

-- La propiedad es

prop_sumaFibsIndiceImpar :: Int -> Property

prop_sumaFibsIndiceImpar n =

n >= 0 ==> sumaFibsIndiceImpar n == fib (2*n)

-- La comprobación es

-- λ> quickCheck prop_sumaFibsIndiceImpar

-- +++ OK, passed 100 tests.

Referencia

Sum of sequence of odd index Fibonacci numbers en ProofWiki.

Pensamiento

El corazón del poeta, tan rico en sonoridades, es casi un insulto a la afonía cordial de la masa.

Antonio Machado

El teorema de Navidad de Fermat

PorJosé A. Alonso 20-diciembre-201927-diciembre-2019

El 25 de diciembre de 1640, en una carta a Mersenne, Fermat demostró la conjetura de Girard: todo primo de la forma 4n+1 puede expresarse de manera única como suma de dos cuadrados. Por eso es conocido como el Teorema de Navidad de Fermat.

Definir las funciones

   representaciones :: Integer -> [(Integer,Integer)]
   primosImparesConRepresentacionUnica :: [Integer]
   primos4nM1 :: [Integer]

representaciones :: Integer -> [(Integer,Integer)]

primosImparesConRepresentacionUnica :: [Integer]

primos4nM1 :: [Integer]

tales que

(representaciones n) es la lista de pares de números naturales (x,y) tales que n = x^2 + y^2 con x <= y. Por ejemplo,

     representaciones  20           ==  [(2,4)]
     representaciones  25           ==  [(0,5),(3,4)]
     representaciones 325           ==  [(1,18),(6,17),(10,15)]
     representaciones 100000147984  ==  [(0,316228)]
     length (representaciones (10^10))    ==  6
     length (representaciones (4*10^12))  ==  7

representaciones 20 == [(2,4)]

representaciones 25 == [(0,5),(3,4)]

representaciones 325 == [(1,18),(6,17),(10,15)]

representaciones 100000147984 == [(0,316228)]

length (representaciones (10^10)) == 6

length (representaciones (4*10^12)) == 7

primosImparesConRepresentacionUnica es la lista de los números primos impares que se pueden escribir exactamente de una manera como suma de cuadrados de pares de números naturales (x,y) con x <= y. Por ejemplo,

     λ> take 20 primosImparesConRepresentacionUnica
     [5,13,17,29,37,41,53,61,73,89,97,101,109,113,137,149,157,173,181,193]

1 2	λ> take 20 primosImparesConRepresentacionUnica [5,13,17,29,37,41,53,61,73,89,97,101,109,113,137,149,157,173,181,193]

primos4nM1 es la lista de los números primos que se pueden escribir como uno más un múltiplo de 4 (es decir, que son congruentes con 1 módulo 4). Por ejemplo,

     λ> take 20 primos4nM1
     [5,13,17,29,37,41,53,61,73,89,97,101,109,113,137,149,157,173,181,193]

1 2	λ> take 20 primos4nM1 [5,13,17,29,37,41,53,61,73,89,97,101,109,113,137,149,157,173,181,193]

El teorema de Navidad de Fermat afirma que un número primo impar p se puede escribir exactamente de una manera como suma de dos cuadrados de números naturales p = x² + y^2 (con x <= y) si, y sólo si, p se puede escribir como uno más un múltiplo de 4 (es decir, que es congruente con 1 módulo 4).

Comprobar con QuickCheck el teorema de Navidad de Fermat; es decir, que para todo número n, los n-ésimos elementos de primosImparesConRepresentacionUnica y de primos4nM1 son iguales.

Soluciones

import Data.Numbers.Primes (primes)
import Test.QuickCheck

-- 1ª definición de representaciones
-- =================================

representaciones :: Integer -> [(Integer,Integer)]
representaciones n =
  [(x,y) | x <- [0..n], y <- [x..n], n == x*x + y*y]

-- 2ª definición de representaciones
-- =================================

representaciones2 :: Integer -> [(Integer,Integer)]
representaciones2 n =
  [(x,raiz z) | x <- [0..raiz (n `div` 2)] 
              , let z = n - x*x
              , esCuadrado z]

-- (esCuadrado x) se verifica si x es un número al cuadrado. Por
-- ejemplo,
--    esCuadrado 25  ==  True
--    esCuadrado 26  ==  False
esCuadrado :: Integer -> Bool
esCuadrado x = x == y * y
  where y = raiz x

-- (raiz x) es la raíz cuadrada entera de x. Por ejemplo,
--    raiz 25  ==  5
--    raiz 24  ==  4
--    raiz 26  ==  5
raiz :: Integer -> Integer 
raiz 0 = 0
raiz 1 = 1
raiz x = aux (0,x)
    where aux (a,b) | d == x    = c
                    | c == a    = a
                    | d < x     = aux (c,b)
                    | otherwise = aux (a,c) 
              where c = (a+b) `div` 2
                    d = c^2

-- 3ª definición de representaciones
-- =================================

representaciones3 :: Integer -> [(Integer,Integer)]
representaciones3 n =
  [(x,raiz3 z) | x <- [0..raiz3 (n `div` 2)] 
               , let z = n - x*x
               , esCuadrado3 z]

-- (esCuadrado x) se verifica si x es un número al cuadrado. Por
-- ejemplo,
--    esCuadrado3 25  ==  True
--    esCuadrado3 26  ==  False
esCuadrado3 :: Integer -> Bool
esCuadrado3 x = x == y * y
  where y = raiz3 x

-- (raiz3 x) es la raíz cuadrada entera de x. Por ejemplo,
--    raiz3 25  ==  5
--    raiz3 24  ==  4
--    raiz3 26  ==  5
raiz3 :: Integer -> Integer
raiz3 x = floor (sqrt (fromIntegral x))
          
-- 4ª definición de representaciones
-- =================================
 
representaciones4 :: Integer -> [(Integer, Integer)]
representaciones4 n = aux 0 (floor (sqrt (fromIntegral n)))
  where aux x y
          | x > y     = [] 
          | otherwise = case compare (x*x + y*y) n of
                          LT -> aux (x + 1) y
                          EQ -> (x, y) : aux (x + 1) (y - 1)
                          GT -> aux x (y - 1)

-- Equivalencia de las definiciones de representaciones
-- ====================================================

-- La propiedad es
prop_representaciones_equiv :: (Positive Integer) -> Bool
prop_representaciones_equiv (Positive n) =
  representaciones  n == representaciones2 n &&
  representaciones2 n == representaciones3 n &&
  representaciones3 n == representaciones4 n

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

-- Comparación de eficiencia de las definiciones de representaciones
-- =================================================================

--    λ> representaciones 3025
--    [(0,55),(33,44)]
--    (2.86 secs, 1,393,133,528 bytes)
--    λ> representaciones2 3025
--    [(0,55),(33,44)]
--    (0.00 secs, 867,944 bytes)
--    λ> representaciones3 3025
--    [(0,55),(33,44)]
--    (0.00 secs, 173,512 bytes)
--    λ> representaciones4 3025
--    [(0,55),(33,44)]
--    (0.00 secs, 423,424 bytes)
--    
--    λ> length (representaciones2 (10^10))
--    6
--    (3.38 secs, 2,188,903,544 bytes)
--    λ> length (representaciones3 (10^10))
--    6
--    (0.10 secs, 62,349,048 bytes)
--    λ> length (representaciones4 (10^10))
--    6
--    (0.11 secs, 48,052,360 bytes)
--
--    λ> length (representaciones3 (4*10^12))
--    7
--    (1.85 secs, 1,222,007,176 bytes)
--    λ> length (representaciones4 (4*10^12))
--    7
--    (1.79 secs, 953,497,480 bytes)

-- Definición de primosImparesConRepresentacionUnica
-- =================================================

primosImparesConRepresentacionUnica :: [Integer]
primosImparesConRepresentacionUnica =
  [x | x <- tail primes
     , length (representaciones4 x) == 1]

-- Definición de primos4nM1
-- ========================

primos4nM1 :: [Integer]
primos4nM1 = [x | x <- primes
                , x `mod` 4 == 1]

-- Teorema de Navidad de Fermat
-- ============================

-- La propiedad es
prop_teoremaDeNavidadDeFermat :: Positive Int -> Bool
prop_teoremaDeNavidadDeFermat (Positive n) =
  primosImparesConRepresentacionUnica !! n == primos4nM1 !! n
             
-- La comprobación es
--    λ> quickCheck prop_teoremaDeNavidadDeFermat
--    +++ OK, passed 100 tests.

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

import Data.Numbers.Primes (primes)

import Test.QuickCheck

-- 1ª definición de representaciones

-- =================================

representaciones :: Integer -> [(Integer,Integer)]

representaciones n =

[(x,y) | x <- [0..n], y <- [x..n], n == x*x + y*y]

-- 2ª definición de representaciones

-- =================================

representaciones2 :: Integer -> [(Integer,Integer)]

representaciones2 n =

[(x,raiz z) | x <- [0..raiz (n `div` 2)]

, let z = n - x*x

, esCuadrado z]

-- (esCuadrado x) se verifica si x es un número al cuadrado. Por

-- ejemplo,

-- esCuadrado 25 == True

-- esCuadrado 26 == False

esCuadrado :: Integer -> Bool

esCuadrado x = x == y * y

where y = raiz x

-- (raiz x) es la raíz cuadrada entera de x. Por ejemplo,

-- raiz 25 == 5

-- raiz 24 == 4

-- raiz 26 == 5

raiz :: Integer -> Integer

raiz 0 = 0

raiz 1 = 1

raiz x = aux (0,x)

where aux (a,b) | d == x = c

| c == a = a

| d < x = aux (c,b)

| otherwise = aux (a,c)

where c = (a+b) `div` 2

d = c^2

-- 3ª definición de representaciones

-- =================================

representaciones3 :: Integer -> [(Integer,Integer)]

representaciones3 n =

[(x,raiz3 z) | x <- [0..raiz3 (n `div` 2)]

, let z = n - x*x

, esCuadrado3 z]

-- (esCuadrado x) se verifica si x es un número al cuadrado. Por

-- ejemplo,

-- esCuadrado3 25 == True

-- esCuadrado3 26 == False

esCuadrado3 :: Integer -> Bool

esCuadrado3 x = x == y * y

where y = raiz3 x

-- (raiz3 x) es la raíz cuadrada entera de x. Por ejemplo,

-- raiz3 25 == 5

-- raiz3 24 == 4

-- raiz3 26 == 5

raiz3 :: Integer -> Integer

raiz3 x = floor (sqrt (fromIntegral x))

-- 4ª definición de representaciones

-- =================================

representaciones4 :: Integer -> [(Integer, Integer)]

representaciones4 n = aux 0 (floor (sqrt (fromIntegral n)))

where aux x y

| x > y = []

| otherwise = case compare (x*x + y*y) n of

LT -> aux (x + 1) y

EQ -> (x, y) : aux (x + 1) (y - 1)

GT -> aux x (y - 1)

-- Equivalencia de las definiciones de representaciones

-- ====================================================

-- La propiedad es

prop_representaciones_equiv :: (Positive Integer) -> Bool

prop_representaciones_equiv (Positive n) =

representaciones n == representaciones2 n &&

representaciones2 n == representaciones3 n &&

representaciones3 n == representaciones4 n

-- La comprobación es

-- λ> quickCheck prop_representaciones_equiv

-- +++ OK, passed 100 tests.

-- Comparación de eficiencia de las definiciones de representaciones

-- =================================================================

-- λ> representaciones 3025

-- [(0,55),(33,44)]

-- (2.86 secs, 1,393,133,528 bytes)

-- λ> representaciones2 3025

-- [(0,55),(33,44)]

-- (0.00 secs, 867,944 bytes)

-- λ> representaciones3 3025

-- [(0,55),(33,44)]

-- (0.00 secs, 173,512 bytes)

-- λ> representaciones4 3025

-- [(0,55),(33,44)]

-- (0.00 secs, 423,424 bytes)

-- λ> length (representaciones2 (10^10))

-- 6

-- (3.38 secs, 2,188,903,544 bytes)

-- λ> length (representaciones3 (10^10))

-- 6

-- (0.10 secs, 62,349,048 bytes)

-- λ> length (representaciones4 (10^10))

-- 6

-- (0.11 secs, 48,052,360 bytes)

-- λ> length (representaciones3 (4*10^12))

-- 7

-- (1.85 secs, 1,222,007,176 bytes)

-- λ> length (representaciones4 (4*10^12))

-- 7

-- (1.79 secs, 953,497,480 bytes)

-- Definición de primosImparesConRepresentacionUnica

-- =================================================

primosImparesConRepresentacionUnica :: [Integer]

primosImparesConRepresentacionUnica =

[x | x <- tail primes

, length (representaciones4 x) == 1]

-- Definición de primos4nM1

-- ========================

primos4nM1 :: [Integer]

primos4nM1 = [x | x <- primes

, x `mod` 4 == 1]

-- Teorema de Navidad de Fermat

-- ============================

-- La propiedad es

prop_teoremaDeNavidadDeFermat :: Positive Int -> Bool

prop_teoremaDeNavidadDeFermat (Positive n) =

primosImparesConRepresentacionUnica !! n == primos4nM1 !! n

-- La comprobación es

-- λ> quickCheck prop_teoremaDeNavidadDeFermat

-- +++ OK, passed 100 tests.

Pensamiento

Dijo Dios: brote la nada
Y alzó su mano derecha,
hasta ocultar su mirada.
Y quedó la nada hecha.

Antonio Machado

Primos o cuadrados de primos

PorJosé A. Alonso 6-diciembre-201917-diciembre-2019

Definir la constante

   primosOcuadradosDePrimos :: [Integer]

1	primosOcuadradosDePrimos :: [Integer]

cuyos elementos son los número primos o cuadrados de primos. Por ejemplo,

   λ> take 20 primosOcuadradosDePrimos
   [2,3,4,5,7,9,11,13,17,19,23,25,29,31,37,41,43,47,49,53]
   λ> primosOcuadradosDePrimos !! (10^6)
   15476729

λ> take 20 primosOcuadradosDePrimos

[2,3,4,5,7,9,11,13,17,19,23,25,29,31,37,41,43,47,49,53]

λ> primosOcuadradosDePrimos !! (10^6)

15476729

Comprobar con QuickCheck que las lista primosOcuadradosDePrimos y unifactorizables (definida en el ejercicio anterior) son iguales.

Soluciones

import Test.QuickCheck
import Data.Numbers.Primes (primeFactors, primes)

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

primosOcuadradosDePrimos :: [Integer]
primosOcuadradosDePrimos =
  filter esPrimoOcuadradoDePrimo [2..]

esPrimoOcuadradoDePrimo :: Integer -> Bool
esPrimoOcuadradoDePrimo n = aux xs
  where xs = primeFactors n
        aux [_]   = True
        aux [x,y] = x == y
        aux _     = False

-- 2ª solución
-- ===========
  
primosOcuadradosDePrimos2 :: [Integer]
primosOcuadradosDePrimos2 = mezcla primes (map (^2) primes)

mezcla :: Ord a => [a] -> [a] -> [a]
mezcla (x:xs) (y:ys) | x < y     = x : mezcla xs (y:ys)
                     | otherwise = y : mezcla (x:xs) ys

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

--    λ> primosOcuadradosDePrimos !! (2*10^4)
--    223589
--    (3.08 secs, 9,829,725,096 bytes)
--    λ> primosOcuadradosDePrimos2 !! (2*10^4)
--    223589
--    (0.04 secs, 73,751,888 bytes)
--
--    λ> primosOcuadradosDePrimos2 !! (5*10^5)
--    7362497
--    (1.29 secs, 3,192,803,040 bytes)

-- Propiedad de equivalencia
-- =========================

-- La propiedad es
prop_equivalencia :: Int -> Property
prop_equivalencia n =
  n >= 0 ==> primosOcuadradosDePrimos2 !! n == unifactorizables !! n

--  unifactorizables es la lísta de los números enteros mayores que 1
--  que se pueden escribir sólo de una forma única como producto de
--  enteros distintos mayores que uno. Por ejemplo,
--     λ> take 20 unifactorizables
--     [2,3,4,5,7,9,11,13,17,19,23,25,29,31,37,41,43,47,49,53]
--     λ> unifactorizables !! 300
--     1873
unifactorizables :: [Integer]
unifactorizables =
  [n | n <- [2..]
     , length (sublistasConProducto n [2..n]) == 1]

-- (sublistasConProducto n xs) es la lista de las sublistas de la
-- lista ordenada estrictamente creciente xs (cuyos elementos son
-- enteros mayores que 1) cuyo producto es el número entero n (con n
-- mayor que 1). Por ejemplo,  
--    λ> sublistasConProducto 72 [2,3,4,5,6,7,9,10,16]
--    [[2,4,9],[3,4,6]]
--    λ> sublistasConProducto 720 [2,3,4,5,6,7,9,10,16]
--    [[2,3,4,5,6],[2,4,9,10],[3,4,6,10],[5,9,16]]
--    λ> sublistasConProducto 2 [4,7]
--    []
sublistasConProducto :: Integer -> [Integer] -> [[Integer]]
sublistasConProducto _ [] = []
sublistasConProducto n (x:xs)
  | x > n     = []
  | x == n    = [[x]]
  | r == 0    = map (x:) (sublistasConProducto q xs)
                ++ sublistasConProducto n xs
  | otherwise = sublistasConProducto n xs
  where (q,r) = quotRem n x

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

import Test.QuickCheck

import Data.Numbers.Primes (primeFactors, primes)

-- 1ª solución

-- ===========

primosOcuadradosDePrimos :: [Integer]

primosOcuadradosDePrimos =

filter esPrimoOcuadradoDePrimo [2..]

esPrimoOcuadradoDePrimo :: Integer -> Bool

esPrimoOcuadradoDePrimo n = aux xs

where xs = primeFactors n

aux [_] = True

aux [x,y] = x == y

aux _ = False

-- 2ª solución

-- ===========

primosOcuadradosDePrimos2 :: [Integer]

primosOcuadradosDePrimos2 = mezcla primes (map (^2) primes)

mezcla :: Ord a => [a] -> [a] -> [a]

mezcla (x:xs) (y:ys) | x < y = x : mezcla xs (y:ys)

| otherwise = y : mezcla (x:xs) ys

-- Comparación de eficiencia

-- =========================

-- λ> primosOcuadradosDePrimos !! (2*10^4)

-- 223589

-- (3.08 secs, 9,829,725,096 bytes)

-- λ> primosOcuadradosDePrimos2 !! (2*10^4)

-- 223589

-- (0.04 secs, 73,751,888 bytes)

-- λ> primosOcuadradosDePrimos2 !! (5*10^5)

-- 7362497

-- (1.29 secs, 3,192,803,040 bytes)

-- Propiedad de equivalencia

-- =========================

-- La propiedad es

prop_equivalencia :: Int -> Property

prop_equivalencia n =

n >= 0 ==> primosOcuadradosDePrimos2 !! n == unifactorizables !! n

-- unifactorizables es la lísta de los números enteros mayores que 1

-- que se pueden escribir sólo de una forma única como producto de

-- enteros distintos mayores que uno. Por ejemplo,

-- λ> take 20 unifactorizables

-- [2,3,4,5,7,9,11,13,17,19,23,25,29,31,37,41,43,47,49,53]

-- λ> unifactorizables !! 300

-- 1873

unifactorizables :: [Integer]

unifactorizables =

[n | n <- [2..]

, length (sublistasConProducto n [2..n]) == 1]

-- (sublistasConProducto n xs) es la lista de las sublistas de la

-- lista ordenada estrictamente creciente xs (cuyos elementos son

-- enteros mayores que 1) cuyo producto es el número entero n (con n

-- mayor que 1). Por ejemplo,

-- λ> sublistasConProducto 72 [2,3,4,5,6,7,9,10,16]

-- [[2,4,9],[3,4,6]]

-- λ> sublistasConProducto 720 [2,3,4,5,6,7,9,10,16]

-- [[2,3,4,5,6],[2,4,9,10],[3,4,6,10],[5,9,16]]

-- λ> sublistasConProducto 2 [4,7]

-- []

sublistasConProducto :: Integer -> [Integer] -> [[Integer]]

sublistasConProducto _ [] = []

sublistasConProducto n (x:xs)

| x > n = []

| x == n = [[x]]

| r == 0 = map (x:) (sublistasConProducto q xs)

++ sublistasConProducto n xs

| otherwise = sublistasConProducto n xs

where (q,r) = quotRem n x

-- La comprobación es

-- λ> quickCheck prop_equivalencia

-- +++ OK, passed 100 tests.

Pensamiento

Despacito y buena letra: el hacer las cosas bien importa más que el hacerlas.

Antonio Machado

Posiciones de conjuntos finitos de naturales

Soluciones

Pensamiento

Infinitud de primos gemelos

Soluciones

Pensamiento

Suma de números de Fibonacci con índice impar

Soluciones

Referencia

Pensamiento

El teorema de Navidad de Fermat

Soluciones

Pensamiento

Primos o cuadrados de primos

Soluciones

Pensamiento

Entradas recientes

Buscar

Correo electrónico