4.2.2.1 Rationale

Before seeing some examples of how a Taylor expansion actually works, we want to understand a little more of the process itself. If we look at Eq. 4.1, we see that the first term in the series is the value of the function at the expansion point a. A Taylor series is nothing other than a function that allows us to understand how much a function changes as we move away from a. The remaining terms of the expansion describe this: the relative change of the function. This first term is important because it defines the offset of the function.

The second expansion term uses the first derivative. It multiplies the value of this derivative with the distance by which x is away from a. The fraction term gives a certain weight to the contribution of the first derivative. Looking at this term, it is easy to understand why it is important to keep x − a very small: the derivative is a linearization of the slope of the function at the point a. The change as we move away from a must be proportional to this slope, which is the origin of this term.

But what do the other terms do? Suppose the function is a rather complicated function for which the derivative is not a constant. In this case, we would make significant mistakes even if x is very close to a. So the obvious solution to this problem is to correct the nonlinearity of the first derivative by finding the derivative to it. This corrects for nonlinearity if this second derivative is a constant. If not, we repeat this process over and over until the corrections finally become unnoticeable. As stated, the fraction terms give a certain weight to these higher derivatives because, obviously the n^th correction is less important than the first correction.

As stated, the more we move away from the point a, the worse the approximation of our function (usually) become. However, in many cases, expanding a Taylor series around a representative point may even be a very good approximation for the whole function.

4.2.2.2 Taylor Series Expansion of the Exponential Function

As a first example, we will derive the Taylor series expansion for the exponential function f (x) = e^x. We will look at the expanded series at the expansion orders n = 1, n = 2, n = 3, and n = 5. We first expand the series around a = 0. According to Eq. 4.1 we find

$\begin{array}{ll}n=1\hfill & f\left(x\right)=1\hfill \\ n=2\hfill & f\left(x\right)=1+x\hfill \\ n=3\hfill & f\left(x\right)=1+x+\frac{1}{2}{x}^2\hfill \\ n=5\hfill & f\left(x\right)=1+x+\frac{1}{2}{x}^2+\frac{1}{6}{x}^3+\frac{1}{24}{x}^4\hfill \end{array}$

The functions are plotted in Fig. 4.2. As you can see, for n = 5, we find a very nice approximation of the exponential function. As long as we evaluate only the function around the expansion point a = 0, we are not making significant mistakes. However, if we wanted to use our approximations in the range of x = −3 (see Fig. 4.3), you can see that we would be making huge mistakes. The reason for this is that expanding the Taylor series around a = 3 leads to significantly different terms. According to Eq. 4.1, we find

$\begin{array}{ll}n=1\hfill & f\left(x\right)={\text{e}}^{-3}\hfill \\ n=2\hfill & f\left(x\right)={\text{e}}^{-3}+{\text{e}}^{-3}\left(x+3\right)\hfill \\ n=3\hfill & f\left(x\right)={\text{e}}^{-3}+{\text{e}}^{-3}\left(x+3\right)+\frac{1}{2}{\text{e}}^{-3}{\left(x+3\right)}^2\hfill \\ n=5\hfill & f\left(x\right)={\text{e}}^{-3}+{\text{e}}^{-3}\left(x+3\right)+\frac{1}{2}{\text{e}}^{-3}{\left(x+3\right)}^2+\frac{1}{6}{\text{e}}^{-3}{\left(x+3\right)}^3+\frac{1}{24}{\text{e}}^{-3}{\left(x+3\right)}^4\hfill \end{array}$

As you can see, these terms are different. When using these equations, we find very good approximations again (see Fig. 4.4).

4.2.2.3 Taylor Series Expansion of the Sine Function

As a final example, we will approximate the sine function using a Taylor series around a = 0. As we will see in a moment, we need higher expansion orders in order to approximate the sine function sufficiently. Here we will use n = 1, n = 3, n = 6, and n = 9, for which we find from Eq. 4.1

$\begin{array}{ll}n=1\hfill & f\left(x\right)=0\hfill \\ n=3\hfill & f\left(x\right)=0+x\hfill \\ n=6\hfill & f\left(x\right)=0+x-\frac{1}{6}{x}^3+\frac{1}{120}{x}^5\hfill \\ n=9\hfill & f\left(x\right)=0+x-\frac{1}{6}{x}^3+\frac{1}{120}{x}^5-\frac{1}{5040}{x}^7\hfill \end{array}$

