Menu Close

Etiqueta: tail

Cálculo de pi mediante el método de Newton

El método de Newton para el cálculo de pi se basa en la relación
Calculo_de_pi_mediante_el_metodo_de_Newton_1
y en el desarrollo del arco seno
Calculo_de_pi_mediante_el_metodo_de_Newton_2
de donde se obtiene la fórmula
Calculo_de_pi_mediante_el_metodo_de_Newton_3

La primeras aproximaciones son

   a(0) = 6*(1/2)                               = 3.0
   a(1) = 6*(1/2+1/(2*3*2^3))                   = 3.125
   a(2) = 6*(1/2+1/(2*3*2^3)+(1*3)/(2*4*5*2^5)) = 3.1390625

Definir las funciones

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

tales que

  • (aproximacionPi n) es la n-ésima aproximación de pi con la fórmula de Newton. Por ejemplo,
     aproximacionPi 0   ==  3.0
     aproximacionPi 1   ==  3.125
     aproximacionPi 2   ==  3.1390625
     aproximacionPi 10  ==  3.1415926468755613
     aproximacionPi 21  ==  3.141592653589793
     pi                 ==  3.141592653589793
  • (grafica xs) dibuja la gráfica de las k-ésimas aproximaciones de pi donde k toma los valores de la lista xs. Por ejemplo, (grafica [1..30]) dibuja
    Calculo_de_pi_mediante_el_metodo_de_Newton_4

Soluciones

import Graphics.Gnuplot.Simple
 
-- 1ª definición
-- =============
 
aproximacionPi :: Int -> Double
aproximacionPi n = 6 * arcsinX
  where arcsinX = 0.5 + sum (take n factoresN)
 
factoresN :: [Double]
factoresN = zipWith (*) (potenciasK 3) fraccionesPI
 
potenciasK :: Double -> [Double]
potenciasK k = (0.5**k)/k : potenciasK (k+2)
 
fraccionesPI :: [Double]
fraccionesPI =
  scanl (*) (1/2) (tail (zipWith (/) [1,3..] [2,4..]))
 
-- 2ª definición
-- =============
 
aproximacionPi2 :: Int -> Double
aproximacionPi2 n = 6 * (serie !! n)
 
serie :: [Double]
serie = scanl1 (+) (zipWith (/)
                            (map fromIntegral numeradores)
                            (map fromIntegral denominadores))
  where numeradores    = 1 : scanl1 (*) [1,3..]
        denominadores  = zipWith (*) denominadores1 denominadores2
        denominadores1 = 2 : scanl1 (*) [2,4..]
        denominadores2 = 1 : [n * 2^n | n <- [3,5..]]
 
grafica :: [Int] -> IO ()
grafica xs = 
    plotList [Key Nothing]
             [(k,aproximacionPi k) | k <- xs]

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

“Mi trabajo siempre trató de unir lo verdadero con lo bello; pero cuando tuve que elegir uno u otro, generalmente elegí lo bello.”

Hermann Weyl.

Medias de dígitos de pi

El fichero Digitos_de_pi.txt contiene el número pi con un millón de decimales; es decir,

   3.1415926535897932384626433832 ... 83996346460422090106105779458151

Definir las funciones

   mediasDigitosDePi        :: IO [Double]
   graficaMediasDigitosDePi :: Int -> IO ()

tales que

  • mediasDigitosDePi es la sucesión cuyo n-ésimo elemento es la media de los n primeros dígitos de pi. Por ejemplo,
     λ> xs <- mediasDigitosDePi
     λ> take 5 xs
     [1.0,2.5,2.0,2.75,4.0]
     λ> take 10 xs
     [1.0,2.5,2.0,2.75,4.0,3.6666666666666665,4.0,4.125,4.0,4.1]
     λ> take 10 <$> mediasDigitosDePi
     [1.0,2.5,2.0,2.75,4.0,3.6666666666666665,4.0,4.125,4.0,4.1]
  • (graficaMediasDigitosDePi n) dibuja la gráfica de los n primeros términos de mediasDigitosDePi. Por ejemplo,
    • (graficaMediasDigitosDePi 20) dibuja
    • (graficaMediasDigitosDePi 200) dibuja
    • (graficaMediasDigitosDePi 2000) dibuja

Soluciones

import Data.Char (digitToInt)
import Data.List (genericLength, inits, tails)
import Graphics.Gnuplot.Simple ( Attribute (Key, PNG)
                               , plotList )
 
-- Definición de mediasDigitosDePi
-- ===============================
 
mediasDigitosDePi :: IO [Double]
mediasDigitosDePi = do
  (_:_:ds) <- readFile "Digitos_de_pi.txt"
  return (medias (digitos ds))
 
-- (digitos cs) es la lista de los digitos de cs. Por ejempplo,
--    digitos "1415926535"  ==  [1,4,1,5,9,2,6,5,3,5]
digitos :: String -> [Int]
digitos = map digitToInt
 
