Series de Fourier (Grupo Eau De Parfum (EDP))

Trabajo realizado por estudiantes
Título	Series de Fourier. Grupo Eau De Parfum (EDP)
Asignatura	EDP
Curso	2023-24
Autores	Lestau Torres, Pablo López Rojo, Celia Muñoz Guijarro, Sofía
Este artículo ha sido escrito por estudiantes como parte de su evaluación en la asignatura

1 Introducción

Para hablar de series de Fourier nos debemos remontar al "problema de la cuerda vibrante". Un problema que fue estudiado por matemáticos del siglo XVIII como d'Alambert, Euler y Daniel Bernouilli entre otros. Este último propuso una solución a este problema mediante la superposición infinita de ondas simples y fue Fourier quien aplicó y perfeccionó estas ideas en 1822 en su libro "Théorie analytique de la Chaleur"; donde sus razonamientos generaron numerosas controversias

En este trabajo abordaremos algunas de ellas y mostraremos la evolución de las conocidas series de Fourier. Pero antes, abordaremos su representación actual empleando el espacio de funciones de cuadrado integrable y espacios de Hilbert.

2 Formulación moderna.

François Marie Charles Fourier

El desarrollo en "serie de Fourier" se realiza sobre funciones de cuadrado integrable, es decir, sobre las funciones que pertenecen al conjunto [math]L^2(Ω;\mathbb{R})[/math], siendo este:

[math]L^2(Ω;\mathbb{R}) = \{ f : Ω \to \mathbb{R} : \int_{Ω} |f(x)|^2 dx \lt \infty \}[/math]

Este conjunto tiene definido el siguiente producto escalar:

[math]\langle f, g \rangle \ \stackrel{\mathrm{def}}{=} \ \frac{1}{2\pi}\int_{-\pi}^{\pi} f(x)\overline{g(x)}\,dx,[/math]

