product

La sucesión de Sylvester

PorJosé A. Alonso 6-mayo-202013-mayo-2020

La sucesión de Sylvester es la sucesión que comienza en 2 y sus restantes términos se obtienen multiplicando los anteriores y sumándole 1.

Definir las funciones

   sylvester        :: Integer -> Integer
   graficaSylvester :: Integer -> Integer -> IO ()

1 2	sylvester :: Integer -> Integer graficaSylvester :: Integer -> Integer -> IO ()

tales que

(sylvester n) es el n-ésimo término de la sucesión de Sylvester. Por ejemplo,

     λ> [sylvester n | n <- [0..7]]
     [2,3,7,43,1807,3263443,10650056950807,113423713055421844361000443]
     λ> length (show (sylvester 25))
     6830085

λ> [sylvester n | n <- [0..7]]

[2,3,7,43,1807,3263443,10650056950807,113423713055421844361000443]

λ> length (show (sylvester 25))

6830085

(graficaSylvester d n) dibuja la gráfica de los d últimos dígitos de los n primeros términos de la sucesión de Sylvester. Por ejemplo,
- (graficaSylvester 3 30) dibuja
- (graficaSylvester 4 30) dibuja
- (graficaSylvester 5 30) dibuja

Nota: Se puede usar programación dinámica para aumentar la eficiencia.

Soluciones

import Data.List               (genericIndex)
import Data.Array              ((!), array)
import Graphics.Gnuplot.Simple (plotList, Attribute (Key, PNG))

-- 1ª solución (por recursión)
-- ===========================

sylvester1 :: Integer -> Integer
sylvester1 0 = 2
sylvester1 n = 1 + product [sylvester1 k | k <- [0..n-1]]

-- 2ª solución (con programación dinámica)
-- =======================================

sylvester2 :: Integer -> Integer
sylvester2 n = v ! n where
  v = array (0,n) [(i,f i) | i <- [0..n]]
  f 0 = 2
  f m = 1 + product [v!k | k <- [0..m-1]]

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

-- Observando que
--    S(n) = 1 + S(0)*S(1)*...*S(n-2)*S(n-1)
--         = 1 + (1 + S(0)*S(1)*...*S(n-2))*S(n-1) - S(n-1)
--         = 1 + S(n-1)*S(n-1) - S(n-1)
--         = 1 + S(n-1)^2 - S(n-1)
-- se obtiene la siguiente definición.
sylvester3 :: Integer -> Integer
sylvester3 0 = 2
sylvester3 n = 1 + x^2 - x
  where x = sylvester3 (n-1)

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

sylvester4 :: Integer -> Integer
sylvester4 n = v ! n where
  v = array (0,n) [(i,f i) | i <- [0..n]]
  f 0 = 2
  f m = 1 + x^2 - x
    where x = v ! (m-1)

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

sylvester5 :: Integer -> Integer
sylvester5 n = sucSylvester5 `genericIndex` n

sucSylvester5 :: [Integer]
sucSylvester5 = iterate (\x -> (x-1)*x+1) 2 

-- La comparación es
--    λ> length (show (sylvester1 23))
--    1707522
--    (6.03 secs, 4,090,415,704 bytes)
--    λ> length (show (sylvester2 23))
--    1707522
--    (0.33 secs, 109,477,296 bytes)
--    λ> length (show (sylvester3 23))
--    1707522
--    (0.35 secs, 109,395,136 bytes)
--    λ> length (show (sylvester4 23))
--    1707522
--    (0.33 secs, 109,402,440 bytes)
--    λ> length (show (sylvester5 23))
--    1707522
--    (0.30 secs, 108,676,256 bytes)

graficaSylvester :: Integer -> Integer -> IO ()
graficaSylvester d n =
  plotList [ Key Nothing
           , PNG ("La_sucesion_de_Sylvester_" ++ show (d,n) ++ ".png")
           ]
           [sylvester5 k `mod` (10^d) | k <- [0..n]]

