Series de Fourier (GRwM)

Trabajo realizado por estudiantes
Título	Series de Fourier. Grupo GRwM
Asignatura	EDP
Curso	2023-24
Autores	Guillermo Gómez Tejedor, Marina Jiménez Barrantes y Rocío Tajuelo Díaz
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 Introducción

En el trabajo que se presenta a continuación vamos a explicar la aproximación de funciones por series trigonométricas. Para ello, vamos a comenzar asentando la base teórica necesaria para entender el funcionamiento de estas aproximaciones. A continuación, vamos a proporcionar algunos ejemplos con el fin de ilustrar lo anterior. Los cálculos y visualizaciones han sido programados con MatLab.

2 Preliminares

En el siglo XIX, Jean-Baptiste Joseph Fourier introdujo una herramienta fundamental en el análisis matemático, las series de Fourier.

En el contexto que trabajaremos, las series de Fourier son series infinitas que convergen puntualmente a funciones periódicas y continuas. Permiten aproximar funciones periódicas mediante una combinación lineal de funciones sinusoidales (senos y cosenos). Serán la herramienta fundamental en nuestro trabajo.

Sin embargo, para entender su importancia y trasfondo, debemos introducir algunos conceptos más generales: los espacios de Hilbert y más tarde adaptarlos al espacio de funciones [math] L^2 [/math] .

Un espacio de Hilbert es un espacio vectorial separable, dotado de un producto escalar, que define una norma, y que es completo con respecto a esa norma. Trataremos espacios de Hilbert de dimensión infinita.

El espacio de funciones [math] L^2 [/math] sobre un un conjunto [math] \Omega [/math], se define como [math] L^2= \{ f: \Omega \longrightarrow \mathbb{R} [/math] tal que [math]\int_{\Omega}|f(x)|^2 \, dx \lt \infty \}[/math] . Además, se define el producto escalar de [math] f, g \in L^2(\Omega) [/math] como [math]\langle f,g\rangle_{L^2(\Omega)}= \int_{\Omega} f(x) \cdot g(x) \, dx [/math]. Gracias a este, podemos definir la norma en [math] L^2 [/math] como [math] \| \cdot \| _ {L^2 (\Omega)} = \sqrt{\langle\cdot , \cdot \rangle_{L^2(\Omega)}} [/math].

El espacio [math] L^2 (\Omega)[/math] con esta norma es un espacio de Hilbert y será sobre el cual trabajaremos a lo largo del documento.

NOTA: Es importante tener en cuenta que esto solo se cumple si consideramos que en [math]L^2[/math] dos funciones son iguales si y solo si son iguales en casi todo punto.

3 Bases de Hilbert

La importancia de los espacios de Hilbert, reside en que poseen propiedades análogas a los espacios vectoriales de dimensión finita, pero en dimensión infinita. Esto permite que estos espacios, que a priori se desvinculan de algunas de las propiedades esenciales de los espacios vectoriales de dimensión finita, sean manipulados con mayor facilidad.

La propiedad esencial de los espacios de Hilbert es que mantienen el concepto de base.

Una familia ortonormal y numerable de elementos de un espacio de Hilbert, cuyo espacio lineal generado es denso en dicho espacio, se llamará base Hilbertiana.

Estas bases son la infraestructura de estos espacios, ya que cada elemento podrá expresarse como una combinación lineal infinita y única de elementos de la base. Estas combinaciones lineales infinitas se denominan series de Fourier.

Y la gran ventaja reside en que el cálculo de los coeficientes de Fourier se obtiene mediante el producto escalar, al igual que en dimensión finita.

4 Base trigonométrica

Si consideramos el espacio [math]L^2([- \pi, \pi ])[/math] de funciones [math]2 \pi [/math] periódicas , una de estas bases es la base trigonométrica [math]\{ 1/2, \cos{(n x)},\sin{(n x)} \}_{n \in \mathbb{N}} ,x\in [-\pi ,\pi][/math].

