Fourierov red

Teorem o konvergenciji

Neka je $f$ po djelovima glatka periodična funkcija s periodom $T=x_1-x_0$ koja na intervalu $[x_0,x_1]$ zadovoljava Dirichletove uvjete:

$f$ je po djelovima neprekidna i njezini prekidi su prve vrste
$f$ je ili monotona ili ima konačno strogih ekstrema

Tada Fourierov red $$ S(x)=\frac{a_0}{2}+\sum_{n=1}^\infty a_n \cos \frac{2 n \pi}{T}x + b_n \sin \frac{2 n \pi}{T}x, $$ gdje je $$ a_0=\frac{2}{T} \int_{x_0}^{x_1} f(x) \, dx, \\ a_n=\frac{2}{T} \int_{x_0}^{x_1} f(x) \cos \frac{2 n \pi}{T}x \,dx,\\ b_n=\frac{2}{T} \int_{x_0}^{x_1} f(x) \sin \frac{2 n \pi}{T}x \,dx, $$ konvergira u svakoj točki $x\in [x_0,x_1]$ i vrijedi $$ S(x)=\begin{cases} f(x), \textrm{ ako je } f \textrm{ neprekinuta u točki } x,\\ \frac{1}{2} [f(x-0)+f(x+0)], \textrm{ inače.} \end{cases} $$

Simboličko računanje

SymPy je paket za simboličko računanje preuzet iz Python-a, Interact je paket za jednostavno manipuliranje parametrima, a Winston je jedan od paketa za crtanje.



In [1]:

    
using SymPy
using Interact
using Winston

Definirajmo simboličku varijablu $x$, simboličku cjelobrojnu (integer) varijablu $n$, i granice intervala $x_0$ i $x_1$.



In [2]:

    
x=Sym("x")
n=symbols("n",integer=true)
x_0=Sym("x_0")
x_1=Sym("x_1")









    Out[2]:





$$x_{1}$$

Računanje koeficijenata



In [3]:

    
a(f,n,x_0,x_1)=2*integrate(f(x)*cos(2*PI*n*x/(x_1-x_0)),(x,x_0,x_1))/(x_1-x_0)
b(f,n,x_0,x_1)=2*integrate(f(x)*sin(2*PI*n*x/(x_1-x_0)),(x,x_0,x_1))/(x_1-x_0)









    Out[3]:





b (generic function with 1 method)

Zadavanje funkcije $f(x)$ i granica intervala $[x_0,x_1]$

Potrebno je koristiti predefiniranu simboličku varijablu PI za razliku od varijabli pi ili $\pi$ kojima je definirana Float64 vrijednost.



In [4]:

    
f(x)=x
x_0=-PI
x_1=PI









    Out[4]:





$$\pi$$

Pogledajmo koeficijente



In [5]:

    
a(f,n,x_0,x_1)









    Out[5]:





$$0$$



In [6]:

    
b(f,n,x_0,x_1)









    Out[6]:





$$\frac{1}{\pi} \begin{cases} 0 & \text{for}\: n = 0 \\- \frac{2 \pi}{n} \left(-1\right)^{n} & \text{otherwise} \end{cases}$$

Računanje i crtanje sume

Za crtanje nam trebaju numeričke granice.



In [7]:

    
x0=N(x_0)
x1=N(x_1)
@manipulate for K=1:10
    S=a(f,0,x_0,x_1)/2+sum([a(f,n,x_0,x_1)*cos(2*PI*n*x/(x_1-x_0))+b(f,n,x_0,x_1)*sin(2*PI*n*x/(x_1-x_0)) for n=1:K])
    g(x)=S(x)
    plot(f,x0,x1)
    oplot(g,x0,x1,"r")
    title("Funkcija i "*string(K)*" članova reda")
end









    











    Out[7]:

Numerička integracija

Fourierove koeficijente čemo izračunati numeričkom integracijom (vidi Numeričko integriranje) koristeći Julia naredbu quadgk.

Ovaj dio je izrađen prema bilježnici lecture-2.ipynb Stevena Johnsona izrađenoj za predmet 18.303.

Definirajmo funkcije sinecoef i coscoef koje numerički računaju koeficijente. Parametar abstol je toleranca numeričke integracije: želimo da je greška mala u odnosu na $\sqrt{\int_{x_0}^{x_1} |f(x)|^2 dx}$.