Fig. 4.5 shows the results of the Taylor series. As we can see, with n = 9, we achieve a decent approximation around a = 0.

4.2.4.1 Infinitesimal Positive Change

Suppose we have expanded the Taylor series around a point a and now want to know how much the function changes as we move from a to a + ∆x. In this case, Eq. 4.1 becomes

$\begin{array}{ll}f\left(a+\Delta x\right)\hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(a+\Delta x-a\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(a+\Delta x-a\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(a+\Delta x-a\right)}^3+\mathrm{...}\hfill \\ \hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(\Delta x\right)}^3+\cdots \hfill \\ \hfill & ={\displaystyle \sum _{n=0}^{\infty }\frac{1}{n!}\frac{d^{n}f\left(a\right)}{d{x}^n}{\left(\Delta x\right)}^n}\hfill \\ \hfill & ={\displaystyle \sum _{n=0}^{n_{\max }}\frac{1}{n!}\frac{d^{n}f\left(a\right)}{d{x}^n}{\left(\Delta x\right)}^n+O_{n_{\max }+1}}\hfill \end{array}$

Note that nothing spectacular has happened here. We now simply have an idea how the function changes as we move the infinitesimal distance ∆x in the positive direction.

4.2.4.2 Infinitesimal Negative Change

In analogy, we now determine the Taylor series as we move from a by the infinitesimally small distance ∆x into the negative direction, i.e., to a − ∆x. In this case Eq. 4.1 becomes

$\begin{array}{ll}f\left(a-\Delta x\right)\hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(a-\Delta x-a\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(a-\Delta x-a\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(a-\Delta x-a\right)}^3+\mathrm{...}\hfill \\ \hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(-\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(-\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(-\Delta x\right)}^3+\cdots \hfill \\ \hfill & ={\displaystyle \sum _{n=0}^{\infty }{\left(-1\right)}^n\frac{1}{n!}\frac{d^{n}f\left(a\right)}{d{x}^n}{\left(\Delta x\right)}^n+O_{n+1}}\hfill \\ \hfill & ={\displaystyle \sum _{n=0}^{n_{\max }}{\left(-1\right)}^n\frac{1}{n!}\frac{d^{n}f\left(a\right)}{d{x}^n}{\left(\Delta x\right)}^n+O_{n_{\max }+1}}\hfill \end{array}$

As you can see, Eq. 4.3 and Eq. 4.4 are almost identical. However, Eq. 4.4 shows alternating negative terms. This is the basic concept that we will use in the next step.

4.2.4.3 Deriving the First Derivative: Central Finite Difference

Now we will simply subtract Eq. 4.3 and Eq. 4.4. In this case, we find

$\begin{array}{ll}f\left(a+\Delta x\right)-f\left(a-\Delta x\right)\hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(\Delta x\right)}^3+\mathrm{...}\hfill \\ \hfill & -\left(f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(-\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(-\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(-\Delta x\right)}^3+\cdots \right)\hfill \\ \hfill & =2\left(\frac{2}{1!}\frac{df}{dx}{\left(\Delta x\right)}^1+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(\Delta x\right)}^3+\cdots \right)\hfill \\ \hfill & =2\left({\displaystyle \sum _{n=0}^{\infty }\frac{1}{\left(2n+1\right)!}\frac{d^{2n+1}f\left(a\right)}{d{x}^{2n+1}}{\left(\Delta x\right)}^{2n+1}+O_{2n+1}}\right)\hfill \end{array}$

In many applications, especially in numerical simulation, Eq. 4.5 is not expanded farther than n = 1, in which case we find

$\begin{array}{ll}f\left(a+\Delta x\right)-f\left(a-\Delta x\right)\hfill & =2\frac{df}{dx}\Delta x\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(a+\Delta x\right)-f\left(a-\Delta x\right)}{2\Delta x}\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(x+\Delta x\right)-f\left(x-\Delta x\right)}{2\Delta x}+O{}_3\hfill \end{array}$

where we have replaced the expansion point a with the variable x. Eq. 4.6 is referred to as the central finite distance approximation. It allows deriving the first derivative of a function if we know the values of the function at the distance x + ∆x and x − ∆x, which is why it is referred to as the central finite distance approximation.

As the formula stands, it does nothing but calculate the slope.

4.2.4.4 Deriving the First Derivative: Forward and Backward Finite Difference

We can derive a second method for putting into relation the first derivative and the values of the function at two given points. For this we consider Eq. 4.3 and expand the series only to n = 1, in which case we find