que lo dota de estructura de espacio de Hilbert [math](L^2(Ω,\mathbb{R});\langle , \rangle_{L^2((Ω,\mathbb{R})})[/math].Además, cumple ser separable, por tanto, admite una base hilbertiana, es decir, ortonormal. En el contexto de un intervalo de la forma [math] [-\pi,\pi][/math] una base trigonométrica típica es:

[math] \mathcal{B} =\left\{ \frac{1}{\sqrt{2\pi}} \right\} \cup \left\{ \frac{1}{\sqrt{\pi}} \cos(nx),\frac{1}{\sqrt{\pi}} \sin(nx) \right\}^{\infty}_{n=1} [/math]

A partir de esta podemos transformar el intervalo a otro de la forma \([a; b]\) mediante el cambio de variable \(h(x) = 2\pi b\) , obteniendo

[math] \mathcal{B} =\left\{ \frac{1}{\sqrt{b-a}} \right\} \cup \left\{ \frac{\sqrt{2}}{\sqrt{b-a}} \cos\left(n\left(\frac{2\pi(x-a)}{b-a}-\pi\right)\right),\frac{\sqrt{2}}{\sqrt{b-a}} \sin\left(n\left(\frac{2\pi(x-a)}{b-a}-\pi\right)\right) \right\}^{\infty}_{n=1} [/math]

La relevancia de estas bases radica en su aplicación en la aproximación de funciones mediante series trigonométricas. Por ejemplo, sobre la base proporcionada, la aproximación de Fourier responderá a

[math]f(t) \approx \sum_{n=1}^{\infty}c_n e_n \quad \text{con} \quad c_n=\langle f, e_n \rangle[/math]

Para comprender esta aproximación, es crucial visualizar las componentes de las bases: las funciones seno, coseno y la constante. Así pues, comenzamos dibujando el elemento constante.

Función [math]f(x)=1/2[/math] para distintos valores de n en [-1,1]

import numpy as np
import matplotlib.pyplot as plt

# Parámetros iniciales
N = 10
x = np.linspace(-1, 1, 400)

# Crear las subfiguras para cada n
for n in range(1, N + 1):
    y = 0.5 * np.ones_like(x)
    
    plt.figure(figsize=(6, 4))
    plt.axhline(0, color='gray')  # Línea horizontal en y=0
    plt.axvline(0, color='gray')  # Línea vertical en x=0
    plt.xlim([-1.2, 1.2])
    plt.ylim([-1.2, 1.2])
    plt.plot(x, y, 'b')
    plt.title(f"$n=${n}", fontsize=10)
    plt.grid(True)

A continuación, mostramos los diez primeros elementos de nuestra base correspondientes a la forma [math]\cos(n\pi x)[/math] :

Función [math]\cos(n\pi x)[/math] para distintos valores de n en [-1,1]

import numpy as np
import matplotlib.pyplot as plt

# Parámetros iniciales
N = 10
x = np.linspace(-1, 1, 400)

# Crear las subfiguras para cada n
for n in range(1, N + 1):
    y = np.cos(np.pi * n * x) 
    
    plt.figure(figsize=(6, 4))
    plt.axhline(0, color='gray')  # Línea horizontal en y=0
    plt.axvline(0, color='gray')  # Línea vertical en x=0
    plt.xlim([-1.2, 1.2])
    plt.ylim([-1.2, 1.2])
    plt.plot(x, y, 'b')
    plt.title(f"$n=${n}", fontsize=10)
    plt.grid(True)

Podemos observar que a medida que aumenta n, la frecuencia de la función coseno aumenta, de hecho, por cada incremento de n, se añade una oscilación completa adicional a la función dentro del intervalo. Por otra parte, el periodo se ve afectado por el factor [math]n\pi [/math] a medida que n incrementa, el periodo en [math][-1,1][/math] disminuye, resultando en más ciclos completos de la función coseno. De forma similar, visualizamos los diez primeros términos para el elemento [math]\sin(n \pi x)[/math] de nuestra base trigonométrica.

Función [math]\sin(n\pi x)[/math] para distintos valores de n en [-1,1]

import numpy as np
import matplotlib.pyplot as plt

# Parámetros iniciales
N = 10
x = np.linspace(-1, 1, 400)

# Crear las subfiguras para cada n
for n in range(1, N + 1):
    y = np.sin(np.pi * n * x) 
    
    plt.figure(figsize=(6, 4))
    plt.axhline(0, color='gray')  # Línea horizontal en y=0
    plt.axvline(0, color='gray')  # Línea vertical en x=0
    plt.xlim([-1.2, 1.2])
    plt.ylim([-1.2, 1.2])
    plt.plot(x, y, 'b')
    plt.title(f"$n=${n}", fontsize=10)
    plt.grid(True)

Análogamente al caso del coseno, es claro observar que a medida que n aumenta, la forma de la onda sigue siento sinusoidal, pero con más ciclos dentro del mismo intervalo, y el número de oscilaciones de la función dentro de este aumenta.

3 Del Debate con Euler al Desafío del Fenómeno de Gibbs

3.1 Análisis de las Discrepancias en la Teoría de las Series de Fourier y su Resolución Moderna

Una vez entendidas las series de Fourier de manera moderna, podemos volver al contexto histórico y comprender las controversias. La primera de ellas viene de la mano de Euler.

Daniel Bernouilli introdujo la solución del problema mencionado como una superposición de ondas más simples, específicamente aquellas de la forma [math] u_n(x, t) = \sin(nx) \cos(nt), \quad \forall n \in \mathbb{N}. [/math] Esto significaba que la solución se podía representar como [math] u(x, t) = \sum_{n=1}^{\infty} a_n \sin(nx) \cos(nt) [/math]. Esta solución fue rechazada por Euler, ya que consideraba matemáticamente inaceptable que cualquier función arbitraria pudiera representarse por medio de una suma trigonométrica. En cierta parte, Euler tenía razón, pues solo se pueden aproximar las funciones [math]f[/math] que se encuentre en el espacio de funciones de cuadrado integrable, como por ejemplo, la función [math]f(x)=x e^{-x}[/math]:

[math] \int_{1}^{3} (x e^{-x})^ 2\, dx = \frac{5(e^4 - 5)}{4 e^6} \approx 0.15368 \lt \infty [/math]

De modo que podemos calcular su serie de Fourier en el intervalo [math][1,3][/math]. Primero, calcularemos la base mediante la cual se va a aproximar. Para ello, utilizaremos la base trigonométrica ortogonal definida previamente con los parámetros [math]a = 1[/math] y [math]b = 3[/math], obteniendo:

[math] \mathcal{B} = \{ e_n \}_{n=1}^{\infty} = \left\{ \frac{1}{\sqrt{2}} \right\} \cup \left\{ \cos\left(n\pi (x-2)\right), \sin\left(n\pi (x-2)\right) \right\}_{n=1}^{\infty} [/math]

Donde tanto la serie de Fourier como sus coeficientes se calculan como se menciona en 2 Formulación moderna. Para realizar estos cálculos se ha implementado el siguiente código de Python, que además grafica la aproximación para diferentes valores de n

Resultado de ejecutar el código. Tres aproximaciones de [math]f(x)[/math] variando el número de términos de la serie de Fourier

Gif de las aproximaciones de f(x) a medida que n aumenta de 1 a 50

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

# Definición de la función f(x)
def f(x):
    return x * np.exp(-x)

# Definición de la base trigonométrica
def base(x, n_terms):
    terms = [1/np.sqrt(2)]
    for n in range(1, n_terms + 1):
        terms.append(np.cos(n * np.pi * (x - 2)))
        terms.append(np.sin(n * np.pi * (x - 2)))
    return terms

# Coeficientes de la serie de Fourier <f(x),en >
def fourier_coefficients(f, n_terms):
    a0, _ = quad(f, 1, 3)/np.sqrt(2)  
    coefficients = [a0]
    for n in range(1, n_terms + 1):
        an, _ = quad(lambda x: f(x) * base(x, n_terms)[2*n-1], 1, 3) #<f(x),cos(nπ(x-2)) >
        bn, _ = quad(lambda x: f(x) * base(x, n_terms)[2*n], 1, 3)  #<f(x),sin(nπ(x-2)) >
        coefficients.append(an)
        coefficients.append(bn)
    return coefficients

# Aproximación de la función utilizando la serie de Fourier
def fourier_series(x, coefficients):
    serie = 0
    for i, coeff in enumerate(coefficients):
        series += coeff * base(x, len(coefficients)//2)[i]
    return series

# Intervalo [1,3]
x_values = np.linspace(1, 3, 1000)

# Aproximaciones con diferentes cantidades de términos
seriefor i, n_terms in enumerate([5, 10, 20]):
    plt.figure()  # Nueva figura para cada iteración
    coefficients = fourier_coefficients(f, n_terms)
    approx_values = [fourier_serie(x, coefficients) for x in x_values]
    plt.plot(x_values, approx_values, label=f'{n_terms} términos', color='blue')

    # Función original
    plt.plot(x_values, f(x_values), label='Función original', color='black', linewidth=1)

    plt.xlabel('$x$')
    plt.ylabel('$f(x)$')
    plt.title(f'Aproximación de f(x) con {n_terms} términos de serie de Fourier')
    plt.legend()
    plt.grid(False)
    plt.show()

Supongamos a continuación que queremos aproximar la función a otro tipo de intervalos, unos de la forma [math][0,L][/math]. Entonces, es interesante estudiar una posible simetría de la misma de modo que podemos extenderla primero al intervalo [math][-L, L][/math]. Esta extensión puede realizarse de dos maneras: de forma par o impar.

Si elegimos extender la función de manera par, esto implica que [math]f(x)[/math] será igual a [math]f(-x)[/math] para [math]x[/math] en el intervalo [math](-L, 0)[/math]. En este caso, la serie de Fourier de la función incluirá únicamente [math]\cos (n\theta)[/math] y un término constante. Entonces, la expresión de la aproximación será: [math] \\ f(x) \approx \frac{a_0}{2} + \sum_{n=1}^{\infty} a_n \cos \left( \frac{n\pi x}{L} \right), \quad a_0 = \frac{2}{L}\int_{0}^{1} f(x) dx ,\quad a_n = \frac{2}{L} \int_{0}^{L} f(x) \cos \left( \frac{n\pi x}{L} \right) dx [/math]

Si por el contrario optamos por la extensión de manera impar, es decir, [math]f(x)=-f(-x)[/math] para [math]x \in (-L, 0)[/math], la serie de Fourier consistirá solo en [math]\sin (n\theta)[/math]. Por tanto,

[math] f(x) \approx \sum_{n=1}^{\infty} b_n \sin \left( \frac{n\pi x}{L} \right), \quad b_n = \frac{2}{L} \int_{0}^{L} f(x) \sin \left( \frac{n\pi x}{L} \right) dx [/math]

Teóricamente, se debe a que para semi intervalos de la forma [math][0,L][/math], los siguientes conjuntos son bases Hilbertianas completas y ortogonales:

[math] \left\{ \frac{1}{2}, \cos \left( \frac{n\pi x}{L} \right) \right\}_{n=1}^{\infty} [/math], [math] \quad \left\{ \sin \left( \frac{n\pi x}{L} \right) \right\}_{n=1}^{\infty} [/math]

Esta simplificación en la composición de términos coseno o seno permite una convergencia más rápida de la serie de Fourier, lo que significa que se puede obtener una buena aproximación de la función original con un número menor de términos, sin tener que hacer uso de la base general para un intervalo [math][a,b][/math].

Es interesante ver un ejemplo. Supongamos la función [math] f(x) = x(1-x) [/math] en el intervalo [math] [0, 1] [/math] La extendemos de forma impar en [math] [-1, 1] [/math], de la siguiente forma

[math] g(x) = \begin{cases} x(x + 1) & \text{si } -1 \leq x \lt 0 \\ x(1 - x) & \text{si } 0 \leq x \leq 1 \end{cases} [/math]

Obteniendo una función continua y diferenciable. Mostraremos la extensión impar mediante este código.

Extensión impar de la función [math] f(x) = x(1-x) [/math]

import numpy as np
import matplotlib.pyplot as plt

# Definir las funciones f1(x) y f2(x)
def f1(x):
    return x * (x + 1)

def f2(x):
    return x * (1 - x)

# Rango de valores de x para f1(x) y f2(x)
x_values_1 = np.linspace(-1, 0, 100)
x_values_2 = np.linspace(0, 1, 100)

# Calcular los valores de las funciones f1(x) y f2(x)
f1_values = f1(x_values_1)
f2_values = f2(x_values_2)

# Crear el gráfico
plt.figure(figsize=(8, 5))

# Graficar f1(x) y f2(x)
plt.plot(x_values_1, f1_values, label='g(x) = -f(-x) = x * (x + 1)', color='blue')
plt.plot(x_values_2, f2_values, label='f(x) = x * (1 - x)', color='red')

# Configuración de ejes y leyenda
plt.xlabel('x')
plt.ylabel('f(x)')
plt.title('Extensión Impar')
plt.legend()

# Mostrar la gráfica
plt.grid(False)
plt.show()

De modo que la suma parcial hasta el término [math] n [/math] de la serie viene dada por la expresión:

[math] f_n(x) = \sum_{k=1}^{n} a_k \sin(k\pi x), [/math]

donde los coeficientes [math] a_k [/math] se calculan como:

[math] a_k = 2 \int_{0}^{1} f(x) \sin(k\pi x) \, dx. [/math]

Para estudiar la serie, es interesante observar cómo esta evoluciona en función del número de términos que tomemos. Por esta razón, se ha implementado un código en Python que calcula la serie de Fourier hasta cierto término [math] n [/math] y la represente junto a la función original.

Serie de Fourier cambiando el número de términos

import numpy as np
import matplotlib.pyplot as plt

# Definir la función f(x) para x en [0, 1]
def f(x):
    return x * (1 - x)

# Coeficientes de la serie de Fourier
def fourier_coefficients(func, num_terms, x):
    coefficients = []
    for n in range(1, num_terms + 1):
        a_n = 2 * np.trapz(func(x) * np.sin(np.pi * n * x), x)
        coefficients.append(a_n)
    return coefficients

# Serie de Fourier a partir de los coeficientes
def fourier_series(x, coefficients):
    series_sum = np.zeros_like(x)
    for n, coeff in enumerate(coefficients, start=1):
        series_sum += coeff * np.sin(np.pi * n * x)
    return series_sum

num_terms_values = [1, 5, 10, 20]  # Valores de num_terms
x = np.linspace(0, 1, 1000)

# Crear las gráficas
fig, axes = plt.subplots(2, 2, figsize=(12, 8))

for i, ax in enumerate(axes.flat):
    num_terms = num_terms_values[i]
    
    # Calcular los coeficientes de Fourier
    coefficients = fourier_coefficients(f, num_terms, x)

    # Calcular la serie de Fourier
    fourier_values = fourier_series(x, coefficients)

    # Graficar la función original (en negro)
    ax.plot(x, f(x), label='Función Original', color='black')

    # Graficar la serie de Fourier (en azul)
    ax.plot(x, fourier_values, label='Serie de Fourier', color='blue')

    # Configuración de los ejes y la leyenda
    ax.set_xlabel('x')
    ax.set_ylabel('f(x)')
    ax.set_title(f'Número de términos = {num_terms}')
    ax.legend()
    ax.grid(False)  # Quitar la cuadrícula de fondo

plt.tight_layout()
plt.show()

Para evaluar adecuadamente la calidad de esta aproximación, es crucial cuantificar su convergencia a la función original. Para ello, calcularemos el error o distancia mediante [math]\sqrt{\int_{0}^{1} |f(x) - f_n(x)|^2 \, dx}[/math] y [math]\sup_{x \in [0,1]} |f(x) - f_n(x)|[/math]), siendo ambas normas de [math] L^2 [/math]. Para mayor comprensión se representará en una gráfica para valores diferentes de n. Por razones visuales es interesante estudiarlo en escala logarítmica.

Errores de la aproximación

# Rango de valores de x para f1(x) y f2(x)
x_values_1 = np.linspace(-1, 0, 100)
x_values_2 = np.linspace(0, 1, 100)
# Función para calcular el error en la norma L2 para un número dado de coeficientes de Fourier
def calculate_l2_norm_error(func, num_terms, x):
    coefficients = fourier_coefficients(func, num_terms, x)
    fourier_series_values = fourier_series(x, coefficients)
    l2_norm_error = np.sqrt(np.trapz(np.abs(func(x) - fourier_series_values)**2, x))
    return l2_norm_error

# Función para calcular el error máximo absoluto para un número dado de coeficientes de Fourier
def calculate_max_absolute_error(func, num_terms, x):
    coefficients = fourier_coefficients(func, num_terms, x)
    fourier_series_values = fourier_series(x, coefficients)
    max_absolute_error = np.max(np.abs(func(x) - fourier_series_values))
    return max_absolute_error


# Rango de valores para el número de coeficientes de Fourier
num_terms_range = np.arange(51) 

# Lista para almacenar los errores en la norma L2 y el error máximo absoluto
l2_norm_errors = []
max_absolute_errors = []

# Calcular los errores para cada número de coeficientes de Fourier
for num_terms in num_terms_range:
    l2_norm_error = calculate_l2_norm_error(f, num_terms, x)
    l2_norm_errors.append(l2_norm_error)
    
    max_absolute_error = calculate_max_absolute_error(f, num_terms, x)
    max_absolute_errors.append(max_absolute_error)


# Crear una figura con dos subtramas
fig, (ax1, ax) = plt.subplots(1, 2, figsize=(10, 5))

# Graficar el error en la norma L2
ax1.plot(num_terms_range, np.log(l2_norm_errors), color='blue')
ax1.set_xlabel('Número de coeficientes de Fourier')
ax1.set_ylabel('log(Error)')
ax1.set_title(r'Error $\left(\int_{0}^{1} |f(x) - f_n(x)|^2 dx\right)^{1/2}$')
ax1.grid(False)

# Graficar el error máximo absoluto
ax.plot(num_terms_range, np.log(max_absolute_errors), color='red')
ax.set_xlabel('Número de coeficientes de Fourier')
ax.set_ylabel('log(Error)')
ax.set_title(r'Error $\sup_{x\in[0,1]} |f(x) - f_n(x)|$')
ax.grid(False)

plt.tight_layout()
plt.show()

Observamos que ambos errores parecen decrecer de forma similar, en función del número de términos que tomemos.

3.2 Limitaciones de las Series de Fourier en Funciones Discontinuas y la Resolución Mediante Sumas de Cesàro

La segunda cuestión abordada por Fourier fue la generalidad en la representación de funciones discontinuas. Para competir con otros métodos, introdujo la técnica de unir los saltos por líneas verticales en sus diagramas en el artículo de 1807, lo que permitió la continuidad visual de las funciones. Sin embargo, esta aproximación no pudo anticipar completamente el 'fenómeno Gibbs', descrito posteriormente por J.W. Gibbs. Este fenómeno reveló que las sumas parciales de las series de Fourier presentan un comportamiento oscilatorio anómalo cerca de puntos de discontinuidad, lo que evidenció una limitación en la teoría de Fourier en términos de convergencia uniforme.

Para estudiar este caso, utilizaremos un ejemplo específico: la función [math] f(x) = \begin{cases} 1 & \text{si } x \leq \frac{1}{2} \\ 0 & \text{si } x \gt \frac{1}{2} \end{cases} [/math], la cual es discontinua en [math] x = \frac{1}{2} [/math]. Sin embargo, haciendo uso de lo mencionado anteriormente, puede extenderse de manera par al intervalo [math] [-1,1] [/math] como [math] g(x) = \begin{cases} 0 & \text{si } x \lt -\frac{1}{2} \\ 1 & \text{si } -\frac{1}{2} \leq x \leq \frac{1}{2} \\ 0 & \text{si } x \gt \frac{1}{2} \end{cases} [/math].

Extensión par de la función [math] f(x) = x(1-x) [/math]

import numpy as np
import matplotlib.pyplot as plt

# Definir la función por partes
def f(x):
    if x < -1/2:
        return 0
    elif -1/2 <= x <= 1/2:
        return 1
    else:
        return 0

# Rango de valores de x
x_values = np.linspace(-1, 1, 1000)

# Calcular los valores de la función f(x)
f_values = [f(x) for x in x_values]

# Crear el gráfico
plt.figure(figsize=(8, 5))

# Graficar f(x) en rojo para x <= 0 y en azul para x > 0
plt.plot(x_values[x_values <= 0], f_values[:len(x_values[x_values <= 0])], color='red', label='Extensión par', linestyle='-')
plt.plot(x_values[x_values > 0], f_values[len(x_values[x_values <= 0]):], color='blue', label='f(x)', linestyle='-')

# Configuración de ejes y leyenda
plt.xlabel('x')
plt.ylabel('g(x)')
plt.title('Extensión par')
plt.legend()

# Mostrar la gráfica
plt.grid(False)
plt.show()

Por tanto, podemos calcular la serie mediante [math]f(x) \approx \frac{a_0}{2} + \sum_{n=1}^{\infty} a_n \cos \left( n\pi x \right)[/math]

[math]a_0 = \frac{2}{1} \int_{0}^{1} f(x) dx, \quad a_n = \frac{2}{1} \int_{0}^{1} f(x) \cos \left( n\pi x \right) dx[/math]

Para estudiar la serie, es interesante observar cómo evoluciona en función del número de términos que tomemos. Por esta razón, se ha implementado un código en Python que calcula la serie de Fourier hasta cierto término [math] n [/math], calculando así los coeficientes.

Serie de Fourier para distinto número de términos

import numpy as np
import matplotlib.pyplot as plt

# Definir la función f(x)
def f(x):
    return np.where(x <= 0.5, 1, 0)

# Función para calcular los coeficientes de Fourier para n términos
def fourier_coefficients(f, n, x):
    coefficients = np.zeros(n + 1)
    for k in range(n + 1):
        coefficients[k] = 2 * np.trapz(f(x) * np.cos(k * np.pi * x), x)
    coefficients[0] = np.trapz(f(x), x)
    return coefficients

# Función para calcular la serie de Fourier para n términos
def fourier_series(x, coefficients):
    series = np.zeros_like(x)
    series += coefficients[0] # Término a_0
    for k in range(1, len(coefficients)):
        series += coefficients[k] * np.cos(k * np.pi * x)
    return series
#GRÁFICAS
# Parámetros
num_terms_values = [5, 10, 20, 50] # Valores de num_terms
x = np.linspace(0, 1, 1000)

# Crear las gráficas
fig, axes = plt.subplots(2, 2, figsize=(12, 8))

for i, ax in enumerate(axes.flat):
    num_terms = num_terms_values[i]

# Calcular los coeficientes de Fourier
    coefficients = fourier_coefficients(f, num_terms, x)

# Calcular la serie de Fourier
    fourier_values = fourier_series(x, coefficients)

# Graficar la función original (en negro)
    ax.plot(x, f(x), label='Función Original', color='black')

# Graficar la serie de Fourier (en azul)
    ax.plot(x, fourier_values, label='Serie de Fourier', color='blue')

# Configuración de los ejes y la leyenda
    ax.set_xlabel('x')
    ax.set_ylabel('f(x)')
    ax.set_title(f'Número de términos = {num_terms}')
    ax.legend()
    ax.grid(False) # Quitar la cuadrícula de fondo

plt.tight_layout()
plt.show()

Para analizar las gráficas resultantes, nos enfocaremos inicialmente en la discontinuidad en [math] x = \frac{1}{2} [/math]. Observamos que la serie de Fourier en dicho punto parece converger a [math] \frac{1}{2} [/math], precisamente el punto medio entre [math] f\left(\frac{1}{2}^-\right) [/math] y [math] f\left(\frac{1}{2}^+\right) [/math]. Este comportamiento es significativo ya que representa cómo la serie de Fourier aborda la discontinuidad en la función original. A este resultado se le conoce como:

[math]\textbf{Teorema (Convergencia puntual de las series de Fourier)}[/math]: Sea [math] f(x) [/math] una función continua a trozos en [math][-L, L][/math], cuya derivada existe y es continua a trozos. Entonces, la serie de Fourier asociada a [math] f(x) [/math] converge en cada punto [math] x \in [-L, L] [/math] y su suma,

[math] f_n(x) = \frac{a_0}{2} + \sum_{n=1}^{\infty} \left( a_n \cos \left( \frac{n \pi x}{L} \right) + b_n \sin \left( \frac{n \pi x}{L} \right) \right), [/math]

verifica las siguientes relaciones:

• [math] f_n(x_0) = f(x_0) [/math] si [math] x_0 \in (-L, L) [/math] y [math] f(x) [/math] es continua en el punto [math] x_0 [/math],
• [math] f_n(x_0) = \frac{f(x_0^+) + f(x_0^-)}{2} [/math] si [math] x_0 \in (-L, L) [/math] y [math] f(x) [/math] es discontinua en el punto [math] x_0 [/math],

donde [math] f(x_0^+) [/math] y [math] f(x_0^-) [/math] denotan los límites laterales de [math] f(x) [/math] en el punto [math] x_0 [/math].

Cuando nos aproximamos a los puntos de discontinuidad, se observa un fenómeno en las series de Fourier: las sumas parciales no se ajustan suavemente al valor medio esperado según el Teorema de Convergencia. Más bien, estas sumas exhiben un comportamiento característico al sobrepasar el valor de la función en los extremos del salto. Este fenómeno peculiar, se identifica como el fenómeno de Gibbs. Curiosamente, incluso al incrementar el número de términos en la suma parcial, la situación persiste sin cambio aparente. Una estrategia para atenuar este efecto es emplear las sumas de Cesàro. Este método implica calcular promedios de las sumas parciales antes de tomar el límite. Se basa en un resultado conocido para sucesiones numéricas que establece que si el límite de la sucesión \(\{b_n\}\) es \(\ell\), entonces el promedio de la sucesión \((b_1 + b_2 + \dots + b_n)/n\) también converge a \(\ell\).

Para una serie \(\sum_{n=1}^{\infty} a_n\), estudiaremos el límite de la sucesión \(\sigma_N = (s_1 + s_2 + \dots + s_N)/N\), donde \(s_k\) es la \(k\)-ésima suma parcial. Este resultado asegura que, si la serie es convergente, \(\sigma_N\) converge al mismo valor que \(s_N\). Aplicando este concepto a nuestras series de Fourier, la fórmula para la suma de Cesàro hasta el término \(N\) sería \(S_N(x) = \frac{1}{N+1} \sum_{n=0}^{N} f_n(x)\).

Para ilustrar esto veremos unas gráficas donde se observan las sumas de Cesàro, variando el número de coeficientes.

Sumas de Cesaro distinto número de términos

import numpy as np
import matplotlib.pyplot as plt

# Definir la función f(x)
def f(x):
    return np.where(x <= 0.5, 1, 0)

# Función para calcular los coeficientes de Fourier para n términos
def fourier_coefficients(f, n, x):
    coefficients = np.zeros(n + 1)
    for k in range(n + 1):
        coefficients[k] = 2 * np.trapz(f(x) * np.cos(k * np.pi * x), x)
    coefficients[0] = np.trapz(f(x), x)
    return coefficients

# Función para calcular la serie de Fourier para n términos
def fourier_series(x, coefficients):
    series = np.zeros_like(x)
    series += coefficients[0] # Término a_0
    for k in range(1, len(coefficients)):
        series += coefficients[k] * np.cos(k * np.pi * x)
    return series

# Función para calcular las sumas de Cesàro
def cesaro_sum(f, x, num_terms):
    cesaro_series = np.zeros_like(x)
    for n in range(num_terms + 1):
        coefficients = fourier_coefficients(f, n, x)
        fourier_partial = fourier_series(x, coefficients)
        cesaro_series += fourier_partial
    return cesaro_series / (num_terms + 1)
#GRÁFICAS
# Parámetros
num_terms_values = [1, 5, 10, 50] # Valores de N
x = np.linspace(0, 1, 1000)

# Crear las gráficas
fig, axes = plt.subplots(2, 2, figsize=(12, 8))

for i, ax in enumerate(axes.flat):
    num_terms = num_terms_values[i]

# Calcular las sumas de Cesàro
    cesaro_values = cesaro_sum(f, x, num_terms)

# Graficar la función original
    ax.plot(x, f(x), label='Función Original', color='black')

# Graficar las sumas de Cesàro
    ax.plot(x, cesaro_values, label=f'Sumas de Cesàro (N={num_terms})', color='blue')

# Configuración de los ejes y la leyenda
    ax.set_xlabel('x')
    ax.set_ylabel('f(x)')
    ax.set_title(f'N={num_terms}')
    ax.legend()
    ax.grid(False)

plt.tight_layout()
plt.show()

Observamos en las gráficas como claramente el fenómeno de Gibbs se ha apaciguado. Por último, de forma análoga a como hicimos en el primer ejemplo, es interesante estudiar el error, para ello hemos calculado los errores de [math] L^2 [/math].

Errores de la aproximación

# Función para calcular el error en la norma L2 para un número dado de coeficientes de Fourier
def calculate_l2_norm_error(func, num_terms, x):
    casaro=cesaro_sum(func,x,num_terms)
    l2_norm_error = np.sqrt(np.trapz(np.abs(func(x) - casaro)**2, x))
    return l2_norm_error

# Función para calcular el error máximo absoluto para un número dado de coeficientes de Fourier
def calculate_max_absolute_error(func, num_terms, x):
    casaro=cesaro_sum(func,x,num_terms)
    max_absolute_error = np.max(np.abs(func(x) - casaro))
    return max_absolute_error

# Parámetros
x = np.linspace(0, 1, 1000)

# Rango de valores para el número de coeficientes de Fourier
num_terms_range = np.arange(101)

# Lista para almacenar los errores en la norma L2 y el error máximo absoluto
l2_norm_errors = []
max_absolute_errors = []

# Calcular los errores para cada número de coeficientes de Fourier
for num_terms in num_terms_range:
l2_norm_error = calculate_l2_norm_error(f, num_terms, x)
l2_norm_errors.append(l2_norm_error)

max_absolute_error = calculate_max_absolute_error(f, num_terms, x)
max_absolute_errors.append(max_absolute_error)


# Crear una figura con dos subtramas
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10, 5))