De forma que, dada una función [math]f \in L^2([- \pi, \pi ])[/math] se tiene que [math] f \approx \frac{a_0}{2}+\sum^{\infty}_{n=1}(a_n\cos{(nx)}+b_n\sin{(nx)})[/math]. Como la base trigonométrica es ortogonal, podemos obtener los coeficientes de la serie de Fourier como

[math] a_0=\frac{1}{\pi} \int^{\pi}_{-\pi} f(x) dx [/math],

[math] a_n=\frac{1}{\pi} \int^{\pi}_{-\pi} f(x) \cos {(n x)}dx [/math],

[math] b_n=\frac{1}{\pi} \int^{\pi}_{-\pi} f(x) \sin {(n x)}dx [/math].

Si consideramos ahora el espacio de funciones [math]L^2([a,b])[/math], al variar el intervalo de [math] [- \pi, \pi ][/math] a [math] [a, b][/math], debemos hacer un cambio de variable de forma que las nuevas funciones de la base tengan periodo igual a la longitud del nuevo intervalo. De este modo, la base trigonométrica pasa a ser [math]\left \{ 1/2, \cos{\left(\frac{2n\pi}{b-a} x\right)},\sin{\left(\frac{2n\pi}{b-a} x \right)} \right\}_{n \in \mathbf{N}} ,x\in [a,b][/math]. Esta base no tiene por qué estar normalizada, sin embargo, que no lo esté solo altera cómo se obtienen los coeficientes de Fourier.

Al igual que antes, sea [math] f \in L^2([a,b])[/math], podemos aproximar la función mediante la serie de Fourier como

[math] f\approx \frac{a_0}{2}+\sum^{\infty}_{n=1}\left(a_n\cos{\left(\frac{2n\pi}{b-a} x\right)}+b_n\sin{\left(\frac{2n\pi}{b-a} x\right)}\right)[/math],

donde:

[math] a_0=\frac{1}{\|1/2\|_{L^2([a,b])}}\int^b_a \frac{f(x)}{2}dx[/math],

[math] a_n=\frac{1}{\left\|\cos{\left(\frac{2n\pi}{b-a} x\right)}\right\|_{L^2([a,b])}}\int^b_a f(x) \cos{\left(\frac{2n\pi}{b-a} x\right)}dx[/math],

[math] b_n=\frac{1}{\left\|\sin{\left(\frac{2n\pi}{b-a} x\right)}\right\|_{L^2([a,b])}}\int^b_a f(x) \sin{\left(\frac{2n\pi}{b-a} x\right)}dx [/math].

En la siguiente gráfica se representa esta base para los 10 primeros términos con [math][a,b]=[-1,1][/math], correspondiendo la recta negra a la función [math]1/2[/math] , y representados en línea continua y punteada [math]\cos{\left(n\pi x\right)}[/math] y [math]\sin{\left(n\pi x\right)}[/math], respectivamente.

f1=@(x,n) cos(n*pi.*x); f2=@(x,n) sin(n*pi.*x);  
xx=linspace(-1,1,1000);  
n = 10;  
colores=["#FF0000"; "#00FF00";"#0000FF"; "#F08080"; "#FFFF00" ; "#A901DB "; "#FF8000"; "#04B404"; "#FF0080"; "#848484"];   

hold on  
yline(1/2,"Color","k",LineWidth=1.3)  
lista=["1/2"];  

for i=1:n  
    plot(xx,f1(xx,i),"Color",colores(i),LineWidth=1.3)  
    plot(xx,f2(xx,i), ":","Color",colores(i),LineWidth=1.3)  
    lista(2*i)="cos("+num2str(i)+"*pi*x)";  
    lista(2*i+1)="sin("+num2str(i)+"*pi*x)";  
end  

