Fourier Series Representation and Animation

(Uses Maple V): plots package, Heaviside, sin, cos, sum, evalf, array, int, animate, proc, lhs, rhs
REF: (1) Monagan, Geddes, Labahn & Vorkoetter, Maple Programming Guide, Springer Verlag, 1996. chap.8, (2) O. Mitchell animation sample in Maple V release 3.

Fourier Series:

Any periodic function, f(t), of period T can be represented by an infinite trigonometric series,

f(t) =

where:

= , is the fundamental frequency in radians/sec, while f is in Hertz.
=nth harmonic frequency,
,,= the formula to calculate these coefficients are presented in the next section.

Practical Fourier Series Representation:

Although the infinite series given above is 'theoretically' valid but it can not be used in real life. For one, it takes an infinite number of harmonic generators to composed a linear fourier series representation of a given signal. What we need is an approximate representation consisting of a finite number of harmonic terms. Like this,

f(t) =

where:
m = mth harmonic required to meet a certain approximation criteria.


and the period T=B-A,

Fourier Series Representation of 'any' function & their plots:
-- and comment on Maple V "inert functions".

Monagan and companies has presented a very sophisticated fourier series code. They use sequence function to generate sequence of frames and superimposed the plot of the original function. User can then observed the progression of the fourier series approximation. In order to accomplished their task they have to use a so-called Maple V inert function. In particular the integration function, Int, employed in their code has a first letter in "capital" I rather than a "lower case", i. Maple treat such function with first letter in capital as an "inert function". Invoking an inert function essentially telling Maple not to evaluate an expression at all, or at least delay it until you specifically say so. A typical usage of inert function is to prevent symbolic evaluation and allow later numerical procedures to be invoked. Monagan's code is so complex that the integration need to be postpone. Here is a sample usage of inert function, Int, and a fully evaluated, int. The left-hand side uses inert integration, Int, while the right-hand side uses a fully evaluated integration, int.

> Int(1/(1+t^3),t) = int(1/(1+t^3),t);
Result,

Notice the left-hand side remains "unevaluated" and the right-hand side is fully integrated.

Given below is a much reduced version of the code. It can take in 'any' function as a function of a single argument of 'any' name. User can specify the range of the argument and indicate the desired 'harmonics'. This code is simple enough that inert integration is not needed here, so we do not use them. Four samples are shown.

> restart;
> fourierREPT:=
>          #func ----------------  any given function,     eg. t^2-6*t+2
>          #trange:name=range ---  range of indept var,    eg. t=0..10
>          #m:posint ------------  positive integer (>0),  eg. 6
> proc( func, trange:name=range, m:posint )
>    local t,A,B,T,n,j,p,q,partsum;
>    A:=lhs(rhs(trange));
>    B:=rhs(rhs(trange));
>    T:=B-A;
>    t:=2*Pi*lhs(trange)/T;
>    partsum:=1/T*evalf(int( func,trange));
>    for n from 1 to m do
>       partsum:=partsum
>          +2/T*evalf(int(func*sin(n*t),trange))*sin(n*t)
>          +2/T*evalf(int(func*cos(n*t),trange))*cos(n*t);
>    od;
> end:
> #############################
> #  4 samples,
> #
> with(plots):
> alias(u=Heaviside):
> f1:=(-t+7)*(u(t)-u(t-7))+7*u(t-7);
> g1:=fourierREPT(f1,t=0..10,8);
> f2:=20*exp(-0.2*t)*sin(2*Pi*t)*(u(t)-u(t-7))+10*u(t-7);
> g2:=fourierREPT(f2,t=0..10,14);
> f3:=20-10*exp(-2*x);
> g3:=fourierREPT(f3,x=0..4,6);
> f4:=t^2-6*t+2;
> g4:=fourierREPT(f4,t=0..10,6);
> plot({f1,g1},t=0..10,title='f1_and_g1');
> plot({f2,g2},t=0..10,title='f2_and_g2');
> plot({f3,g3},x=0..4,title='f3_and_g3');
> plot({f4,g4},t=0..10,title='f4_and_g4');

The functions, their fourier series representation and their plots.


Fourier Series Animation:

This sample written by Ocie Mitchell in 1993, illustrate Fourier Series animation convergence to the given function. The function f:=x->1-x^2 on the right side of the arrow maybe replaced by any function define in the interval (0,1). N the number of fourier terms can also be change by the user.

> with(plots):
> af:=(t,k)->Heaviside(t-k)*(t-k)-Heaviside(t-1-k)*(t-1-k);
> N:=20;
> f:=x->1-x^2;
> bf:=array(1..N);
> for k to N do bf[k]:=evalf(2*int(sin(k*Pi*x)*f(x),x=0..1)) od:
> ts:=(x,t)->sum(af(t,n)*bf[n]*sin(n*Pi*x),n=1..N);
> animate(ts(x,t),x=0..1,t=0..N,frames=3*N,numpoints=5*N);
> animate({f(x),ts(x,t)},x=0..1,t=0..N,frames=3*N,numpoints=5*N);
> 

When executed it gave the following result,