Ir al contenido principal

Polinomios de Bell

B₀(x) = 1 (polinomio unidad)
Bₙ(x) = x·[Bₙ(x) + Bₙ'(x)]

Por ejemplo,

B₀(x) = 1                     = 1
B₁(x) = x·(1+0)               = x
B₂(x) = x·(x+1)               = x²+x
B₃(x) = x·(x²+x+2x+1)         = x³+3x²+x
B₄(x) = x·(x³+3x²+x+3x²+6x+1) = x⁴+6x³+7x²+x

Definir la función

polBell :: Integer -> Polinomio Integer

tal que (polBell n) es el polinomio de Bell de grado n. Por ejemplo,

polBell 4                    ==  x^4 + 6*x^3 + 7*x^2 + 1*x
coeficiente 2 (polBell 4)    ==  7
coeficiente 2 (polBell 30)   ==  536870911
coeficiente 1 (polBell 1000) == 1
length (show (coeficiente 9 (polBell 1000)))  ==  949

Notas: Se usa la librería I1M.PolOperaciones que se encuentra aquí y se describe aquí. Además, en el último ejemplo se usa la función coeficiente tal que (coeficiente k p) es el coeficiente del término de grado k en el polinomio p definida por

coeficiente :: Num a => Int -> Polinomio a -> a
coeficiente k p | k == n                 = coefLider p
                | k > grado (restoPol p) = 0
                | otherwise              = coeficiente k (restoPol p)
                where n = grado p

Soluciones

import Data.List (genericIndex)
import I1M.PolOperaciones

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

polBell1 :: Integer -> Polinomio Integer
polBell1 0 = polUnidad
polBell1 n = multPol (consPol 1 1 polCero) (sumaPol p (derivada p))
  where p = polBell1 (n-1)

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

polBell2 :: Integer -> Polinomio Integer
polBell2 n = sucPolinomiosBell `genericIndex` n

sucPolinomiosBell :: [Polinomio Integer]
sucPolinomiosBell = iterate f polUnidad
  where f p = multPol (consPol 1 1 polCero) (sumaPol p (derivada p))

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

-- La comparación es
--    λ> length (show (coeficiente 9 (polBell1 1000)))
--    949
--    (4.52 secs, 1002959992 bytes)
--    λ> length (show (coeficiente 9 (polBell2 1000)))
--    949
--    (4.01 secs, 1002475712 bytes)

-- (coeficiente k p) es el coeficiente del término de grado k en el
-- polinomio p.
coeficiente :: Num a => Int -> Polinomio a -> a
coeficiente k p | k == n                 = coefLider p
                | k > grado (restoPol p) = 0
                | otherwise              = coeficiente k (restoPol p)
                where n = grado p