In [8]:

    
Tf=x1-x0
sinecoef(f, m, x0, x1, Tf) = 2 * quadgk(x -> f(x) * sin(2*m*π*x/Tf)/Tf, x0,x1, abstol=1e-8 * sqrt(quadgk(x->abs2(f(x)),x0,x1)[1]))[1]
coscoef(f, m, x0, x1, Tf) = 2 * quadgk(x -> f(x) * cos(2*m*π*x/Tf)/Tf, x0,x1, abstol=1e-8 * sqrt(quadgk(x->abs2(f(x)),x0,x1)[1]))[1]

# i druga funkcija koja računa na vektoru prirodnih brojeva
sinecoef(f, M::AbstractVector,x0,x1,Tf) = Float64[sinecoef(f,m,x0,x1,Tf) for m in M]
coscoef(f, M::AbstractVector,x0,x1,Tf) = Float64[coscoef(f,m,x0,x1,Tf) for m in M]









    Out[8]:





coscoef (generic function with 2 methods)



In [9]:

    
# First, define a function to evaluate N terms of the sine series, given the coefficients b
function fouriersum(a, b, x, T)
    f = a[1]/2
    for n = 1:length(b)
        f += a[n+1]* cos(2*n*π*x/T) + b[n] * sin(2*n*π*x/T)
    end
    return f
end
fouriersum(a, b, X::AbstractVector, T) = Float64[fouriersum(a, b,x,T) for x in X]









    Out[9]:





fouriersum (generic function with 2 methods)



In [11]:

    
# Uzmimo više točaka (inače "nestane" Gibbsov efekt!!)
X=linspace(x0,x1,1000)
@manipulate for n=1:1:99
    a = coscoef(f, 0:n, x0, x1, Tf)
    b = sinecoef(f, 1:n, x0, x1, Tf)
    plot(f,x0,x1,"b",X, fouriersum(a, b, X, Tf), "r-")
    title("Funkcija i "*string(n)*" članova reda")
end









    











    Out[11]:

Primjeri



In [13]:

    
f1(x)=x^2
x_0=-1
x_1=1









    Out[13]:





1



In [15]:

    
a(f1,n,x_0,x_1)









    Out[15]:





$$\begin{cases} \frac{2}{3} & \text{for}\: n = 0 \\\frac{4 \left(-1\right)^{n}}{\pi^{2} n^{2}} & \text{otherwise} \end{cases}$$



In [16]:

    
x0=N(x_0)
x1=N(x_1)
@manipulate for K=1:10
    S=a(f1,0,x_0,x_1)/2+sum([a(f,n,x_0,x_1)*cos(2*PI*n*x/(x_1-x_0))+b(f1,n,x_0,x_1)*sin(2*PI*n*x/(x_1-x_0)) for n=1:K])
    g(x)=S(x)
    plot(f1,x0,x1)
    oplot(g,x0,x1,"r")
    title("Funkcija i "*string(K)*" članova reda")
end









    











    Out[16]:



In [19]:

    
Tf=x1-x0
X=linspace(x0,x1,1000)
@manipulate for n=1:1:99
    a = coscoef(f1, 0:n, x0, x1, Tf)
    b = sinecoef(f1, 1:n, x0, x1, Tf)
    plot(f1,x0,x1,"b",X, fouriersum(a, b, X, Tf), "r-")
    title("Funkcija i "*string(n)*" članova reda")
end









    











    Out[19]:



In [20]:

    
# Funkcija zadana po djelovima
f2(x)=x+abs(x)
x_0=-1
x_1=1









    Out[20]:





1

Na ovu funkciju, zadanu po djelovima, ne mogu se direktno primijeniti simboličke formule, već bi ih trebalo "ručno" prilagoditi na odgovarajuči interval. Numeričko računanje radi. Funkcija se još može zadati i kao

f(x)=sympy.Piecewise((0,x<=0), (x,x>0))

ili

f(x)=x*sympy.Heaviside(x)



In [22]:

    
# a(f2,n,x_0,x_1)



In [25]:

    
x0=N(x_0)
x1=N(x_1)
Tf=x1-x0
X=linspace(x0,x1,1000)
@manipulate for n=1:2:99
    a = coscoef(f2, 0:n, x0, x1, Tf)
    b = sinecoef(f2, 1:n, x0, x1, Tf)
    plot(f2,x0,x1,"b",X, fouriersum(a, b, X, Tf), "r-")
    title("Funkcija i "*string(n)*" članova reda")
end









    











    Out[25]:



In [ ]: