Menu Close

Etiqueta: Polinomios

Valores de polinomios representados con vectores

Los polinomios se pueden representar mediante vectores usando la librería Data.Array. En primer lugar, se define el tipo de los polinomios (con coeficientes de tipo a) mediante

   type Polinomio a = Array Int a

Como ejemplos, definimos el polinomio

   ej_pol1 :: Array Int Int
   ej_pol1 = array (0,4) [(0,6),(1,2),(2,-5),(3,0),(4,7)]

que representa a 6 + 2x – 5x^2 + 7x^4 y el polinomio

   ej_pol2 :: Array Int Double
   ej_pol2 = array (0,4) [(0,6.5),(1,2),(2,-5.2),(3,0),(4,7)]

que representa a 6.5 + 2x – 5.2x^2 + 7x^4

Definir la función

   valor :: Num a => Polinomio a -> a -> a

tal que (valor p b) es el valor del polinomio p en el punto b. Por ejemplo,

   valor ej_pol1 0  ==  6
   valor ej_pol1 1  ==  10
   valor ej_pol1 2  ==  102
   valor ej_pol2 0  ==  6.5
   valor ej_pol2 1  ==  10.3
   valor ej_pol2 3  ==  532.7
   length (show (valor (listArray (0,5*10^5) (repeat 1)) 2)) == 150516

Soluciones

