2.10 Some Fourier Series

```#
#
#
#                                           WORKSHEET # 21
#
#                                           Some Fourier Series
#
#
> restart;
--------------------------------------------------------------------------------
# For a function f(x) defined on the interval [0,2*Pi] the Fourier series
# is a linear combination of the periodic functions cos(n*x), sin(n*x)
# and a constant term.
--------------------------------------------------------------------------------
> f(x)=(1/2)*a[0]+Sum(a[n]*cos(n*x)+b[n]*sin(n*x),n=1..infinity);

/infinity                                \
| -----                                  |
|  \                                     |
f(x) = 1/2 a[0] + |   )     (a[n] cos(n x) + b[n] sin(n x))|
|  /                                     |
| -----                                  |
\ n = 1                                  /
--------------------------------------------------------------------------------
# The coefficients in this formula are given by the following formulas.
# The formula for a[n] is valid for n=0,1,... and that for b[n] is valid for
# n=1,2,...
--------------------------------------------------------------------------------
> a[n]=(1/Pi)*Int(f(x)*cos(n*x),x=0..2*Pi);

2 Pi
/
|
|    f(x) cos(n x) dx
|
/
0
a[n] = ----------------------
Pi
--------------------------------------------------------------------------------
> b[n]=(1/Pi)*Int(f(x)*sin(n*x),x=0..2*Pi);

2 Pi
/
|
|    f(x) sin(n x) dx
|
/
0
b[n] = ----------------------
Pi
--------------------------------------------------------------------------------
# If we assume that f(x) is actually represented by the infinite series
# above and we can integrate term-by-term in this series then the
# formulas for the coefficients follow from the following orthogonal
# properties of sin(x) and cos(x):
#
--------------------------------------------------------------------------------
> weknow:=cos(2*m*Pi)=1,sin(2*m*Pi)=0,cos(2*n*Pi)=1,sin(2*n*Pi)=0;\

weknow := cos(2 m Pi) = 1, sin(2 m Pi) = 0, cos(2 n Pi) = 1, sin(2 n Pi) = 0
--------------------------------------------------------------------------------
> subs(weknow,Int(sin(m*x)*sin(m*x),x=0..2*Pi)=int(sin(m*x)*sin(m*x),x=0..2*Pi));

2 Pi
/
|            2
|    sin(m x)  dx = Pi
|
/
0
--------------------------------------------------------------------------------
> subs(weknow,Int(cos(m*x)*cos(m*x),x=0..2*Pi)=int(cos(m*x)*cos(m*x),x=0..2*Pi));

2 Pi
/
|            2
|    cos(m x)  dx = Pi
|
/
0
--------------------------------------------------------------------------------
> Int(sin(m*x)*cos(n*x),x=0..2*Pi)=0;

2 Pi
/
|
|    sin(m x) cos(n x) dx = 0
|
/
0
--------------------------------------------------------------------------------
# All other integrals involving sin(m*x) and cos(n*x) are zero.  That is
# if m does not equal n then:
--------------------------------------------------------------------------------
> Int(sin(m*x)*sin(n*x),x=0..2*Pi)=0,Int(cos(m*x)*cos(n*x),x=0..2*Pi)=0;

2 Pi                            2 Pi
/                               /
|                               |
|    sin(m x) sin(n x) dx = 0,  |    cos(m x) cos(n x) dx = 0
|                               |
/                               /
0                               0
--------------------------------------------------------------------------------
# Now we will compute the Fourier series for some simple functions.
--------------------------------------------------------------------------------
> f:=x*(2*Pi-x);

f := x (2 Pi - x)
--------------------------------------------------------------------------------
> F:=plot(f,x=0..2*Pi):
--------------------------------------------------------------------------------
> with(plots):
--------------------------------------------------------------------------------
> display(F,title=`f(x)=x*(2*Pi-x)`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> B:=int(f*sin(n*x),x=0..2*Pi);

n sin(2 n Pi) Pi + cos(2 n Pi)     2
B := - 2 ------------------------------ + ----
3                   3
n                   n
--------------------------------------------------------------------------------
> subs(weknow,B);

0
--------------------------------------------------------------------------------
# Thus all the sine coefficients are zero.
--------------------------------------------------------------------------------
> A:=int(f*cos(n*x),x=0..2*Pi);

- n cos(2 n Pi) Pi + sin(2 n Pi)      Pi
A := 2 -------------------------------- - 2 ----
3                      2
n                      n
--------------------------------------------------------------------------------
> A:=subs(weknow,A);

Pi
A := - 4 ----
2
n
--------------------------------------------------------------------------------
> a0:=int(f,x=0..2*Pi)/(2*Pi);

2
a0 := 2/3 Pi
--------------------------------------------------------------------------------
> a:=unapply(A/Pi,n);

4
a := n -> - ----
2
n
--------------------------------------------------------------------------------
> a(5);

-4/25
--------------------------------------------------------------------------------
> S:=n->a0+sum(a(k)*cos(k*x),k=1..n);

/  n                \
|-----              |
| \                 |
S := n -> a0 + |  )   a(k) cos(k x)|
| /                 |
|-----              |
\k = 1              /
--------------------------------------------------------------------------------
> S(1);

2
2/3 Pi  - 4 cos(x)
--------------------------------------------------------------------------------
> S(2);

2
2/3 Pi  - 4 cos(x) - cos(2 x)
--------------------------------------------------------------------------------
> S(3);

2
2/3 Pi  - 4 cos(x) - cos(2 x) - 4/9 cos(3 x)
--------------------------------------------------------------------------------
> sun.rise;

sunrise
--------------------------------------------------------------------------------
> subscript.10;

subscript10
--------------------------------------------------------------------------------
> for j from 1 to 5 do S(j) od;j:='j':

2
2/3 Pi  - 4 cos(x)

2
2/3 Pi  - 4 cos(x) - cos(2 x)

2
2/3 Pi  - 4 cos(x) - cos(2 x) - 4/9 cos(3 x)

2
2/3 Pi  - 4 cos(x) - cos(2 x) - 4/9 cos(3 x) - 1/4 cos(4 x)

2
2/3 Pi  - 4 cos(x) - cos(2 x) - 4/9 cos(3 x) - 1/4 cos(4 x) - 4/25 cos(5 x)
--------------------------------------------------------------------------------
> for j from 1 to 10 do P.j:=plot(S(j),x=0..2*Pi) od:
--------------------------------------------------------------------------------
> F:=plot(f,x=0..2*Pi,style=POINT):
--------------------------------------------------------------------------------
> display(F,title=`f(x)=x*(2Pi-x)`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> display([F,P1],title=`first partial sum`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> display([F,P2],title=`second partial sum`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> display([F,P3],title=`third partial sum`);\
display([F,P4],title=`fourth partial sum`);\
display([F,P3],title=`fifth partial sum`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> for j from 1 to 20 do print(j,evalf(subs(x=0,S(j)-f))) od;j:='j':

1, 2.579736270

2, 1.579736270

3, 1.135291826

4, .885291826

5, .725291826

6, .614180714

7, .532548061

8, .470048061

9, .420665345

10, .380665345

11, .347607494

12, .319829716

13, .296161077

14, .275752914

15, .257975136

16, .242350136

17, .228509306

18, .216163627

19, .205083294

20, .195083294
--------------------------------------------------------------------------------
> for j from 1 to 20 do print(j,evalf(subs(x=Pi,S(j)-f))) od;j:='j':

1, .710131866

2, -.289868134

3, .1545763104

4, -.0954236896

5, .0645763104

6, -.0465348007

7, .03509785236

8, -.02740214764

9, .02198056841

10, -.01801943159

11, .01503841965

12, -.01273935813

13, .01092928092

14, -.00947888235

15, .00829889543

16, -.00732610457

17, .00651472588

18, -.00583095313

19, .00524937928

20, -.00475062072
--------------------------------------------------------------------------------
> fourier:=a0+Sum(a(k)*cos(k*x),k=1..infinity);

/infinity             \
| -----               |
2   |  \          cos(k x)|
fourier := 2/3 Pi  + |   )     - 4 --------|
|  /              2   |
| -----          k    |
\ k = 1               /
--------------------------------------------------------------------------------
> subs(x=0,fourier);

/infinity           \
| -----             |
2   |  \          cos(0)|
2/3 Pi  + |   )     - 4 ------|
|  /             2  |
| -----         k   |
\ k = 1             /
--------------------------------------------------------------------------------
> eval(");

/infinity       \
| -----         |
2   |  \          4 |
2/3 Pi  + |   )     - ----|
|  /          2 |
| -----      k  |
\ k = 1         /
--------------------------------------------------------------------------------
> evalf(");\

-8
.3*10
--------------------------------------------------------------------------------
> Sum(1/k^2,k=1..infinity)=sum(1/k^2,k=1..infinity);

infinity
-----
\        1          2
)     ---- = 1/6 Pi
/        2
-----    k
k = 1
--------------------------------------------------------------------------------
> mod2Pi:=x->x-2*Pi*floor(x/(2*Pi));

x
mod2Pi := x -> x - 2 Pi floor(1/2 ----)
Pi
--------------------------------------------------------------------------------
> fperiodic:=x->mod2Pi(x)*(2*Pi-mod2Pi(x));

fperiodic := x -> mod2Pi(x) (2 Pi - mod2Pi(x))
--------------------------------------------------------------------------------
> Fper:=plot(fperiodic(x),x=-4*Pi..4*Pi,style=POINT,numpoints=200):
--------------------------------------------------------------------------------
> display(Fper,title=`periodic extension of f(x)`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
>  for j from 1 to 5 do Pper.j:=plot(S(j),x=-4*Pi..4*Pi,numpoints=200) od:
--------------------------------------------------------------------------------
> display([Fper,Pper1],title=`first partial sum`);display([Fper,Pper2],title=`second partial sum`);\
display([Fper,Pper3],title=`third partial sum`);
--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
#

** Maple V Graphics **

--------------------------------------------------------------------------------
> fperiodic(Pi);

2
Pi
--------------------------------------------------------------------------------
> approx:=a0+sum(a(k)*cos(k*Pi),k=1..infinity);

/infinity              \
| -----                |
2   |  \          cos(k Pi)|
approx := 2/3 Pi  + |   )     - 4 ---------|
|  /               2   |
| -----           k    |
\ k = 1                /
--------------------------------------------------------------------------------
> wealsoknow:=cos(k*Pi)=(-1)^k;

k
wealsoknow := cos(k Pi) = (-1)
--------------------------------------------------------------------------------
> subs(wealsoknow,approx);

/infinity          \
| -----           k|
2   |  \          (-1) |
2/3 Pi  + |   )     - 4 -----|
|  /             2 |
| -----         k  |
\ k = 1            /
--------------------------------------------------------------------------------
> eval(");

2
2/3 Pi  + 4 hypergeom([1, 1, 1], [2, 2], -1)
--------------------------------------------------------------------------------
> ?hypergeom
--------------------------------------------------------------------------------
> evalf("");
Error, (in evalf/Hypergeom) Too many iterations without convergence

--------------------------------------------------------------------------------
> odd:=sum(1/(2*k-1)^2,k=1..  infinity);

2
odd := 1/8 Pi
--------------------------------------------------------------------------------
> even:=sum(1/(2*k)^2,k=1..  infinity);

2
even := 1/24 Pi
--------------------------------------------------------------------------------
> odd-even;

2
1/12 Pi
--------------------------------------------------------------------------------
> a0+4*";

2
Pi
--------------------------------------------------------------------------------

```