Series de Fourier. Yan, Otelo, Miguel

Trabajo realizado por estudiantes
Título	Series de Fourier. Grupo 6-A
Asignatura	EDP
Curso	2023-24
Autores	Yan Wang Otelo Gallego Ayala Miguel Cazorla Pedraza
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

A lo largo de este artículo trataremos de entender la aproximación de funciones por series trigonométricas, más concretamente, por Series de Fourier. Para ello aproximaremos varias funciones con diferentes propiedades, para así llegar a comprender mejor como funcionan estas aproximaciones.

1 Series de Fourier

No tendría sentido intentar entender la aproximación de funciones por series trigonométricas, sin si quiera saber que son estás series. Así, comenzamos con una breve explicación teórica de las Series de Fourier.

Tomamos el espacio de Hilbert [math] L^2([-\pi,\pi]) [/math], con producto escalar asociado [math] ( f,g ) _{L^2} = \int_{-\pi}^{\pi} f \cdot g \hspace{5px} dx [/math], norma asociada a dicho producto escalar [math] || f ||_{L^2} = \sqrt{\langle f,g \rangle _{L^2}} [/math] y la base ortogonal [math] \{\frac{1}{2},cos(nx),sen(nx)\}_{n\in \mathbb{N}} [/math] compuesta por funciones [math] 2\pi [/math] periódicas (omitiremos la comprobación de que la base es ortogonal y de que el espacio es de Hilbert). Por ser un sistema ortogonal, podemos aproximar cualquier función [math] f(x) \in L^2([-\pi,\pi]) [/math] mediante la serie

[math] f(x)∼ \frac{a_0}{2} + \sum_{n=1}^{N} [a_nsen(nx) + b_nsen(nx)] [/math]

donde los coeficientes [math] a_0, a_n [/math] y [math] b_n [/math] reciben el nombre de Coeficientes de Fourier y vienen dados (por ser la base antes mencionada una base ortognal) por las expresiones

[math] a_0=\frac{\langle f , 1/2 \rangle _{L^2}}{|| \frac{1}{2}||^{2}_{L^2}}, \hspace{30px} a_n=\frac{\langle f , sen(nx) \rangle _{L^2}}{||sen(nx)||^{2}_{L^2}} \hspace{15px} \text{y} \hspace{15px} b_n=\frac{\langle f , cos(nx) \rangle _{L^2}}{||cos(nx)||^{2}_{L^2}}. [/math]

Para finalizar con esta introducción teórica, veamos rápidamente como cambian los elementos de nuestra base trigonométrica dependiendo del intervalo de definición de las funciones de [math] L^2 [/math]. Si el intervalo es la forma [math] [-T,T] [/math], nótese que mediante el cambio de variable [math] y=\frac{xT}{\pi} [/math] en la base anterior obtenemos [math] \{\frac{1}{2},cos(\frac{ny\pi}{T}),sen(\frac{ny\pi}{T})\}_{n\in \mathbb{N}} [/math], siendo esta una base ortogonal de [math] L^2([-T,T]) [/math]. Esta idea se puede generalizar a intervalos de la forma [math] [a,b] [/math] para cualesquiera [math] a [/math] y [math] b [/math] reales. Bastaría tomar la base trigonométrica del espacio [math] L^2([-\frac{b-a}{2},\frac{b-a}{2}]) [/math] , puesto que este intervalo tiene la misma longitud que [math] [a,b] [/math], y los elementos de dicha base son funciones periódicas cuyo periodo es precisamente la longitud [math] b-a [/math].

2 La base trigonométrica

Diez primeros términos de la base trigonométrica

Antes de comenzar con la aproximación de funciones, lo primero que haremos será dibujar los diez primeros términos de nuestra base trigonométrica, [math] \{\frac{1}{2},cos(nx\pi),sen(nx\pi)\}_{n\in \mathbb{N}} [/math] , con el fin de hacernos una idea de la forma geométrica de los elementos de nuestra base. Para dibujarlos utilizamos el siguiente código de matlab:

x = linspace(-1, 1, 1000); % Generar valores de x en el intervalo [-1, 1]
n = 1:10; % Números de términos

% Graficar los términos en una sola imagen con subgráficas
figure;

for i = 1:length(n)
    subplot(5, 2, i); % Crear una subgráfica
    plot(x, 1/2*ones(size(x)), 'r--', x, cos(n(i)*pi*x), 'g', x, sin(n(i)*pi*x), 'b'); % Graficar los términos
    title(['Términos para n = ', num2str(n(i))]); % Título de la subgráfica
    xlabel('x'); % Etiqueta del eje x
    ylabel(['y_', num2str(n(i))]); % Etiqueta del eje y
    legend({'1/2', ['cos(', num2str(n(i)), '\pi x)'], ['sin(', num2str(n(i)), '\pi x)']}); % Leyenda
