Take the 2-minute tour ×

I'd like to plot some partial sums for a Fourier Series problem, but I am not sure if the output I am getting is correct. I want to be able to plot the partial sums and the function on the same graph. Here is something that I attempted.

 s[n_, x_] :=  8/4 + 3/(9 \[Pi]) Sum[(6 (-1)^k)/(k \[Pi]) 
                Cos[(k \[Pi] x)/2] + (16 (-1)^k + 13)/(\[Pi] k) 
                Sin[(k \[Pi] x)/2],{k, 0, n}]

 partialsums = Table[s[n, x], {n, 1, 5}];

 f[x_] = Piecewise[{{-x^3-2x,-2 < x < 0},{-1+x,0<= x <= 2}}]

 Plot[Evaluate[partialsums], {x, -4 Pi, 4 Pi}]

Any ideas about the best method to tackle something like this?

share|improve this question
There's a k in the denominator of the sum (actually, two) and your index for k goes from 0 to n, so the sum is always indeterminate. –  David Skulsky Mar 2 '12 at 18:54
@DavidSkulsky: Should it be {k,1,n}. I actually only want to plot a certain amount of partial sums. Does this number go where the n is located? –  night owl Mar 2 '12 at 18:58
Once you fix the sum, you should be able to plot what you want by just doing Plot[Append[partialsums,f[x]],{x,-4 Pi, 4 Pi}]. –  David Skulsky Mar 2 '12 at 18:59
add comment

3 Answers

up vote 11 down vote accepted

In the definition of s you're summing from k==0. Since the summand has a term 1/k this gives a divide-by-zero error when calculating the partial sums. The sum should in fact start from k==1 (the zeroth coefficient is taken care of by the constant term in front of the sum). The first few approximations then look like

s[n_, x_] := 
 8/4 + 3/(9 \[Pi]) Sum[(6 (-1)^k)/(k \[Pi]) Cos[(k \[Pi] x)/
        2] + (16 (-1)^k + 13)/(\[Pi] k) Sin[(k \[Pi] x)/2], {k, 1, n}]

Plot[Evaluate[partialsums], {x, -4, 4}]

Mathematica graphics

To compare this with f we plot the s[5,x] and f in the same plot:

Plot[{s[5, x], f[x]}, {x, -2, 2}]

Mathematica graphics

which doesn't seem right to me, so I suspect you made a mistake somewhere in calculating the coefficients.

Calculating coefficients by hand

We could use FourierSeries to calculate the partial sums, but this is very slow, and doesn't give produce the general equation for the coefficients. Therefore it's better to calculate the coefficients by hand. The $n$-th coefficient can be calculated according to

coeff[0] = 1/4 Integrate[f[x], {x, -2, 2}];
coeff[n_] = 1/4 Integrate[f[x] Exp[I Pi n x/2], {x, -2, 2}]
(1/(2 n^4 \[Pi]^4))E^(-I n \[Pi]) (-48 + 48 E^(I n \[Pi]) - 
 48 I n \[Pi] + 28 n^2 \[Pi]^2 - 6 E^(I n \[Pi]) n^2 \[Pi]^2 + 
 2 E^(2 I n \[Pi]) n^2 \[Pi]^2 + 12 I n^3 \[Pi]^3 - 
 I E^(I n \[Pi]) n^3 \[Pi]^3 - I E^(2 I n \[Pi]) n^3 \[Pi]^3)

Then the partial sums are given by

series[m_, x_] := Sum[Exp[-I Pi n x/2] coeff[n], {n, -m, m}]

Plotting the first few approximations:

Plot[Evaluate[Table[series[j, x], {j, 0, 5}]], {x, -6, 6}]

Mathematica graphics

To see how this compares with the original function f:

Plot[Evaluate[{series[5, x], f[Mod[x, 4, -2]]}], {x, -4, 4}]

Mathematica graphics

which looks a lot better than the before.

