Ejercicio

Ecuación diofántica con primos (OME1995 P4)

PorJosé A. Alonso 4-mayo-202111-mayo-2021

El enunciado del problema 4 de la OME (Olimpiada Matemática Española) del 1995 es

Siendo p un número primo, halla las soluciones enteras de la ecuación:

p.(x + y) = x.y

Definir la función

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

1	soluciones :: Integer -> [(Integer, Integer)]

tal que (soluciones p) es la lista de los pares de enteros (x,y) tales que p.(x + y) = x.y. Por ejemplo,

   λ> take 3 (soluciones 7)
   [(0,0),(14,14),(-42,6)]
   λ> sum [x+y | (x,y) <- take 6 (soluciones (primes !! (10^6)))]
   185830404

λ> take 3 (soluciones 7)

[(0,0),(14,14),(-42,6)]

λ> sum [x+y | (x,y) <- take 6 (soluciones (primes !! (10^6)))]

185830404

Soluciones

import Data.List (sort)
import Data.Numbers.Primes (primes)

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

soluciones :: Integer -> [(Integer, Integer)]
soluciones p =
  [(x,y) | (x,y) <- pares,
           p * (x + y) == x * y]

-- pares es la lista de los pares de números enteros ordenados por la
-- suma de los valores absolutos de sus componentes. Por ejemplo,
--    λ> take 38 pares
--    [(0,0),(-1,0),(0,-1),(0,1),(1,0),(-2,0),(-1,-1),(-1,1),(0,-2),
--     (0,2),(1,-1),(1,1),(2,0),(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),
--     (0,-3),(0,3),(1,-2),(1,2),(2,-1),(2,1),(3,0),(-4,0),(-3,-1),
--     (-3,1),(-2,-2),(-2,2),(-1,-3),(-1,3),(0,-4),(0,4),(1,-3),(1,3),
--     (2,-2),(2,2)]
pares :: [(Integer,Integer)]
pares = concatMap paresEnterosSuma [0..]

-- (paresEnterosSuma n) es la lista de pares de enteros (x,y) tales que la suma
-- de los valores absolutos de x e y es igual a n. Por ejemplo,
--    λ> paresEnterosSuma 3
--    [(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),(0,-3),(0,3),(1,-2),(1,2),
--     (2,-1),(2,1),(3,0)]
paresEnterosSuma :: Integer -> [(Integer,Integer)]
paresEnterosSuma n = concatMap aux [-n..n]
  where aux k | m == 0    = [(k,0)]
              | otherwise = [(k,-m),(k,m)]
          where m = n - abs k

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

-- Observando los siguientes cálculos:
--    take 6 (soluciones 2)  == [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]
--    take 6 (soluciones 3)  == [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]
--    take 6 (soluciones 5)  == [(0,0),(10,10),(-20,4),(4,-20),(6,30),(30,6)]
--    take 6 (soluciones 7)  == [(0,0),(14,14),(-42,6),(6,-42),(8,56),(56,8)]
--    take 6 (soluciones 11) == [(0,0),(22,22),(-110,10),(10,-110),(12,132),(132,12)]

soluciones2 :: Integer -> [(Integer, Integer)]
soluciones2 2 = [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]
soluciones2 3 = [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]
soluciones2 p =
  [(0,0), (2*p,2*p), (p*(1-p),p-1), (p-1,p*(1-p)), (p+1,p*(p+1)), (p*(p+1),p+1)]

-- Comparación de equivalencia
-- ===========================

-- La propiedad es
prop_equivalencia :: Bool
prop_equivalencia =
  and [ take 6 (soluciones p) == soluciones2 p
      | p <- takeWhile (<= 37) primes]

-- La comprobación es
--    λ> prop_equivalencia
--    True

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

-- La comparación es
--    λ> take 6 (soluciones 37)
--    [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]
--    (3.88 secs, 2,521,196,352 bytes)
--    λ> take 6 (soluciones2 37)
--    [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]
--    (0.01 secs, 144,192 bytes)

import Data.List (sort)

import Data.Numbers.Primes (primes)

-- 1ª solución

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

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

soluciones p =

[(x,y) | (x,y) <- pares,

p * (x + y) == x * y]

-- pares es la lista de los pares de números enteros ordenados por la

-- suma de los valores absolutos de sus componentes. Por ejemplo,

-- λ> take 38 pares

-- [(0,0),(-1,0),(0,-1),(0,1),(1,0),(-2,0),(-1,-1),(-1,1),(0,-2),

-- (0,2),(1,-1),(1,1),(2,0),(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),

-- (0,-3),(0,3),(1,-2),(1,2),(2,-1),(2,1),(3,0),(-4,0),(-3,-1),

-- (-3,1),(-2,-2),(-2,2),(-1,-3),(-1,3),(0,-4),(0,4),(1,-3),(1,3),

-- (2,-2),(2,2)]

pares :: [(Integer,Integer)]

pares = concatMap paresEnterosSuma [0..]

-- (paresEnterosSuma n) es la lista de pares de enteros (x,y) tales que la suma

-- de los valores absolutos de x e y es igual a n. Por ejemplo,

-- λ> paresEnterosSuma 3

-- [(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),(0,-3),(0,3),(1,-2),(1,2),

-- (2,-1),(2,1),(3,0)]

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

paresEnterosSuma n = concatMap aux [-n..n]

where aux k | m == 0 = [(k,0)]

| otherwise = [(k,-m),(k,m)]

where m = n - abs k

-- 2ª solución

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

-- Observando los siguientes cálculos:

-- take 6 (soluciones 2) == [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]

-- take 6 (soluciones 3) == [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]

-- take 6 (soluciones 5) == [(0,0),(10,10),(-20,4),(4,-20),(6,30),(30,6)]

-- take 6 (soluciones 7) == [(0,0),(14,14),(-42,6),(6,-42),(8,56),(56,8)]

-- take 6 (soluciones 11) == [(0,0),(22,22),(-110,10),(10,-110),(12,132),(132,12)]

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

soluciones2 2 = [(0,0),(-2,1),(1,-2),(4,4),(3,6),(6,3)]

soluciones2 3 = [(0,0),(-6,2),(2,-6),(6,6),(4,12),(12,4)]

soluciones2 p =

[(0,0), (2*p,2*p), (p*(1-p),p-1), (p-1,p*(1-p)), (p+1,p*(p+1)), (p*(p+1),p+1)]

-- Comparación de equivalencia

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

-- La propiedad es

prop_equivalencia :: Bool

prop_equivalencia =

and [ take 6 (soluciones p) == soluciones2 p

| p <- takeWhile (<= 37) primes]

-- La comprobación es

-- λ> prop_equivalencia

-- True

-- Comparación de eficiencia

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

-- La comparación es

-- λ> take 6 (soluciones 37)

-- [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]

-- (3.88 secs, 2,521,196,352 bytes)

-- λ> take 6 (soluciones2 37)

-- [(0,0),(74,74),(-1332,36),(36,-1332),(38,1406),(1406,38)]

-- (0.01 secs, 144,192 bytes)

Nuevas soluciones

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

Ordenación de los pares de enteros

PorJosé A. Alonso 3-mayo-202110-mayo-2021

Los pares de números enteros se pueden ordenar por la suma de los valores absolutos de sus componentes. Los primeros pares en dicha ordenación son

   (0,0),
   (-1,0),(0,-1),(0,1),(1,0),
   (-2,0),(-1,-1),(-1,1),(0,-2),(0,2),(1,-1),(1,1),(2,0), ...

(0,0),

(-1,0),(0,-1),(0,1),(1,0),

(-2,0),(-1,-1),(-1,1),(0,-2),(0,2),(1,-1),(1,1),(2,0), ...

Definir la lista

   pares :: [(Integer, Integer)]

1	pares :: [(Integer, Integer)]