import Data.List (genericIndex)

import Data.Array ((!), array)

import Graphics.Gnuplot.Simple (plotList, Attribute (Key, PNG))

-- 1ª solución (por recursión)

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

sylvester1 :: Integer -> Integer

sylvester1 0 = 2

sylvester1 n = 1 + product [sylvester1 k | k <- [0..n-1]]

-- 2ª solución (con programación dinámica)

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

sylvester2 :: Integer -> Integer

sylvester2 n = v ! n where

v = array (0,n) [(i,f i) | i <- [0..n]]

f 0 = 2

f m = 1 + product [v!k | k <- [0..m-1]]

-- 3ª solución

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

-- Observando que

-- S(n) = 1 + S(0)*S(1)*...*S(n-2)*S(n-1)

-- = 1 + (1 + S(0)*S(1)*...*S(n-2))*S(n-1) - S(n-1)

-- = 1 + S(n-1)*S(n-1) - S(n-1)

-- = 1 + S(n-1)^2 - S(n-1)

-- se obtiene la siguiente definición.

sylvester3 :: Integer -> Integer

sylvester3 0 = 2

sylvester3 n = 1 + x^2 - x

where x = sylvester3 (n-1)

-- 4ª solución

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

sylvester4 :: Integer -> Integer

sylvester4 n = v ! n where

v = array (0,n) [(i,f i) | i <- [0..n]]

f 0 = 2

f m = 1 + x^2 - x

where x = v ! (m-1)

-- 5ª solución

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

sylvester5 :: Integer -> Integer

sylvester5 n = sucSylvester5 `genericIndex` n

sucSylvester5 :: [Integer]

sucSylvester5 = iterate (\x -> (x-1)*x+1) 2

-- La comparación es

-- λ> length (show (sylvester1 23))

-- 1707522

-- (6.03 secs, 4,090,415,704 bytes)

-- λ> length (show (sylvester2 23))

-- 1707522

-- (0.33 secs, 109,477,296 bytes)

-- λ> length (show (sylvester3 23))

-- 1707522

-- (0.35 secs, 109,395,136 bytes)

-- λ> length (show (sylvester4 23))

-- 1707522

-- (0.33 secs, 109,402,440 bytes)

-- λ> length (show (sylvester5 23))

-- 1707522

-- (0.30 secs, 108,676,256 bytes)

graficaSylvester :: Integer -> Integer -> IO ()

graficaSylvester d n =

plotList [ Key Nothing

, PNG ("La_sucesion_de_Sylvester_" ++ show (d,n) ++ ".png")

]

[sylvester5 k `mod` (10^d) | k <- [0..n]]

Cálculo de pi mediante la variante de Euler de la serie armónica

PorJosé A. Alonso 9-abril-202026-marzo-2022

En el artículo El desarrollo más bello de Pi como suma infinita, Miguel Ángel Morales comenta el desarrollo de pi publicado por Leonhard Euler en su libro «Introductio in Analysis Infinitorum» (1748).

El desarrollo es el siguiente

y se obtiene a partir de la serie armónica

modificando sólo el signo de algunos términos según el siguiente criterio:

Dejamos un + cuando el denominador de la fracción sea un 2 o un primo de la forma 4m-1.
Cambiamos a – si el denominador de la fracción es un primo de la forma 4m+1.
Si el número es compuesto ponemos el signo que quede al multiplicar los signos correspondientes a cada factor.

Por ejemplo,

la de denominador 3 = 4×1-1 lleva un +,
la de denominador 5 = 4×1+1 lleva un -,
la de denominador 13 = 4×3+1 lleva un -,
la de denominador 6 = 2×3 lleva un + (porque los dos llevan un +),
la de denominador 10 = 2×5 lleva un – (porque el 2 lleva un + y el 5 lleva un -) y
la de denominador 50 = 5x5x2 lleva un + (un – por el primer 5, otro – por el segundo 5 y un + por el 2).

