Introduction
This is an introduction to an approximation technique called Fourier approximation. Given a particular function \(y = f(x)\), a Fourier approximation of \(f(x)\) is defined in the following way.
First, for each \(k = 0, 1, 2, \ldots\) we compute the following numbers (“coefficients”)
\(\displaystyle a_k = \frac1{\pi}\int_{-\pi}^{\pi} f(x)\cos(k x)\ dx,\) and \(\displaystyle b_k = \frac1{\pi}\int_{-\pi}^{\pi} f(x)\sin(k x)\ dx.\)
Notice that, with \(k = 0\), we get that \(\frac12 a_0\) is the average value of \(f(x)\) on the interval \([-\pi, \pi]\) since \(\cos(0)=1\). Also, \(b_0 = 0\) since \(\sin(0) = 0\). Now the \(n^{th}\) Fourier approximation is defined as
\(\displaystyle f_n(x) = \frac{a_0}{2} + \sum_{k=1}^n \left(a_k\cos(k x) + b_k\sin(k x)\right). \)
For example, the function \(f_2(x)\) will be \(\frac{a_0}{2} + a_1\cos(x) + b_1\sin(x) + a_2\cos(2x) + b_2\sin(2x)\).
Constructing the Fourier Approximation in SageMath
Let’s begin by constructing Fourier approximations for an example function, \(f(x) = x^2\). Once we’ve evaluated a line with the definition f(x) = x^2
, we would store the coefficient \(a_2\), for example, in SageMath by setting a variable equal to
1/pi*integral( f(x)*cos(2*x), x, -pi, pi )
\(\left(1\right)\)
As you will see below, we assign an entire list
of these coefficients to the variable name a
. The coefficient \(a_2\) will be in position 2 of the list, and we can get that coefficient from the list by typing a[2]
.
Avoiding a pitfall. When \(f(x)\) is a polynomial of \(x\) or \(e^x\) then computing the exact integrals for Fourier coefficients (using anti-derivatives) is relatively quick. For other functions, like the one in Exercise 4, it is more difficult. In order to sidestep the issue, in this lab we’ll use approximate integration to compute Fourier coefficients. We replace the command from \((1)\) with the command below (the [0]
makes the output be the approximate integral value only, ignoring the error bound).
numerical_integral(1/pi*f(x)*cos(2*x), -pi, pi)[0]
Exercise 1. For the function \(f(x) = x^2\), find the Fourier coefficients \(a_0, a_1\), and \(b_1\), displaying your output in your Sage worksheet.
Say we want to find the first 20 coefficients, \(a_k\) and \(b_k\) for \(k\) between \(0\) and \(20\). Then writing down one line at a time like this will be quite a pain! There is a nice succinct way to make a list of all the coefficients that you want. Say that you want to get these coefficients for \(k\) between \(0\) and \(6\). Then you can evaluate the following lines.
a=[ numerical_integral(1/pi*f(x)*cos(k*x), -pi, pi)[0] for k in range(7) ] b=[ numerical_integral(1/pi*f(x)*sin(k*x), -pi, pi)[0] for k in range(7) ]
(Make sure to scroll sideways with the line above so that you get the entire line.)
These lines make two lists of numbers, one called a
and one called b
. There are 7 numbers in each list; for example, list b
has the numbers \(b_0, b_1, b_2, b_3, b_4, b_5, b_6\), in that order. After SageMath computes these, you can look at \(b_3\) by typing b[3]
. (Note where the 7 appears in range(7)
. This means that values of k
go through the list [0, 1, 2, 3, 4, 5, 6]
.)
Once we have these numbers, if we want to define the \(5^{th}\) Fourier approximation, we can do so with
f5(x) = a[0]/2 + sum( [a[k]*cos(k*x) + b[k]*sin(k*x) for k in range(1,6)] )
What is being done in this SageMath definition of f5(x)
? Start by looking at it, and comparing to the definition of \(f_n(x)\) in the introduction.
There is a list of functions that equal coefficient a[k]
times cos(k*x)
plus b[k]*sin(k*x)
, one for each k
; the range
command in that list forces k
to start at 1 and go up until 5 (instead of starting at 0). Then I put the list inside of the command sum( )
to add the functions together.1 This is added to a[0]/2
.
After computing the Fourier approximation f5(x)
for the function \(f(x) = x^2\), we can plot it together with f(x)
over the interval \([-\pi, \pi]\). This lets us see how close its values are to those of \(f(x)\). Displaying a plot of the functions together, \(f_5(x)\) in red and the original \(f(x)\) in blue, can be done as follows.
plot( [f5(x), f(x)], -pi, pi, color=['red', 'blue'] )
Exercises 2 – 4. For each of the following functions, compute the \(20^{th}\) Fourier approximation \(f_{20}(x)\). Then plot the functions \(f(x)\) and \(f_{20}(x)\) in the same plot over the interval \([-\pi, \pi]\).
-
-
- \(f(x) = x^3.\)
- \(f(x) = e^x.\)
- \(f(x) = \dfrac{x^2-x}{1+x^4}.\)
-
Interval of Approximation
When plotting the functions above, you may have noticed that, over \(x\) values near the endpoints of \([-\pi, \pi]\), the \(y\)-values of the approximation were not very close to the \(y\)-values of \(f(x)\). It gets worse if you compare the two functions outside of the interval \([-\pi, \pi]\) (the interval used for integration). Can this be fixed?
Consider a function \(f(x)\) defined on an interval \([-c, c]\). Recall that multiplying the input by a constant affects the graph by stretching it horizontally. In other words, if \(x_1 = 2x_0\) for a number \(x_0\) in the interval \([-c, c]\), then the value \(y_1 = f(x_1)\) is the same as \(f(2x_0)\).
If we can approximate the function \(f(2x)\) well on the interval \([-\pi, \pi]\), then we should be able to use that to get a good approximation to \(f(x)\) on the interval \([-2\pi, 2\pi]\).
Setup for Exercise 5. Let \(g(x) = e^x\) and let \(f(x) = g(2x)\). The function \(f_{50}(x)\), a Fourier approximation of \(f(x)\), approximates \(g(2x)\) on the interval \([-\pi, \pi]\). So, on the interval \([-2\pi, 2\pi]\), the function \(f_{50}(x/2)\) should be a good approximation to \(g(x)\).
Exercise 5. Compute the Fourier approximation \(f_{50}(x)\) of the function \(g(2x)\), as explained above. Then, display a plot of \(e^x\) along with the good approximation to it, on the interval \([-2\pi, 2\pi]\).
- You can always add up a list of things this way in SageMath (or in Python). ↩