cuyos elementos son los pares de números enteros con la ordenación anterior. Por ejemplo,

   λ> take 38 pares
   [(0,0),(-1,0),(0,-1),(0,1),(1,0),(-2,0),(-1,-1),(-1,1),(0,-2),
    (0,2),(1,-1),(1,1),(2,0),(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),
    (0,-3),(0,3),(1,-2),(1,2),(2,-1),(2,1),(3,0),(-4,0),(-3,-1),
    (-3,1),(-2,-2),(-2,2),(-1,-3),(-1,3),(0,-4),(0,4),(1,-3),(1,3),
    (2,-2),(2,2)]

λ> take 38 pares

[(0,0),(-1,0),(0,-1),(0,1),(1,0),(-2,0),(-1,-1),(-1,1),(0,-2),

(0,2),(1,-1),(1,1),(2,0),(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),

(0,-3),(0,3),(1,-2),(1,2),(2,-1),(2,1),(3,0),(-4,0),(-3,-1),

(-3,1),(-2,-2),(-2,2),(-1,-3),(-1,3),(0,-4),(0,4),(1,-3),(1,3),

(2,-2),(2,2)]

Soluciones

import Data.List (sort)
import Test.QuickCheck (Property, (==>), quickCheck)

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

pares :: [(Integer, Integer)]
pares =
  concat [sort (concat [aux (x,n-x) | x <- [0..n]])
         | n <- [0..]]
  where
    aux (0,0) = [(0,0)]
    aux (0,y) = [(0,y), (0,-y)]
    aux (x,0) = [(x,0), (-x,0)]
    aux (x,y) = [(x,y), (-x,y), (x,-y), (-x,-y)]

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

pares2 :: [(Integer,Integer)]
pares2 = concatMap paresEnterosSuma [0..]

-- (paresEnterosSuma n) es la lista de pares de enteros (x,y) tales que la suma
-- de los valores absolutos de x e y es igual a n. Por ejemplo,
--    λ> paresEnterosSuma 3
--    [(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),(0,-3),(0,3),(1,-2),(1,2),
--     (2,-1),(2,1),(3,0)]
paresEnterosSuma :: Integer -> [(Integer,Integer)]
paresEnterosSuma n = concatMap aux [-n..n]
  where aux k | m == 0    = [(k,0)]
              | otherwise = [(k,-m),(k,m)]
          where m = n - abs k

-- Comprobación de equivalencia
-- ============================

-- La propiedad es
prop_pares :: Int -> Property
prop_pares n =
  n >= 0 ==>
  pares !! n == pares2 !! n

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

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

-- La comparación es
--    λ> pares !! (10^7)
--    (304,-1932)
--    (9.15 secs, 11,031,475,520 bytes)
--    λ> pares2 !! (10^7)
--    (304,-1932)
--    (2.84 secs, 3,401,450,320 bytes)

import Data.List (sort)

import Test.QuickCheck (Property, (==>), quickCheck)

-- 1ª solución

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

pares :: [(Integer, Integer)]

pares =

concat [sort (concat [aux (x,n-x) | x <- [0..n]])

| n <- [0..]]

where

aux (0,0) = [(0,0)]

aux (0,y) = [(0,y), (0,-y)]

aux (x,0) = [(x,0), (-x,0)]

aux (x,y) = [(x,y), (-x,y), (x,-y), (-x,-y)]

-- 2ª solución

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

pares2 :: [(Integer,Integer)]

pares2 = concatMap paresEnterosSuma [0..]

-- (paresEnterosSuma n) es la lista de pares de enteros (x,y) tales que la suma

-- de los valores absolutos de x e y es igual a n. Por ejemplo,

-- λ> paresEnterosSuma 3

-- [(-3,0),(-2,-1),(-2,1),(-1,-2),(-1,2),(0,-3),(0,3),(1,-2),(1,2),

-- (2,-1),(2,1),(3,0)]

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

paresEnterosSuma n = concatMap aux [-n..n]

where aux k | m == 0 = [(k,0)]

| otherwise = [(k,-m),(k,m)]

where m = n - abs k

-- Comprobación de equivalencia

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

-- La propiedad es

prop_pares :: Int -> Property

prop_pares n =

n >= 0 ==>

pares !! n == pares2 !! n

-- La comprobación es

-- λ> quickCheck prop_pares

-- +++ OK, passed 100 tests.

-- Comparación de eficiencia

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

-- La comparación es

-- λ> pares !! (10^7)

-- (304,-1932)

-- (9.15 secs, 11,031,475,520 bytes)

-- λ> pares2 !! (10^7)

-- (304,-1932)

-- (2.84 secs, 3,401,450,320 bytes)

Nuevas soluciones

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

Números potencias perfectas de la suma de sus dígitos

PorJosé A. Alonso 30-abril-20217-mayo-2021

El número 2401 es una potencia de la suma de sus dígitos, ya que dicha suma es 7 y 7^4 = 2401.

Definir la lista

   potenciaSumaDigitos :: [Integer]

1	potenciaSumaDigitos :: [Integer]

cuyos elementos son los números que son potencias de las sumas de sus dígitos. Por ejemplo,

   λ> take 17 potenciaSumaDigitos
   [0,1,2,3,4,5,6,7,8,9,81,512,2401,4913,5832,17576,19683]

1 2	λ> take 17 potenciaSumaDigitos [0,1,2,3,4,5,6,7,8,9,81,512,2401,4913,5832,17576,19683]

Soluciones

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

potenciaSumaDigitos :: [Integer]
potenciaSumaDigitos = 0 : filter esPotenciaSumaDigitos [0..]

-- (esPotenciaSumaDigitos n) se verifica si n es una potencia de la suma
-- de sus dígitos. Por ejemplo,
--    esPotenciaSumaDigitos 2401  ==  True
--    esPotenciaSumaDigitos 2402  ==  False
esPotenciaSumaDigitos :: Integer -> Bool
esPotenciaSumaDigitos n =
  or [n == x^k | k <- [1..n]]
  where x = sumaDigitos n

-- (sumaDigitos n) es la suma de los dígitos de n. Por ejemplo,
--    sumaDigitos 2021  ==  5
sumaDigitos :: Integer -> Integer
sumaDigitos = sum . digitos

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,
--    digitos 2021  ==  [2,0,2,1]
digitos :: Integer -> [Integer]
digitos x = [read [c] | c <- show x]

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

potenciaSumaDigitos2 :: [Integer]
potenciaSumaDigitos2 = 0 : 1 : filter esPotenciaSumaDigitos2 [2..]

-- (esPotenciaSumaDigitos2 n) se verifica si n es una potencia de la suma
-- de sus dígitos. Por ejemplo,
--    esPotenciaSumaDigitos2 2401  ==  True
--    esPotenciaSumaDigitos2 2402  ==  False
esPotenciaSumaDigitos2 :: Integer -> Bool
esPotenciaSumaDigitos2 n =
  n == x^k
  where x = sumaDigitos n
        k = round (logBase (fromIntegral x) (fromIntegral n))

-- 1ª solución

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

potenciaSumaDigitos :: [Integer]

potenciaSumaDigitos = 0 : filter esPotenciaSumaDigitos [0..]

-- (esPotenciaSumaDigitos n) se verifica si n es una potencia de la suma

-- de sus dígitos. Por ejemplo,

-- esPotenciaSumaDigitos 2401 == True

-- esPotenciaSumaDigitos 2402 == False

esPotenciaSumaDigitos :: Integer -> Bool

esPotenciaSumaDigitos n =

or [n == x^k | k <- [1..n]]

where x = sumaDigitos n

-- (sumaDigitos n) es la suma de los dígitos de n. Por ejemplo,

-- sumaDigitos 2021 == 5

sumaDigitos :: Integer -> Integer

sumaDigitos = sum . digitos

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,

-- digitos 2021 == [2,0,2,1]

digitos :: Integer -> [Integer]

digitos x = [read [c] | c <- show x]

-- 2ª solución

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

potenciaSumaDigitos2 :: [Integer]

