Menu Close

Categoría: Medio

Operaciones con series de potencias

José A. Alonso---5 Comentarios

Una serie de potencias es una serie de la forma

   a₀ + a₁x + a₂x² + a₃x³ + ...

Las series de potencias se pueden representar mediante listas infinitas. Por ejemplo, la serie de la función exponencial es

   e^x = 1 + x + x²/2! + x³/3! + ...

y se puede representar por [1, 1, 1/2, 1/6, 1/24, 1/120, …]

Las operaciones con series se pueden ver como una generalización de las de los polinomios.

En lo que sigue, usaremos el tipo (Serie a) para representar las series de potencias con coeficientes en a y su definición es

   type Serie a = [a]

Definir las siguientes funciones

   opuesta      :: Num a => Serie a -> Serie a
   suma         :: Num a => Serie a -> Serie a -> Serie a
   resta        :: Num a => Serie a -> Serie a -> Serie a
   producto     :: Num a => Serie a -> Serie a -> Serie a
   cociente     :: Fractional a => Serie a -> Serie a -> Serie a
   derivada     :: (Num a, Enum a) => Serie a -> Serie a
   integral     :: (Fractional a, Enum a) => Serie a -> Serie a
   expx         :: Serie Rational

tales que

  • (opuesta xs) es la opuesta de la serie xs. Por ejemplo,
     λ> take 7 (opuesta [-6,-4..])
     [6,4,2,0,-2,-4,-6]
  • (suma xs ys) es la suma de las series xs e ys. Por ejemplo,
     λ> take 7 (suma [1,3..] [2,4..])
     [3,7,11,15,19,23,27]
  • (resta xs ys) es la resta de las series xs es ys. Por ejemplo,
     λ> take 7 (resta [3,5..] [2,4..])
     [1,1,1,1,1,1,1]
     λ> take 7 (resta ([3,7,11,15,19,23,27] ++ repeat 0) [1,3..])
     [2,4,6,8,10,12,14]
  • (producto xs ys) es el producto de las series xs e ys. Por ejemplo,
     λ> take 7 (producto [3,5..] [2,4..])
     [6,22,52,100,170,266,392]
  • (cociente xs ys) es el cociente de las series xs e ys. Por ejemplo, 
     λ> take 7 (cociente ([6,22,52,100,170,266,392] ++ repeat 0) [3,5..])
     [2.0,4.0,6.0,8.0,10.0,12.0,14.0]
  • (derivada xs) es la derivada de la serie xs. Por ejemplo,
     λ> take 7 (derivada [2,4..])
     [4,12,24,40,60,84,112]
  • (integral xs) es la integral de la serie xs. Por ejemplo, 
     λ> take 7 (integral ([4,12,24,40,60,84,112] ++ repeat 0))
     [0.0,4.0,6.0,8.0,10.0,12.0,14.0]
  • expx es la serie de la función exponencial. Por ejemplo,
     λ> take 8 expx
     [1 % 1,1 % 1,1 % 2,1 % 6,1 % 24,1 % 120,1 % 720,1 % 5040]
     λ> take 8 (derivada expx)
     [1 % 1,1 % 1,1 % 2,1 % 6,1 % 24,1 % 120,1 % 720,1 % 5040]
     λ> take 8 (integral expx)
     [0 % 1,1 % 1,1 % 2,1 % 6,1 % 24,1 % 120,1 % 720,1 % 5040]

Soluciones

type Serie a = [a] 
 
opuesta :: Num a => Serie a -> Serie a
opuesta = map negate
 
suma :: Num a => Serie a -> Serie a -> Serie a
suma = zipWith (+)
 
resta :: Num a => Serie a -> Serie a -> Serie a
resta xs ys = suma xs (opuesta ys)
 
producto :: Num a => Serie a -> Serie a -> Serie a
producto (x:xs) zs@(y:ys) = 
    x*y : suma (producto xs zs) (map (x*) ys)
 
cociente :: Fractional a => Serie a -> Serie a -> Serie a
cociente (x:xs) (y:ys) = zs 
    where zs = x/y : map (/y) (resta xs (producto zs ys))  
 