-- (medias xs) es la lista de las medias de los segmentos iniciales de
-- xs. Por ejemplo,
--    λ> medias [1,4,1,5,9,2,6,5,3,5]
--    [1.0,2.5,2.0,2.75,4.0,3.6666666666666665,4.0,4.125,4.0,4.1]
medias :: [Int] -> [Double]
medias xs = map media (tail (inits xs))
 
-- (media xs) es la media aritmética de xs. Por ejemplo,
--    media [1,4,1,5,9]  ==  4.0
media :: [Int] -> Double
media xs = fromIntegral (sum xs) / genericLength xs
 
-- Definición de graficaMediasDigitosDePi
-- ======================================
 
graficaMediasDigitosDePi :: Int -> IO ()
graficaMediasDigitosDePi n = do
  xs <- mediasDigitosDePi
  plotList [ Key Nothing
           , PNG ("Medias_de_digitos_de_pi_" ++ show n ++ ".png")
           ]
           (take n xs)

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

Es el mejor de los buenos
quien sabe que en esta vida
todo es cuestión de medida:
un poco más, algo menos.

Antonio Machado

Primero no consecutivo

Definir la función

   primeroNoConsecutivo :: (Eq a,Enum a) => [a] -> Maybe a

tal que (primeroNoConsecutivo xs) es el primer elemento de la lista xs que no es igual al siguiente de su elemento anterior en xs o Nothing si tal elemento no existe. Por ejemplo

   primeroNoConsecutivo [1,2,3,4,6,7,8] == Just 6
   primeroNoConsecutivo "bcdfg"         == Just 'f'
   primeroNoConsecutivo "bcdef"         == Nothing

Soluciones

import Data.Maybe (listToMaybe)
 
-- 1ª solución
primeroNoConsecutivo1 :: (Eq a, Enum a) => [a] -> Maybe a
primeroNoConsecutivo1 xs
  | null ys   = Nothing
  | otherwise = Just (head ys)
  where ys = [y | (z,y) <- zip xs (tail xs), y /= succ z]
 
-- 2ª solución
primeroNoConsecutivo2 :: (Eq a, Enum a) => [a] -> Maybe a
primeroNoConsecutivo2 xs = 
  listToMaybe [y | (z,y) <- zip xs (tail xs), y /= succ z]
 
-- 3ª solución
primeroNoConsecutivo3 :: (Eq a,Enum a) => [a] -> Maybe a
primeroNoConsecutivo3 (x:y:zs)
  | succ x /= y = Just y 
  | otherwise   = primeroNoConsecutivo3 (y:zs)
primeroNoConsecutivo3 _ = Nothing
 
-- 4ª solución
primeroNoConsecutivo :: (Eq a,Enum a) => [a] -> Maybe a
primeroNoConsecutivo [] = Nothing
primeroNoConsecutivo (x:ys) = aux x ys
  where aux _ [] = Nothing
        aux x' (z:zs) | z == succ x' = aux z zs
                      | otherwise    = Just z

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

“La única enseñanza que un profesor puede dar, en mi opinión, es la de pensar delante de sus alumnos.”

Henri Lebesgue.

Producto de Fibonaccis consecutivos

Los números de Fibonacci son los números F(n) de la siguiente sucesión

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

que comienza con 0 y 1 y los siguientes términos son las sumas de los dos anteriores.

Un número x es el producto de dos números de Fibonacci consecutivos si existe un n tal que

   F(n) * F(n+1) = x

y su prueba es (F(n),F(n+1),True). Por ejemplo, 714 es el producto de dos números de Fibonacci consecutivos ya que

F(8) = 21, F(9) = 34 y 714 = 21 * 34.

Su prueba es (21, 34, True).

Un número x no es el producto de dos números de Fibonacci consecutivos si no existe un n tal que

   F(n) * F(n+1) = x

y su prueba es (F(m),F(m+1),False) donde m es el menor número tal que

   F(m) * F(m+1) > x

Por ejemplo, 800 no es el producto de dos números de Fibonacci consecutivos, ya que

 F(8) = 21, F(9) = 34, F(10) = 55 y 21 * 34 < 800 < 34 * 55.

Su prueba es (34, 55, False),

Definir la función

   productoFib :: Integer -> (Integer, Integer, Bool)

tal que (productoFib x) es la prueba de que es, o no es, el producto de dos números de Fibonacci consecutivos. Por ejemplo,

   productoFib 714  == (21,  34, True)
   productoFib 800  == (34,  55, False)
   productoFib 4895 == (55,  89, True)
   productoFib 5895 == (89, 144, False)

Soluciones

-- 1ª solución
-- ===========
 