potenciaSumaDigitos2 = 0 : 1 : filter esPotenciaSumaDigitos2 [2..]

-- (esPotenciaSumaDigitos2 n) se verifica si n es una potencia de la suma

-- de sus dígitos. Por ejemplo,

-- esPotenciaSumaDigitos2 2401 == True

-- esPotenciaSumaDigitos2 2402 == False

esPotenciaSumaDigitos2 :: Integer -> Bool

esPotenciaSumaDigitos2 n =

n == x^k

where x = sumaDigitos n

k = round (logBase (fromIntegral x) (fromIntegral n))

Nuevas soluciones

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

Órbita con raíz entera (OME1997 P4)

PorJosé A. Alonso 29-abril-20216-mayo-2021

El enunciado del problema 4 de la OME (Olimpiada Matemática Española) del 1997 es

Sea p un número primo. Determinar todos los enteros k tales que sqrt(k² – k*p) es natural.

Definir las funciones

   orbita        :: Integer -> [Integer]
   orbitaDePrimo :: Integer -> [Integer]

1 2	orbita :: Integer -> [Integer] orbitaDePrimo :: Integer -> [Integer]

tales que

(orbita n) es la lista de todos los enteros k tales que sqrt(k² – k*n) es natural. Por ejemplo,

     take 4  (orbita 6)   == [0,-2,6,8]
     take 5  (orbita 36)  == [0,-12,36,48,-64]
     take 6  (orbita 9)   == [0,-3,9,12,-16,25]
     take 8  (orbita 27)  == [0,-9,27,36,-48,75,-169,196]
     take 10 (orbita 111) == [0,-37,111,148,-289,400,-972,1083,-3025,3136]

take 4 (orbita 6) == [0,-2,6,8]

take 5 (orbita 36) == [0,-12,36,48,-64]

take 6 (orbita 9) == [0,-3,9,12,-16,25]

take 8 (orbita 27) == [0,-9,27,36,-48,75,-169,196]

take 10 (orbita 111) == [0,-37,111,148,-289,400,-972,1083,-3025,3136]

(orbitaDePrimo p) es la lista de todos los enteros k tales que sqrt(k² – k*p) es natural, suponiendo que p es un número primo. Por ejemplo,

     orbitaDePrimo 5                  == [0,-4,5,9]
     orbitaDePrimo (primes !! (10^6)) == [0,15485867,-59953011442489,59953026928356]

1 2	orbitaDePrimo 5 == [0,-4,5,9] orbitaDePrimo (primes !! (10^6)) == [0,15485867,-59953011442489,59953026928356]

Soluciones

import Data.Numbers.Primes (primes)

orbita :: Integer -> [Integer]
orbita n =
  [k | k <- enteros,
       k^2 - k*n >= 0,
       esCuadrado (k^2 - k*n)]

-- entero es la lista de los números enteros. Por ejemplo,
--    λ> take 20 enteros
--    [0,-1,1,-2,2,-3,3,-4,4,-5,5,-6,6,-7,7,-8,8,-9,9,-10]
enteros :: [Integer]
enteros = 0 : concat [[-n,n] | n <- [1..]]

-- (esCuadrado x) se verifica si x es un cuadrado perfecto. Por
-- ejemplo,
--    esCuadrado 16  ==  True
--    esCuadrado 27  ==  False
esCuadrado :: Integer -> Bool
esCuadrado x =
  (raizEntera x)^2 == x

-- (raizEntera x) es el mayor entero cuyo cuadrado es menor o igual que
-- x. Por ejemplo,
--    raizEntera 16  ==  4
--    raizEntera 27  ==  5
raizEntera :: Integer -> Integer
raizEntera x = aux (1,x)
    where aux (a,b) | d == x    = c
                    | c == a    = c
                    | d < x     = aux (c,b)
                    | otherwise = aux (a,c)
              where c = (a+b) `div` 2
                    d = c^2

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

orbitaDePrimo1 :: Integer -> [Integer]
orbitaDePrimo1 2 = take 2 (orbita 2)
orbitaDePrimo1 p = take 4 (orbita p)

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

-- Basada en los siguientes cálculos
--    orbitaDePrimo1 2  == [0,2]
--    orbitaDePrimo1 3  == [0,-1,3,4]
--    orbitaDePrimo1 5  == [0,-4,5,9]
--    orbitaDePrimo1 7  == [0, 7,  -9, 16]
--    orbitaDePrimo1 11 == [0,11, -25, 36]
--    orbitaDePrimo1 13 == [0,13, -36, 49]
--    orbitaDePrimo1 17 == [0,17, -64, 81]
--    orbitaDePrimo1 19 == [0,19, -81,100]
--    orbitaDePrimo1 23 == [0,23,-121,144]

orbitaDePrimo2 :: Integer -> [Integer]
orbitaDePrimo2 p
  | p == 2    = [0,2]
  | p <= 5    = [0, -((p-1) `div` 2)^2, p, ((p+1) `div` 2)^2]
  | otherwise = [0, p, -((p-1) `div` 2)^2, ((p+1) `div` 2)^2]

-- Comprobación de equivalencia
-- ============================

-- La comprobación es
--    λ> and [orbitaDePrimo1 n == orbitaDePrimo2 n | n <- take 30 primes]
--    True

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

-- La comparación es
--    λ> orbitaDePrimo1 (primes !! 100)
--    [0,547,-74529,75076]
--    (4.94 secs, 4,471,368,256 bytes)
--    λ> orbitaDePrimo2 (primes !! 100)
--    [0,547,-74529,75076]
--    (0.01 secs, 302,096 bytes)

import Data.Numbers.Primes (primes)

orbita :: Integer -> [Integer]

orbita n =

[k | k <- enteros,

k^2 - k*n >= 0,

esCuadrado (k^2 - k*n)]

-- entero es la lista de los números enteros. Por ejemplo,

-- λ> take 20 enteros

-- [0,-1,1,-2,2,-3,3,-4,4,-5,5,-6,6,-7,7,-8,8,-9,9,-10]

enteros :: [Integer]

enteros = 0 : concat [[-n,n] | n <- [1..]]

-- (esCuadrado x) se verifica si x es un cuadrado perfecto. Por

-- ejemplo,

-- esCuadrado 16 == True

-- esCuadrado 27 == False

esCuadrado :: Integer -> Bool

esCuadrado x =

(raizEntera x)^2 == x

-- (raizEntera x) es el mayor entero cuyo cuadrado es menor o igual que

-- x. Por ejemplo,

-- raizEntera 16 == 4

-- raizEntera 27 == 5

raizEntera :: Integer -> Integer

raizEntera x = aux (1,x)

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

| c == a = c

| d < x = aux (c,b)

| otherwise = aux (a,c)

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

d = c^2

-- 1ª definición de orbitaDePrimo

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

orbitaDePrimo1 :: Integer -> [Integer]

orbitaDePrimo1 2 = take 2 (orbita 2)

orbitaDePrimo1 p = take 4 (orbita p)

-- 2ª definición de orbitaDePrimo

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

-- Basada en los siguientes cálculos

-- orbitaDePrimo1 2 == [0,2]

-- orbitaDePrimo1 3 == [0,-1,3,4]

-- orbitaDePrimo1 5 == [0,-4,5,9]

-- orbitaDePrimo1 7 == [0, 7, -9, 16]

-- orbitaDePrimo1 11 == [0,11, -25, 36]

-- orbitaDePrimo1 13 == [0,13, -36, 49]

-- orbitaDePrimo1 17 == [0,17, -64, 81]

-- orbitaDePrimo1 19 == [0,19, -81,100]

-- orbitaDePrimo1 23 == [0,23,-121,144]

orbitaDePrimo2 :: Integer -> [Integer]

orbitaDePrimo2 p

| p == 2 = [0,2]

| p <= 5 = [0, -((p-1) `div` 2)^2, p, ((p+1) `div` 2)^2]