# Graficar el error en la norma L2
ax1.plot(num_terms_range, l2_norm_errors, label='Error en la norma L2', color='blue')
ax1.set_xlabel('Número de coeficientes de Fourier')
ax1.set_ylabel('Error')
ax1.set_title(r'Error $\left(\int_{0}^{1} |f(x) - f_n(x)|^2 dx\right)^{1/2}$')
ax1.grid(False)

# Graficar el error máximo absoluto
ax2.plot(num_terms_range, max_absolute_errors, label='Error máximo absoluto', color='red')

ax2.set_xlabel('Número de coeficientes de Fourier')
ax2.set_ylabel('Error')
ax2.set_title(r'Error $\sup_{x\in[0,1]} |f(x) - f_n(x)|$')
ax2.grid(False)

plt.tight_layout()
plt.show()

4 Nuevas bases: base compleja

El estudio sobre las series de Fourier y las bases trigonométricas continúa y aparece la conceptualización de las bases trigonométricas complejas. Esto no fue obra directa de Fourier, sino de matemáticos que le sucedieron. Inspirados por la fórmula de Euler, [math]e^{ix} = \cos(x) + i\sin(x)[/math], que demuestra cómo las funciones exponenciales complejas pueden expresarse a través de sus componentes real e imaginario (correspondientes a las funciones coseno y seno, respectivamente). Estos matemáticos adoptaron las exponenciales complejas para simplificar y extender la formulación de las series de Fourier. Por lo que definieron la base trigonométrica compleja :[math]\{e^{inx}\}_{n \in\mathbb{Z}}[/math] en el intervalo [[math]-\pi,\pi[/math]]