Definir las funciones

  aproximacionPi :: Int -> Double
  grafica        :: Int -> IO ()

1 2	aproximacionPi :: Int -> Double grafica :: Int -> IO ()

tales que

(aproximacionPi n) es la aproximación de pi obtenida sumando los n primeros términos de la serie de Euler. Por ejemplo.

     aproximacionPi 1        ==  1.0
     aproximacionPi 10       ==  2.3289682539682537
     aproximacionPi 100      ==  2.934318000847734
     aproximacionPi 1000     ==  3.0603246224585128
     aproximacionPi 10000    ==  3.1105295744825403
     aproximacionPi 100000   ==  3.134308801935256
     aproximacionPi 1000000  ==  3.1395057903490806

aproximacionPi 1 == 1.0

aproximacionPi 10 == 2.3289682539682537

aproximacionPi 100 == 2.934318000847734

aproximacionPi 1000 == 3.0603246224585128

aproximacionPi 10000 == 3.1105295744825403

aproximacionPi 100000 == 3.134308801935256

aproximacionPi 1000000 == 3.1395057903490806

(grafica n) dibuja la gráfica de las aproximaciones de pi usando k sumando donde k toma los valores de la lista [100,110..n]. Por ejemplo, al evaluar (grafica 4000) se obtiene

Soluciones

import Data.Numbers.Primes
import Graphics.Gnuplot.Simple

-- 1ª definición
-- =============

aproximacionPi :: Int -> Double
aproximacionPi n =
  sum [1 / fromIntegral (k * signo k) | k <- [1..n]] 

signoPrimo :: Int -> Int
signoPrimo 2 = 1
signoPrimo p | p `mod` 4 == 3 = 1
             | otherwise      = -1

signo :: Int -> Int
signo n | isPrime n = signoPrimo n
        | otherwise = product (map signoPrimo (primeFactors n))

-- 2ª definición
-- =============

aproximacionPi2 :: Int -> Double
aproximacionPi2 n = serieEuler !! (n-1)

serieEuler :: [Double]
serieEuler =
  scanl1 (+) [1 / fromIntegral (n * signo n) | n <- [1..]]

-- Definición de grafica
-- =====================

grafica :: Int -> IO ()
grafica n = 
    plotList [Key Nothing]
             [(k,aproximacionPi2 k) | k <- [100,110..n]]

import Data.Numbers.Primes

import Graphics.Gnuplot.Simple

-- 1ª definición

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

aproximacionPi :: Int -> Double

aproximacionPi n =

sum [1 / fromIntegral (k * signo k) | k <- [1..n]]

signoPrimo :: Int -> Int

signoPrimo 2 = 1

signoPrimo p | p `mod` 4 == 3 = 1

| otherwise = -1

signo :: Int -> Int

signo n | isPrime n = signoPrimo n

| otherwise = product (map signoPrimo (primeFactors n))

-- 2ª definición

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

aproximacionPi2 :: Int -> Double

aproximacionPi2 n = serieEuler !! (n-1)

serieEuler :: [Double]

serieEuler =

scanl1 (+) [1 / fromIntegral (n * signo n) | n <- [1..]]

-- Definición de grafica

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

grafica :: Int -> IO ()

grafica n =

plotList [Key Nothing]

[(k,aproximacionPi2 k) | k <- [100,110..n]]

Otras soluciones

Se pueden escribir otras soluciones en los comentarios.
El código se debe escribir entre una línea con <pre lang="haskell"> y otra con </pre>

Productos de dos y tres números consecutivos

PorJosé A. Alonso 3-abril-202010-abril-2020

Definir la función

   productos :: Integer -> Integer -> [[Integer]]

1	productos :: Integer -> Integer -> [[Integer]]

tal que (productos n x) es las listas de n elementos consecutivos cuyo producto es x. Por ejemplo,

   productos 2 6     ==  [[2,3]]
   productos 3 6     ==  [[1,2,3]]
   productos 4 1680  ==  [[5,6,7,8]]
   productos 2 5     ==  []