| otherwise = [0, p, -((p-1) `div` 2)^2, ((p+1) `div` 2)^2]

-- Comprobación de equivalencia

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

-- La comprobación es

-- λ> and [orbitaDePrimo1 n == orbitaDePrimo2 n | n <- take 30 primes]

-- True

-- Comparación de eficiencia

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

-- La comparación es

-- λ> orbitaDePrimo1 (primes !! 100)

-- [0,547,-74529,75076]

-- (4.94 secs, 4,471,368,256 bytes)

-- λ> orbitaDePrimo2 (primes !! 100)

-- [0,547,-74529,75076]

-- (0.01 secs, 302,096 bytes)

Nuevas soluciones

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

Números iguales a potencias de las sumas de sus cifras (OME1999 P2)

PorJosé A. Alonso 28-abril-20215-mayo-2021

El enunciado del problema 2 de la OME (Olimpiada Matemática Española) del 1998 es

Hallar todos los números naturales de 4 cifras, escritos en base 10, que sean iguales al cubo de la suma de sus cifras.

Definir la función

   especiales :: Integer -> Integer -> [Integer]

1	especiales :: Integer -> Integer -> [Integer]

tal que (especiales a b) es la lista de los números de a cifras que son iguales la suma de sus cifras elevada a b. Por ejemplo,

   especiales 5 3  ==  [17576,19683]
   especiales 6 4  ==  [234256,390625,614656]

1 2	especiales 5 3 == [17576,19683] especiales 6 4 == [234256,390625,614656]

Usando la función anterior, calcular las soluciones del problema de la Olimpiada.

Soluciones

import Data.Char (digitToInt)

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

especiales :: Integer -> Integer -> [Integer]
especiales a b =
  [n | n <- [10^(a-1)..10^a-1],
       n == (sumaDigitos n)^b]

-- (sumaDigitos n) es la suma de los dígitos de n. Por ejemplo,
--    sumaDigitos 2021  ==  5
sumaDigitos :: Integer -> Integer
sumaDigitos = sum . digitos

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,
--    digitos 2021  ==  [2,0,2,1]
digitos :: Integer -> [Integer]
digitos x = [read [c] | c <- show x]

-- Cálculo de la solución del problema de la Olimpiada:
--    λ> especiales 4 3
--    [4913,5832]

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