lista(2)="cos(pi*x)"; lista(3)="sin(pi*x)"; legend(lista)  
title('Base trigonométrica'); xlabel('eje X'); ylabel('eje Y')

5 Aproximación de funciones impares

Siguiendo con nuestro objetivo de ilustrar la aproximación de funciones por series trigonométricas, vamos a considerar la función [math]f(x)=x(1-x)[/math] en el intervalo [math][0,1][/math].

En primer lugar obtenemos la extensión impar de [math]f[/math] a [math][-1,1][/math]:

[math] g(x)= \left \{ \begin{array}{ll} f(x), & \quad x\in (0,1], \\ 0, & \quad x=0, \\ -f(-x), & \quad x \in [-1,0). \end{array} \right. [/math]

Calculamos ahora la serie de Fourier para [math]g[/math] en [math][-1,1][/math]. Debemos tener en cuenta que la función [math]f[/math] es continua y cumple que el intervalo [math][0,1][/math] se puede dividir en un conjunto de subintervalos finitos en los cuales la función es monónotona, es decir, [math]f[/math] cumple la condición de Dirichlet. Gracias a esto, podemos asegurar que [math]\forall x_{0} \in [0,1][/math] la serie de Fourier va a converger puntualmente a [math]f(x_{0})[/math].

Observamos que

[math] a_0=\frac{1}{\|1/2\|_{L^2([a,b])}}\int^b_a \frac{f(x)}{2}\, dx=0,[/math]

[math] a_k=\frac{1}{\left\|\cos{\left(\frac{2n\pi}{b-a} x\right)}\right\|_{L^2([a,b])}} \int_{1}^{-1} \cos{(k\pi x)} f(x) \, dx [/math].

Esto se debe a que la integral de funciones impares sobre intervalos de la forma [math][-T,T][/math] es nula. Con ello como la multiplicación de una función par con una impar resulta en otra función impar, al calcular los coeficientes de [math]\cos{(k\pi x)}[/math] y [math]\frac{1}{2}[/math] obtenemos valores nulos.

Por tanto, se tiene que [math]g(x)= \sum_{k=1}^{\infty} b_k \sin{(k\pi x)}[/math].

Teniendo esto en cuenta, vamos a emplear las funciones impares de la base trigonométrica en [math][-1,1][/math], es decir [math]\{ \sin{(k\pi x)} \}_{k\in \mathbb{N}}[/math].

Definimos [math] f_n (x)[/math] como la suma de los primeros [math]n [/math] términos de la serie de Fourier, es decir, [math]f_n(x)= \sum_{k=1}^{n} b_k \sin{(k\pi x)}[/math], con [math]b_k = 2\int_{0}^{1} f(x)\sin{(k\pi x)} \,dx [/math].

A continuación, se representa gráficamente [math] f(x)[/math] y [math]f_n (x)[/math] para [math]n=1,5,10[/math] y el error obtenido en las normas [math] L^2[/math] y uniforme en función del número de términos de la serie. Los coeficientes de Fourier se obtienen numéricamente resolviendo las integrales usando la fórmula del trapecio con una división bastante fina [math](10^{-3})[/math] con la función [math]trapz[/math] que proporciona MatLab.

Representación gráfica de las sumas parciales de la serie de Fourier de f.

Errores de las sumas parciales de la serie de Fourier de f con las normas L^2 y del supremo en función de n.

f=@(x)(x.*(1-x));  %Definimos la función 
inicio=0; fin=1; %límites del intervalo 
g=@(x)(f(-x).*(sign(x)-1)/2+f(x).*(sign(x)+1)/2);  %Extensión impar de f  
inicio_impar=-fin; %límite a la izq 
div=10^-3; 

%Coeficientes de Fourier 
N2=100; N=[1;5;10;N2]; longitud = (fin-inicio_impar); 
a=zeros(max(N),1); b=zeros(max(N),1); 
X=inicio_impar:div:fin; 
G=g(X); 
a_0=trapz(X,G/2); %Coeficiente a_0 
F=a_0*ones(max(N),length(X)); 
NL2=zeros(N2,1); 
Nsup=zeros(N2,1); 

for i=1:max(N) 
    cc=cos(i*2*pi*X/longitud); 
    ss=sin(i*2*pi*X/longitud); 
    %Coeficientes a_k y b_k 
    a(i)=trapz(X,G.*cc); 
    b(i)=trapz(X,G.*ss);   
    F(i:max(N),:)=F(i:max(N),:)+a(i)*ones(max(N)-i+1,1)*cc + b(i)*ones(max(N)-i+1,1)*ss; %valores de f_n 
    NL2(i) = sqrt(trapz(X,abs(G-F(i,:)).^2)/2); %Error norma L^2 
    Nsup(i) = max(abs(G-F(i,:))); %Error norma supremo 
end 

figure(1) 
name="Extensión impar de f";
for j=1:length(N)-1
    subplot(1,length(N)-1,j)
    plot(X,G,"k","LineWidth",2)
    subtitle("f_{"+num2str(N(j))+"}")
    hold on
    plot(X,real(F(N(j),:)),"Color","r","LineWidth",1)
    legend(name,"f_{"+num2str(N(j))+"}");
    xlabel("eje x"); ylabel("eje y");
end

%Función error respecto de n 
figure(2) 
semilogy(1:N2,NL2,"Color","g","LineWidth",1.5) 
hold on  
semilogy(1:N2,Nsup,"Color","m","LineWidth",1.5) 
title("Errores de la base trigonométrica escala logarítmica");  
legend("Norma de L^2","Norma del supremo"); xlabel("N"); ylabel("Error")

6 Aproximación de funciones discontinuas y pares

A continuación vamos a estudiar la aproximación de una función discontinua. En este caso, vamos a trabajar con una función con una discontinuidad de salto finito.

En los puntos de continuidad de la función tendremos que la serie de Fourier va a converger al punto [math]f(x_{0}) [/math]. Sin embargo, si [math]x_{0} [/math] es un punto de discontinuidad, la serie converge en [math]x_{0} [/math] a [math]\frac{f(x_{0}^-) + f(x_{0}^+)}{2 }[/math].

Al intentar aproximar de forma continua una discontinuidad de salto, se obtienen pendientes extremas, que resultan en grandes oscilaciones en torno a esos puntos de discontinuidad. A este comportamiento oscilatorio de la serie de Fourier se le conoce como fenómeno de Gibbs.

Para observar este fenómeno, vamos a considerar la función discontinua [math]f(x) =1_{x \leq 1/2}(x)[/math]. Lo primero que observamos es que esta función se puede extender a [math][-1,1][/math] de forma par.

Por un razonamiento similar al de la extensión impar, la serie de Fourier en este caso tiene únicamente los términos de la base trigonométrica que son funciones pares, es decir, [math]1/2[/math] y [math]\cos{(n \pi x)}[/math].

Al aproximar esta función obtenemos los siguientes resultados:

f=@(x)(x<=1/2);  %Definimos la función escalón en el 1/2
inicio=0; fin=1; %límites del intervalo

%Extensión par
g=@(x)(f(-x).*(x<0) + f(x).*(x>=0));  %Extensión par de f 
inicio_par=-fin; %límite a la izquierda
div=10^-3;

%Coeficientes de Fourier
N2=100; N=[1;5;10;N2]; longitud = (fin-inicio_par);
a=zeros(max(N),1); b=zeros(max(N),1);
X=inicio_par:div:fin;
G=g(X);
a_0=trapz(X,G/2);
F=a_0*ones(max(N),length(X));
S=a_0*transpose(2:max(N)+1)*ones(1,length(X));
NL2=zeros(N2,1); Nsup=zeros(N2,1); NL2S=zeros(N2,1); NsupS=zeros(N2,1);

for i=1:max(N)
    cc=cos(i*2*pi*X/longitud);
    ss=sin(i*2*pi*X/longitud);
    %Coefs
    a(i)=trapz(X,G.*cc);
    b(i)=trapz(X,G.*ss);  
    %f_n
    F(i:max(N),:)=F(i:max(N),:)+a(i)*ones(max(N)-i+1,1)*cc + b(i)*ones(max(N)-i+1,1)*ss;
    %Errores de f_n
    NL2(i) = sqrt( trapz(X, abs( G-F(i,:) ).^2)/2 );
    Nsup(i) = max( abs( G-F(i,:) ) );
    %S_n
    sat=transpose(1:max(N)-i+1);
    S(i:max(N),:)=S(i:max(N),:)+a(i)*sat*cc + b(i)*sat*ss;
    S(i,:)=S(i,:)/(i+1);
    %Errores de S_n
    NL2S(i) = sqrt( trapz(X, abs( G-S(i,:) ).^2)/2 );
    NsupS(i) = max( abs( G-S(i,:) ) );
end

figure(1)
colors=["g","b","m"];
plot(X,G,"k","LineWidth",1.5)
name=["Extensión par de f"];
hold on
for j=1:length(N)-1
    plot(X,real(F(N(j),:)),"Color",colors(j),"LineWidth",1.5)
    name(j+1)="f_{"+num2str(N(j))+"}";
end
legend(name); xlabel("eje x"); ylabel("eje y")
title("Aproximaciones con la Base Trigonométrica")

%Representamos
figure(2)
plot(1:N2,NL2,"Color","g","LineWidth",1.5)
hold on 
plot(1:N2,Nsup,"Color","m","LineWidth",1.5)
title("Errores de la base trigonométrica")
legend("Norma de L^2","Norma del supremo"); xlabel("N"); ylabel("Error")

%%Sumas de Cesàro
figure(3)
title("Sumas de Cesàro")
colors=["g","b","m"];
plot(X,G,"k","LineWidth",1.5)
name=["Extensión par de f"];
hold on
for j=1:length(N)-1
    plot(X,real(S(N(j),:)),"Color",colors(j),"LineWidth",1.5)
    name(j+1)="S_{"+num2str(N(j))+"}";
end
legend(name); xlabel("eje x"); ylabel("eje y")

figure(4)
plot(1:N2,NL2S,"Color","r","LineWidth",1.5)
hold on 
plot(1:N2,NsupS,"Color","b","LineWidth",1.5)
title("Errores de las sumas de Cesàro")
legend("Norma de L^2","Norma del supremo"); xlabel("N"); ylabel("Error")

7 Cambio de intervalo

Otro caso para analizar es realizar un cambio de intervalo. Para ello vamos a considerar el intervalo [math][1,3][/math]. Dado que la longitud del intervalo es 2, la base trigonométrica asociada es [math]\left\{ 1/2, \cos{(n\pi x)},\sin{(n\pi x)}\right\}_{n \in \mathbb{N}} [/math] . Sea [math]f(x)=xe^{-x}[/math]. Aproximamos la función de la misma forma que en los casos anteriores.

Gráficamente se puede observar la representación de [math]f(x)[/math] y [math]f_n(x)[/math] con [math]n=5,10,20[/math].

f=@(x)(x.*(exp(-x)));  %Definimos la función
inicio=1; fin=3; %límites del intervalo
div=10^-3;

%Coeficientes de Fourier
N2=0; N=[5;10;20;N2]; longitud = (fin-inicio);
a=zeros(max(N),1); b=zeros(max(N),1);
X=inicio:div:fin;
G=f(X);
a_0=trapz(X,G/2);
F=a_0*ones(max(N),length(X));

for i=1:max(N)
    cc=cos(i*2*pi*X/longitud);
    ss=sin(i*2*pi*X/longitud);
    %Coefs
    a(i)=trapz(X,G.*cc);
    b(i)=trapz(X,G.*ss);  
    %f_n
    F(i:max(N),:)=F(i:max(N),:)+a(i)*ones(max(N)-i+1,1)*cc + b(i)*ones(max(N)-i+1,1)*ss;
end

figure(1)
colors=["g","b","m"];
plot(X,G,"k","LineWidth",1.5)
name=["Extensión impar de f"];
hold on
for j=1:length(N)-1
    plot(X,real(F(N(j),:)),"Color",colors(j),"LineWidth",1.5)
    name(j+1)="f_{"+num2str(N(j))+"}";
end
legend(name); xlabel("eje x"); ylabel("eje y")
title("Aproximaciones con la Base Trigonométrica")

8 Base trigonométrica compleja

Por último, vamos a considerar la base trigonométrica compleja, que viene dada por [math] \left\{ e^{\frac{2\pi inx}{b-a}}\right\} _{n \in \mathbb{Z} }[/math] en [math][a, b] [/math].

Cabe destacar que, al igual que ocurría en el caso real, es una base ortonormal si el intervalo sobre el que está definido tiene longitud 1.

Nótese que esta definición de base extiende la dada en los reales, pues podemos escribir todos los elementos de la base real como combinación de los elementos de la compleja:

[math] \sin{(n \pi x)} = \frac{e^{i n x}- e^{-i n x}}{2i} [/math] ,

[math] cos{(n \pi x)} = \frac{e^{i n x}+ e^{-i n x}}{2} [/math],

[math] \frac{1}{2}=\frac{1}{2}e^{0i} [/math].

En el espacio de funciones [math] L^2 ([a,b]) [/math] que toman valores complejos, el producto escalar se define como [math] \langle f,g \rangle_{L^2}= \int_{a} ^{b} f(x) \overline{g(x)} \, dx \in \mathbb{C} [/math] . En este caso los coeficientes de Fourier son números complejos, incluso aunque la función tome valores reales.

Veamos ahora un ejemplo al emplear esta base. Para ello, vamos a considerar la función [math] f(x) = 4x(1/2-x)^2[/math] y aproximarla. Comenzamos calculando la base trigonométrica adaptada al intervalo [math] [0,1] [/math] que será [math] \{ e^{i 2\pi nx}\} _{n \in \mathbb{Z} }[/math].

En la siguiente gráfica se representa la aproximación considerando los primeros 5, 10 y 20 términos de la serie.

f=@(x)4*x.*(1/2-x).^2;  %Definimos la función
inicio=0; fin=1; %límites del intervalo
div=10^-3;

%Coeficientes de Fourier
N2=100; N=[5;10;20;N2]; longitud = (fin-inicio);
a=zeros(max(N),1); b=zeros(max(N),1);
X=inicio:div:fin;
G=f(X);
a_0=trapz(X,G.*exp(2i*0*X/longitud));
F=a_0*ones(max(N),length(X));

for j=1:max(N)
    eplus=exp(j*2i*pi*X/longitud);
    emin=exp(-j*2i*pi*X/longitud);
    %Coefs
    a(j)=trapz(X,G.*eplus);
    b(j)=trapz(X,G.*emin);  
    %f_n
    F(j:max(N),:)=F(j:max(N),:)+a(j)*ones(max(N)-j+1,1)*eplus + b(j)*ones(max(N)-j+1,1)*emin;
end

figure(1)
colors=["g","b","m"];
plot(X,G,"k","LineWidth",1.5)
name=["f"];
hold on
for j=1:length(N)-1
    plot(X,flip(F(N(j),:)),"Color",colors(j),"LineWidth",1.5)
    name(j+1)="f_{"+num2str(N(j))+"}";
end
legend(name); xlabel("eje x"); ylabel("eje y")
title("Aproximaciones con la Base Trigonométrica Compleja")