productos 2 6 == [[2,3]]

productos 3 6 == [[1,2,3]]

productos 4 1680 == [[5,6,7,8]]

productos 2 5 == []

Comprobar con QuickCheck que si n > 0 y x > 0, entonces

   productos n (product [x..x+n-1]) == [[x..x+n-1]]

1	productos n (product [x..x+n-1]) == [[x..x+n-1]]

Usando productos, definir la función

   productosDe2y3consecutivos :: [Integer]

1	productosDe2y3consecutivos :: [Integer]

cuyos elementos son los números naturales (no nulos) que pueden expresarse simultáneamente como producto de dos y tres números consecutivos. Por ejemplo,

   head productosDe2y3consecutivos  ==  6

1	head productosDe2y3consecutivos == 6

Nota. Según demostró Mordell en 1962, productosDe2y3consecutivos sólo tiene dos elementos.

Soluciones

import Test.QuickCheck

-- 1ª definición
productos1 :: Integer -> Integer -> [[Integer]]
productos1 n x =
  [[y..y+n-1] | y <- [1..x]
              , product [y..y+n-1] == x]

-- 2ª definición
productos2 :: Integer -> Integer -> [[Integer]]
productos2 n x =
  [[z..z+n-1] | z <- [1..y]
              , product [z..z+n-1] == x]
  where y = head (filter (\y -> y^n >= x) [2..])

productos :: Integer -> Integer -> [[Integer]]
productos = productos2

prop_productos n x =
  n > 0 && x > 0 ==> productos n (product [x..x+n-1]) == [[x..x+n-1]]

-- La comprobación es
--    λ> quickCheck prop_productos
--    +++ OK, passed 100 tests.
--    (0.10 secs, 26409644 bytes)

productosDe2y3consecutivos :: [Integer]
productosDe2y3consecutivos =
  [x | x <- [1..] 
     , (not . null) (productos 2 x)
     , (not . null) (productos 3 x)]

-- El cálculo es
--    λ> take 2 productosDe2y3consecutivos
--    [6,210]
--    λ> productos 2 210
--    [[14,15]]
--    λ> productos 3 210
--    [[5,6,7]]

import Test.QuickCheck

-- 1ª definición

productos1 :: Integer -> Integer -> [[Integer]]

productos1 n x =

[[y..y+n-1] | y <- [1..x]

, product [y..y+n-1] == x]

-- 2ª definición

productos2 :: Integer -> Integer -> [[Integer]]

productos2 n x =

[[z..z+n-1] | z <- [1..y]

, product [z..z+n-1] == x]

where y = head (filter (\y -> y^n >= x) [2..])

productos :: Integer -> Integer -> [[Integer]]

productos = productos2

prop_productos n x =

n > 0 && x > 0 ==> productos n (product [x..x+n-1]) == [[x..x+n-1]]

-- La comprobación es

-- λ> quickCheck prop_productos

-- +++ OK, passed 100 tests.

-- (0.10 secs, 26409644 bytes)

productosDe2y3consecutivos :: [Integer]

productosDe2y3consecutivos =

[x | x <- [1..]

, (not . null) (productos 2 x)

, (not . null) (productos 3 x)]

-- El cálculo es

-- λ> take 2 productosDe2y3consecutivos

-- [6,210]

-- λ> productos 2 210

-- [[14,15]]

-- λ> productos 3 210

-- [[5,6,7]]

Otras soluciones

Se pueden escribir otras soluciones en los comentarios.
El código se debe escribir entre una línea con <pre lang="haskell"> y otra con </pre>

Pensamiento

«El verdadero viaje de descubrimiento no consiste en buscar nuevos paisajes sino en tener nuevos ojos.»

Marcel Proust.

Números de Bell

PorJosé A. Alonso 5-febrero-202012-febrero-2020