especiales2 :: Integer -> Integer -> [Integer]
especiales2 c e = [z | z <- map (^e) [a..min (9*c) b], z == (f z)^e]
  where [c', e'] = map fromIntegral [c,e]
        [a , b ] = [ceiling (10**((c'+p)/e')) | p <- [-1,0]]
        f = sum . map (toInteger . digitToInt) . show

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

-- La comparación es
--    λ> especiales 6 4
--    [234256,390625,614656]
--    (7.41 secs, 21,808,919,208 bytes)
--    λ> especiales2 6 4
--    [234256,390625,614656]
--    (0.01 secs, 168,072 bytes)

import Data.Char (digitToInt)

-- 1ª solución

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

especiales :: Integer -> Integer -> [Integer]

especiales a b =

[n | n <- [10^(a-1)..10^a-1],

n == (sumaDigitos n)^b]

-- (sumaDigitos n) es la suma de los dígitos de n. Por ejemplo,

-- sumaDigitos 2021 == 5

sumaDigitos :: Integer -> Integer

sumaDigitos = sum . digitos

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,

-- digitos 2021 == [2,0,2,1]

digitos :: Integer -> [Integer]

digitos x = [read [c] | c <- show x]

-- Cálculo de la solución del problema de la Olimpiada:

-- λ> especiales 4 3

-- [4913,5832]

-- 2ª solución

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

especiales2 :: Integer -> Integer -> [Integer]

especiales2 c e = [z | z <- map (^e) [a..min (9*c) b], z == (f z)^e]

where [c', e'] = map fromIntegral [c,e]

[a , b ] = [ceiling (10**((c'+p)/e')) | p <- [-1,0]]

f = sum . map (toInteger . digitToInt) . show

-- Comparación de eficiencia

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

-- La comparación es

-- λ> especiales 6 4

-- [234256,390625,614656]

-- (7.41 secs, 21,808,919,208 bytes)

-- λ> especiales2 6 4

-- [234256,390625,614656]

-- (0.01 secs, 168,072 bytes)

Nuevas soluciones

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

Máximos de una función recursiva (OME2002 P3)

PorJosé A. Alonso 27-abril-20214-mayo-2021

El enunciado del problema 5 de la OME (Olimpiada Matemática Española) del 2002 es

La función g se define sobre los números naturales y satisface las condiciones:

g(1) = 1

g(2n) = g(n)

g(2n + 1) = g(2n) + 1

Sea n un número natural tal que 1 ≤ n ≤ 2002. Calcula el valor máximo M de g(n). Calcula también cuántos valores de n satisfacen g(n) = M.

Los valores de la función g para n de 1 a 30 son

   1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4

1	1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4

Definir la función

   maximoG :: Integer -> Integer

1	maximoG :: Integer -> Integer

tal que (maximoG m) es el máximo de los valores de g(n) para n en {1, 2,…, m}. Por ejemplo,

   maximoG 30           ==  4
   maximoG (10^(10^5))  ==  332192

1 2	maximoG 30 == 4 maximoG (10^(10^5)) == 332192

Usando la función maximoG, calcular los valores pedidos en el problema.

Soluciones

import Data.List (genericLength, genericTake, group)
import Test.QuickCheck (Property, (==>), quickCheck)

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

maximoG :: Integer -> Integer
maximoG m = maximum (genericTake m sucesionG)

-- sucesionG es la lista de los valores de g. Por ejemplo,
--    λ> take 30 sucesionG
--    [1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4]
sucesionG :: [Integer]
sucesionG = map g [1..]

-- (g n) es el valor de g(n).
g :: Integer -> Integer
g 1 = 1
g n | even n    = g (n `div` 2)
    | otherwise = g (n `div` 2) + 1

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

-- Observando los siguientes cálculos
--   λ> map maximoG [1..40]
--   [1,1,2,2,2,2,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,5,5,5,5,5,5,5,5,5,5]
--   λ> take 10 (map length (group (map maximoG [1..])))
--   [2,4,8,16,32,64,128,256,512,1024]

maximoG2 :: Integer -> Integer
maximoG2 m = head [x | x <- [1..], m + 1 < 2^x] - 1

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

maximoG3 :: Integer -> Integer
maximoG3 m = genericLength (takeWhile (<m+2) potenciasDeDos) - 1

-- potenciasDeDos es la lista de las potencias de dos. Por ejemplo,
--    take 10 potenciasDeDos  ==  [1,2,4,8,16,32,64,128,256,512]
potenciasDeDos :: [Integer]
potenciasDeDos = iterate (*2) 1

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

maximoG4 :: Integer -> Integer
maximoG4 m = floor (logBase 2 (fromIntegral (m+1)))

-- Comprobación de equivalencia
-- ============================

-- La propiedad es
prop_maximoG :: Integer -> Property
prop_maximoG m =
  m > 0 ==>
  all (== (maximoG m))
      [ maximoG2 m,
        maximoG3 m,
        maximoG4 m]

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

-- Nota: Aunque QuickCheck no ha encontrado ningún contraejemplo, la
-- definición maximoG4 que usa números decimales falla para números muy
-- grandes. Por ejemplo,
--    λ> maximoG4 (10^308)
--    1023
--    λ> maximoG4 (10^309)
--    1797693134862315907729305190789024733617976978942306572734300...
--    λ> maximoG3 (10^308)
--    1023
--    λ> maximoG3 (10^309)
--    1026

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

-- La comparación es
--    λ> maximoG (10^6)
--    19
--    (8.80 secs, 6,015,772,760 bytes)
--    λ> maximoG2 (10^6)
--    19
--    (0.02 secs, 106,752 bytes)
--    λ> maximoG3 (10^6)
--    19
--    (0.02 secs, 102,592 bytes)
--    λ> maximoG4 (10^6)
--    19
--    (0.01 secs, 98,464 bytes)
--
--    λ> maximoG2 (10^308)
--    1023
--    (0.03 secs, 3,266,096 bytes)
--    λ> maximoG3 (10^308)
--    1023
--    (0.02 secs, 266,856 bytes)
--    λ> maximoG4 (10^308)
--    1023
--    (0.02 secs, 103,176 bytes)
--
--    λ> maximoG2 (10^16789)
--    55771
--    (1.70 secs, 1,201,022,600 bytes)
--    λ> maximoG3 (10^16789)
--    55771
--    (0.24 secs, 214,872,864 bytes)

-- Cálculos del problema
-- =====================

-- El máximo se calcula como sigue:
--    λ> maximum (take 2002 sucesionG)
--    10
-- Por tanto, el máximo es 10.

-- Los valores de n tales que g(n) es el máximo se calcula como sigue:
--    λ> [n | n <- [1..2002], g n == 10]
--    [1023,1535,1791,1919,1983]

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

import Data.List (genericLength, genericTake, group)

import Test.QuickCheck (Property, (==>), quickCheck)

-- 1ª solución

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

maximoG :: Integer -> Integer

maximoG m = maximum (genericTake m sucesionG)

-- sucesionG es la lista de los valores de g. Por ejemplo,

-- λ> take 30 sucesionG

-- [1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,1,2,2,3,2,3,3,4,2,3,3,4,3,4,4]

sucesionG :: [Integer]

sucesionG = map g [1..]

-- (g n) es el valor de g(n).

g :: Integer -> Integer

g 1 = 1

g n | even n = g (n `div` 2)

| otherwise = g (n `div` 2) + 1

-- 2ª solución

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

-- Observando los siguientes cálculos

-- λ> map maximoG [1..40]

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

-- λ> take 10 (map length (group (map maximoG [1..])))

-- [2,4,8,16,32,64,128,256,512,1024]

maximoG2 :: Integer -> Integer

maximoG2 m = head [x | x <- [1..], m + 1 < 2^x] - 1

-- 3ª solución

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

maximoG3 :: Integer -> Integer

maximoG3 m = genericLength (takeWhile (<m+2) potenciasDeDos) - 1

-- potenciasDeDos es la lista de las potencias de dos. Por ejemplo,

-- take 10 potenciasDeDos == [1,2,4,8,16,32,64,128,256,512]

potenciasDeDos :: [Integer]

potenciasDeDos = iterate (*2) 1

-- 4ª solución

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

maximoG4 :: Integer -> Integer

maximoG4 m = floor (logBase 2 (fromIntegral (m+1)))

-- Comprobación de equivalencia

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

-- La propiedad es

prop_maximoG :: Integer -> Property

prop_maximoG m =

m > 0 ==>

all (== (maximoG m))

[ maximoG2 m,

maximoG3 m,

maximoG4 m]

-- La comprobación es

-- λ> quickCheck prop_maximoG

-- +++ OK, passed 100 tests.

-- Nota: Aunque QuickCheck no ha encontrado ningún contraejemplo, la

-- definición maximoG4 que usa números decimales falla para números muy

-- grandes. Por ejemplo,

-- λ> maximoG4 (10^308)

-- 1023

-- λ> maximoG4 (10^309)

-- 1797693134862315907729305190789024733617976978942306572734300...

-- λ> maximoG3 (10^308)

-- 1023

-- λ> maximoG3 (10^309)

-- 1026

-- Comparación de eficiencia

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

-- La comparación es

-- λ> maximoG (10^6)

-- 19

-- (8.80 secs, 6,015,772,760 bytes)

-- λ> maximoG2 (10^6)

-- 19

-- (0.02 secs, 106,752 bytes)

-- λ> maximoG3 (10^6)

-- 19

-- (0.02 secs, 102,592 bytes)

-- λ> maximoG4 (10^6)

-- 19

-- (0.01 secs, 98,464 bytes)

-- λ> maximoG2 (10^308)

-- 1023

-- (0.03 secs, 3,266,096 bytes)

-- λ> maximoG3 (10^308)

-- 1023

-- (0.02 secs, 266,856 bytes)

-- λ> maximoG4 (10^308)

-- 1023

-- (0.02 secs, 103,176 bytes)

-- λ> maximoG2 (10^16789)

-- 55771

-- (1.70 secs, 1,201,022,600 bytes)

-- λ> maximoG3 (10^16789)

-- 55771

-- (0.24 secs, 214,872,864 bytes)

-- Cálculos del problema

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

-- El máximo se calcula como sigue:

-- λ> maximum (take 2002 sucesionG)

-- 10

-- Por tanto, el máximo es 10.

-- Los valores de n tales que g(n) es el máximo se calcula como sigue:

-- λ> [n | n <- [1..2002], g n == 10]

-- [1023,1535,1791,1919,1983]

Nuevas soluciones

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

Productos de cuatro consecutivos (OME2006 P5)

PorJosé A. Alonso 26-abril-20213-mayo-2021

El enunciado del problema 5 de la OME (Olimpiada Matemática Española) del 2006 es

Probar que el producto de cuatro naturales consecutivos no puede ser ni cuadrado ni cubo perfecto.

Definir la lista

   productos :: [Integer]

1	productos :: [Integer]

cuyos elementos son los productos de cuatro enteros positivos consecutivos. Por ejemplo,

   λ> take 12 productos
   [24,120,360,840,1680,3024,5040,7920,11880,17160,24024,32760]

1 2	λ> take 12 productos [24,120,360,840,1680,3024,5040,7920,11880,17160,24024,32760]

Comprobar con QuickCheck que los elementos de la lista productos no son ni cuadrados ni cubos perfectos.

Soluciones

import Test.QuickCheck

productos :: [Integer]
productos =
  [product [n..n+3] | n <- [1..]]

-- La propiedad es
prop_productos :: Int -> Property
prop_productos n =
  n >= 0 ==>
  not (esCuadrado x) && not (esCubo x)
  where x = productos !! n

-- (esCuadrado x) se verifica si x es un cuadrado perfecto. Por
-- ejemplo,
--    esCuadrado 16  ==  True
--    esCuadrado 27  ==  False
esCuadrado :: Integer -> Bool
esCuadrado x = (raizEntera 2 x)^2 == x

-- (esCubo x) se verifica si x es un cubo perfecto. Por ejemplo,
--    esCubo 8  ==  True
--    esCubo 9  ==  False
esCubo :: Integer -> Bool
esCubo x = (raizEntera 3 x)^3 == x

-- (raizEntera n x) es el mayor entero cuya potencia n-ésima es menor o
-- igual que x. Por ejemplo,
--    raizEntera 2 15  ==  3
--    raizEntera 2 16  ==  4
--    raizEntera 2 17  ==  4
raizEntera :: Int -> Integer -> Integer
raizEntera n x = aux (1,x)
    where aux (a,b) | d == x    = c
                    | c == a    = c
                    | d < x     = aux (c,b)
                    | otherwise = aux (a,c)
              where c = (a+b) `div` 2
                    d = c^n

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

import Test.QuickCheck

productos :: [Integer]

productos =

[product [n..n+3] | n <- [1..]]

-- La propiedad es

prop_productos :: Int -> Property

prop_productos n =

n >= 0 ==>

not (esCuadrado x) && not (esCubo x)

where x = productos !! n

-- (esCuadrado x) se verifica si x es un cuadrado perfecto. Por

-- ejemplo,

-- esCuadrado 16 == True

-- esCuadrado 27 == False

esCuadrado :: Integer -> Bool

esCuadrado x = (raizEntera 2 x)^2 == x

-- (esCubo x) se verifica si x es un cubo perfecto. Por ejemplo,

-- esCubo 8 == True

-- esCubo 9 == False

esCubo :: Integer -> Bool

esCubo x = (raizEntera 3 x)^3 == x

-- (raizEntera n x) es el mayor entero cuya potencia n-ésima es menor o

-- igual que x. Por ejemplo,

-- raizEntera 2 15 == 3

-- raizEntera 2 16 == 4

-- raizEntera 2 17 == 4

raizEntera :: Int -> Integer -> Integer

raizEntera n x = aux (1,x)

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

| c == a = c

| d < x = aux (c,b)

| otherwise = aux (a,c)

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

d = c^n

-- La comprobación es

-- λ> quickCheck prop_productos

-- +++ OK, passed 100 tests.

Nuevas soluciones

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

Mayores potencias de 5 que dividen a la sucesión (OME2014 P4)

PorJosé A. Alonso 23-abril-202130-abril-2021

El enunciado del problema 4 de la OME (Olimpiada Matemática Española) del 2014 es

Sea {x(n)} la sucesión de enteros positivos definida por x(1) = 2 y x(n+1) = 2*x(n)^3+x(n) para todo n >= 1. Determinar la mayor potencia de 5 que divide al número x(2014)^2+1.

Definir la función

   mayorExponente :: Integer -> Integer

1	mayorExponente :: Integer -> Integer

tal que (mayorExponente n) es el mayor m tal que 5^m divide al número x(n)^2+1.

Soluciones

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

mayorExponente :: Integer -> Integer
mayorExponente n =
  mayorExponenteDe5QueDivide ((termino n)^2 + 1)

-- (termino n) es el término n-ésimo de la sucesión. Por ejemplo,
--    termino 2  ==  18
--    termino 3  ==  11682
--    termino 4  ==  3188464624818
termino :: Integer -> Integer
termino 1 = 2
termino n = 2*x^3+x
  where x = termino (n-1)

-- (mayorExponenteDe5QueDivide n) es el mayor exponente de 5 que divide a n. Por
-- ejemplo,
--    mayorExponenteDe5QueDivide 50  ==  2
mayorExponenteDe5QueDivide :: Integer -> Integer
mayorExponenteDe5QueDivide n
  | n `mod` 5 /= 0 = 0
  | otherwise      = 1 + mayorExponenteDe5QueDivide (n `div` 5)

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

-- Con el siguente cálculo
--    λ> [mayorExponente n | n <- [1..5]]
--    [1,2,3,4,5]
-- se observa que el mayor exponente de 5 que divide a x(n)^2+1 es n.

mayorExponente2 :: Integer -> Integer
mayorExponente2 = id

-- 1ª solución

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

mayorExponente :: Integer -> Integer

mayorExponente n =

mayorExponenteDe5QueDivide ((termino n)^2 + 1)

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

-- termino 2 == 18

-- termino 3 == 11682

-- termino 4 == 3188464624818

termino :: Integer -> Integer

termino 1 = 2

termino n = 2*x^3+x

where x = termino (n-1)

-- (mayorExponenteDe5QueDivide n) es el mayor exponente de 5 que divide a n. Por

-- ejemplo,

-- mayorExponenteDe5QueDivide 50 == 2

mayorExponenteDe5QueDivide :: Integer -> Integer

mayorExponenteDe5QueDivide n

| n `mod` 5 /= 0 = 0

| otherwise = 1 + mayorExponenteDe5QueDivide (n `div` 5)

-- 2ª solución

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

-- Con el siguente cálculo

-- λ> [mayorExponente n | n <- [1..5]]

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

-- se observa que el mayor exponente de 5 que divide a x(n)^2+1 es n.

mayorExponente2 :: Integer -> Integer

mayorExponente2 = id

Nuevas soluciones

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

Permutaciones divisibles (OME2016 P5)

PorJosé A. Alonso 22-abril-202129-abril-2021

El enunciado del problema 5 de la OME (Olimpiada Matemática Española) del 2016 es

De entre todas las permutaciones (a(1), a(2),…, a(n)) del conjunto {1, 2,…, n},(n ≥ 1 entero), se consideran las que cumplen que 2(a(1) + a(2) +···+ a(m)) es divisible por m, para cada m = 1, 2,…, n. Calcular el número total de estas permutaciones.

Llamaremos permutaciones divisibles a las que cumplen la propiedad anterior. Por ejemplo, [2,3,4,1] es una permutación divisible de {1,2,3,4} ya que es una permutación del conjunto y se cumplen las condiciones:

2*2 = 4 es divisible por 1,
2*(2+3) = 10 es divisible por 2
2*(2+3+4) = 18 es divisible por 3.
2*(2+3+4+1) = 20 es divisible por 4.

Definir las siguientes funciones:

   permutacionesDivisibles  :: Integer -> [[Integer]]
   nPermutacionesDivisibles :: Integer -> Integer

1 2	permutacionesDivisibles :: Integer -> [[Integer]] nPermutacionesDivisibles :: Integer -> Integer

tales que

(permutacionesDivisibles n) es la lista de las permutaciones divisibles de {1,2,…,n}. Por ejemplo,

     λ> permutacionesDivisibles 2
     [[1,2],[2,1]]
     λ> permutacionesDivisibles 3
     [[1,2,3],[1,3,2],[2,1,3],[2,3,1],[3,1,2],[3,2,1]]
     λ> permutacionesDivisibles 4
     [[1,2,3,4],[1,3,2,4],[2,1,3,4],[2,3,1,4],[2,3,4,1],[2,4,3,1],
      [3,1,2,4],[3,2,1,4],[3,2,4,1],[3,4,2,1],[4,2,3,1],[4,3,2,1]]
     λ> length (permutacionesDivisibles 20)
     786432

λ> permutacionesDivisibles 2

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

λ> permutacionesDivisibles 3

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

λ> permutacionesDivisibles 4

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

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

λ> length (permutacionesDivisibles 20)

786432

(nPermutacionesDivisibles n) es el número de permutaciones divisibles de {1,2,…,n}. Por ejemplo,

     nPermutacionesDivisibles 4  ==  12
     nPermutacionesDivisibles (10^8) `mod` (10^9)  ==  340332032

1 2	nPermutacionesDivisibles 4 == 12 nPermutacionesDivisibles (10^8) `mod` (10^9) == 340332032

Soluciones

import Data.List (genericLength, permutations, sort)

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

permutacionesDivisibles :: Integer -> [[Integer]]
permutacionesDivisibles n =
  sort (filter aux (permutations [1..n]))
  where
    aux xs = and [(2*x) `mod` n == 0 | (x,n) <- zip (sumasParciales xs) [1..]]

-- (sumasParciales xs) es la lista de las suas parciales de xs. Por ejemplo,
--    sumasParciales [1..9]  ==  [1,3,6,10,15,21,28,36,45]
sumasParciales :: [Integer] -> [Integer]
sumasParciales = scanl1 (+)

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

permutacionesDivisibles2 :: Integer -> [[Integer]]
permutacionesDivisibles2 = sort . aux
  where aux n
          | n < 4     = permutations [1..n]
          | otherwise = [xs ++ [n] | xs <- xss] ++
                        [map (+1) xs ++ [1] | xs <- xss]
          where xss = aux (n-1)

-- Comprobación de equivalencia
-- ============================

-- La comprobación para los n primeros valores es
--    λ> and [permutacionesDivisibles2 n == permutacionesDivisibles n | n <- [1..9]]
--    True

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

-- La comparación es
--    λ> length (permutacionesDivisibles 10)
--    768
--    (8.42 secs, 10,420,281,776 bytes)
--    λ> length (permutacionesDivisibles2 10)
--    768
--    (0.02 secs, 2,250,152 bytes)
--
--    λ> length (permutacionesDivisibles2 20)
--    786432
--    (14.33 secs, 4,458,839,736 bytes)

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

--    nPermutacionesDivisibles 5  ==  24
nPermutacionesDivisibles :: Integer -> Integer
nPermutacionesDivisibles =
  genericLength . permutacionesDivisibles

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

nPermutacionesDivisibles2 :: Integer -> Integer
nPermutacionesDivisibles2 =
  genericLength . permutacionesDivisibles2

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

-- Observando los siguientes cálculos:
--    λ> [nPermutacionesDivisibles2 n | n <- [1..12]]
--    [1,2,6,12,24,48,96,192,384,768,1536,3072]
--    λ> 1 : 2 : [3*2^(n-2) | n <- [3..12]]
--    [1,2,6,12,24,48,96,192,384,768,1536,3072]

nPermutacionesDivisibles3 :: Integer -> Integer
nPermutacionesDivisibles3 n
  | n < 3     = n
  | otherwise = 3*2^(n-2)

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

-- La comparación es
--    λ> nPermutacionesDivisibles 10
--    768
--    (7.95 secs, 10,704,422,592 bytes)
--    λ> nPermutacionesDivisibles2 10
--    768
--    (0.01 secs, 2,269,872 bytes)
--    λ> nPermutacionesDivisibles3 10
--    768
--    (0.01 secs, 102,560 bytes)
--
--    λ> nPermutacionesDivisibles2 20
--    786432
--    (13.00 secs, 4,477,717,968 bytes)
--    λ> nPermutacionesDivisibles3 20
--    786432
--    (0.01 secs, 106,848 bytes)

import Data.List (genericLength, permutations, sort)

-- 1ª definición de permutacionesDivisibles

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

permutacionesDivisibles :: Integer -> [[Integer]]

permutacionesDivisibles n =

sort (filter aux (permutations [1..n]))

where

aux xs = and [(2*x) `mod` n == 0 | (x,n) <- zip (sumasParciales xs) [1..]]

-- (sumasParciales xs) es la lista de las suas parciales de xs. Por ejemplo,

-- sumasParciales [1..9] == [1,3,6,10,15,21,28,36,45]

sumasParciales :: [Integer] -> [Integer]

sumasParciales = scanl1 (+)

-- 2ª definición de permutacionesDivisibles

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

permutacionesDivisibles2 :: Integer -> [[Integer]]

permutacionesDivisibles2 = sort . aux

where aux n

| n < 4 = permutations [1..n]

| otherwise = [xs ++ [n] | xs <- xss] ++

[map (+1) xs ++ [1] | xs <- xss]

where xss = aux (n-1)

-- Comprobación de equivalencia

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

-- La comprobación para los n primeros valores es

-- λ> and [permutacionesDivisibles2 n == permutacionesDivisibles n | n <- [1..9]]

-- True

-- Comparación de eficiencia

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

-- La comparación es

-- λ> length (permutacionesDivisibles 10)

-- 768

-- (8.42 secs, 10,420,281,776 bytes)

-- λ> length (permutacionesDivisibles2 10)

-- 768

-- (0.02 secs, 2,250,152 bytes)

-- λ> length (permutacionesDivisibles2 20)

-- 786432

-- (14.33 secs, 4,458,839,736 bytes)

-- 1ª definición de nPermutacionesDivisibles

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

-- nPermutacionesDivisibles 5 == 24

nPermutacionesDivisibles :: Integer -> Integer

nPermutacionesDivisibles =

genericLength . permutacionesDivisibles

-- 2ª definición de nPermutacionesDivisibles

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

nPermutacionesDivisibles2 :: Integer -> Integer

nPermutacionesDivisibles2 =

genericLength . permutacionesDivisibles2

-- 3ª definición de nPermutacionesDivisibles

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

-- Observando los siguientes cálculos:

-- λ> [nPermutacionesDivisibles2 n | n <- [1..12]]

-- [1,2,6,12,24,48,96,192,384,768,1536,3072]

-- λ> 1 : 2 : [3*2^(n-2) | n <- [3..12]]

-- [1,2,6,12,24,48,96,192,384,768,1536,3072]

nPermutacionesDivisibles3 :: Integer -> Integer

nPermutacionesDivisibles3 n

| n < 3 = n

| otherwise = 3*2^(n-2)

-- Comparación de eficiencia

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

-- La comparación es

-- λ> nPermutacionesDivisibles 10

-- 768

-- (7.95 secs, 10,704,422,592 bytes)

-- λ> nPermutacionesDivisibles2 10

-- 768

-- (0.01 secs, 2,269,872 bytes)

-- λ> nPermutacionesDivisibles3 10

-- 768

-- (0.01 secs, 102,560 bytes)

-- λ> nPermutacionesDivisibles2 20

-- 786432

-- (13.00 secs, 4,477,717,968 bytes)

-- λ> nPermutacionesDivisibles3 20

-- 786432

-- (0.01 secs, 106,848 bytes)

Nuevas soluciones

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

Números fibonaccianos

PorJosé A. Alonso 21-abril-202128-abril-2021

El enunciado del segundo problema de este mes de la RSME es el siguiente:

Un número de al menos tres cifras se denomina fibonacciano si sus cifras, a partir de la tercera, son iguales a la suma de las dos cifras anteriores. Por ejemplo, 5279 es un número fibonacciano, pues su tercera cifra, 7, es suma de las dos anteriores (5+2) y su cuarta cifra, 9, también (2+7).

Te daremos el problema por válido si respondes bien a estas dos cuestiones:
a) ¿cuántas cifras como máximo puede tener un número fibonacciano?
b) ¿cuántos números fibonaccianos hay?