end

Estas representaciones muestran algo que ya podíamos intuir por la forma de los elementos de nuestra base; a medida que el valor de [math] n [/math] aumenta, el periodo de las funciones [math] cos(nx\pi) [/math] y [math] sen(nx\pi) [/math] disminuye.

3 Aproximación de función continua

Extensión impar de la función f(x)=x(1-x)

Una vez representados algunos de los términos de nuestra base, usaremos estos para aproximar mediante una serie la función [math] f(x)=x(1-x) [/math] en el intervalo [math] [0,1] [/math]. Para facilitarnos la tarea, primero extenderemos de forma impar la función [math] f(x) [/math] al intervalo [math] [-1,1] [/math]. Nótese que los elementos de la forma [math]\frac{1}{2} [/math] y [math] cos(nx\pi) [/math] son impares, y por la definición antes explicada de los coeficientes de Fourier se tiene que en el dominio [math] [-1,1] [/math] dichos coeficientes serán nulos, por lo que la serie estará formada únicamente por términos de la forma [math] b_{n}sen(nx\pi) [/math].

% Definir la función
f = @(x) x.*(1-x);

% Definir el intervalo
x = linspace(0, 1, 100); % Generar 100 puntos entre 0 y 1
x2 = linspace(-1, 0, 100);

% Evaluar la función en los puntos del intervalo
y = f(x);
y1=-f(x);
% Graficar la función
hold on
plot(x2, y1, 'b', 'LineWidth', 2);
plot(x,y,'b','LineWidth', 2)
xlabel('x');
ylabel('Extensión impar def(x) ');
title('Gráfico de f(x) = x(1-x) extendida');
grid on;

Representación gráfica de [math] f(x) [/math] y de las aproximaciones [math] f_{1}(x), f_{5}(x) [/math] y [math] f_{10}(x) [/math]

Errores en las normas [math] L^2 [/math] y uniforme en función de [math] n [/math]

Una vez aclarado el por qué de la extensión impar, utilizamos el método del trapecio para calcular numéricamente los coeficientes de Fourier y representamos gráficamente tanto la función [math] f(x) [/math] como las aproximaciones [math] f_{1}(x), f_{5}(x) [/math] y [math] f_{10}(x) [/math]. Además, calculamos el error en las normas [math] L^2 [/math] y uniforme en función del número [math] n [/math] de términos de la serie, y graficamos también los errores obtenidos.

Todo esto queda escrito en el siguiente código de Matlab:

f1 = @(x)(1-x).*x;

coeficientes = zeros(1, 10);
N=1000;
a=0; b=1;
h=(b-a)/(N-1);
u=a:h:b; 
 w=ones(N,1);
w(1)=1/2; w(N)=1/2;

for k = 1:10 
    f=(f1(u).*sin(k.*pi.*u))';
    coeficientes(1,k) = h*w'*f;
end

coeficientes = 2*coeficientes;

y = zeros(10,1000);
y(1,:) = coeficientes(1,1).*sin(pi*u);
for k = 2:10
    y(k,:) = y(k-1,:) + coeficientes(1,k).*sin(k*pi*u);
end

figure(1)
hold on
plot(u,f1(u), "black")
plot(u,y(1,:), "g")
plot(u,y(5,:), "blue")
plot(u,y(10,:), 'red')
hold off
title("Aproximación de f por series de Fourier")
legend("f(x) =(1-x)x", "n =1", "n=5", "n=10")

m = 10;
errores_L2 = zeros(1,m); errores_uniforme = zeros(1,m);
for n=1:m
    integrando = (abs(f1(u)-y(n,:))).^2;
    g = integrando';
    int = h*w'*g;
    errores_L2(1,n) = sqrt(int);
    errores_uniforme(1,n) = max(abs(f1(u)-y(n,:)));
end
figure(2)
hold on
v = 1:1:m;
plot(v,errores_L2)
plot(v,errores_uniforme)
hold off
title("Gráfica de errores")
legend("error en norma L^2", "error en norma uniforme")

En la primera imagen podemos ver como ya a partir de [math] n=5 [/math] la gráfica de [math] f(x) [/math] queda prácticamente solapada por la aproximación [math] f_n(x) [/math]. Esto se confirma en el segundo gráfico, donde vemos que los errores para [math] n=3 [/math] son del orden de [math] 10^{-10} [/math], esto es, que las funciones [math] f_n(x) [/math] aproximan muy correctamente la función [math] f(x) [/math].

5 Cambio de intervalo y aproximación de [math] f(x) = xe^{−x} [/math]