Una partición de un conjunto A es un conjunto de subconjuntos no vacíos de A, disjuntos dos a dos y cuya unión es A. Por ejemplo, el conjunto {1, 2, 3} tiene exactamente 5 particiones:

   {{1}, {2}, {3}}
   {{1,2}, {3}}
   {{1,3}, {2}}
   {{1}, {2,3}}
   {{1,2,3}}

{{1}, {2}, {3}}

{{1,2}, {3}}

{{1,3}, {2}}

{{1}, {2,3}}

El n-ésimo número de Bell, B(n), es el número de particiones de un conjunto de n elementos. Por lo visto anteriormentem B(3) = 5.

Definir las funciones

   particiones :: [a] -> [[[a]]]
   bell :: Integer -> Integer

1 2	particiones :: [a] -> [[[a]]] bell :: Integer -> Integer

tales que

(particiones xs) es el conjunto de las particiones de xs. Por ejemplo,

     λ> particiones [1,2]
     [[[1,2]],[[1],[2]]]
     λ> particiones [1,2,3]
     [[[1,2,3]],[[1],[2,3]],[[1,2],[3]],[[2],[1,3]],[[1],[2],[3]]]
     λ> particiones "abcd"
     [["abcd"],["a","bcd"],["ab","cd"],["b","acd"],["abc","d"],["bc","ad"],
      ["ac","bd"],["c","abd"],["a","b","cd"],["a","bc","d"],["a","c","bd"],
      ["ab","c","d"],["b","ac","d"],["b","c","ad"],["a","b","c","d"]]

λ> particiones [1,2]

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

λ> particiones [1,2,3]

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

λ> particiones "abcd"

[["abcd"],["a","bcd"],["ab","cd"],["b","acd"],["abc","d"],["bc","ad"],

["ac","bd"],["c","abd"],["a","b","cd"],["a","bc","d"],["a","c","bd"],

["ab","c","d"],["b","ac","d"],["b","c","ad"],["a","b","c","d"]]

(bell n) es el n-ésimo número de Bell. Por ejemplo,

     λ> bell 3
     5
     λ> map bell [0..10]
     [1,1,2,5,15,52,203,877,4140,21147,115975]

λ> bell 3

λ> map bell [0..10]

[1,1,2,5,15,52,203,877,4140,21147,115975]

Comprobar con QuickCheck que (bell n) es equivalente a la función B(n) definida por

$B(0) = 1$
$B(n) = \displaystyle \sum_{k=0}^{n-1} \binom{n-1}{k} B(k)$

Soluciones

import Data.List (genericLength)
import Test.QuickCheck

-- Definición de particiones
-- =========================

particiones :: [a] -> [[[a]]]
particiones [] = [[]]
particiones (x:xs) =
  concat [([x] : yss) : inserta x yss | yss <- ysss]
  where ysss = particiones xs

-- (inserta x yss) es la lista obtenida insertando x en cada uno de los
-- elementos de yss. Por ejemplo, 
--    λ> inserta 1 [[2,3],[4],[5,6,7]]
--    [[[1,2,3],[4],[5,6,7]],[[2,3],[1,4],[5,6,7]],[[2,3],[4],[1,5,6,7]]]
inserta :: a -> [[a]] -> [[[a]]]
inserta _ []       = []
inserta x (ys:yss) = ((x:ys):yss) : [ys : zs | zs <- inserta x yss] 

-- Definición de Bell
-- ==================

bell :: Integer -> Integer
bell n = genericLength (particiones [1..n])

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

prop_Bell :: Integer -> Property
prop_Bell n =
  n >= 0 ==> bell n == b n

b :: Integer -> Integer
b 0 = 1
b n = sum [comb (n-1) k * b k | k <- [0..n-1]]

comb :: Integer -> Integer -> Integer
comb n k = product [n-k+1..n] `div` product [1..k]

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

import Data.List (genericLength)

import Test.QuickCheck