import Data.List (foldl')
import Data.Array (Array, (!), array, assocs, bounds, elems, listArray)
import Test.QuickCheck
 
type Polinomio a = Array Int a
 
ej_pol1 :: Array Int Int
ej_pol1 = array (0,4) [(1,2),(2,-5),(4,7),(0,6),(3,0)]
 
ej_pol2 :: Array Int Double
ej_pol2 = array (0,4) [(1,2),(2,-5.2),(4,7),(0,6.5),(3,0)]
 
-- 1ª solución
-- ===========
 
valor1 :: Num a => Polinomio a -> a -> a
valor1 p b = sum [(p!i)*b^i | i <- [0..n]]
  where (_,n) = bounds p
 
-- 2ª solución
-- ===========
 
valor2 :: Num a => Polinomio a -> a -> a
valor2 p b = sum [(p!i)*b^i | i <- [0..length p - 1]]
 
-- 3ª solución
-- ===========
 
valor3 :: Num a => Polinomio a -> a -> a
valor3 p b = sum [v*b^i | (i,v) <- assocs p]
 
-- 4ª solución
-- ===========
 
valor4 :: Num a => Polinomio a -> a -> a
valor4 = valorLista4 . elems
 
valorLista4 :: Num a => [a] -> a -> a
valorLista4 xs b =
  sum [(xs !! i) * b^i | i <- [0..length xs - 1]]
 
-- 5ª solución
-- ===========
 
valor5 :: Num a => Polinomio a -> a -> a
valor5 = valorLista5 . elems
 
valorLista5 :: Num a => [a] -> a -> a
valorLista5 []     _ = 0
valorLista5 (x:xs) b = x + b * valorLista5 xs b
 
-- 6ª solución
-- ===========
 
valor6 :: Num a => Polinomio a -> a -> a
valor6 = valorLista6 . elems
 
valorLista6 :: Num a => [a] -> a -> a
valorLista6 xs b = aux xs
  where aux []     = 0
        aux (y:ys) = y + b * aux ys
 
-- 7ª solución
-- ===========
 
valor7 :: Num a => Polinomio a -> a -> a
valor7 = valorLista7 . elems
 
valorLista7 :: Num a => [a] -> a -> a
valorLista7 xs b = foldr (\y r -> y + b * r) 0 xs
 
-- 8ª solución
-- ===========
 
valor8 :: Num a => Polinomio a -> a -> a
valor8 = valorLista8 . elems
 
valorLista8 :: Num a => [a] -> a -> a
valorLista8 xs b = aux 0 (reverse xs)
  where aux r []     = r
        aux r (y:ys) = aux (y + r * b) ys
 
-- 9ª solución
-- ===========
 
valor9 :: Num a => Polinomio a -> a -> a
valor9 = valorLista9 . elems
 
valorLista9 :: Num a => [a] -> a -> a
valorLista9 xs b = aux 0 (reverse xs)
  where aux = foldl (\ r y -> y + r * b)
 
-- 10ª solución
-- ============
 
valor10 :: Num a => Polinomio a -> a -> a
valor10 p b =
  foldl (\ r y -> y + r * b) 0 (reverse (elems p))
 
-- 11ª solución
-- ============
 
valor11 :: Num a => Polinomio a -> a -> a
valor11 p b =
  foldl' (\ r y -> y + r * b) 0 (reverse (elems p))
 
-- 12ª solución
-- ============
 
valor12 :: Num a => Polinomio a -> a -> a
valor12 p b =
  sum (zipWith (*) (elems p) (iterate (* b) 1))
 
-- 13ª solución
-- ============
 
valor13 :: Num a => Polinomio a -> a -> a
valor13 p b =
  foldl' (+) 0 (zipWith (*) (elems p) (iterate (* b) 1))
 
-- Equivalencia de las definiciones
-- ================================
 
-- La propiedad es
prop_valor :: [Integer] -> Integer -> Bool
prop_valor xs b =
  all (== valor1 p b)
      [f p b | f <- [valor2,
                     valor3,
                     valor4,
                     valor5,
                     valor6,
                     valor7,
                     valor8,
                     valor9,
                     valor10,
                     valor11,
                     valor12,
                     valor13]]
  where p = listArray (0, length xs - 1) xs
 
-- La comprobación es
--    λ> quickCheck prop_valor
--    +++ OK, passed 100 tests.
 
-- Comparación de eficiencia
-- =========================
 
-- La comparación es
--    λ> length (show (valor1 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (7.62 secs, 2,953,933,864 bytes)
--    λ> length (show (valor2 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (8.26 secs, 2,953,933,264 bytes)
--    λ> length (show (valor3 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (7.49 secs, 2,954,733,184 bytes)
--    λ> length (show (valor4 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (84.80 secs, 2,956,333,712 bytes)
--    λ> length (show (valor5 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.34 secs, 1,307,347,416 bytes)
--    λ> length (show (valor6 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.26 secs, 1,308,114,752 bytes)
--    λ> length (show (valor7 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.21 secs, 1,296,843,456 bytes)
--    λ> length (show (valor8 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.28 secs, 1,309,591,744 bytes)
--    λ> length (show (valor9 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.27 secs, 1,299,191,672 bytes)
--    λ> length (show (valor10 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (1.30 secs, 1,299,191,432 bytes)
--    λ> length (show (valor11 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (0.23 secs, 1,287,654,752 bytes)
--    λ> length (show (valor12 (listArray (0,10^5) (repeat 1)) 2))
--    30104
--    (0.75 secs, 1,309,506,968 bytes)

El código se encuentra en GitHub.

La elaboración de las soluciones se encuentran en el siguiente vídeo

Polinomios de Fibonacci

La sucesión de polinomios de Fibonacci se define por

   p(0) = 0
   p(1) = 1
   p(n) = x*p(n-1) + p(n-2)

Los primeros términos de la sucesión son

   p(2) = x
   p(3) = x^2 + 1
   p(4) = x^3 + 2*x
   p(5) = x^4 + 3*x^2 + 1

Definir la lista

   sucPolFib :: [Polinomio Integer]

tal que sus elementos son los polinomios de Fibonacci. Por ejemplo,

   λ> take 7 sucPolFib
   [0,1,1*x,x^2 + 1,x^3 + 2*x,x^4 + 3*x^2 + 1,x^5 + 4*x^3 + 3*x]
   λ> sum (map grado (take 3000 sucPolFib2))
   4495501

Comprobar con QuickCheck que el valor del n-ésimo término de sucPolFib para x=1 es el n-ésimo término de la sucesión de Fibonacci 0, 1, 1, 2, 3, 5, 8, …

Nota. Limitar la búsqueda a ejemplos pequeños usando

   quickCheckWith (stdArgs {maxSize=5}) prop_polFib

Soluciones

import Data.List (genericIndex)
import I1M.PolOperaciones
import Test.QuickCheck
 
-- 1ª solución
-- ===========
 
sucPolFib :: [Polinomio Integer]
sucPolFib = [polFibR n | n <- [0..]]
 
polFibR :: Integer -> Polinomio Integer
polFibR 0 = polCero
polFibR 1 = polUnidad
polFibR n = 
  sumaPol (multPol (consPol 1 1 polCero) (polFibR (n-1)))
          (polFibR (n-2))
 
-- 2ª definición (dinámica)
-- ========================
 
sucPolFib2 :: [Polinomio Integer]
sucPolFib2 = 
  polCero : polUnidad : zipWith f (tail sucPolFib2) sucPolFib2
  where f p = sumaPol (multPol (consPol 1 1 polCero) p)
 
-- La propiedad es
prop_polFib :: Integer -> Property
prop_polFib n = 
    n >= 0 ==> valor (polFib n) 1 == fib n
    where polFib n = sucPolFib2 `genericIndex` n
          fib n    = fibs `genericIndex` n
 
fibs :: [Integer]
fibs = 0 : 1 : zipWith (+) fibs (tail fibs)
 
-- La comprobación es
--    ghci> quickCheckWith (stdArgs {maxSize=5}) prop_polFib
--    +++ OK, passed 100 tests.

Polinomio digital

Definir la función

   polinomioDigital :: Int -> Polinomio Int

tal que (polinomioDigital n) es el polinomio cuyos coeficientes son los dígitos de n. Por ejemplo,

   λ> polinomioDigital 5703 
   5*x^3 + 7*x^2 + 3

Nota: Este ejercicio debe realizarse usando únicamente las funciones de la librería I1M.Pol que se encuentra aquí y se describe aquí.

Soluciones

import I1M.Pol
 
polinomioDigital :: Int -> Polinomio Int
polinomioDigital = creaPolDensa . digitos
 
-- (digitos n) es la lista de las dígitos de n. Por ejemplo,        
--    dígitos 142857  ==  [1,4,2,8,5,7]
digitos :: Int -> [Int]
digitos n = [read [x]| x <- show n]
 
-- (creaPolDensa xs) es el polinomio cuya representación densa es
-- xs. Por ejemplo, 
--    creaPolDensa [7,0,0,4,0,3]  ==  7*x^5 + 4*x^2 + 3
creaPolDensa :: (Num a, Eq a) => [a] -> Polinomio a
creaPolDensa []     = polCero
creaPolDensa (x:xs) = consPol (length xs) x (creaPolDensa xs)

Polinomios de Fibonacci

La sucesión de polinomios de Fibonacci se define por

   p(0) = 0
   p(1) = 1
   p(n) = x*p(n-1) + p(n-2)

Los primeros términos de la sucesión son

   p(2) = x
   p(3) = x^2 + 1
   p(4) = x^3 + 2*x
   p(5) = x^4 + 3*x^2 + 1

Definir la lista

   sucPolFib :: [Polinomio Integer]

tal que sus elementos son los polinomios de Fibonacci. Por ejemplo,

   λ> take 7 sucPolFib
   [0,1,1*x,x^2 + 1,x^3 + 2*x,x^4 + 3*x^2 + 1,x^5 + 4*x^3 + 3*x]
   λ> sum (map grado (take 3000 sucPolFib2))
   4495501

Comprobar con QuickCheck que el valor del n-ésimo término de sucPolFib para x=1 es el n-ésimo término de la sucesión de Fibonacci 0, 1, 1, 2, 3, 5, 8, …

Nota. Limitar la búsqueda a ejemplos pequeños usando

   quickCheckWith (stdArgs {maxSize=5}) prop_polFib

Soluciones

Polinomio cromático de un grafo

El polinomio cromático de un grafo calcula el número de maneras en las cuales puede ser coloreado el grafo usando un número de colores dado, de forma que dos vértices adyacentes no tengan el mismo color.

En el caso del grafo completo de n vértices, su polinomio cromático es

   P(n,x) = x(x-1)(x-2) ... (x-(n-1))

Por ejemplo,

   P(3,x) = x(x-1)(x-2)      = x^3 - 3*x^2 + 2*x
   P(4,x) = x(x-1)(x-2)(x-3) = x^4 - 6*x^3 + 11*x^2 - 6*x

Lo que significa que P(4)(x) es el número de formas de colorear el grafo completo de 4 vértices con x colores. Por tanto,

   P(4,2) =  0 (no se puede colorear con 2 colores)
   P(4,4) = 24 (hay 24 formas de colorearlo con 4 colores)

Definir la función

   polGC:: Int -> Polinomio Int

tal que (polGC n) es el polinomio cromático del grafo completo de n vértices. Por ejemplo,

   polGC 4  ==  x^4 + -6*x^3 + 11*x^2 + -6*x
   polGC 5  ==  x^5 + -10*x^4 + 35*x^3 + -50*x^2 + 24*x

Comprobar con QuickCheck que si el número de colores (x) coincide con el número de vértices del grafo (n), el número de maneras de colorear el grafo es n!.

Nota 1. Al hacer la comprobación limitar el tamaño de las pruebas como se indica a continuación

   ghci> quickCheckWith (stdArgs {maxSize=7}) prop_polGC
   +++ OK, passed 100 tests.

Nota 2: Este ejercicio debe realizarse usando únicamente las funciones de la librería de polinomios (I1M.PolOperaciones) que se describe aquí y se encuentra aquí.

Soluciones

import I1M.PolOperaciones
import Test.QuickCheck
 
-- 1ª solución
-- ===========
 
polGC :: Int -> Polinomio Int
polGC 0 = consPol 0 1 polCero
polGC n = polGC (n-1) `multPol` consPol 1 1 (consPol 0 (-n+1) polCero)
 
-- 2ª solución
-- ===========
 
polGC2 :: Int -> Polinomio Int
polGC2 n = multLista (map polMon [0..n-1])
 
-- (polMon n) es el monomio x-n. Por ejemplo,
--    polMon 3  ==  1*x + -3
polMon:: Int -> Polinomio Int
polMon n = consPol 1 1 (consPol 0 (-n) polCero)
 
-- (multLista ps) es el producto de la lista de polinomios ps.
multLista :: [Polinomio Int] -> Polinomio Int
multLista []     = polUnidad
multLista (p:ps) = multPol p (multLista ps)
 
-- 3ª solución (por plegado)
-- =========================
 
polGC3 :: Int -> Polinomio Int
polGC3 n = foldl multPol polUnidad
           [consPol 1 1 (consPol 0 (-i) polCero) | i <- [0..n-1]]
 
-- Comprobación
-- ============
 
-- La propiedad es
prop_polGC :: Int -> Property
prop_polGC n = 
    n > 0 ==> valor (polGC n) n == product [1..n]
 
-- La comprobación es
--    ghci> quickCheckWith (stdArgs {maxSize=7}) prop_polGC
--    +++ OK, passed 100 tests.
--    (0.04 secs, 7785800 bytes)

Polinomios pares

Un polinomio de coeficientes enteros se dirá par si todos sus coeficientes son números pares. Por ejemplo, el polinomio 2x³ – 4x² + 8 es par y el x² + 2x + 10 no lo es.

Definir el predicado

   parPol :: Integral a => Polinomio a -> Bool

tal que (parPol p) se verifica si p es un polinomio par. Por ejemplo,

   ghci> parPol (consPol 3 2 (consPol 2 (-4) (consPol 0 8 polCero)))
   True
   ghci> parPol (consPol 2 1 (consPol 1 2 (consPol 0 10 polCero)))
   False

Comprobar con QuickCheck que la suma de un polinomio con él mismo es un polinomio par.

Nota: Este ejercicio debe realizarse usando la librería I1M.Pol que se encuentra aquí y se describe aquí.

Soluciones

import I1M.PolOperaciones
import Test.QuickCheck
 
parPol :: Integral a => Polinomio a -> Bool
parPol p = esPolCero p || (even (coefLider p) && parPol (restoPol p))
 
-- La propiedad es
prop_parPol :: Integral a => Polinomio a -> Bool
prop_parPol p =
    parPol (sumaPol p p)
 
-- La comprobación es 
--    ghci> quickCheck prop_parPol
--    +++ OK, passed 100 tests.

Cociente entero de polinomios

El cociente entero de un polinomio P(x) por un monomio axⁿ es el polinomio que se obtiene a partir de los términos de P(x) con un grado mayor o igual que n, realizando la división entera entre sus coeficientes y el coeficiente del monomio divisor y restando el valor de n al de sus grados. Por ejemplo,

  • El cociente entero de 4x⁴ + 6x³ + 7x² + 5x + 2 por el monomio 3x² se obtiene a partir de los términos 4x⁴ + 6x³ + 7x² realizando la división entera entre sus coeficientes y el número 3 y restando 2 a sus grados. De esta forma se obtiene x² + 2x + 2
  • El cociente entero de 6x⁵ + 2x⁴ + 8x³ + 5x² + 8x + 4 por el monomio 4x³ se obtiene a partir de los términos 6x⁵ + 2x⁴ + 8x³ realizando la división entera entre sus coeficientes y el número 4 y restando 3 a sus grados. De esta forma se obtiene x² + 2

Definir la función

   cocienteEntero :: Polinomio Int -> Int -> Int -> Polinomio Int

tal que (cocienteEntero p a n) es el cociente entero del polinomio p por el monomio de grado n y coeficiente a. Por ejemplo,

   ghci> let listaApol xs = foldr (\(n,b) -> consPol n b) polCero xs
   ghci> cocienteEntero (listaApol [(4,4),(3,6),(2,7),(1,5),(0,2)]) 3 2
   x^2 + 2*x + 2
   ghci> cocienteEntero (listaApol [(5,6),(4,2),(3,8),(2,5),(1,8),(0,4)]) 4 3
   x^2 + 2

Nota: Este ejercicio debe realizarse usando únicamente las funciones de la librería I1M.Pol que se encuentra aquí y se describe aquí.

Soluciones

import I1M.Pol
 
cocienteEntero :: Polinomio Int -> Int -> Int -> Polinomio Int
cocienteEntero p a n
    | grado p < n = polCero
    | otherwise   = consPol (grado p - n) (coefLider p `div` a)
                            (cocienteEntero (restoPol p) a n)