>

Approximation by Trigonometric Polynomials.

Ingnore this stuff:

 > restart; with(plots): with(LinearAlgebra):

Warning, the name changecoords has been redefined

 > ip := proc(f,g) int(f*g, x=-Pi..Pi);

 > end;

 > nrm2 := proc(f) ip(f,f); end;

 > nrm := proc(f) sqrt(nrm2(f)); end;

 > append := proc(x, llist) [op(llist), x]; end;

We're working in the space of continuse functions on the interval .  The inner product is

 > ip(f(x),g(x));

n and m will denote postive integers

 > assume(n,posint);

 > assume(m, posint);

 >

Consider the functions 1, cos(nx), sin(nx).  Here's an illustration of what changing n does.

 > animate(plot,[sin(n*x), x=-Pi..Pi], n=1..10);

 >

Compute the inner product

 > ip(sin(n*x), sin(m*x));

L:et's see why this is true:

 > Int(sin(n*x)*sin(m*x), x) = int(sin(n*x)*sin(m*x), x);

 > rh := rhs(%);

 > Int(sin(n*x)*sin(m*x), x=-Pi..Pi) = subs(x=Pi, rh)-subs(x=-Pi, rh);

 > value(rhs(%));

 > ip(cos(n*x), cos(m*x));

 > ip(cos(n*x), sin(m*x));

 > ip(1, cos(n*x));

 > ip(1, sin(n*x));

Thus, these functions are othogonal.  To make them ON, we normalize them.

 > nrm(1);

 > u1 := 1/nrm(1);

 > nrm(cos(n*x));

 > nrm(sin(n*x));

 > basis := [u1];

Here is an ON set of functions:

 > for j from 1 to 5 do basis := append(cos(j*x)/sqrt(Pi), basis); basis := append(sin(j*x)/sqrt(Pi), basis); od:

 > basis;

 > for i from 1 to nops(basis) do

 > for j from i to nops(basis) do print(i,j, ip(basis[i], basis[j])); od;

 > od;

 > proj := proc(f, k)

 > end;

proj(f, k) computes the projection of f onto the space spanned by the first k basis elements.

Consider the following function:

 > myfun := expand(x*(x-2)*(x+1));

 > plot1 := plot(myfun, x=-Pi..Pi): plot1;

 > pr[1] := proj(myfun, 1);

 > pr[2] := proj(myfun, 2);

 > pr[3] := proj(myfun, 3);

 > for j from 4 to nops(basis) do pr[j]:= proj(myfun, j); od;

Plot these projections to see how well they approximate the function.

 > for j from 1 to nops(basis) do plt[j]:= plot(pr[j], x=-Pi..Pi, color=blue); od:

 > for j from 1 to nops(basis) do plt[j] := display([plot1,plt[j]]); od:

 > display([plot1, seq(plt[j], j=1..nops(basis))], insequence=true);

Try it with more basis elements!

 > basis := [u1];

 > for j from 1 to 50 do basis := append(cos(j*x)/sqrt(Pi), basis); basis := append(sin(j*x)/sqrt(Pi), basis); od:

 > pr[1] := ip(myfun, basis[1])*basis[1];

 > for j from 2 to nops(basis) do pr[j] := pr[j-1]+ip(myfun, basis[j])*basis[j]; od:

 > for j from 1 to nops(basis) do plt[j]:= plot(pr[j], x=-Pi..Pi, color=blue); od:

 > for j from 1 to nops(basis) do plt[j] := display([plot1,plt[j]]); od:

 > display([plot1, seq(plt[j], j=1..nops(basis))], insequence=true);

Evaluate the norm of the difference between f and the projections:

 > for j from 1 to nops(basis) do

 > nn2 := evalf(Int((myfun-pr[j])^2,x=-Pi..Pi)); nn := sqrt(nn2); norms[j] := nn; od:

 >

 > for j from 1 to nops(basis) do print(j, norms[j]); od;

 >