-- Definición de particiones

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

particiones :: [a] -> [[[a]]]

particiones [] = [[]]

particiones (x:xs) =

concat [([x] : yss) : inserta x yss | yss <- ysss]

where ysss = particiones xs

-- (inserta x yss) es la lista obtenida insertando x en cada uno de los

-- elementos de yss. Por ejemplo,

-- λ> inserta 1 [[2,3],[4],[5,6,7]]

-- [[[1,2,3],[4],[5,6,7]],[[2,3],[1,4],[5,6,7]],[[2,3],[4],[1,5,6,7]]]

inserta :: a -> [[a]] -> [[[a]]]

inserta _ [] = []

inserta x (ys:yss) = ((x:ys):yss) : [ys : zs | zs <- inserta x yss]

-- Definición de Bell

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

bell :: Integer -> Integer

bell n = genericLength (particiones [1..n])

-- Propiedad

-- =========

prop_Bell :: Integer -> Property

prop_Bell n =

n >= 0 ==> bell n == b n

b :: Integer -> Integer

b 0 = 1

b n = sum [comb (n-1) k * b k | k <- [0..n-1]]

comb :: Integer -> Integer -> Integer

comb n k = product [n-k+1..n] `div` product [1..k]

-- La comprobación es

-- λ> quickCheckWith (stdArgs {maxSize=10}) prop_Bell

-- +++ OK, passed 100 tests.

Otras soluciones

Se pueden escribir otras soluciones en los comentarios.
El código se debe escribir entre una línea con <pre lang=»haskell»> y otra con </pre>

Pensamiento

«Cambiemos nuestra actitud tradicional en la construcción de programas. En lugar de imaginar que nuestra tarea principal es indicarle a una computadora lo que debe hacer, concentrémonos más bien en explicarle a los seres humanos lo que queremos que haga una computadora.»

Donald Knuth.

Acotación del primorial

PorJosé A. Alonso 31-enero-20207-febrero-2020

El primorial de un número natural n es el producto de todos los números primos menores o iguales a n. Por ejemplo, el primorial de 5 es 30 porque el producto de los primos menores o iguales que 5 es

   2 * 3 * 5 = 30

1	2 * 3 * 5 = 30

La propiedad de Erdös de acotación de los primoriales afirma que

Para todo número natural n, su primorial es menor o igual que 4ⁿ.

Definir las funciones

   primorial :: Integer -> Integer
   primoriales :: [Integer]

1 2	primorial :: Integer -> Integer primoriales :: [Integer]

tales que

(primorial n) es el primorial de n. Por ejemplo,

     primorial 3  ==  6
     primorial 5  ==  30
     primorial 8  ==  210

primorial 3 == 6

primorial 5 == 30

primorial 8 == 210

primoriales es la sucesión de los primoriales. Por ejemplo,

   λ> take 15 primoriales
   [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]

1 2	λ> take 15 primoriales [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]

Comprobar con QuickCheck la propiedad de Erdös de acotación de los primoriales.

Soluciones

import Data.Numbers.Primes
import Test.QuickCheck

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

primorial :: Integer -> Integer
primorial n = product (takeWhile (<= n) primes)

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

primorial2 :: Integer -> Integer
primorial2 0 = 1
primorial2 n | gcd n x == 1 = n*x
             | otherwise    = x
  where x = primorial2 (n-1)

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

--    λ> length (show (primorial (5*10^5)))
--    216852
--    (1.65 secs, 2,472,977,584 bytes)
--    λ> length (show (primorial2 (5*10^5)))
--    216852
--    (3.56 secs, 2,719,162,272 bytes)

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

--    λ> take 15 primoriales
--    [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]
primoriales :: [Integer]
primoriales = map primorial [0..]

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

--    λ> take 15 primoriales2
--    [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]
primoriales2 :: [Integer]
primoriales2 = map primorial2 [0..]

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

