|
|
| Línea 1: |
Línea 1: |
| − | import numpy as np
| |
| − | import matplotlib.pyplot as plt
| |
| − | from scipy.integrate import trapezoid
| |
| | | | |
| − | def base_trigonometrica(n):
| |
| − | #Puntos donde se grafica
| |
| − | x = np.linspace(-1,1,1000)
| |
| − | plt.figure(figsize = (10,6))
| |
| − |
| |
| − | #Primer término de la serie
| |
| − | plt.plot(x, np.full_like(x, 0.5), label='1/2 (n=0)', linewidth=2, color='black')
| |
| − |
| |
| − | #Términos trigonométricos
| |
| − | for n in range(1,n):
| |
| − | plt.plot(x, np.cos(n*np.pi*x), color = 'r', linestyle='-')
| |
| − | plt.plot(x, np.sin(n*np.pi*x), color = 'b', linestyle='-')
| |
| − |
| |
| − | plt.title('Primeros términos de la base trigonométrica {1/2, cos(nπx), sin(nπx)}')
| |
| − | plt.xlabel('x')
| |
| − | plt.ylabel('f(x)')
| |
| − | plt.axhline(0, color='black', linewidth=0.5)
| |
| − | plt.axvline(0, color='black', linewidth=0.5)
| |
| − | plt.grid(True, linestyle=':', alpha=0.6)
| |
| − | plt.legend(['1/2','cos','sin'], loc = 'upper right', bbox_to_anchor = (1,1))
| |
| − | plt.tight_layout()
| |
| − |
| |
| − | # Mostrar resultado
| |
| − | plt.show()
| |
| − |
| |
| − |
| |
| − | def serie_fourier_definitiva(f,L,n):
| |
| − |
| |
| − | #Puntos donde se grafica y sus valores exactos
| |
| − | x = np.linspace(-L/2,L/2,1000)
| |
| − | f_ev = f(x)
| |
| − |
| |
| − | def serie_fourier(funcion,x_val,n_terms):
| |
| − | fn = np.zeros_like(x_val)
| |
| − | f_eval = funcion(x_val)
| |
| − |
| |
| − | c0 = trapezoid(2/L * f_eval,x_val)
| |
| − | fn = c0/2
| |
| − | for k in range(1,n_terms+1):
| |
| − | #Coeficientes
| |
| − | cn = trapezoid(2/L * f_eval * np.cos(2*k*np.pi*x_val/L),x_val)
| |
| − | dn = trapezoid(2/L * f_eval * np.sin(2*k*np.pi*x_val/L),x_val)
| |
| − |
| |
| − | #Sumar términos
| |
| − | fn += cn* np.cos(2*k*np.pi*x_val/L) + dn * np.sin(2*k*np.pi*x_val/L)
| |
| − |
| |
| − | return fn
| |
| − |
| |
| − | #Gráfica
| |
| − | plt.figure(figsize = (10,6))
| |
| − | plt.plot(x,f_ev,'k--', label='f(x) original', linewidth=2)
| |
| − |
| |
| − | for i in n:
| |
| − | fn_ev = serie_fourier(f,x,i)
| |
| − | plt.plot(x,fn_ev, label=f'n = {i}')
| |
| − |
| |
| − | plt.title('Aproximación de f(x) mediante serie de Fourier')
| |
| − | plt.legend()
| |
| − | plt.grid(True, alpha=0.3)
| |
| − | plt.show()
| |
| − |
| |
| − | #Cálculo de errores
| |
| − | n_error = []
| |
| − | for i in range(1,100):
| |
| − | n_error.append(i)
| |
| − |
| |
| − | error_l2 = []
| |
| − | error_inf = []
| |
| − |
| |
| − | for i in n_error:
| |
| − | f_n = serie_fourier(f,x,i)
| |
| − | dif = f_ev - f_n
| |
| − |
| |
| − | error_l2.append(np.sqrt(trapezoid(dif**2,x)))
| |
| − | error_inf.append(round(float(np.max(np.abs(dif))),6))
| |
| − |
| |
| − | #Gráfica
| |
| − | plt.figure(figsize = (10,6))
| |
| − | plt.plot(n_error,error_l2, linewidth=2)
| |
| − | plt.title('Error $L^2$')
| |
| − | plt.legend()
| |
| − | plt.grid(True, alpha=0.3)
| |
| − | plt.show()
| |
| − |
| |
| − | def comparar_fourier_cesaro(func, L, n_terms, nombre_funcion="Función"):
| |
| − | """
| |
| − | Calcula y grafica la aproximación de Fourier estándar vs Sumas de Cesàro.
| |
| − |
| |
| − | func: Función original a aproximar.
| |
| − | L: Longitud del intervalo [-L/2, L/2].
| |
| − | n_terms: Número máximo de coeficientes (n).
| |
| − |
| |
| − | """
| |
| − |
| |
| − | def get_coefficients(n):
| |
| − | # Integración por método del trapecio como se indica en el documento
| |
| − | x_int = np.linspace(-L/2, L/2, 2000)
| |
| − | y_int = func(x_int)
| |
| − |
| |
| − | # Coeficiente c0
| |
| − | c0 = (2/L) * trapezoid(y_int, x_int)
| |
| − |
| |
| − | an = []
| |
| − | bn = []
| |
| − | for i in range(1, n + 1):
| |
| − | cos_term = np.cos(2 * np.pi * i * x_int / L)
| |
| − | sin_term = np.sin(2 * np.pi * i * x_int / L)
| |
| − | an.append((2/L) * trapezoid(y_int * cos_term, x_int))
| |
| − | bn.append((2/L) * trapezoid(y_int * sin_term, x_int))
| |
| − | return c0, an, bn
| |
| − |
| |
| − | # Obtener coeficientes
| |
| − | c0, an, bn = get_coefficients(n_terms)
| |
| − |
| |
| − | # Construir sumas parciales S_k para Cesàro
| |
| − | x_range = np.linspace(-L/2,L/2,1000)
| |
| − | sumas_parciales = []
| |
| − | current_sum = np.full_like(x_range, c0 / 2)
| |
| − | sumas_parciales.append(current_sum.copy())
| |
| − |
| |
| − | for i in range(n_terms):
| |
| − | term = an[i] * np.cos(2 * np.pi * (i+1) * x_range / L) + \
| |
| − | bn[i] * np.sin(2 * np.pi * (i+1) * x_range / L)
| |
| − | current_sum += term
| |
| − | sumas_parciales.append(current_sum.copy())
| |
| − |
| |
| − | # Calcular Suma de Cesàro (Promedio de las sumas parciales)
| |
| − | sigma_n = np.mean(sumas_parciales, axis=0)
| |
| − |
| |
| − | # Gráfica
| |
| − | plt.figure(figsize=(12, 6))
| |
| − | plt.plot(x_range, func(x_range), 'k--', label=f"Original: {nombre_funcion}", alpha=0.6)
| |
| − | plt.plot(x_range, current_sum, label=f"Fourier Estándar (n={n_terms})", color='red', alpha=0.5)
| |
| − | plt.plot(x_range, sigma_n, label=f"Suma de Cesàro (n={n_terms})", color='blue', linewidth=2)
| |
| − |
| |
| − | plt.title(f"Aproximación de {nombre_funcion}: Fourier vs Cesàro")
| |
| − | plt.xlabel("x")
| |
| − | plt.ylabel("f(x)")
| |
| − | plt.grid(True, linestyle='--', alpha=0.7)
| |
| − | plt.show()
| |
| − |
| |
| − |
| |
| − |
| |
| − | #Lo que se ejecuta---------------------------------------------------------
| |
| − | base_trigonometrica(6)
| |
| − |
| |
| − | def f(x):
| |
| − | lista = []
| |
| − | for i in x:
| |
| − | if i < 0:
| |
| − | lista.append(0)
| |
| − | else:
| |
| − | lista.append(1)
| |
| − | return np.array(lista)
| |
| − |
| |
| − | def MW(x):
| |
| − | n = 200
| |
| − | a = 1/2
| |
| − | b = 13
| |
| − |
| |
| − | suma = 0
| |
| − | for i in range(n):
| |
| − | suma += a**i * np.cos(b**i * np.pi*x)
| |
| − |
| |
| − | return suma
| |
| − |
| |
| − | def g(x):
| |
| − | return x**2
| |
| − |
| |
| − | serie_fourier_definitiva(f,10,[1,5,10,100])
| |
| − | comparar_fourier_cesaro(f, 10, 50, "Función Discontinua (Escalón)")
| |
| − | comparar_fourier_cesaro(MW, 10, 100, "Monstruo de Weierstrass")
| |