Para estudiar y comprender mejor este nuevo concepto, vamos a calcular la base trigonométrica adaptada al intervalo [0,1]. Esta nueva base es de la forma [math]\{ e^{2\pi i n x} \}_{n \in \mathbb{Z}}[/math]. Y por tanto los coeficientes y la serie de Fourier responden a

[math]c_n = \frac{1}{2\pi} \int_{0}^{1} f(x) \cdot e^{-2\pi i n x} \, dx[/math]

[math]f(x) \approx \sum_{n=-\infty}^{\infty} c_n e^{2\pi i n x}[/math]

Tras esta definición aproximaremos la función [math]f(x) = 4x\left(\frac{1}{2} - x\right)^2 [/math] con los primeros 5, 10 y 20 términos de la serie.

Aproximación de la función [math]f(x) = 4x\left(\frac{1}{2} - x\right)^2 [/math]

import numpy as np
import matplotlib.pyplot as plt

# Definir la función f(x)
def f(x):
    return 4 * x * (1/2 - x)**2

# Coeficientes de Fourier complejos c_n
def fourier_coefficients_complex(f, n, x):
    x = np.linspace(0, 1, 1000)
    c_n = np.trapz(f(x) * np.exp(-2 * np.pi * 1j * n * x), x)
    return c_n