derivada :: (Num a, Enum a) => Serie a -> Serie a
derivada (_:xs) = zipWith (*) xs [1..]
 
integral :: (Fractional a, Enum a) => Serie a -> Serie a
integral xs = 0 : zipWith (/) xs [1..]
 
expx :: Serie Rational
expx = map (1/) (map fromIntegral factoriales)
 
-- factoriales es la lista de los factoriales. Por ejemplo, 
--    take 7 factoriales  ==  [1,1,2,6,24,120,720]
factoriales :: [Integer]
factoriales = 1 : scanl1 (*) [1..]
Seguir leyendo →

Las sucesiones de Loomis

José A. Alonso---4 Comentarios

La sucesión de Loomis generada por un número entero positivo x es la sucesión cuyos términos se definen por

  • f(0) es x
  • f(n) es la suma de f(n-1) y el producto de los dígitos no nulos de f(n-1)

Los primeros términos de las primeras sucesiones de Loomis son

  • Generada por 1: 1, 2, 4, 8, 16, 22, 26, 38, 62, 74, 102, 104, 108, 116, 122, …
  • Generada por 2: 2, 4, 8, 16, 22, 26, 38, 62, 74, 102, 104, 108, 116, 122, 126, …
  • Generada por 3: 3, 6, 12, 14, 18, 26, 38, 62, 74, 102, 104, 108, 116, 122, 126, …
  • Generada por 4: 4, 8, 16, 22, 26, 38, 62, 74, 102, 104, 108, 116, 122, 126, 138, …
  • Generada por 5: 5, 10, 11, 12, 14, 18, 26, 38, 62, 74, 102, 104, 108, 116, 122, …

Se observa que a partir de un término todas coinciden con la generada por 1. Dicho término se llama el punto de convergencia. Por ejemplo,

  • la generada por 2 converge a 2
  • la generada por 3 converge a 26
  • la generada por 4 converge a 4
  • la generada por 5 converge a 26

Definir las siguientes funciones

   sucLoomis           :: Integer -> [Integer]
   convergencia        :: Integer -> Integer
   graficaConvergencia :: [Integer] -> IO ()