$\begin{array}{ll}f\left(a+\Delta x\right)\hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}\Delta x\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(a+\Delta x\right)-f\left(a\right)}{\Delta x}\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(x+\Delta x\right)-f\left(x\right)}{\Delta x}+O{}_2\hfill \end{array}$

where we have again replaced x for a as the independent variable. Eq. 4.7 is referred to as the forward finite distance approximation. It allows deriving the first derivative of a function if we know the values of the function at x and x + ∆x (and v.v.).

We can perform the same operation starting from Eq. 4.4, in which case we find (again expanding only to n = 1)

$\begin{array}{ll}f\left(a-\Delta x\right)\hfill & =f\left(a\right)-\frac{1}{1!}\frac{df}{dx}\Delta x\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(a\right)-f\left(a-\Delta x\right)}{\Delta x}\hfill \\ \frac{df}{dx}\hfill & =\frac{f\left(x\right)-f\left(x-\Delta x\right)}{\Delta x}+O{}_2\hfill \end{array}$

where we have again replaced x for a as the independent variable. Eq. 4.8 is referred to as the backward finite distance approximation. It allows deriving the first derivative of a function if we know the values of the function at x and x − ∆x (and v.v.).

4.2.4.5 Deriving the Second Derivative

In the next step we will add Eq. 4.3 and Eq. 4.4. Subtraction resulted in the central finite distance approximation of the first derivative. Interestingly, something quite different happens in the case of the addition. Here we find

$\begin{array}{ll}f\left(a+\Delta x\right)-f\left(a-\Delta x\right)\hfill & =f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(\Delta x\right)}^3+\mathrm{...}\hfill \\ \hfill & +f\left(a\right)+\frac{1}{1!}\frac{df}{dx}{\left(-\Delta x\right)}^1+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(-\Delta x\right)}^2+\frac{1}{3!}\frac{d^{3}f}{d{x}^3}{\left(-\Delta x\right)}^3+\cdots \hfill \\ \hfill & =2\left(f\left(a\right)+\frac{1}{2!}\frac{d^{2}f}{d{x}^2}{\left(\Delta x\right)}^2+\frac{1}{4!}\frac{d^{4}f\left(a\right)}{d{x}^4}{\left(\Delta x\right)}^4+\cdots \right)\hfill \\ \hfill & =2\left(f\left(a\right)+{\displaystyle \sum _{n=0}^{\infty }\frac{1}{\left(2n\right)!}\frac{d^{2n}f\left(a\right)}{d{x}^{2n}}{\left(\Delta x\right)}^{2n}+O_{2n}}\right)\hfill \end{array}$

Again, we will expand the series only to the first term (n = 1, in this case as the series does not start at 0), in which case we find

$\begin{array}{ll}f\left(a+\Delta x\right)+f\left(a-\Delta x\right)\hfill & =2f\left(a\right)+\frac{d^{2}f}{d{x}^2}{\left(\Delta x\right)}^2\hfill \\ \frac{d^{2}f}{d{x}^2}\hfill & =\frac{f\left(a+\Delta x\right)+f\left(a-\Delta x\right)-2f\left(a\right)}{{\left(\Delta x\right)}^2}\hfill \\ \frac{d^{2}f}{d{x}^2}\hfill & ={\frac{f\left(x+\Delta x\right)+f\left(x-\Delta x\right)-2f\left(x\right)}{{\left(\Delta x\right)}^2}}_3\hfill \end{array}$

where we have again replaced a by x as the independent variable. In this case we have found a correlation between the second derivative and the values of the function at x, x + ∆x, and x − ∆x. If we know these values, we can calculate the value of the second derivative. Equally, if we know the value of the second derivative at x, we can calculate the values of the function at x, x + ∆x, and x − ∆x.

4.2.4.6 Implications

We have now derived four very important equations that allow us to derive the first derivative (see Eq. 4.6, Eq. 4.7, and Eq. 4.8) and the second derivative (see Eq. 4.10) of a function at the location x if we know the values of the function at x + ∆x and x − ∆x. More importantly, however, these functions also enable us to do the reverse process: if we know the values of the derivatives at the location x, we are able to derive the change in the actual function between x + ∆x and x − ∆x. If we know the value of the function at a boundary (usually because a boundary condition is given), we can calculate the changes from this absolute value if we know one of the derivatives.