productoFib :: Integer -> (Integer, Integer, Bool)
productoFib n
  | c == n    = (a,b,True)
  | otherwise = (a,b,False)
  where (a,b,c) = head (dropWhile (\(x,y,z) -> z < n) productos) 
 
-- fibs es la sucesión de números de Fibonacci. Por ejemplo,
--    take 14 fibs  ==  [0,1,1,2,3,5,8,13,21,34,55,89,144,233]
fibs :: [Integer]
fibs = 0 : 1 : zipWith (+) fibs (tail fibs)
 
-- productos es la lista de las ternas (a,b,c) tales que a y b son dos
-- números de Fibonacci consecutivos y c es su producto. Por ejemplo,
--    λ> take 7 productos
--    [(0,1,0),(1,1,1),(1,2,2),(2,3,6),(3,5,15),(5,8,40),(8,13,104)]
productos :: [(Integer,Integer,Integer)]
productos = [(x,y,x*y) | (x,y) <- zip fibs (tail fibs)] 
 
-- 2ª solución
-- ===========
 
productoFib2 :: Integer -> (Integer, Integer, Bool)
productoFib2 n = aux 0 1 n
  where
    aux a b c
        | a * b >= c = (a, b, a * b == c)
        | otherwise  = aux b (a + b) c
 
-- 3ª solución
-- ===========
 
productoFib3 :: Integer -> (Integer, Integer, Bool)
productoFib3 x = aux 0 1
  where
    aux a b | a * b >= x = (a, b, x == a * b)
            | otherwise  = aux b (a + b)
 
-- Comparación de eficiencia
-- =========================
 
-- La comparación es
--    λ> let (x,_,_) = productoFib (10^20000) in length (show x)
--    10000
--    (1.15 secs, 323,396,360 bytes)
--    λ> let (x,_,_) = productoFib2 (10^20000) in length (show x)
--    10000
--    (1.10 secs, 317,268,672 bytes)
--    λ> let (x,_,_) = productoFib3 (10^20000) in length (show x)
--    10000
--    (1.08 secs, 314,972,440 bytes)

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 placer que obtenemos de la música proviene de contar, pero contando inconscientemente. La música no es más que aritmética inconsciente.”

Gottfried Wilhelm Leibniz.

La conjetura de Mertens

Un número entero n es libre de cuadrados si no existe un número primo p tal que p² divide a n; es decir, los factores primos de n son todos distintos.

La función de Möbius μ(n) está definida para todos los enteros positivos como sigue:

  • μ(n) = 1 si n es libre de cuadrados y tiene un número par de factores primos.
  • μ(n) = -1 si n es libre de cuadrados y tiene un número impar de factores primos.
  • μ(n) = 0 si n no es libre de cuadrados.

Sus primeros valores son 1, -1, -1, 0, -1, 1, -1, 0, 0, 1, …

La función de Mertens M(n) está definida para todos los enteros positivos como la suma de μ(k) para 1 ≤ k ≤ n. Sus primeros valores son 1, 0, -1, -1, -2, -1, -2, -2, …

La conjetura de Mertens afirma que

Para todo entero x mayor que 1, el valor absoluto de la función de Mertens en x es menor que la raíz cuadrada de x.

La conjetura fue planteada por Franz Mertens en 1897. Riele Odlyzko, demostraronen 1985 que la conjetura de Mertens deja de ser cierta más o menos a partir de 10^{10^{64}}, cifra que luego de algunos refinamientos se redujo a 10^{10^{40}}.

Definir las funciones

   mobius :: Integer -> Integer
   mertens :: Integer -> Integer
   graficaMertens :: Integer -> IO ()

tales que

  • (mobius n) es el valor de la función de Möbius en n. Por ejemplo,
     mobius 6   ==  1
     mobius 30  ==  -1
     mobius 12  ==  0
  • (mertens n) es el valor de la función de Mertens en n. Por ejemplo,
     mertens 1     ==  1
     mertens 2     ==  0
     mertens 3     ==  -1
     mertens 5     ==  -2
     mertens 661   ==  -11
     mertens 1403  ==  11
  • (graficaMertens n) dibuja la gráfica de la función de Mertens, la raíz cuadrada y el opuestos de la raíz cuadrada para los n primeros n enteros positivos. Por ejemplo, (graficaMertens 1000) dibuja

Comprobar con QuickCheck la conjetura de Mertens.

Nota: El ejercicio está basado en La conjetura de Merterns y su relación con un número tan raro como extremada y colosalmente grande publicado por @Alvy la semana pasada en Microsiervos.

Soluciones

import Data.Numbers.Primes (primeFactors)
import Test.QuickCheck
import Graphics.Gnuplot.Simple
 
mobius :: Integer -> Integer
mobius n | tieneRepetidos xs = 0
         | otherwise         = (-1)^(length xs)
  where xs = primeFactors n
 
