Polinomios de Chebyshev

 

http://www.sc.ehu.es/sbweb/fisica3/especial/chebyshev/chebyshev.html

Polinomios de Chebyshev

En esta página, estudiamos las sorprendentes propiedades de los polinomios de Chebyshev de primera especie T_n(x), los más importantes. Finalizando con la definición y algunas propiedades de los polinomios de Chebyshev de segunda especie U_n(x)

Definimos un conjunto de polinomios T_n(x)=cos(nθ) de grado n, donde cosθ=x en el intervalo -1≤x≤1.

Para determinar la forma de estos polinomios recordaremos la fórmula cos(a+b)=cos(a)cos(b)-sin(a)sin(b)

$\begin{array}{c}\{\begin{array}{c}\cos \left((n+1)\theta \right)=\cos \left(n\theta \right)\cos \theta -\sin \left(n\theta \right)\sin \theta \\ \cos \left((n-1)\theta \right)=\cos \left(n\theta \right)\cos \theta +\sin \left(n\theta \right)\sin \theta \end{array}\\ \cos \left((n+1)\theta \right)+\cos \left((n-1)\theta \right)=2\cos \left(n\theta \right)\cos \theta \\ T_{n+1}\left(x\right)+T_{n-1}\left(x\right)=2x{T}_n\left(x\right)\end{array}$

Los primeros polinomios son:

T₀(x)=cos(0)=1
T₁(x)=cosθ=x;

El resto de polinomios se obtiene empleando la relación de recurrencia:

T_n+1(x)=2xT_n(x)-T_n-1(x)

T₂(x)=cos(2θ)=2x²-1
T₃(x)=cos(3θ)=2x(2x²-1)-x=4x³-3x
T₄(x)= cos(4θ)=2x(4x³-3x)-( 2x²-1)=8x⁴-8x²+1
T₅(x)= cos(5θ)=2x(8x⁴-8x²+1)-(4x³-3x)=16x⁵-20x³+5x
T₆(x)= .....

Creamos una función recursiva que nos devuelve el valor del polinomio T_n(x) de grado n para valores dados de x.

function y=chebyshev(n,x) %n=1,2,3...
    if n==0
        y=ones(1,length(x));
    elseif n==1
        y=x;
    else
        y=2*x.*chebyshev(n-1,x)-chebyshev(n-2,x);
    end
end

Otra forma equivalente de crear una función que nos devuela el valor del polinomio T_n(x) de grado n para valores dados de x es la siguiente:

function T=chebyshev(n,x) %n=0,1,2,3,...
    m=length(x);
    T=zeros(1,m);
    T1=ones(1,m);
    T2=x;
    if n==1
        T=T2;
    else
        for j=2:n
            T=2*x.*T2-T1;
            T1=T2;
            T2=T;
        end
    end
end

Mediante el siguiente script, representamos gráficamente los polinomios T₁(x), T₂(x), T₃(x) y T₄(x)

hold on
for n=1:4
    fplot(@(x) chebyshev(n,x),[-1,1], 'displayName',num2str(n))
end
title('Polinomios de Chebyshev')
legend('-DynamicLegend','location','southeast')
xlabel('x')
ylabel('y')
grid on
hold off 

Se puede obtener los polinomios de Chebyshev T₀(x), T₁(x), T₂(x), T₃(x), T₄(x) para valores de x mediante el siguiente código.

n=5;    
x=-1:0.01:1; %vector de valores de x
m=length(x);
T=zeros(m,n);
T(:,1)=ones(1,m);  %T0(x)
T(:,2)=x;          %T1(x)
hold on
for j=3:n           %de T2(x) a T4(x)
    T(:,j)=2*x'.*T(:,j-1)-T(:,j-2);
end

Creamos una función que nos devuelve el vector de los coeficientes, [a₁,a₂,...a_n,a_n+1], (de mayor a menor grado) del polinomio T_n(x).

function p=chebypoly(n) %n es grado 0,1,2...
    p=zeros(n+1,1);
    p1=[1];
    p2=[1 0];
    if n==0
        p=p1;
     elseif n==1
       p=p2;
    else
       for i=2:n
            p=2*[p2 0]-[0 0 p1];
            p1=p2;
            p2=p;
        end   
    end       
end

La función poly2sym de Math Symbolic convierte los coeficientes del polinomio en la expresión algebraica a₁xⁿ+a₂x^n-1+...a_nx+a_n+1. La función ezplot representa gráficamente los polinomios

syms x;
hold on
for n=1:4
     Tk=poly2sym(chebypoly(n));
     ezplot(Tk,[-1,1]);
end
hold off
title('Polinomios de Chebyshev')
xlabel('x')
ylabel('y')
grid on

Esta misma reepresentación se puede obtener utilizando la función chebyshevT de MATLAB

hold on
fplot(chebyshevT(1:4,x),[-1,1]);
hold off
title('Polinomios de Chebyshev')
xlabel('x')
ylabel('y')
grid on

Propiedades de los polinomios de Chebyshev

Los polinomios de Chebyshev están definidos en el intervalo [-1,1].

Simetría

Los polinomios de índice par son funciones pares f(x)=f(-x) y los polinomios de índice impar son funciones impares f(x)=-f(-x)

Raíces

Obtenemos las raíces del polinomio T_n(x) de tres formas distintas:

• Mediante la definición del polinomio T_n(x)

Las raíces del polinomio T_n(x), se obtienen igualando cos(nθ)=0

$\begin{array}{c}n\theta =(2k-1)\frac{\pi }{2}k=\mathrm{1,2,3....}n\\ \theta =(2k-1)\frac{\pi }{2n}k=\mathrm{1,2,3....}n\\ x_k=\cos \left(\frac{(2k-1)\pi }{2n}\right)k=\mathrm{1,2....}n\end{array}$

Observamos que

$\begin{array}{c}{x}_1=\cos \left(\frac{\pi }{2n}\right)\\ x_n=\cos \left(\frac{\left(2n-1\right)\pi }{2n}\right)=\cos \left(\pi -\frac{\pi }{2n}\right)=-x_1\\ x_2=\cos \left(\frac{3\pi }{2n}\right)\\ x_{n-1}=\cos \left(\frac{\left(2(n-1)-1\right)\pi }{2n}\right)=\cos \left(\pi -\frac{3\pi }{2n}\right)=-x_2\\ x_k=-x_{n+1-k}\end{array}$

Calculamos las raíces del polinomio T₅(x)

>> k=1:5;
>> cos((2*k-1)*pi/(2*5))
ans =    0.9511    0.5878    0.0000   -0.5878   -0.9511

• Mediante la función roots de MATLAB, si disponemos del vector de los coeficientes del polinomio, que calculamos mediante la función chebypoly

>> roots(chebypoly(5))
ans =
         0
   -0.9511
   -0.5878
    0.9511
    0.5878

• Mediante la función solve de Math Symbolic que nos proporciona las raíces exactas del polinomio

>> T5=poly2sym(chebypoly(5))
T5 =16*x^5 - 20*x^3 + 5*x
>> r=solve(T5)
r =
                        0
  (5^(1/2)/8 + 5/8)^(1/2)
  (5/8 - 5^(1/2)/8)^(1/2)
 -(5^(1/2)/8 + 5/8)^(1/2)
 -(5/8 - 5^(1/2)/8)^(1/2)
>> double(r)
ans =
         0
    0.9511
    0.5878
   -0.9511
   -0.5878

Representamos las raíces del polinomio de grado 5, T₅(x), mediante puntos de color rojo a lo largo del eje X

n=5;
k=1:n;
ang=(2*k-1)*pi/(2*n);
xk=cos(ang);
yk=sin(ang);
hold on
plot(cosd(1:180),sind(1:180),'b')
for i=1:length(xk)
    line([0,xk(i)],[0,yk(i)])
    hg=line([xk(i),xk(i)],[0,yk(i)]);
    set(hg,'linestyle','--');
end
line([-1,1],[0,0])
plot(xk,0,'ro','markersize',4,'markerfacecolor','r')
hold off
axis equal
title('Raíces')

Máximos y mínimos

Los extremos (máximo o mínimo) del polinomio T_n(x)=cos(nθ) donde cosθ=x son 1 y -1, respectivamente.

$\frac{d{T}_n\left(x\right)}{dx}=-n\sin \left(n\theta \right)\frac{d\theta }{dx}=\frac{n\sin \left(n\theta \right)}{\sin \theta }$

se producen cuando sin(nθ)=0, es decir nθ=kπ, k=1,2…n-1. En las posiciones.

$x_k^{\prime }=\cos \left(\frac{k\pi }{n}\right)k=\mathrm{1,2...}n-1$

Se excluye k=0 y k=n por que

$\underset{\theta \to \mathrm{0,}\pi }{\lim }\frac{\sin \left(n\theta \right)}{\sin \theta }=\mathrm{n}$

>> k=1:4;
>> cos(k*pi/5)
ans = 0.8090    0.3090   -0.3090   -0.8090

>> T5=poly2sym(chebypoly(5)) 
T5 =16*x^5 - 20*x^3 + 5*x
>> dT5=diff(T5) %calcula la derivada
dT5 =80*x^4 - 60*x^2 + 5
>> r=solve(dT5)
r =
   5^(1/2)/4 + 1/4
   5^(1/2)/4 - 1/4
   1/4 - 5^(1/2)/4
 - 5^(1/2)/4 - 1/4
>> y=subs(T5,x,r); %valor del polinomio T5 para cada una de las raíces r.
>> simplify(y) 
ans =
 -1
  1
 -1
  1
>> double(r)
ans =
    0.8090
    0.3090
   -0.3090
   -0.8090

Relaciones entre dos polinomiosT_i(x) y T_j(x),

$\underset{-1}{\overset{1}{\int }}\frac{T_i\left(x\right)T_j\left(x\right)}{\sqrt{1-x^2}}dx=\{\begin{array}{c}0i\ne j\\ \pi i=j=0\\ \frac{\pi }{2}i=j\ne 0\end{array}$

>> syms x;
>> T5=poly2sym(chebypoly(5))
T5 =16*x^5 - 20*x^3 + 5*x
>> T4=poly2sym(chebypoly(4))
T4 =8*x^4 - 8*x^2 + 1
>> f=T4*T5/sqrt(1-x^2);
>> int(f,x,-1,1)
ans =0
 
>> f=T5*T5/sqrt(1-x^2);
>> int(f,x,-1,1)
ans =pi/2
 
>> f=1/sqrt(1-x^2);
>> int(f,x,-1,1)
ans =pi

Relación de ortogonalidad

Como hemos visto, las raíces del polinomio de grado n, T_n(x) son

$x_k=\cos \left(\frac{(2k-1)\pi }{2n}\right)k=\mathrm{1,2....}n$

Sean dos polinomios, T_i(x) y T_j(x), i,j<n vamos a comprobar que se cumple la siguiente relación que es importante para la aproximación de funciones f(x)

$\sum _{k=1}^n{T}_i\left(x_k\right)T_j\left(x_k\right)=\{\begin{array}{c}0i\ne j\\ ni=j=0\\ \frac{n}{2}i=j\ne 0\end{array}$

Esta es una relación análoga a la existente en las series de Fourier

n=7;
k=1:n;
xk=cos((2*k-1)*pi/(2*n));
y5=chebyshev(5,xk);
y4=chebyshev(4,xk);
suma=sum(y5.*y4)
suma=sum(y5.*y5)
suma=sum(y4.*y4)

En la ventana de comandos obtenemos el siguiente resultado

suma = -2.4980e-015
suma =    3.5000
suma =    3.5000

De otra forma, mediante la función chebypoly

>> n=7;
>> k=1:n;
>> xk=cos((2*k-1)*pi/(2*n));
>> T5=chebypoly(5);
>> T4=chebypoly(4);
>> suma=sum(polyval(T5,xk).*polyval(T4,xk))
suma = -1.2768e-015
>> suma=sum(polyval(T5,xk).*polyval(T5,xk))
suma =    3.5000
>> suma=sum(polyval(T4,xk).*polyval(T4,xk))
suma =    3.5000

Utilizando Math Symbolic obtenemos valores exactos

>> syms x;
>> T5=poly2sym(chebypoly(5));
>> T4=poly2sym(chebypoly(4));
>> n=sym('7');  %raices del polinomio T7
>> k=1:n;
>> xk=cos((2*k-1)*pi/(2*n));
>> suma=sum(subs(T5,x,xk).*subs(T4,x,xk))
suma =0 
>> suma=sum(subs(T5,x,xk).*subs(T5,x,xk));
>> simplify(suma)
ans =cos(pi/7) - cos((2*pi)/7) + cos((3*pi)/7) + 3
>> double(ans)
ans =    3.5000  %que es n/2
>> suma=sum(subs(T4,x,xk).*subs(T4,x,xk));
>> simplify(suma)
ans =cos((2*pi)/7) - cos(pi/7) - cos((3*pi)/7) + 4
>> double(ans)
ans =    3.5000 %que es n/2

Polinomios de Chebyshev de segunda especie

Definimos un conjunto de polinomios

$U_n\left(x\right)=\frac{\sin \left((n+1)\theta \right)}{\sin \theta }$

grado n, donde cosθ=x en el intervalo -1≤x≤1.

Para determinar la forma de estos polinomios recordaremos la fórmula sin(a+b)=sin(a)cos(b)+sin(b)cos(a)

$\begin{array}{c}\sin \left((n+2)\theta \right)=\sin \left((n+1)\theta \right)\cos \theta +\cos \left((n+1)\theta \right)\sin \theta =\\ \sin \left((n+1)\theta \right)\cos \theta +\left\{\cos \left(n\theta \right)\cos \theta -\sin \left(n\theta \right)\sin \theta \right\}\sin \theta =\\ \sin \left((n+1)\theta \right)\cos \theta +\cos \left(n\theta \right)\cos \theta \sin \theta -\sin \left(n\theta \right){\sin }^2\theta =\\ \sin \left((n+1)\theta \right)\cos \theta +\cos \left(n\theta \right)\cos \theta \sin \theta +\sin \left(n\theta \right){\cos }^2\theta -\sin \left(n\theta \right)\end{array}$

Agrupando términos

$\begin{array}{c}\sin \left((n+2)\theta \right)+\sin \left(n\theta \right)=\sin \left((n+1)\theta \right)\cos \theta +\left\{\sin \left(n\theta \right)\cos \theta +\cos \left(n\theta \right)\sin \theta \right\}\cos \theta =\\ =2\cos \theta \cdot \sin \left((n+1)\theta \right)\end{array}$

Dividiendo entre sinθ

$\begin{array}{c}\frac{\sin \left((n+2)\theta \right)}{\sin \theta }+\frac{\sin \left(n\theta \right)}{\sin \theta }=2\cos \theta \frac{\sin \left((n+1)\theta \right)}{\sin \theta }\\ U_{n+1}\left(x\right)+U_{n-1}\left(x\right)=2x{U}_n\left(x\right)\end{array}$

Esta es la relación de recurrencia, para utilizarla hemos de determinar lo dos primeros polinomios, U₀(x) y U₁(x)

Los primeros polinomios son:

$\begin{array}{c}{U}_0\left(x\right)=\frac{\sin \left(\theta \right)}{\sin \theta }=1\\ U_1\left(x\right)=\frac{\sin \left(2\theta \right)}{\sin \theta }=2\cos \theta =2x\end{array}$

U₂(x)=2x·U₁(x)-U₀(x)=4x²-1

Alternativamente, podríamos calcular U₂(x) a partir de su definición

$\begin{array}{c}{U}_2\left(x\right)=\frac{\sin \left(3\theta \right)}{\sin \theta }=\frac{\sin \left(2\theta \right)\cos \theta +\cos \left(2\theta \right)\sin \theta }{\sin \theta }=\frac{2\sin \theta {\cos }^2\theta +{\cos }^2\theta \sin \theta -{\sin }^3\theta }{\sin \theta }=\\ \frac{2\sin \theta {\cos }^2\theta +{\cos }^2\theta \sin \theta -(1-{\cos }^2\theta )\sin \theta }{\sin \theta }=4x^2-1\end{array}$

Creamos una función recursiva que nos devuelve el valor del polinomio U_n(x) de grado n para valores dados de x.

function y=chebyshev_2(n,x) %n=0,1,2,3...
    if n==0
        y=ones(1,length(x));
    elseif n==1
        y=2*x;
    else
        y=2*x.*chebyshev_2(n-1,x)-chebyshev_2(n-2,x);
    end
end

Mediante el siguiente script, representamos gráficamente los polinomios U₀(x), U₁(x), U₂(x), U₃(x) y U₄(x)

hold on
for n=0:4
    fplot(@(x) chebyshev_2(n,x),[-1,1], 'displayName',num2str(n))
end
title('Polinomios de Chebyshev')
legend('-DynamicLegend','location','southeast')
xlabel('x')
ylabel('y')
grid on
hold off 

Una representación gráfica similar se obtiene con la función chebyshevU de MATLAB

syms x;
hold on
fplot(chebyshevU(0:4,x),[-1,1]);
hold off
title('Polinomios de Chebyshev')
xlabel('x')
ylabel('y')
grid on

En la figura, vemos la simetría

$U_n(-x)={\left(-1\right)}^n{U}_n\left(x\right)$

Cuando n es un número par, 0,2,4... la función es simétrica respecto del eje Y. Cuando n es un número impar, 1,3,5...la función es antisimétrica

Cuando x=cosθ=1, θ=0, 2π,...,

$U_n\left(1\right)=\underset{\theta \to 0}{\lim }\frac{\sin \left((n+1)\theta \right)}{\sin \theta }=n+1$

>> syms n x;
>> limit(sin((n+1)*x)/sin(x),x,0)
ans =n + 1

Teniendo en cuenta la simetría de las funciones U_n(x)

$U_n(-1)={\left(-1\right)}^n\left(n+1\right)$

Creamos una función que nos devuelve el vector de los coeficientes, [a₁,a₂,...a_n,a_n+1], (de mayor a menor grado) del polinomio U_n(x).

function p=chebypoly_2(n) %n es grado 1,2...
    p=zeros(n,1);
    p0=[0];
    p1=[1];
    for i=2:n+1
        p=2*[p1,0]-[0,p0];
        p0=[0,p1];
        p1=p;
    end
end

La función poly2sym de Math Symbolic convierte los coeficientes del polinomio en la expresión algebraica a₁xⁿ+a₂x^n-1+...a_nx+a_n+1. La función ezplot representa gráficamente los polinomios

hold on
for n=1:4
     Un=poly2sym(chebypoly_2(n));
     ezplot(Un,[-1,1]);
end
hold off
ylim([-4,5])
title('Polinomios de Chebyshev')
xlabel('x')
ylabel('y')
grid on

Raíces

Obtenemos las raíces del polinomio U_n(x) de tres formas distintas:

• Mediante su definición

$\begin{array}{c}{U}_n\left(x\right)=\frac{\sin \left((n+1)\theta \right)}{\sin \theta }=0\\ {\theta }_m=\frac{m\pi }{n+1}m=\mathrm{1,2...}n\\ x_m=\cos {\theta }_m\end{array}$

n=4;
m=1:n;
th=m*pi/(n+1);
x=cos(th);
disp(x)

    0.8090    0.3090   -0.3090   -0.8090

• Mediante la función roots de MATLAB, si disponemos del vector de los coeficientes del polinomio, que calculamos mediante la función chebypoly_2

n=4;
p=chebypoly_2(n);
x=roots(p);
disp(x')

   -0.8090    0.8090   -0.3090    0.3090

• Mediante la función solve de Math Symbolic que nos proporciona las raíces exactas del polinomio

n=4;
p=poly2sym(chebypoly_2(n));
x=solve(p);
disp(x)

 - 5^(1/2)/4 - 1/4
   1/4 - 5^(1/2)/4
   5^(1/2)/4 - 1/4
   5^(1/2)/4 + 1/4

$-\frac{1}{4}\left(\sqrt{5}+1\right)-\frac{1}{4}\left(\sqrt{5}-1\right)\frac{1}{4}\left(\sqrt{5}-1\right)\frac{1}{4}\left(\sqrt{5}+1\right)$

Relaciones entre dos polinomios U_i(x) y U_j(x),

$\underset{-1}{\overset{1}{\int }}{U}_i\left(x\right)\cdot U_j\left(x\right)\sqrt{1-x^2}dx=\{\begin{array}{c}0i\ne j\\ \frac{\pi }{2}i=j\end{array}$

>> syms x;
>>  U5=poly2sym(chebypoly_2(5))
U5 =32*x^5 - 32*x^3 + 6*x
>>  U4=poly2sym(chebypoly_2(4))
U4 =16*x^4 - 12*x^2 + 1

>>  f=U4*U5*sqrt(1-x^2);
>>  int(f,x,-1,1)
ans =0
 
>>   f=U5*U5*sqrt(1-x^2);
>>  int(f,x,-1,1)
ans =pi/2 

Como hemos visto, las raíces del polinomio de grado n, U_n(x) son

$x_m=\cos \left(\frac{m\pi }{n+1}\right)m=\mathrm{1,2...}n$

Sean dos polinomios, U_i(x) y U_j(x), i,j<n vamos a comprobar que se cumple la siguiente relación

$\sum _{m=1}^n{U}_i\left(x_m\right)U_j\left(x_m\right)\left(1-x_m^2\right)=\{\begin{array}{c}0i\ne j\\ \frac{n+1}{2}i=j\end{array}$

n=7;
m=1:n;
xm=cos(m*pi/(n+1)); %raíces
y5=chebyshev_2(5,xm);
y4=chebyshev_2(4,xm);
suma=sum((y5.*y4).*(1-xm.^2))
suma=sum((y5.*y5).*(1-xm.^2))
suma=sum((y4.*y4).*(1-xm.^2))

suma =  -8.8818e-16
suma =    4.0000
suma =    4.0000