As we will see, in many problems in practical physics and especially fluid mechanics, we often know something about a function’s second derivative because the most important PDE in fluid mechanics, the Navier-Stokes equation, is a second-order PDE (see Eq. 11.40) of the velocity. In many cases, we want to derive the velocity, which is why finding equations that relate the second derivative of a function and the function itself is very important.

It is important to keep one thing in mind here. As we have seen, a Taylor series expansion only does a good job approximating a function around the value it was expanded around. Therefore, we will need to make sure that the values for ∆x we are working with are very small. Obviously, the higher the local changes in the functions become, the less precise a Taylor series approximation that takes into account only the linear terms (as do Eq. 4.6, Eq. 4.7, Eq. 4.8, and Eq. 4.10) becomes.

4.3.2.1 Rationale

Again, before seeing some examples of the Fourier series, we want to understand how the Fourier series expansion actually works. First, the Fourier series uses three types of coefficients: a₀, a_n, and b_n. The coefficient a₀ is obtained if n = 0 is applied to the term a_n cos (n x) + b_n sin (n x). In this case the sine function is 0 and the cosine function 1, in which case the term a₀ is obtained from Eq. 4.20. This is why this constant is given with a₀ and not b₀. If we look at Eq. 4.17, we see that the function is thus approximated with sine and cosine functions of increasing frequency. For n = 1, we obtain the functions cos x and sin x, for n = 2, we obtain cos ∗ (2 x) and sin ( 2 x), and so on.

The cosine and the sine functions are “weighted” with the constants a_n and b_n, respectively. They determine how well a given sine or cosine function at the given frequency will approximate the function. If you look carefully at Eq. 4.20 and Eq. 4.21, you will see why. These functions effectively multiply f (x), the function we want to approximate, with the sine or cosine function of the given frequency. Now, if this specific sine or cosine approximates the function well, they will overlap, which is why the total surface area of the product of these two functions should be a high value that is either positive or negative at most points of the domain (see Fig. 4.6a). If we take the integral of such a function, we expect a large value. As this integral gives rise to a_n and b_n, we expect these weighting factors to be high (see Eq. 4.20 and Eq. 4.21). On the other hand, if the sine or cosine function does a poor job of approximating the function, the product of the two will either be essentially 0 or feature positive and negative regions. If we take the integral of such a function, we expect a significantly smaller integral and thus a smaller weighting factor (see Fig. 4.6b).

4.3.2.2 Fourier Series Expansion of the Square Wave Function

We have already briefly introduced the square wave function as a good example of how the Fourier series expansion works (see Fig. 4.6). In this section, we perform the Fourier series expansion of the square wave function, which is defined as

$f_{\text{square}}\left(x\right)=\{\begin{array}{ll}1\hfill & x<\pi \hfill \\ -1\hfill & \text{otherwise}\hfill \end{array}$

Now let us exemplary expand the Fourier series for this function up to n = 1. From Eq. 4.19 we find

$\begin{array}{l}{a}_0=\frac{1}{\pi }{\int }_0^{2\pi }f\left(x\right)dx\hfill \\ =\frac{1}{\pi }\left({\int }_0^{\pi }1dx+{\int }_{\pi }^{2\pi }\left(-1\right)dx\right)\hfill \\ =\frac{1}{\pi }\left({\int }_0^{\pi }1dx-{\int }_{\pi }^{2\pi }\left(1\right)dx\right)\hfill \\ =0\hfill \end{array}$

Furthermore, we find from Eq. 4.20 for n = 1

$\begin{array}{l}{a}_1=\frac{1}{\pi }{\int }_0^{2\pi }f\left(x\right)\cos xdx\hfill \\ =\frac{1}{\pi }\left({\int }_0^{\pi }1\cos xdx+{\int }_{\pi }^{2\pi }\left(-1\right)\cos xdx\right)\hfill \\ =\frac{1}{\pi }\left(\left(\sin \pi -\sin \theta \right)+\left(\sin 2\pi -\sin \pi \right)\right)\hfill \\ =0\hfill \end{array}$

which tells us that a₁ = 0. In fact, we find that a_2n+1 = 0 for all n, which leaves only the terms b_n as relevant terms. From Eq. 4.21 we find for n = 1

$\begin{array}{l}{b}_1=\frac{1}{\pi }{\int }_0^{2\pi }f\left(x\right)\sin xdx\hfill \\ =\frac{1}{\pi }\left({\int }_0^{\pi }1\sin xdx+{\int }_{\pi }^{2\pi }\left(1\right)\sin xdx\right)\hfill \\ =\frac{1}{\pi }\left(\left(\cos \pi -\cos \theta \right)+\left(\cos 2\pi -\cos \pi \right)\right)\hfill \\ =\frac{4}{\pi }\hfill \end{array}$

So in total we find for the Fourier expansion of the square function

$f_{\text{square}}\left(x\right)=\frac{4}{\pi }\sin x$

Higher order expansions can be calculated likewise. For these we find

$\begin{array}{ll}n=1\hfill & f_{\text{square}}\left(x\right)=\frac{4\sin x}{\pi }\hfill \\ n=3\hfill & f_{\text{square}}\left(x\right)=\frac{4\sin x}{\pi }+\frac{4}{3}\cdot \frac{4\sin \left(3x\right)}{\pi }\hfill \\ n=5\hfill & f_{\text{square}}\left(x\right)=\frac{4\sin x}{\pi }+\frac{4}{3}\cdot \frac{4\sin \left(3x\right)}{\pi }+\frac{4}{5}\cdot \frac{4\sin \left(5x\right)}{\pi }\hfill \\ n=7\hfill & f_{\text{square}}\left(x\right)=\frac{4\sin x}{\pi }+\frac{4}{3}\cdot \frac{4\sin \left(3x\right)}{\pi }+\frac{4}{5}\cdot \frac{4\sin \left(5x\right)}{\pi }+\frac{4}{7}\cdot \frac{4\sin \left(7x\right)}{\pi }\hfill \\ n=9\hfill & f_{\text{square}}\left(x\right)=\frac{4\sin x}{\pi }+\frac{4}{3}\cdot \frac{4\sin \left(3x\right)}{\pi }+\frac{4}{5}\cdot \frac{4\sin \left(5x\right)}{\pi }+\frac{4}{7}\cdot \frac{4\sin \left(7x\right)}{\pi }+\frac{4}{9}\cdot \frac{4\sin \left(9x\right)}{\pi }\hfill \end{array}$

The expansion functions together with the square wave function are shown in Fig. 4.7. As you can see, an expansion to n = 1 gives only a very vague approximation of the square wave function. If more and more sine functions of higher orders are added, the approximation becomes gradually better until for n = 9 one can already see the square wave function quite clearly.

Gibbs Phenomenon. As you can see in Fig. 4.7, the Fourier transform of a piecewise continuous function tends to form periodic oscillations at locations of discontinuities. This phenomenon is named after Willard Gibbs¹ who described it in 1899 [5].

4.3.2.3 Fourier Series Expansion of a Triangular-Like Function

As a next example, we will look at a somewhat odd function that vaguely resembles a triangular function. The function is given by

$f_{\text{triang}}\left(x\right)=\{\begin{array}{ll}1\hfill & x<\frac{x}{\pi }\hfill \\ -1\hfill & x\ge 1-\frac{x}{\pi }\hfill \end{array}$

We can repeat the indicated process to determine the coefficients a₀, a_n, and b_n. The Fourier series expansions to a given expansion order n can then be found to be

$\begin{array}{ll}n=1\hfill & f_{\text{triang}}\left(x\right)=\frac{4\cos x}{{\pi }^2}+\frac{2\sin x}{\pi }\hfill \\ n=3\hfill & f_{\text{triang}}\left(x\right)=\frac{4\cos x}{{\pi }^2}+\frac{2\sin x}{\pi }-\frac{4\cos \left(3x\right)}{9\pi }+\frac{2\sin \left(3x\right)}{3\pi }\hfill \\ n=5\hfill & f_{\text{triang}}\left(x\right)=\frac{4\cos x}{\pi }+\frac{2\sin x}{\pi }-\frac{4\cos \left(3x\right)}{9{\pi }^2}+\frac{2\sin \left(3x\right)}{3\pi }-\frac{4\cos \left(5x\right)}{25{\pi }^2}+\frac{2\sin \left(5x\right)}{5\pi }\hfill \\ n=7\hfill & f_{\text{triang}}\left(x\right)=\frac{4\cos x}{{\pi }^2}+\frac{2\sin x}{\pi }-\frac{4\cos \left(3x\right)}{9{\pi }^2}+\frac{2\sin \left(3x\right)}{3\pi }-\frac{4\cos \left(5x\right)}{25{\pi }^2}+\frac{2\sin \left(5x\right)}{5\pi }-\frac{4\cos \left(7x\right)}{49{\pi }^2}+\frac{2\sin \left(7x\right)}{7\pi }\hfill \\ n=9\hfill & f_{\text{triang}}\left(x\right)=\frac{4\cos x}{{\pi }^2}+\frac{2\sin x}{\pi }-\frac{4\cos \left(3x\right)}{9{\pi }^2}+\frac{2\sin \left(3x\right)}{3\pi }-\frac{4\cos \left(5x\right)}{25{\pi }^2}+\frac{2\sin \left(5x\right)}{5\pi }-\frac{4\cos \left(7x\right)}{49{\pi }^2}+\frac{2\sin \left(7x\right)}{7\pi }-\frac{4\cos \left(9x\right)}{81{\pi }^2}+\frac{2\sin \left(9x\right)}{9\pi }\hfill \end{array}$

The expansion of the Fourier series is shown in Fig. 4.8. As you can see, at expansion orders around or above n = 9, we can approximate the function to a very high degree.

4.3.2.4 Discontinuities

When you carefully look at both cases we just discussed, you will notice one important point: both functions have discontinuities, i.e., sharp steps at which the function cannot be differentiated. We have briefly touched upon differentiability when discussing the Taylor series for which we stated that the function must be continuous and the derivatives must be defined in the domain of interest, otherwise, the Taylor expansion will fail. This is not the case for the Fourier series. In fact, if a function has a discontinuity at x₀ and the values of the function when approaching “from the left side” (which is often written as f (x → −x₀)) and “from the right side” (which is then written as f (x → + x₀)) are not equal, the Fourier series will deliver the arithmetic mean of the two values. This is interesting and can be clearly seen, e.g., when looking at Fig. 4.7. The discontinuity here is at x = π, and the two limit values are 1 and −1, in which case the Fourier series expansion has the value 0.

This is a very important property of the Fourier series. This series is stable even if the function does something very “funky” within the domain in which we are investigating it.

4.3.2.5 Maple Listing for the Fourier Expansion of f (x)

Just as we did for the Taylor series, we will now convert the equations we found into a Maple worksheet that will allow us to create Fourier series. The listing is shown in listing 4.2.

1  restart:read "Core.txt":
2
3  #makes a Fourier expansion
4  expandfourier:=proc(f,variable,expansionorder) local i,store;
5   return unapply(
6     1/Pi*int(f,variable=0..2*Pi)
7     + sum(
8         1/Pi*int(f*cos(i*variable),variable=0..2*Pi)*cos(i*variable)
9       + 1/Pi*int(f*sin(i*variable),variable=0..2*Pi)*sin(i*variable)
10    ,i=1..n),variable);
11 end proc:
12
13 #make a square wave function
14 f_square:=x->piecewise(x<Pi,1,x>=Pi,-1):
15
16 #make the Fourier expansion
17 expandfourier(f_square(x),x,1);
18 expandfourier(f_square(x),x,3);
19 expandfourier(f_square(x),x,5);
20 expandfourier(f_square(x),x,7);
21 expandfourier(f_square(x),x,9);
22
23 #and plot us some functions
24 quickplot([expandfourier(f_square(x),x,1),x->f_square(x)],0..2*Pi,-2..2,"x","y",["n=1","square function"]);
25 quickplot([expandfourier(f_square(x),x,3),x->f_square(x)],0..2*Pi,-2..2,"x","y",["n=3","square function"]);
26 quickplot([expandfourier(f_square(x),x,5),x->f_square(x)],0..2*Pi,-2..2,"x","y",["n=5","square function"]);
27 quickplot([expandfourier(f_square(x),x,9),x->f_square(x)],0..2*Pi,-2..2,"x","y",["n=9","square function"]);
28
29 #next we look at a quasi-triangle function
30 f_triangle:=x->piecewise(x<Pi,x/Pi,x>=Pi,1-x/Pi):
31
32 #make the Fourier expansion
33 expandfourier(f_triangle(x),x,1);
34 expandfourier(f_triangle(x),x,3);
35 expandfourier(f_triangle(x),x,5);
36 expandfourier(f_triangle(x),x,7);
37 expandfourier(f_triangle(x),x,9);
38
39 #and plot us some functions
40 quickplot([expandfourier(f_triangle(x),x,1),x->f_triangle(x)],0..2*Pi,-2..2,"x","y",["n=1","triang. function"]);
41 quickplot([expandfourier(f_triangle(x),x,3),x->f_triangle(x)],0..2*Pi,-2..2,"x","y",["n=3","triang. function"]);
42 quickplot([expandfourier(f_triangle(x),x,5),x->f_triangle(x)],0..2*Pi,-2..2,"x","y",["n=5","triang. function"]);
43 quickplot([expandfourier(f_triangle(x),x,9),x->f_triangle(x)],0..2*Pi,-2..2,"x","y",["n=9","triang. function"]);

We will briefly walk through the listing. In line 1 we restart the Maple server and include Core.txt. Next we define the function that performs the Fourier series expansion (see line 4). This function takes three arguments: the function f to expand, the independent variable variable, and the expansion order expansionorder. This function first calculates the coefficients a₀ (see line 6), a_n (see line 8) and b_n (see line 9) and form their products with the sine and cosine functions, respectively, before summing them up to the given expansion order (see line 7). The sum is then turned into a function of the independent variable (see line 5) and returned.

In the following, we use this function. First, we define the square wave function Eq. 4.22 in line 14 and then find the Fourier series with the expansion order n = 1, n = 3, n = 5, n = 7, and n = 9 (see line 17). The output of these code lines has already been presented. Next we plot these functions (see line 24), which produces the plots shown in Fig. 4.7.

In the next step, we define the triangular-like function given by Eq. 4.23 (see line 30). Again, we perform the Fourier series expansion (see line 33) for which we already presented the results. These functions are then plotted (see line 40), which produces the plots shown in Fig. 4.8.

As you can see, it is not very difficult to use these functions. Feel free to modify them and evaluate different functions.

4.3.3.1 Period Length

When working with periodic functions, we come across the angular frequency ${\omega }_o=2\pi v=\frac{2\pi }{T}$, which is nothing other than the reciprocal of the period length T . The angular frequency converts the physical length T or the frequency ν into a function of π, which is required for periodic functions.

In order to convert our physical length x_max into an expression of π, we will do a similar thing. We will replace the independent variable x by the independent variable $\tilde{x}=\frac{2\pi }{x_{\text{max}}}x$. If we now look at our function, we find the values $f\left(\tilde{x}=0\right)=0$ and $f\left(\tilde{x}=x_{\text{max}}\right)=\frac{2\pi }{x_{\text{max}}}=2\pi$, both of which would become 0 if the sine function were applied to them. By using this scheme, we are able to convert our function into an expression that can be approximated using a sine function.

Please note that the integration boundaries can be chosen more or less arbitrarily as long as they span a full period, i.e., a full interval. The interval can be shifted and offset in order to match the function best. This is due to the fact that we will make the function periodic. Thus it does not matter at which point of the interval we start, as long as we integrate along a full interval shifting the boundaries will amount to the same integral.

4.3.3.2 Fourier Expansion

Using our nonperiodic functions, we rewrite the coefficients for the Fourier expansion slightly resulting in

$f\left(x\right)=\frac{a_0}{2}+{\displaystyle \sum _{n=1}^{\infty }\left(a_n\cos \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)+b_n\sin \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)\right)}$