tieneRepetidos :: [Integer] -> Bool
tieneRepetidos xs =
  or [x == y | (x,y) <- zip xs (tail xs)]
 
mertens :: Integer -> Integer
mertens n = sum (map mobius [1..n])
 
-- Definición de graficaMertens
-- ============================
 
graficaMertens :: Integer -> IO ()
graficaMertens n = do
  plotLists [ Key Nothing
            , Title "Conjetura de Mertens"
            , PNG "La_conjetura_de_Mertens.png"
            ]
            [ [mertens k | k <- [1..n]]
            , raices
            , map negate raices
            ]
 
  where
    raices = [ceiling (sqrt k) | k <- [1..fromIntegral n]]
 
-- Conjetura de Mertens
-- ====================
 
-- La conjetura es
conjeturaDeMertens :: Integer -> Property
conjeturaDeMertens n =
  n > 1
  ==>
  abs (mertens n) < ceiling (sqrt n')
  where n' = fromIntegral n
 
-- La comprobación es
--    λ> quickCheck conjeturaDeMertens
--    +++ 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

“El control de la complejidad es la esencia de la programación informática.”

Brian Kernighan.

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

Infinitud de primos gemelos

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

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)

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.

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

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

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)

Por ejemplo,

   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)

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.

Referencia

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

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]

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
  • 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]
  • 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]

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.

Pensamiento

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

Antonio Machado

Mayor producto de n dígitos consecutivos de un número

Definir la función

   mayorProducto :: Int -> Integer -> Integer

tal que (mayorProducto n x) es el mayor producto de n dígitos consecutivos del número x (suponiendo que x tiene al menos n dígitos). Por ejemplo,

   mayorProducto 2 325                  ==  10
   mayorProducto 5 11111                ==  1
   mayorProducto 5 113111               ==  3
   mayorProducto 5 110111               ==  0
   mayorProducto 5 10151112             ==  10
   mayorProducto 5 101511124            ==  10
   mayorProducto 5 (product [1..1000])  ==  41472

Nota: Este ejercicio está basado en el problema 8 del Proyecto Euler

Soluciones

import Data.List (inits, tails)
import Data.Char (digitToInt)
 
-- 1ª solución
-- ===========
 
mayorProducto :: Int -> Integer -> Integer
mayorProducto n x =
  maximum [product xs | xs <- segmentos n (digitos x)]
 
-- (digitos x) es la lista de las digitos del número x. Por ejemplo, 
--    digitos 325  ==  [3,2,5]
digitos :: Integer -> [Integer]
digitos x = map (toInteger . digitToInt) (show x)
 
-- (segmentos n xs) es la lista de los segmentos de longitud n de la
-- lista xs. Por ejemplo,
--    segmentos 2 [3,5,4,6]  ==  [[3,5],[5,4],[4,6]]
segmentos :: Int -> [Integer] -> [[Integer]]
segmentos n xs = take (length xs - n + 1) (map (take n) (tails xs))
 
-- 2ª solución
-- ===========
 
mayorProducto2 :: Int -> Integer -> Integer
mayorProducto2 n x = maximum (aux ns)
    where ns     = [read [d] | d <- show x]
          aux xs | length xs < n = []
                 | otherwise     = product (take n xs) : aux (tail xs)
 
-- 3ª solución
-- ===========
 
mayorProducto3 :: Int -> Integer -> Integer
mayorProducto3 n = maximum
                 . map (product . take n)
                 . filter ((>=n) . length) 
                 . tails
                 . digitos
 
-- 4ª solución
-- ===========
 
mayorProducto4 :: Int -> Integer -> Integer
mayorProducto4 n = maximum  
                 . map (product . map (fromIntegral . digitToInt)) 
                 . filter ((==n) . length) 
                 . concatMap inits
                 . tails 
                 . show
 
-- Comparación de eficiencia
-- =========================
 
--    λ> mayorProducto 5 (product [1..500])
--    28224
--    (0.01 secs, 1,645,256 bytes)
--    λ> mayorProducto2 5 (product [1..500])
--    28224
--    (0.03 secs, 5,848,416 bytes)
--    λ> mayorProducto3 5 (product [1..500])
--    28224
--    (0.03 secs, 1,510,640 bytes)
--    λ> mayorProducto4 5 (product [1..500])
--    28224
--    (1.85 secs, 10,932,551,216 bytes)
--    
--    λ> mayorProducto 5 (product [1..7000])
--    46656
--    (0.10 secs, 68,590,808 bytes)
--    λ> mayorProducto2 5 (product [1..7000])
--    46656
--    (1.63 secs, 157,031,432 bytes)
--    λ> mayorProducto3 5 (product [1..7000])
--    46656
--    (1.55 secs, 65,727,176 bytes)

Pensamiento

“El control de la complejidad es la esencia de la programación.” ~ B.W. Kernigan