# Aproximación de la función f(x) con los primeros n términos de la serie de Fourier compleja
def fourier_series_complex(f, n, x):
    series_sum = 0
    for i in range(-n, n + 1):
        c_n = fourier_coefficients_complex(f, i, x)
        series_sum += c_n * np.exp(2 * np.pi * 1j * i * x)
    return series_sum

# Parámetros
num_terms_values = [5, 10, 20, 50]  # Valores de num_terms
x = np.linspace(0, 1, 1000)

# Crear las gráficas
fig, axes = plt.subplots(2, 2, figsize=(12, 8))

for i, ax in enumerate(axes.flat):
    num_terms = num_terms_values[i]
    
    # Calcular los coeficientes de Fourier
    coefficients = fourier_coefficients_complex(f, num_terms, x)

    # Calcular la serie de Fourier
    fourier_values = fourier_series_complex(f, num_terms, x)

    # Graficar la función original (en negro)
    ax.plot(x, f(x), label='Función Original', color='black')

    # Graficar la serie de Fourier (en azul)
    ax.plot(x, fourier_values, label='Serie de Fourier', color='blue')

    # Configuración de los ejes y la leyenda
    ax.set_xlabel('x')
    ax.set_ylabel('f(x)')
    ax.set_title(f'Número de términos = {num_terms}')
    ax.legend()
    ax.grid(False)  # Quitar la cuadrícula de fondo

plt.tight_layout()
plt.show()