$=\frac{a_0}{2}+{\displaystyle \sum _{n=1}^{n_{\max }}\left(a_n\cos \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)+b_n\sin \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)\right)}+O_{n,\max }$

$a_0=\frac{2\pi }{x_{\max }-x_{\min }}{\int }_{x_{\min }}^{x_{\max }}f\left(x\right)dx$

$a_n=\frac{2}{x_{\max }-x_{\min }}{\int }_{x_{\min }}^{x_{\max }}f\left(x\right)\cos \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)dx$

$b_n=\frac{2}{x_{\max }-x_{\min }}{\int }_{x_{\min }}^{x_{\max }}f\left(x\right)\sin \left(n\frac{2\pi }{x_{\max }-x_{\min }}x\right)dx$

4.3.3.3 Fourier Expansion of a Constant to a Sine Series in One Dimension

A Constant Function. We will now derive a very important Fourier series that we will need later, i.e., the Fourier series of a constant. The function we want to expand is defined as

$f_{\text{constant}}\left(x\right)=\{\begin{array}{ll}{c}_0\hfill & 0\le x\le a\hfill \\ -c_0\hfill & \text{otherwise}\hfill \end{array}$

which we convert to a function of π as detailed in section 4.3.3

$f_{\text{constant}}\left(\pi \frac{x}{a}\right)=\{\begin{array}{ll}{c}_0\hfill & 0\le x<a\hfill \\ -c_0\hfill & \text{otherwise}\hfill \end{array}$

As you can see, Eq. 4.30 looks very much like the square function. Please note that we are only interested in the function in the interval 0 ≤ x ≤ a. In most cases we want the function to have the value 0 at the boundaries, i.e., at x = 0 and x = a. We define the function’s value to be a constant c₀ in this interval. We set the function to the negative value outside this interval as displayed in Fig. 4.9, which will make the terms a₀ and a_n disappear, as we will see in a moment. We also use the interval 0 . . . 2a instead of only 0 . . . a, which results in the term $\frac{\pi }{a}$ instead of $\frac{2\pi }{a}$.

Fourier Series of the Constant in One Dimension. We will now derive the terms required for the Fourier according to Eq. 4.26, Eq. 4.27, and Eq. 4.28. In this case we find

$\begin{array}{l}{a}_0=\frac{2}{2a}{\int }_0^{2a}f\left(x\right)dx\hfill \\ =\frac{1}{a}\left(c_{0}a-c_{0}a\right)\hfill \\ =0\hfill \end{array}$

$\begin{array}{l}{a}_{\mathrm{n}}=\frac{2}{2a}{\int }_0^{2a}f\left(x\right)\cos \left(n\pi \frac{x}{a}\right)dx\hfill \\ =\frac{c_0}{a}{\int }_0^{2a}\left(\cos \left(n\pi \frac{x}{a}\right)dx\right)\hfill \\ =\frac{c_0}{a}\frac{a}{n\pi }\left({\left[\sin \left(n\pi \frac{x}{a}\right)\right]}_0^a-{\left[\sin \left(\frac{n\pi }{a}\right)\right]}_0^{2a}\right)\hfill \\ =\frac{c_0}{n\pi }\left(\left(\sin \left(n\pi \right)-\sin 0\right)-\left(\sin \left(n2\pi \right)-\sin \pi \right)\right)\hfill \\ =0\hfill \end{array}$

$\begin{array}{l}{b}_n=\frac{2}{2a}{\int }_0^{2a}f\left(x\right)\sin \left(nx\right)dx\hfill \\ =\frac{c_0}{a}{\int }_0^{2a}\left(\sin \left(n\pi \frac{x}{a}\right)dx\right)\hfill \\ =\frac{c_0}{a}\frac{a}{n\pi }\left({\left[-\cos \left(n\pi \frac{x}{a}\right)\right]}_0^a+{\left[\cos \left(n\pi \frac{x}{a}\right)\right]}_0^{2a}\right)\hfill \\ =\frac{c_0}{n\pi }\left(-\left(\cos \left(n\pi \right)-\cos 0\right)+\left(\cos \left(2\pi n\right)-\cos \left(n\pi \right)\right)\right)\hfill \\ =\frac{c_0}{n\pi }\left(-\cos \left(n\pi \right)+1+1-\cos \left(n\pi \right)\right)\hfill \\ =\frac{c_0}{n\pi }\left(2-2\cos \left(n\pi \right)\right)\hfill \end{array}$

We will simplify the expression for b_n a little more taking into account that, for all n = 0, 2, 4 . . ., cos (nπ) = 1 and therefore b_n = 0 because the term in the bracket is 0. On the other hand, for all n = 1, 3, 4 . . ., the expression cos (nπ) = −1, in which case the term in the bracket becomes 4. Therefore our Fourier series must be evaluated only for all odd values of n. From this we find the Fourier series to Eq. 4.30 to be

$\begin{array}{ll}{f}_{\text{constant}}\left(\pi \frac{x}{a}\right)\hfill & ={\displaystyle \sum _{n=1,3,5\cdots }^{\infty }\frac{4}{n\pi }{c}_0\sin }\left(n\pi \frac{x}{a}\right)\hfill \\ \hfill & =c_0{\displaystyle \sum _{n=0}^{\infty }\frac{4}{\left(2n+1\right)\pi }\sin }\left(\left(2n+1\right)\pi \frac{x}{a}\right)\hfill \end{array}$

Eq. 4.34 is the function we are looking for. It expresses a constant in one dimension as a Fourier series. It is plotted together with the constant function Eq. 4.34 in Fig. 4.10. As you can see, the function does a pretty good job at approximating the constant function.

4.3.3.4 Fourier Expansion of a Constant to a Sine Series in Two Dimensions

In the next step, we will do an expansion of a constant in two dimensions. This function is defined as