Edit: Real coefficients

Here, coeff[n] are the coefficients for the Fourier series in exponential form, but these can be easily converted to the coefficients for the $\cos$ and $\sin$ series, a_n and b_n, by doing something like

a[0] = coeff[0];
a[n_] = Simplify[ComplexExpand[coeff[n] + coeff[-n]]];
b[n_] = Simplify[ComplexExpand[I (coeff[n] - coeff[-n])]];
share|improve this answer
Thanks, This is a very interesting way or technique of doing this. I also like that about how to covert back and forth between the exponential form and the trigonometric forms. That will be very useful. :). Very neat and organized approach. –  night owl Mar 4 '12 at 9:37
add comment

You can use FourierCoefficient to pre-calculate the Fourier coefficients to arbitrary degree and then use the result very effectively.

There may be some issues with zero-th degree, therefore I excluded this using Piecewise. Here's the main code block:

f[x_] := Piecewise[{
    {-x^3 - 2 x, -2 < x < 0},
    {-1 + x, 0 <= x <= 2}},
    Module[{x, fp},
    (* Set parameters so that the integration runs
       from -2 to 2 *)
    fp = {0, -Pi/2};
    fc = Piecewise[{
        {FourierCoefficient[f[x], x, 0, FourierParameters -> fp], # == 0},
        {ComplexExpand@FourierCoefficient[f[x], x, #, FourierParameters -> fp], True}
    }] &;
    fc = Evaluate /@ fc;

The output (fc) of this is something very unpleasant to look at; however, there's nothing Fourier-related left there, all the hard math has been done already, and only a bunch of elementary functions remain. fc is now a function of one argument that gives you the $n$-th Fourier coefficient of f in no time.

(* Calculate the first 2001 Fourier coefficients *)
AbsoluteTiming[Table[fc[n], {n, -1000, 1000}]] // First
0.777059 seconds

To convert this back to the function, you have to do the partial sum with your hands, for example here are the second, fourth and eightieth partial sums:

myPartialSums = Table[
    (* The Pi/2 compensates for the custom FourierParameters, see
       documentation of FourierSeries/FourierParameters
       under "more info" *)
    Re[Sum[fc[k] Exp[-Pi/2 I k t], {k, -n, n}]],
    {n, {2, 4, 80}}
A very large output has been generated,
but we're luckily not interested in it
anyway but would rather plot it
    {f[t]} ~Join~ myPartialSums,
    {t, -10, 10},
    PlotRange -> All, 
    Evaluated -> True, 
    PlotStyle -> {Thick, Automatic, Automatic, Automatic}

plot of Fourier series partial sums

share|improve this answer
add comment

If the summation index goes from 1 to n, as you suggested (?) above, and you make the mod to the Plot function I mentioned



f[x_] = Piecewise[{{-x^3 - 2 x, -2 < x < 0}, {-1 + x, 0 <= x <= 2}}]

Plot[Evaluate@Tooltip[Append[partialsums, f[x]]], {x, -4 Pi, 4 Pi}]

Then you get the following figure:

enter image description here

share|improve this answer
Your plot is running into a known issue with mma which I can't seem to find a reference to. Since Plot has the Attribute HoldAll, it is treating Append[partialsums, f[x]] as a single function, and hence colors all the lines the same color. To get around this either use Evaluate@Append[...] or Plot[#, ...]& @ Append[...]. –  rcollyer Mar 2 '12 at 19:21
Interesting-I was wondering about that! Thanks! –  David Skulsky Mar 2 '12 at 19:33
Thats nice David. I was wondering what was going on with the colors. I was fooling around with it to get something like this and then I seen that rcollyer figured it out. Thank You, that will be a useful option to remember when doing these. Thank You! :) –  night owl Mar 2 '12 at 20:18
@nightowl and David, I found the discussion on this. –  rcollyer Mar 5 '12 at 2:21
add comment

Your Answer