En la definición de fibonacciano la suma de las cifras tiene que menor que 10, pero podemos generalizarlo sustituyendo 10 por número n. Dichos números de llaman fibonaccianos generalizados acotados por n. Por ejemplo, 571219315081 es un fibonacciano generalizado acotado por 100 ya que la sucesión de sus dígitos es 5, 7, 12 (= 5+7), 19 (= 7+12), 31 (= 12+19) 50 (=19+31) y 81 (=31+50).

Definir las funciones

   esFibonacciano :: Integer -> Bool
   fibonaccianos  :: [Integer]
   fibonaccianosG :: Integer -> [Integer]

esFibonacciano :: Integer -> Bool

fibonaccianos :: [Integer]

fibonaccianosG :: Integer -> [Integer]

tales que

(esFibonacciano n) se verifica si n es un número fibonacciano. Por ejemplo,

     esFibonacciano 5279    ==  True
     esFibonacciano 527916  ==  False

1 2	esFibonacciano 5279 == True esFibonacciano 527916 == False

fibonaccianos es la lista de los números fibonaccianos. Por ejemplo,

     λ> take 60 fibonaccianos
     [101,112,123,134,145,156,167,178,189,202,213,224,235,246,257,268,
      279,303,314,325,336,347,358,369,404,415,426,437,448,459,505,516,
      527,538,549,606,617,628,639,707,718,729,808,819,909,1011,1123,
      1235,1347,1459,2022,2134,2246,2358,3033,3145,3257,3369,4044,4156]

