Berechnen Sie die Fourier-Reihe mit dem Trigonometrie Ansatz

Question

Ich versuche, die Fourier-Reihen-Funktion gemäß den folgenden Formeln zu implementieren:

... wo ...

... und ...

Hier ist meine Herangehensweise an das Problem:

import numpy as np 
import pylab as py 

# Define "x" range. 
x = np.linspace(0, 10, 1000) 

# Define "T", i.e functions' period. 
T = 2 
L = T/2 

# "f(x)" function definition. 
def f(x): 
    return np.sin(np.pi * 1000 * x) 

# "a" coefficient calculation. 
def a(n, L, accuracy = 1000): 
    a, b = -L, L 
    dx = (b - a)/accuracy 
    integration = 0 
    for i in np.linspace(a, b, accuracy): 
     x = a + i * dx 
     integration += f(x) * np.cos((n * np.pi * x)/L) 
    integration *= dx 
    return (1/L) * integration 

# "b" coefficient calculation. 
def b(n, L, accuracy = 1000): 
    a, b = -L, L 
    dx = (b - a)/accuracy 
    integration = 0 
    for i in np.linspace(a, b, accuracy): 
     x = a + i * dx 
     integration += f(x) * np.sin((n * np.pi * x)/L) 
    integration *= dx 
    return (1/L) * integration 

# Fourier series. 
def Sf(x, L, n = 10): 
    a0 = a(0, L) 
    sum = 0 
    for i in np.arange(1, n + 1): 
     sum += ((a(i, L) * np.cos(n * np.pi * x)) + (b(i, L) * np.sin(n * np.pi * x))) 
    return (a0/2) + sum  

# x axis. 
py.plot(x, np.zeros(np.size(x)), color = 'black') 

# y axis. 
py.plot(np.zeros(np.size(x)), x, color = 'black') 

# Original signal. 
py.plot(x, f(x), linewidth = 1.5, label = 'Signal') 

# Approximation signal (Fourier series coefficients). 
py.plot(x, Sf(x, L), color = 'red', linewidth = 1.5, label = 'Fourier series') 

# Specify x and y axes limits. 
py.xlim([0, 10]) 
py.ylim([-2, 2]) 

py.legend(loc = 'upper right', fontsize = '10') 

py.show() 

... und hier ist das, was ich nach dem Ergebnis Plotten:

ich gelesen habe die How to calculate a Fourier series in Numpy? und ich habe diesen Ansatz bereits implementiert. Es funktioniert gut, aber es verwenden die Expotential-Methode, wo ich auf Trigonometrie-Funktionen und die rechteckige Methode im Fall der Berechnung der Integraionen für a_{n} und b_{n} Koeffizienten konzentrieren möchte.

Vielen Dank im Voraus.

UPDATE (GELÖST)

Schließlich ist hier ein funktionierendes Beispiel des Codes. Allerdings werde ich mehr Zeit darauf verwenden, also wenn es etwas gibt, das verbessert werden kann, wird es getan werden.

from __future__ import division 
import numpy as np 
import pylab as py 

# Define "x" range. 
x = np.linspace(0, 10, 1000) 

# Define "T", i.e functions' period. 
T = 2 
L = T/2 

# "f(x)" function definition. 
def f(x): 
    return np.sin((np.pi) * x) + np.sin((2 * np.pi) * x) + np.sin((5 * np.pi) * x) 

# "a" coefficient calculation. 
def a(n, L, accuracy = 1000): 
    a, b = -L, L 
    dx = (b - a)/accuracy 
    integration = 0 
    for x in np.linspace(a, b, accuracy): 
     integration += f(x) * np.cos((n * np.pi * x)/L) 
    integration *= dx 
    return (1/L) * integration 

# "b" coefficient calculation. 
def b(n, L, accuracy = 1000): 
    a, b = -L, L 
    dx = (b - a)/accuracy 
    integration = 0 
    for x in np.linspace(a, b, accuracy): 
     integration += f(x) * np.sin((n * np.pi * x)/L) 
    integration *= dx 
    return (1/L) * integration 

# Fourier series. 
def Sf(x, L, n = 10): 
    a0 = a(0, L) 
    sum = np.zeros(np.size(x)) 
    for i in np.arange(1, n + 1): 
     sum += ((a(i, L) * np.cos((i * np.pi * x)/L)) + (b(i, L) * np.sin((i * np.pi * x)/L))) 
    return (a0/2) + sum 

# x axis. 
py.plot(x, np.zeros(np.size(x)), color = 'black') 

# y axis. 
py.plot(np.zeros(np.size(x)), x, color = 'black') 

# Original signal. 
py.plot(x, f(x), linewidth = 1.5, label = 'Signal') 

# Approximation signal (Fourier series coefficients). 
py.plot(x, Sf(x, L), '.', color = 'red', linewidth = 1.5, label = 'Fourier series') 

# Specify x and y axes limits. 
py.xlim([0, 5]) 
py.ylim([-2.2, 2.2]) 

py.legend(loc = 'upper right', fontsize = '10') 

py.show() 

Answer 1

Betrachten Sie Ihren Code in einer anderen Art und Weise, Block für Block zu entwickeln. Sie sollten überrascht sein, wenn ein Code wie dieser beim ersten Versuch funktionieren würde. Debugging ist eine Option, wie @ tom10 sagte. Die andere Möglichkeit besteht darin, den Code Schritt für Schritt im Interpreter zu entwickeln, besser noch mit ipython.

Oben erwarten Sie, dass b_1000 nicht Null ist, da der Eingang f(x) eine Sinuskurve mit einem 1000 darin ist. Sie erwarten auch, dass alle anderen Koeffizienten Null sind, oder?

Dann sollten Sie nur auf die Funktion b(n, L, accuracy = 1000) konzentrieren. Wenn man es anschaut, laufen 3 Dinge falsch. Hier sind einige Hinweise.

• die Multiplikation von dx ist innerhalb der Schleife. Dessen sicher?
• in der Schleife, i soll eine ganze Zahl sein, oder? Ist es wirklich eine ganze Zahl? durch Prototyping oder Debugging würden Sie dies
• vorsichtig sein, wenn Sie (1/L) oder einen ähnlichen Ausdruck schreiben. Wenn Sie Python2.7 verwenden, tun Sie wahrscheinlich falsch. Wenn nicht, verwenden Sie mindestens from __future__ import division an der Spitze Ihrer Quelle. Lesen Sie this PEP, wenn Sie nicht wissen, worüber ich spreche.

Wenn Sie diese 3 Punkte adressieren, wird b() funktionieren. Dann denke an a in ähnlicher Weise.