En este apartado pretendemos aproximar la función [math] f(x) = xe^{−x} [/math] definida en el intervalo [math] [1,3] [/math]. Aplicando lo explicado en la introducción teórica a inicios de este artículo, la base trigonométrica de [math] L^2([1,3]) [/math] por periodicidad y por la longitud de los intervalos será la misma que la de [math] L^2([-1,1]) [/math], esto es, [math] \{\frac{1}{2},cos(nx\pi),sen(nx\pi)\}_{n\in \mathbb{N}} [/math].

'Aproximación de [math] f(x) = xe^{−x} [/math] con [math] n= 5 [/math], [math] 10 [/math] y [math] 20 [/math]'

A continuación, de la misma forma que se hizo en apartados anteriores, usamos dicha base para aproximar la función [math] f(x) = xe^{−x} [/math]. Esto queda representado mediante el siguiente código de Matlab:

f =@(x) x.*exp(-x);
m = 20;
coeficientes1 = zeros(1, m); coeficientes2 = zeros(1, m);
N=1000;                         %Number of points
a=1; b=3;                     %Extremes of the interval
h=(b-a)/(N);
u=a:h:b;                       %coordinates of the partition
 w=ones(N+1,1);                 %weights vector
w(1)=1/2; w(N+1)=1/2;

coef0 = h*w'*((f(u))');
for k = 1:m 
    f1=(f(u).*sin(k.*pi.*u))';
    f2=(f(u).*cos(k.*pi.*u))'; %function
    coeficientes1(1,k) = h*w'*f1;              
    coeficientes2(1,k) = h*w'*f2;
end

y = zeros(m,N+1);
y(1,:) = coef0/2+ coeficientes1(1,1).*sin(pi*u) + coeficientes2(1,1).*cos(pi*u);
for k = 2:m
    y(k,:) = y(k-1,:) + coeficientes1(1,k).*sin(k*pi*u) + coeficientes2(1,k).*cos(k*pi*u);
end


figure(1)
hold on
plot(u,f(u),'black')
plot(u,y(5,:), "blue")
plot(u,y(10,:), 'red')
plot(u,y(m,:), 'yellow')
xlabel('x');
ylabel('f(x), f_{5}(x), f_{10}(x), f_{20}(x)');
title('Aproximación de f(x) = xe−x con n= 5, 10 y 20');
hold off
legend("f", "f_{5}", "f_{10}", "f_{20}")

En este gráfico podemos ver como a medida que aumenta la [math] n [/math] la aproximación de la serie trigonométrica se acerca más a la función [math] f(x)[/math]. Además, vemos como de nuevo

6 Base trigonométrica compleja

En este apartado nos piden trabajar con la base trigonométrica compleja [math]\{e^{inx}\}_{n \in \mathbb{Z}}[/math] en [math][-\pi, \pi][/math]. En concreto, buscamos aproximar la función [math]f(x) = 4x(1/2-x)^2[/math] en el intervalo [math][0,1][/math]. Por lo que la base que usaremos será [math]\{e^{2\pi inx}\}_{n \in \mathbb{Z}}[/math]

Debemos tener en cuenta, como indica el enunciado, que en el espacio [math]L^2(-\pi, \pi)[/math] de funciones el producto escalar es:

[math] ( f,g ) _{L^2} = \int_{-\pi}^{\pi} f(x) \overline{g(x)} \hspace{5px} dx [/math]

A continuación se muestra el código empleado para aproximar la función por series con los primeros 5, 10 y 20 términos:

'Aproximación de f(x) = xe−x con n= 5, 10 y 20'

f =@(x) 4*x.*(1/2-x).^2; 
m = 20;

coeficientes1 = zeros(1, 2*m+1); %vector de coefs
N=1000;                          %regla del trapecio para la integral
a=0; b=1;                     
h=(b-a)/(N);
u=a:h:b;                       
w=ones(N+1,1);                 
w(1)=1/2; w(N+1)=1/2;

for k = 1:m
    fun1=(conj(f(u)).*(exp(-2*1i*pi*k.*u)))';      %coefs para n<0
    coeficientes1(1,2*k) = h*w'*fun1;
    fun2=(conj(f(u)).*(exp(2*1i*pi*k.*u)))';       %coefs para n<0     
    coeficientes1(1,2*k+1) = h*w'*fun2;
end
fun3=(f(u))';                                      %coeficiente 0
coeficientes1(1,1) = h*w'*fun3;

y = zeros(m+1,N+1);                                %calculamos las sumas
y(1,:) = coeficientes1(1,1);
for k = 1:m
    y(k+1,:) = y(k,:) + coeficientes1(1,2*k).*(exp(-2*1i*k*pi.*u)) + coeficientes1(1,2*k+1).*(exp(2*1i*k*pi.*u));
end


figure(1)                                          %rep gráfica
hold on
plot(u,f(u),'black')
plot(u,y(5,:), "blue")
plot(u,y(10,:), 'red')
plot(u,y(m,:), 'yellow')
hold off
legend("f", "n=5", "n=10", "n ="+ num2str(m))