λ> take 60 fibonaccianos

[101,112,123,134,145,156,167,178,189,202,213,224,235,246,257,268,

279,303,314,325,336,347,358,369,404,415,426,437,448,459,505,516,

527,538,549,606,617,628,639,707,718,729,808,819,909,1011,1123,

1235,1347,1459,2022,2134,2246,2358,3033,3145,3257,3369,4044,4156]

(fibonaccianosG n) es la lista de los números fibonaccianos generalizados acotados por n. Por ejemplo,

     λ> take 60 (fibonaccianosG 100)
     [101,112,123,134,145,156,167,178,189,202,213,224,235,246,257,268,
      279,303,314,325,336,347,358,369,404,415,426,437,448,459,505,516,
      527,538,549,606,617,628,639,707,718,729,808,819,909,1011,1123,
      1235,1347,1459,1910,2022,2134,2246,2358,2810,2911,3033,3145,3257]
     λ> take 12 (drop 60 (fibonaccianosG 10))
     [4268,5055,5167,5279,6066,6178,7077,7189,8088,9099,10112,11235]
     λ> take 12 (drop 60 (fibonaccianosG 100))
     [3369,3710,3811,3912,4044,4156,4268,4610,4711,4812,4913,5055]
     λ> length (fibonaccianosG (10^40))
     16888
     λ> length (show (last (fibonaccianosG (10^40))))
     3943