tales que

  • (sucLoomis x) es la sucesión de Loomis generada por x. Por ejemplo,
     λ> take 15 (sucLoomis 1)
     [1,2,4,8,16,22,26,38,62,74,102,104,108,116,122]
     λ> take 15 (sucLoomis 2)
     [2,4,8,16,22,26,38,62,74,102,104,108,116,122,126]
     λ> take 15 (sucLoomis 3)
     [3,6,12,14,18,26,38,62,74,102,104,108,116,122,126]
     λ> take 15 (sucLoomis 4)
     [4,8,16,22,26,38,62,74,102,104,108,116,122,126,138]
     λ> take 15 (sucLoomis 5)
     [5,10,11,12,14,18,26,38,62,74,102,104,108,116,122]
     λ> take 15 (sucLoomis 20)
     [20,22,26,38,62,74,102,104,108,116,122,126,138,162,174]
     λ> take 15 (sucLoomis 100)
     [100,101,102,104,108,116,122,126,138,162,174,202,206,218,234]
     λ> sucLoomis 1 !! (2*10^5)
     235180736652
  • (convergencia x) es el término de convergencia de la sucesioń de Loomis generada por x xon la geerada por 1. Por ejemplo, 
     convergencia  2      ==  2
     convergencia  3      ==  26
     convergencia  4      ==  4
     convergencia 17      ==  38
     convergencia 19      ==  102
     convergencia 43      ==  162
     convergencia 27      ==  202
     convergencia 58      ==  474
     convergencia 63      ==  150056
     convergencia 81      ==  150056
     convergencia 89      ==  150056
     convergencia (10^12) ==  1000101125092
  • (graficaConvergencia xs) dibuja la gráfica de los términos de convergencia de las sucesiones de Loomis generadas por los elementos de xs. Por ejemplo, (graficaConvergencia ([1..50]) dibuja
    Las_sucesiones_de_Loomis_1
    y graficaConvergencia ([1..148] \ [63,81,89,137]) dibuja
    Las_sucesiones_de_Loomis_2

Soluciones

import Data.List               ((\\))
import Data.Char               (digitToInt)
import Graphics.Gnuplot.Simple (plotList, Attribute (Key, Title, XRange, PNG))
 
-- 1ª definición de sucLoomis
-- ==========================
 
sucLoomis :: Integer -> [Integer]
sucLoomis x = map (loomis x) [0..]
 
loomis :: Integer -> Integer -> Integer
loomis x 0 = x
loomis x n = y + productoDigitosNoNulos y
  where y = loomis x (n-1)
 
productoDigitosNoNulos :: Integer -> Integer
productoDigitosNoNulos = product . digitosNoNulos
 
digitosNoNulos :: Integer -> [Integer]
digitosNoNulos x =
  [read [c] | c <- show x, c /= '0']
 
-- 2ª definición de sucLoomis
-- ==========================
 
sucLoomis2 :: Integer -> [Integer]
sucLoomis2 = iterate siguienteLoomis 
 
siguienteLoomis :: Integer -> Integer
siguienteLoomis y = y + productoDigitosNoNulos y
 
-- 3ª definición de sucLoomis
-- ==========================
 
sucLoomis3 :: Integer -> [Integer]
sucLoomis3 =
  iterate ((+) <*> product .
           map (toInteger . digitToInt) .
           filter (/= '0') . show)
 
-- Comparación de eficiencia
-- =========================
 
--    λ> sucLoomis 1 !! 30000
--    6571272766
--    (2.45 secs, 987,955,944 bytes)
--    λ> sucLoomis2 1 !! 30000
--    6571272766
--    (2.26 secs, 979,543,328 bytes)
--    λ> sucLoomis3 1 !! 30000
--    6571272766
--    (0.31 secs, 88,323,832 bytes)
 
-- 1ª definición de convergencia
-- =============================
 
convergencia1 :: Integer -> Integer
convergencia1 x =
  head (dropWhile noEnSucLoomisDe1 (sucLoomis x))
 
noEnSucLoomisDe1 :: Integer -> Bool
noEnSucLoomisDe1 x = not (pertenece x sucLoomisDe1)
 
sucLoomisDe1 :: [Integer]
sucLoomisDe1 = sucLoomis 1
 
pertenece :: Integer -> [Integer] -> Bool
pertenece x ys =
  x == head (dropWhile (<x) ys)
 
-- 2ª definición de convergencia
-- =============================
 
convergencia2 :: Integer -> Integer
convergencia2 = aux (sucLoomis3 1) . sucLoomis3
 where aux as@(x:xs) bs@(y:ys) | x == y    = x
                               | x < y     = aux xs bs
                               | otherwise = aux as ys
 
-- 3ª definición de convergencia
-- =============================
 
convergencia3 :: Integer -> Integer
convergencia3 = head . interseccion (sucLoomis3 1) . sucLoomis3
 
-- (interseccion xs ys) es la intersección entre las listas ordenadas xs
-- e ys. Por ejemplo,
--    λ> take 10 (interseccion (sucLoomis3 1) (sucLoomis3 2))
--    [2,4,8,16,22,26,38,62,74,102]
interseccion :: Ord a => [a] -> [a] -> [a]
interseccion = aux
  where aux as@(x:xs) bs@(y:ys) = case compare x y of
                                    LT ->     aux xs bs
                                    EQ -> x : aux xs ys
                                    GT ->     aux as ys
        aux _         _         = []                           
 
-- 4ª definición de convergencia
-- =============================
 
convergencia4 :: Integer -> Integer
convergencia4 x = perteneceA (sucLoomis3 x) 1
  where perteneceA (y:ys) n | y == c    = y
                            | otherwise = perteneceA ys c
          where c = head $ dropWhile (< y) $ sucLoomis3 n
 
-- Comparación de eficiencia
-- =========================
 
--    λ> convergencia1 (10^4)
--    150056
--    (2.94 secs, 1,260,809,808 bytes)
--    λ> convergencia2 (10^4)
--    150056
--    (0.03 secs, 700,240 bytes)
--    λ> convergencia3 (10^4)
--    150056
--    (0.03 secs, 1,165,496 bytes)
--    λ> convergencia4 (10^4)
--    150056
--    (0.02 secs, 1,119,648 bytes)
--    
--    λ> convergencia2 (10^12)
--    1000101125092
--    (1.81 secs, 714,901,080 bytes)
--    λ> convergencia3 (10^12)
--    1000101125092
--    (1.92 secs, 744,932,184 bytes)
--    λ> convergencia4 (10^12)
--    1000101125092
--    (1.82 secs, 941,053,328 bytes)
 
-- Definición de graficaConvergencia
-- ==================================
 
graficaConvergencia :: [Integer] -> IO ()
graficaConvergencia xs =
  plotList [ Key Nothing
           , Title "Convergencia de sucesiones de Loomis"
           , XRange (fromIntegral (minimum xs),fromIntegral (maximum xs))
           , PNG "Las_sucesiones_de_Loomis_2.png"
           ]
           [(x,convergencia2 x) | x <- xs]
Seguir leyendo →

Operaciones con polinomios como diccionarios

José A. Alonso---3 Comentarios

Los polinomios se pueden representar mediante diccionarios con los exponentes como claves y los coeficientes como valores.

El tipo de los polinomios con coeficientes de tipo a se define por

   type Polinomio a = M.Map Int a

Dos ejemplos de polinomios (que usaremos en los ejemplos) son

   3 + 7x - 5x^3
   4 + 5x^3 + x^5

se definen por

  ejPol1, ejPol2 :: Polinomio Int
  ejPol1 = M.fromList [(0,3),(1,7),(3,-5)]
  ejPol2 = M.fromList [(0,4),(3,5),(5,1)]

Definir las funciones

   sumaPol :: (Num a, Eq a) => Polinomio a -> Polinomio a -> Polinomio a
   multPorTerm :: Num a => (Int,a) -> Polinomio a -> Polinomio a
   multPol :: (Eq a, Num a) => Polinomio a -> Polinomio a -> Polinomio a

tales que

  • (sumaPol p q) es la suma de los polinomios p y q. Por ejemplo,
     λ> sumaPol ejPol1 ejPol2
     fromList [(0,7),(1,7),(5,1)]
     λ> sumaPol ejPol1 ejPol1
     fromList [(0,6),(1,14),(3,-10)]
  • (multPorTerm (n,a) p) es el producto del término ax^n por p. Por ejemplo, 
     λ> multPorTerm (2,3) (M.fromList [(0,4),(2,1)])
     fromList [(2,12),(4,3)]
  • (multPol p q) es el producto de los polinomios p y q. Por ejemplo, 
     λ> multPol ejPol1 ejPol2
     fromList [(0,12),(1,28),(3,-5),(4,35),(5,3),(6,-18),(8,-5)]
     λ> multPol ejPol1 ejPol1
     fromList [(0,9),(1,42),(2,49),(3,-30),(4,-70),(6,25)]
     λ> multPol ejPol2 ejPol2
     fromList [(0,16),(3,40),(5,8),(6,25),(8,10),(10,1)]

Soluciones

import qualified Data.Map as M
 
type Polinomio a = M.Map Int a 
 
ejPol1, ejPol2 :: Polinomio Int
ejPol1 = M.fromList [(0,3),(1,7),(3,-5)]
ejPol2 = M.fromList [(0,4),(3,5),(5,1)]
 
sumaPol :: (Num a, Eq a) => Polinomio a -> Polinomio a -> Polinomio a
sumaPol p q = 
    M.filter (/=0) (M.unionWith (+) p q)
 
multPorTerm :: Num a => (Int,a) -> Polinomio a -> Polinomio a
multPorTerm (n,a) p =
    M.map (*a) (M.mapKeys (+n) p)
 
multPol :: (Eq a, Num a) => Polinomio a -> Polinomio a -> Polinomio a
multPol p q
    | M.null p  = M.empty
    | otherwise = sumaPol (multPorTerm t q) (multPol r q)
    where (t,r) = M.deleteFindMin p
Seguir leyendo →

Cadenas de divisores

José A. Alonso---4 Comentarios

Una cadena de divisores de un número n es una lista donde cada elemento es un divisor de su siguiente elemento en la lista. Por ejemplo, las cadenas de divisores de 12 son [2,4,12], [2,6,12], [2,12], [3,6,12], [3,12], [4,12], [6,12] y [12].

Definir la función

   cadenasDivisores :: Int -> [[Int]]

tal que (cadenasDivisores n) es la lista de las cadenas de divisores de n. Por ejemplo,

   λ> cadenasDivisores 12
   [[2,4,12],[2,6,12],[2,12],[3,6,12],[3,12],[4,12],[6,12],[12]]
   λ> length (cadenaDivisores 48)
   48
   λ> length (cadenaDivisores 120)
   132

Soluciones

import Data.List (sort)
import Data.Numbers.Primes (isPrime)
 
-- 1ª definición
-- =============
 
cadenasDivisores :: Int -> [[Int]]
cadenasDivisores n = sort (extiendeLista [[n]])
    where extiendeLista []           = []
          extiendeLista ((1:xs):yss) = xs : extiendeLista yss
          extiendeLista ((x:xs):yss) =
              extiendeLista ([y:x:xs | y <- divisores x] ++ yss)
 
-- (divisores x) es la lista decreciente de los divisores de x distintos
-- de x. Por ejemplo,
--    divisores 12  ==  [6,4,3,2,1]
divisores :: Int -> [Int]
divisores x = 
    [y | y <- [a,a-1..1], x `mod` y == 0]
    where a = x `div` 2
 
-- 2ª definición
-- =============
 
cadenasDivisores2 :: Int -> [[Int]]
cadenasDivisores2 = sort . aux
    where aux 1 = [[]]
          aux n = [xs ++ [n] | xs <- concatMap aux (divisores n)]
 
-- 3ª definición
-- =============
 
cadenasDivisores3 :: Int -> [[Int]]
cadenasDivisores3 = sort . map reverse . aux
    where aux 1 = [[]]
          aux n = map (n:) (concatMap aux (divisores3 n))
 
-- (divisores3 x) es la lista creciente de los divisores de x distintos
-- de x. Por ejemplo,
--    divisores3 12  ==  [1,2,3,4,6]
divisores3 :: Int -> [Int]
divisores3 x = 
    [y | y <- [1..a], x `mod` y == 0]
    where a = x `div` 2
 
-- 1ª definición de nCadenasDivisores
-- ==================================
 
nCadenasDivisores1 :: Int -> Int
nCadenasDivisores1 = length . cadenasDivisores
 
-- 2ª definición de nCadenasDivisores
-- ==================================
 
nCadenasDivisores2 :: Int -> Int
nCadenasDivisores2 1 = 1
nCadenasDivisores2 n = 
    sum [nCadenasDivisores2 x | x <- divisores n]
Seguir leyendo →

Sucesión fractal

José A. Alonso---5 Comentarios

La sucesión fractal

   0, 0, 1, 0, 2, 1, 3, 0, 4, 2, 5, 1, 6, 3, 7, 0, 8, 4, 9, 2, 
   10, 5, 11, 1, 12, 6, 13, 3, 14, 7, 15, ...

está construida de la siguiente forma:

  • los términos pares forman la sucesión de los números naturales 
     0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, ...
  • los términos impares forman la misma sucesión original
     0, 0, 1, 0, 2, 1, 3, 0, 4, 2, 5, 1, 6, 3, 7, ...

Definir las funciones

   sucFractal     :: [Integer]
   sumaSucFractal :: Integer -> Integer

tales que

  • sucFractal es la lista de los términos de la sucesión fractal. Por ejemplo, 
     take 20 sucFractal   == [0,0,1,0,2,1,3,0,4,2,5,1,6,3,7,0,8,4,9,2]
     sucFractal !! 30     == 15
     sucFractal !! (10^7) == 5000000
  • (sumaSucFractal n) es la suma de los n primeros términos de la sucesión fractal. Por ejemplo,
     sumaSucFractal 10      == 13
     sumaSucFractal (10^5)  == 1666617368
     sumaSucFractal (10^10) == 16666666661668691669
     sumaSucFractal (10^15) == 166666666666666166673722792954
     sumaSucFractal (10^20) == 1666666666666666666616666684103392376198
     length (show (sumaSucFractal (10^15000))) == 30000
     sumaSucFractal (10^15000) `mod` (10^9)    == 455972157

Soluciones

 
-- 1ª definición de sucFractal
-- ===========================
 
sucFractal1 :: [Integer]
sucFractal1 = 
  map termino [0..]
 
-- (termino n) es el término n de la secuencia anterior. Por ejemplo,
--   termino 0            ==  0
--   termino 1            ==  0
--   map termino [0..10]  ==  [0,0,1,0,2,1,3,0,4,2,5]
termino :: Integer -> Integer
termino 0 = 0
termino n 
  | even n    = n `div` 2
  | otherwise = termino (n `div` 2)
 
-- 2ª definición de sucFractal
-- ===========================
 
sucFractal2 :: [Integer]
sucFractal2 =
  0 : 0 : mezcla [1..] (tail sucFractal2)
 
-- (mezcla xs ys) es la lista obtenida intercalando las listas infinitas
-- xs e ys. Por ejemplo,
--    take 10 (mezcla [0,2..] [0,-2..])  ==  [0,0,2,-2,4,-4,6,-6,8,-8]
mezcla :: [Integer] -> [Integer] -> [Integer]
mezcla (x:xs) (y:ys) =
  x : y : mezcla xs ys
 
-- Comparación de eficiencia de definiciones de sucFractal
-- =======================================================
 
--    λ> sum (take (10^6) sucFractal1)
--    166666169612
--    (5.56 secs, 842,863,264 bytes)
--    λ> sum (take (10^6) sucFractal2)
--    166666169612
--    (1.81 secs, 306,262,616 bytes)
 
-- En lo que sigue usaremos la 2ª definición
sucFractal :: [Integer]
sucFractal = sucFractal2
 
-- 1ª definición de sumaSucFractal
-- ===============================
 
sumaSucFractal1 :: Integer -> Integer
sumaSucFractal1 n =
  sum (map termino [0..n-1])
 
-- 2ª definición de sumaSucFractal
-- ===============================
 
sumaSucFractal2 :: Integer -> Integer
sumaSucFractal2 n =
  sum (take (fromIntegral n) sucFractal)
 
-- 3ª definición de sumaSucFractal
-- ===============================
 
sumaSucFractal3 :: Integer -> Integer
sumaSucFractal3 0 = 0
sumaSucFractal3 1 = 0
sumaSucFractal3 n
  | even n    = sumaN (n `div` 2) + sumaSucFractal3 (n `div` 2)
  | otherwise = sumaN ((n+1) `div` 2) + sumaSucFractal3 (n `div` 2)
  where sumaN n = (n*(n-1)) `div` 2
 
-- Comparación de eficiencia de definiciones de sumaSucFractal
-- ===========================================================
 
--    λ> sumaSucFractal1 (10^6)
--    166666169612
--    (5.25 secs, 810,622,504 bytes)
--    λ> sumaSucFractal2 (10^6)
--    166666169612
--    (1.72 secs, 286,444,048 bytes)
--    λ> sumaSucFractal3 (10^6)
--    166666169612
--    (0.01 secs, 0 bytes)
--    
--    λ> sumaSucFractal2 (10^7)
--    16666661685034
--    (17.49 secs, 3,021,580,920 bytes)
--    λ> sumaSucFractal3 (10^7)
--    16666661685034
--    (0.01 secs, 0 bytes)
Seguir leyendo →