I. TELECOMUNICACION. ALGEBRA LINEAL Sesi´on de laboratorio 1. Polinomios y algoritmo de Horner Esta primera sesi´on intenta familiarizar al estudiante con la representaci´ on y manejo de los polinomios en MATLAB.
1
Representaci´ on de polinomios
Los polinomios en MATLAB vienen representados por vectores, cuyas componentes guardan los coeficientes del polinomio en potencias de x y en orden decreciente. As´ı, la expresi´on >> p = [1
3
− 2];
>> tiene m´as de un sentido en MATLAB: puede representar el vector (fila) p = (1, 3, −2) (el punto y coma ; impide que la expresi´on del vector aparezca expl´ıcitamente en pantalla). pero tambi´en representa el polinomio p(x) = x2 + 3x − 2. El n´ umero de componentes de p determina el grado del polinomio, de forma que grado de p = n´ umero de componentes de p −1. Varios comandos de MATLAB van asociados a la manipulaci´on con polinomios. Los m´as elementales se refieren a las operaciones de producto y divisi´on, a su evaluaci´ on y al c´alculo aproximado de sus ra´ıces.
2
Operaciones de producto y divisi´ on de polinomios
El comando para realizar el producto de dos polinomios es conv. Por ejemplo, los vectores >> p = [1 3
− 2]; q = [1
− 1];
>> representan, respectivamente, los polinomios p(x) = x2 + 3x − 2, q(x) = x − 1. El vector >> r = conv(p, q) >> r = 1 2
−5
2
representa el producto r(x) = p(x)q(x) = x3 + 2x2 − 5x + 2. La instrucci´on >> help
conv
proporciona una breve explicaci´on del comando, sus variantes, as´ı como comandos relacionados. Uno de ellos es el que implementa la divisi´on, deconv. Muestra como resultado los polinomios
cociente y resto de la divisi´on de los polinomios de entrada. >> b = [1
3
− 2]; a = [1
− 1];
>> [q, r] = deconv(b, a) >> q = 1
4
r= 0
0 2
As´ı, la divisi´on de b(x) = x2 + 3x − 2 entre a(x) = x − 1 proporciona el cociente q(x) = x + 4 y el resto r(x) = 2, de modo que b(x) = a(x)q(x) + r(x).
3
Evaluaci´ on y determinaci´ on de ra´ıces
Algunas operaciones asociadas a polinomios no pueden realizarse exactamente y MATLAB implementa algoritmos num´ericos de aproximaci´ on. Tal es el caso de la determinaci´on de las ra´ıces. El comando roots proporciona una aproximaci´ on a los ceros de un polinomio. Por ejemplo, las instrucciones >> p = [1
3
− 2];
>> roots(p) ans = −3.5616 0.5616 muestran aproximaciones a las ra´ıces del polinomio p(x) = x2 + 3x − 2 que, por la f´ormula de resoluci´on de ecuaciones de segundo grado, son x± =
√ ´ 1³ −3 ± 17 . 2
Cuando el resultado de una instrucci´on ejecutada no lleva asociado expl´ıcitamente una variable (como en el uso de roots en el ejemplo anterior) MATLAB asigna una por defecto, llamada ans (de la palabra inglesa answer). Por otro lado, MATLAB siempre realiza las operaciones con la misma precisi´on, si bien la representaci´ on de los n´ umeros puede variar. Esto viene controlado por el comando format. As´ı, por ejemplo, la sucesi´on >> f ormat long >> p = [1
3
e
− 2]; roots(p)
ans = −3.561552812808830e + 00 0.5615528128088303e − 01
genera una representaci´on de las ra´ıces con m´as decimales y en notaci´on cient´ıfica. La instrucci´on help format se˜ nala las posibles variantes del comando. La evaluaci´on de polinomios viene implementada por el comando polyval. Las instrucciones >> p = [1
3
− 2]; c = 1;
>> pc = polyval(p, c) pc = 2 dan como resultado la evaluaci´on del polinomio p(x) = x2 + 3x − 2 en el punto c = 1, que se guarda en la variable pc = p(c) = 2. La aplicaci´on m´as elemental del comando polyval utiliza el algoritmo de Horner. Como ilustraci´on general de la construcci´on de un algoritmo y su programaci´on en MATLAB, vamos a implementar nuestro propio programa de evaluaci´ on de un polinomio. La forma m´as habitual de programar un algoritmo comienza elaborando el llamado pseudoc´odigo. Se trata de un conjunto de instrucciones que muestra claramente los pasos de que consta el proceso, de manera que su traslaci´on a cualquier lenguaje de programaci´on sea m´as o menos directa. En nuestro caso, y en su versi´ on m´as sencilla, para evaluar un polinomio P (x) = a0 + a1 x + · · · + aN xN , en un punto c, puede escribirse P (c) = a0 + c(a1 + c(a2 + · · · + c(aN −2 + c(aN −1 + caN )) · · ·)), en una serie de multiplicaciones encajadas. El algoritmo, desarrollado aparentemente de forma independiente por Ruffini (1804), Budan (1807 y 1813) y Horner (1819), se basa en la evaluaci´ on sucesiva de los par´entesis, a partir del m´as interno. Si ponemos bN = aN y calculamos bk = ak + cbk+1 ,
k = N − 1, . . . , 0,
entonces b0 = P (c). Un pseudoc´odigo para este algoritmo podr´ıa ser el siguiente: Pseudoc´ odigo del algoritmo de Horner Dados a0 , . . . , aN y c bN = aN Desde k = N − 1 hasta 0 bk = ak + cbk+1 P (c) = b0 La traslaci´on al lenguaje de MATLAB tendr´ıa la forma siguiente:
function pc=horner(a,c) % Algoritmo de Horner para evaluar un polinomio % de la forma % a(1) + a(2)x + a(3)x2 + ...a(m)xm−1 % El vector a guarda los coeficientes % c es el punto de evaluaci´ on m = length(a); b(m) = a(m); for k = m − 1 : −1 : 1 b(k) = a(j) + c ∗ b(k + 1); end pc = b(1); Es necesario realizar algunos comentarios: (i) Todo programa que represente un algoritmo en MATLAB debe encabezarse con una instrucci´on del tipo function [a1 , . . . , aN ] = nombre (b1 , . . . , bM ) donde a1 , . . . , aN son los datos de salida y b1 , . . . , bM son las variables de entrada (ii) Las l´ıneas precedidas de % son de comentario. MATLAB las ignora a efectos de computaci´on y son las que aparecen al teclear en la ventana de comandos >> help nombre Son muy u ´tiles para que el usuario reconozca, en pocas palabras, lo que realiza el programa. (iii) En el interior de un fichero de funci´on (el programa) se pueden declarar variables, realizar operaciones, crear condicionales, ciclos, llamar a funciones ya creadas, bien de MATLAB o del usuario, etc. As´ı, el comando length calcula la longitud de un vector. En nuestro caso, hay que tener en cuenta que el vector representa los coeficientes del polinomio y que MATLAB no admite vectores con componentes no positivas. De esta manera a = [a(1), . . . , a(m)] representa el polinomio a(1) + a(2)x + a(3)x2 + ...a(m)xm−1 y el pseudoc´odigo debe adaptarse a ello. (iv) El ciclo for k = m − 1 : −1 : 1 b(k) = a(j) + c ∗ b(k + 1); end implementa en MATLAB el cuerpo central del algoritmo. Recu´erdese lo mencionado en el apartado (iii), que influye en los valores l´ımite (en comp´araci´ on con el pseudoc´odigo) y obs´ervese c´omo se implementa un ciclo decreciente en MATLAB. As´ımismo, la u ´ltima instrucci´on del programa asigna el resultado del algoritmo a la variable de salida.
4
Algunos ejercicios
1. Realiza el producto y la divisi´on de los polinomios siguientes en MATLAB: (i) p(x) = x3 + x2 + x + 1, (ii) p(x) = x3 − 1,
q(x) = x2 + 2.
q(x) = x − 1.
(iii) p(x) = x4 + x2 + 2,
q(x) = x3 + x − 1.
2. Calcula una aproximaci´on a las ra´ıces de los polinomios siguientes (expres´andolas en dos formatos) y eval´ ualos en el punto indicado. (i) p(x) = x3 + x2 + x + 1, (ii) p(x) = x5 + x3 x, (iii) p(x) = x4 + x2 + 2,
c = −1.
c = i. c = −2.
3. Compara el psudoc´odigo del algoritmo de Horner con la siguiente variante Pseudoc´ odigo del algoritmo de Horner. Versi´ on 2 Dados a0 , . . . , aN y c, Desde k = N − 1 hasta 0 ak ← ak + cak+1 P (c) = a0 Elabora el correspondiente programa MATLAB y ejec´ utalo con p(x) = x3 + x2 + x + 1 y c = −1.