--    λ> take 15 primoriales3
--    [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]
primoriales3 :: [Integer]
primoriales3 = scanl1 f [1..]
  where f x n | gcd n x == 1 = n*x
              | otherwise    = x

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

--    λ> minimum (take 5000 primoriales)
--    1
--    (1.56 secs, 4,857,760,464 bytes)
--    λ> minimum (take 5000 primoriales2)
--    1
--    (9.39 secs, 10,942,848,240 bytes)
--    λ> minimum (take 5000 primoriales3)
--    1
--    (0.01 secs, 5,575,024 bytes)
--    
--    λ> minimum (take 6000 primoriales)
--    1
--    (2.22 secs, 7,013,937,248 bytes)
--    λ> minimum (take 6000 primoriales3)
--    1
--    (0.01 secs, 6,737,328 bytes)

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

prop_primorial :: Integer -> Property
prop_primorial n =
  n >= 0 ==> primorial n <= 4^n

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

import Data.Numbers.Primes

import Test.QuickCheck

-- 1ª definición de primorial

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

primorial :: Integer -> Integer

primorial n = product (takeWhile (<= n) primes)

-- 2ª definición de primorial

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

primorial2 :: Integer -> Integer

primorial2 0 = 1

primorial2 n | gcd n x == 1 = n*x

| otherwise = x

where x = primorial2 (n-1)

-- Comparación de eficiencia

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

-- λ> length (show (primorial (5*10^5)))

-- 216852

-- (1.65 secs, 2,472,977,584 bytes)

-- λ> length (show (primorial2 (5*10^5)))

-- 216852

-- (3.56 secs, 2,719,162,272 bytes)

-- 1ª definición de primoriales

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

-- λ> take 15 primoriales

-- [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]

primoriales :: [Integer]

primoriales = map primorial [0..]

-- 2ª definición de primoriales

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

-- λ> take 15 primoriales2

-- [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]

primoriales2 :: [Integer]

primoriales2 = map primorial2 [0..]

-- 3ª definición de primoriales

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

-- λ> take 15 primoriales3

-- [1,1,2,6,6,30,30,210,210,210,210,2310,2310,30030,30030]

primoriales3 :: [Integer]

primoriales3 = scanl1 f [1..]

where f x n | gcd n x == 1 = n*x

| otherwise = x

-- Comparación de eficiencia

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

-- λ> minimum (take 5000 primoriales)

-- 1

-- (1.56 secs, 4,857,760,464 bytes)

-- λ> minimum (take 5000 primoriales2)

-- 1

-- (9.39 secs, 10,942,848,240 bytes)

-- λ> minimum (take 5000 primoriales3)

-- 1

-- (0.01 secs, 5,575,024 bytes)

-- λ> minimum (take 6000 primoriales)

-- 1

-- (2.22 secs, 7,013,937,248 bytes)

-- λ> minimum (take 6000 primoriales3)

-- 1

-- (0.01 secs, 6,737,328 bytes)

-- Propiedad

-- =========

prop_primorial :: Integer -> Property

prop_primorial n =

n >= 0 ==> primorial n <= 4^n

-- La comprobación es

-- λ> quickCheck prop_primorial

-- +++ OK, passed 100 tests.

Otras soluciones

Se pueden escribir otras soluciones en los comentarios.
El código se debe escribir entre una línea con <pre lang=»haskell»> y otra con </pre>

Pensamiento

«Las matemáticas son la reina de las ciencias y la teoría de los números es la reina de las matemáticas.»

Carl Friedrich Gauss.

La sucesión de Sylvester

Soluciones

Cálculo de pi mediante la variante de Euler de la serie armónica

Soluciones

Otras soluciones

Productos de dos y tres números consecutivos

Soluciones

Otras soluciones

Pensamiento

Números de Bell

Soluciones

Otras soluciones

Pensamiento

Acotación del primorial

Soluciones

Otras soluciones

Pensamiento

Entradas recientes

Buscar

Correo electrónico