λ> take 60 (fibonaccianosG 100)

[101,112,123,134,145,156,167,178,189,202,213,224,235,246,257,268,

279,303,314,325,336,347,358,369,404,415,426,437,448,459,505,516,

527,538,549,606,617,628,639,707,718,729,808,819,909,1011,1123,

1235,1347,1459,1910,2022,2134,2246,2358,2810,2911,3033,3145,3257]

λ> take 12 (drop 60 (fibonaccianosG 10))

[4268,5055,5167,5279,6066,6178,7077,7189,8088,9099,10112,11235]

λ> take 12 (drop 60 (fibonaccianosG 100))

[3369,3710,3811,3912,4044,4156,4268,4610,4711,4812,4913,5055]

λ> length (fibonaccianosG (10^40))

16888

λ> length (show (last (fibonaccianosG (10^40))))

3943

Usando las funciones anteriores, calcular cuántas cifras como máximo puede tener un número fibonacciano y cuántos números fibonaccianos hay.

Soluciones

import Data.List (inits, isPrefixOf, sort)

-- Definición de esFibonacciano
-- ============================

esFibonacciano :: Integer -> Bool
esFibonacciano n =
  n > 99 && ds `isPrefixOf` drop 2 fs
  where (a:b:ds) = digitos n
        fs = a : b : zipWith (+) fs (tail fs)

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,
--    digitos 325  ==  [3,2,5]
digitos :: Integer -> [Int]
digitos n =
  [read [c] | c <- show n]

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

fibonaccianos :: [Integer]
fibonaccianos =
  filter esFibonacciano [100..10112358]

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

fibonaccianos2 :: [Integer]
fibonaccianos2 =
  sort (concat [aux a b | a <- [1..9], b <- [0..9]])
  where
    aux a b = [digitosAnumero zs | zs <- drop 3 (inits (takeWhile (<10) (fibs a b)))]

-- (fibs a b) es la sucesión de números de Fibonacci cuyos dos primeros
-- términos son a y b. Por ejemplo,
--    take 9 (fibs 2 4)  ==  [2,4,6,10,16,26,42,68,110]
--    take 9 (fibs 3 6)  ==  [3,6,9,15,24,39,63,102,165]
fibs :: Integer -> Integer -> [Integer]
fibs a b = fs
  where fs = a : b : zipWith (+) fs (tail fs)

-- (digitosAnumero xs) es el número cuya lista de dígitos es xs. Por
-- ejemplo,
--    digitosAnumero [8,1,6,4,9]  ==  81649
digitosAnumero :: [Integer] -> Integer
digitosAnumero xs =
  read (concatMap show xs)

-- fibonaccianosG
-- ==============

fibonaccianosG :: Integer -> [Integer]
fibonaccianosG n =
  sort (concat [fibonaccianosGAux n a b | a <- [1..9], b <- [0..9]])

fibonaccianosGAux :: Integer -> Integer -> Integer -> [Integer]
fibonaccianosGAux n a b =
  [digitosAnumero zs | zs <- drop 3 (inits (takeWhile (<n) (fibs a b)))]

-- Comprobación de la generalización
--    λ> fibonaccianos2 == fibonaccianosG 10
--    True

--    λ> length (fibonaccianosG (10^40))
--    16888
--    (2.41 secs, 10,546,873,616 bytes)
--    λ> length (show (last (fibonaccianosG (10^40))))
--    3943
--    (2.41 secs, 10,547,101,856 bytes)

-- Cálculos
-- ========

-- El cálculo del máximo número de cifras de los números fibonaccianos
-- es
--      λ> maximum [length (show n) | n <- fibonaccianos2]
--      8
--      λ> length (show (last fibonaccianos2))
--      8

-- El cálculo de la cantidad de números fibonaccianos es
--      λ> length fibonaccianos2
--      84

import Data.List (inits, isPrefixOf, sort)

-- Definición de esFibonacciano

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

esFibonacciano :: Integer -> Bool

esFibonacciano n =

n > 99 && ds `isPrefixOf` drop 2 fs

where (a:b:ds) = digitos n

fs = a : b : zipWith (+) fs (tail fs)

-- (digitos n) es la lista de los dígitos de n. Por ejemplo,

-- digitos 325 == [3,2,5]

digitos :: Integer -> [Int]

digitos n =

[read [c] | c <- show n]

-- 1ª definición de fibonaccianos

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

fibonaccianos :: [Integer]

fibonaccianos =

filter esFibonacciano [100..10112358]

-- 2ª definición de fibonaccianos

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

fibonaccianos2 :: [Integer]

fibonaccianos2 =

sort (concat [aux a b | a <- [1..9], b <- [0..9]])

where

aux a b = [digitosAnumero zs | zs <- drop 3 (inits (takeWhile (<10) (fibs a b)))]

-- (fibs a b) es la sucesión de números de Fibonacci cuyos dos primeros

-- términos son a y b. Por ejemplo,

-- take 9 (fibs 2 4) == [2,4,6,10,16,26,42,68,110]

-- take 9 (fibs 3 6) == [3,6,9,15,24,39,63,102,165]

fibs :: Integer -> Integer -> [Integer]

fibs a b = fs

where fs = a : b : zipWith (+) fs (tail fs)

-- (digitosAnumero xs) es el número cuya lista de dígitos es xs. Por

-- ejemplo,

-- digitosAnumero [8,1,6,4,9] == 81649

digitosAnumero :: [Integer] -> Integer

digitosAnumero xs =

read (concatMap show xs)

-- fibonaccianosG

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

fibonaccianosG :: Integer -> [Integer]

fibonaccianosG n =

sort (concat [fibonaccianosGAux n a b | a <- [1..9], b <- [0..9]])

fibonaccianosGAux :: Integer -> Integer -> Integer -> [Integer]

fibonaccianosGAux n a b =

[digitosAnumero zs | zs <- drop 3 (inits (takeWhile (<n) (fibs a b)))]

-- Comprobación de la generalización

-- λ> fibonaccianos2 == fibonaccianosG 10

-- True

-- λ> length (fibonaccianosG (10^40))

-- 16888

-- (2.41 secs, 10,546,873,616 bytes)

-- λ> length (show (last (fibonaccianosG (10^40))))

-- 3943

-- (2.41 secs, 10,547,101,856 bytes)

-- Cálculos

-- ========

-- El cálculo del máximo número de cifras de los números fibonaccianos

-- es

-- λ> maximum [length (show n) | n <- fibonaccianos2]

-- 8

-- λ> length (show (last fibonaccianos2))

-- 8

-- El cálculo de la cantidad de números fibonaccianos es

-- λ> length fibonaccianos2

-- 84

Nuevas soluciones

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

Ecuación diofántica con primos (OME1995 P4)

Soluciones

Nuevas soluciones

Ordenación de los pares de enteros

Soluciones

Nuevas soluciones

Números potencias perfectas de la suma de sus dígitos

Soluciones

Nuevas soluciones

Órbita con raíz entera (OME1997 P4)

Soluciones

Nuevas soluciones

Números iguales a potencias de las sumas de sus cifras (OME1999 P2)

Soluciones

Nuevas soluciones

Máximos de una función recursiva (OME2002 P3)

Soluciones

Nuevas soluciones

Productos de cuatro consecutivos (OME2006 P5)

Soluciones

Nuevas soluciones

Mayores potencias de 5 que dividen a la sucesión (OME2014 P4)

Soluciones

Nuevas soluciones

Permutaciones divisibles (OME2016 P5)

Soluciones

Nuevas soluciones

Números fibonaccianos

Soluciones

Nuevas soluciones

Entradas recientes

Buscar

Correo electrónico