Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to Mathematica tutorial for the first course APMA0330
Return to Mathematica tutorial for the second course APMA0340
Return to the main page for the first course APMA0330
Return to the main page for the second course APMA0340
Return to Part V of the course APMA0340
Introduction to Linear Algebra with Mathematica

Glossary

Preface

The Gibbs phenomenon was first noticed and analyzed by the English mathematician Henry Wilbraham (1825--1883) in 1848, and rediscovered by an American scientist J. Willard Gibbs (1839--1903) 50 years later. The term "Gibbs phenomenon" was introduced by the American mathematician Maxime Bôcher in 1906. The history of this discovery can be found in an article by Hewitt&Hewitt. This phenomenon assures that the Fourier series of a piecewise continuously differentiable periodic function oscillates near the jump of discontinuity by overshooting and undershooting it by about 9%. This lack of improvement in the approximations near the discontinuity manifested in the continual presence of the overshoot or undershoot is called the Gibbs phenomenon. Henry Wilbraham received a BA in 1846 and an MA in 1849 from Cambridge.

Gibbs Phenomenon

Josiah Willard Gibbs

Let \( F_N (x) \) be the finite Fourier sum for the periodic function f(x) with N+1 terms:

\[ F_N (x) = \frac{a_0}{2} + \sum_{k=1}^N \left( a_k \cos \frac{k\pi x}{\ell} + b_k \sin \frac{k\pi x}{\ell} \right) , \]
where the Fourier coefficients 𝑎k and bk were defined previously. The Gibbs phenomenon is the peculiar manner in which the Fourier series of a piecewise continuously differentiable periodic function behaves at a jump of discontinuity. The partial Fourier sums ripple near every point of discontinuity in an amount proportional to the finite jump. So finite Fourier sums (as well as other eigenfunction series) provide an overshoot/undershoot (or "ringing") occurring at discontinuities of finite length. Mathematically, it reflects the difference between uniform convergence and pointwise convergence. The Nth partial sum of the Fourier series has large oscillations near the jump, which increases/decreases the maximum/minimum of the partial sum of the function itself. The overshoot does not die out as N increases, but instead approaches a finite limit. If a function f(x) has a discontinuity at the point \( x_0 \) of amount \( c = f(x_0 +0) - f(x_0 -0) , \) then finite Fourier sums experience overshoot and undershoot in a neighborhood of the amount c w where w = 0.0894899… is the Wilbraham constant and c is the value of jump at that point. These undershoots or overshoots cannot be eliminated by increasing the number of terms in the finite Fourier partial sums; actually these finite sums swap the range
\[ \frac{2}{\pi}\,\mbox{Si}(\pi ) \times \left\vert f(x_0 +0) - f(x_0 -0) \right\vert = 1.1789797444721675\ldots \times \left\vert f(x_0 +0) - f(x_0 -0) \right\vert \]
when \( N\mapsto \infty . \) Here Si is the sine integral
\[ \mbox{Si}(x) = \int_0^{x} \frac{\sin t}{t}\,{\text d}t , \]
a special function that can be evaluated by Mathematica with a one line built-in command
N[SinIntegral[Pi]]
1.851937051982466
Division by π/2
N[Pi/2]
1.5707963267948966
yields
N[SinIntegral[Pi]/Pi*2]
1.1789797444721672
So the overshoot/undershoot of a finite Fourier series for a function having a unit jump of discontinuity will exceed the actual value of the function at that point by the Wilbraham constant
\[ w = \frac{1}{2} \left( \frac{2}{\pi}\,\mbox{Si}(\pi ) -1 \right) = \frac{1}{2} \left( \frac{2}{\pi}\,\int_0^{\pi} \frac{\sin t}{t}\,{\text d}t -1 \right) \approx 0.08948987223608362 , \]
which is evaluated with 20 decimal places by Mathematica as
w = N[(SinIntegral[Pi]/Pi*2 - 1)/2, 20]
0.089489872236083635116

If the function is a periodic continuation on all real axis, we do not observe the Gibbs phenomenon. The considered examples show that if a function is a piecewise continuous with finite jumps, its Fourier coefficients decay at a rate of 1/k; if a continuous function has corners (not differentiable), then its Fourier coefficients decay at a rate of 1/k2. The more derivatives exist, the faster Fourier coefficients decay. The general case provides the following lemma.

Riemann--Lebesgue lemma: If a periodic function f of period T is absolutely integrable (so \( \int_{T} \left\vert f(x) \right\vert {\text d}x < \infty ), \) then its Fourier coefficients tend to zero.

Derivation of Wilbraham constant (Click to see the proof)

The shape of the Dirichlet kernel causes two distinct problems at a discontinuity:

As N → ∞, the overshoot does not disappear; it simply compresses closer to the discontinuity because the peak location is 𝑥 = π/N. The Gibbs phenomenon occurs because the Dirichlet kernel contains negative side lobes and fails to absorb oscillations under integration. Increasing N narrows the width of the ringing ripples but fails to reduce the maximum height of the overshoot.

Example 1: Consider a piecewise continuous function (which is actually piecewise constant):

\[ f(x) = \begin{cases} 1 , & \quad \mbox{if $x\in (-2,-1)$, } \\ 0 , & \quad \mbox{if $x\in (-1,0)$, } \\ 2, & \quad \mbox{for }x\in (0,2) . \end{cases} \]
We calculate its Fourier coefficients with Mathematica:

x[t_]= Piecewise[{{1,-2<t<-1},{0,-1<t<0},{2,0<t<2}}]
a0=Integrate[1,{x,-2,-1}]/2 + Integrate[1,{x,0,2}];
an = Integrate[Cos[n*Pi*x/2], {x, -2, -1}] /2 + Integrate[Cos[n*Pi*x/2], {x, 0, 2}];
bn = Integrate[Sin[n*Pi*x/2], {x, -2, -1}] /2 + Integrate[Sin[n*Pi*x/2], {x, 0, 2}];
fourier[m_] = 5/4 + Sum[an*Cos[n*Pi*x/2] + bn*Sin[n*Pi*x/2], {n, 1, m}];
Plot[{fourier[50], fourier[10]}, {x, -2.1, 2.1}, PlotStyle -> {{Thick, Blue}, {Thick, Orange}}, Ticks -> {{-2, -1, 1, 2}, {2.18, 1.09, 0.91, -0.18}}]
\[ a_0 = \frac{1}{2} \,\int_{-2}^{-1} {\text d}x + \int_0^2 {\text d}x = \frac{5}{2} , \quad a_n = - \frac{1}{n\pi}\,\sin\frac{n\pi}{2} , \quad b_n = \frac{1}{n\pi} \left[ (-1)^n -2 + \cos \frac{n\pi}{2} \right] . \]
                 
       Fourier approximation with n=50 vs n=10 terms.            Fourier approximation with n=100 vs n=10 terms.

In the graph shown, the orange curve was calculated with n=10 terms while the blue curve was calculated with n=50 terms. As a result, when the number of terms in the Fourier series is increased, the overshoot and undershoot cannot be removed, but the “ripples” will move closer and closer (horizontally) to the point of discontinuity. There are four points of discontinuity on the interval [-2,2]; at three of them (x = -2, -1, 2) the given function experiences jump of discontinuity of 1 and at the point x = 0 it has jump of value 2.

At point x = -2, Fourier series has undershoot value of \( 1 - 0.0894899 \approx 0.91051 \) because 0.0894899 ≈ 0.1789797444721675 /2 = 2*Si(π)/π - 1. Correspondingly, at point x = -1, we have overshoot 1 + 0.0894899 ≈ 1.0894899… . At point x = 2, the overshoot will be about 2.0894899 and undershoot will be about 0.91051. However, at point x = 0, the jump of discontinuity is 2; therefore, we expect Gibbs phenomenon at this point to be

Changing the FourierParameters setting allows control over the limits of integration on the coefficients, and therefore, the base frequency of the series. Some work will go into calculating what values for parameters will give the proper limits of integration. There is ordinarily no reason to change the first parameter from its default setting, which is 1.

It has the following Fourier series:

curve5 = FourierTrigSeries[x[t],t,5,FourierParameters->{1,Pi/2}]
Out[33]= 3/4 - Cos[(\[Pi] t)/2]/\[Pi] + Cos[(3 \[Pi] t)/2]/(3 \[Pi]) - Cos[(5 \[Pi] t)/2]/(5 \[Pi]) + Sin[(\[Pi] t)/2]/\[Pi] + Sin[\[Pi] t]/\[Pi] + Sin[(3 \[Pi] t)/2]/(3 \[Pi]) + Sin[(5 \[Pi] t)/2]/(5 \[Pi])
When plotted, it looks like this:
Plot[{x[t],curve5},{t,-2,2}, PlotStyle -> Thick, Ticks-> {{-2,-1,0,1,2},{-0.08,1.08}}]
   ■
End of Example 1

Example 2: Consider the continuous function on the interval (−π,π):

\[ f(x) = \begin{cases} 0, & \ \mbox{ if } -\pi < x \le 0, \\ x, & \ \mbox{ if } 0 \le x < \pi . \end{cases} \]
It has a Fourier series expansion on [-π,π]:
\[ f(x) \sim \frac{\pi}{4} - \frac{2}{\pi} \sum_{k\ge 0} \frac{1}{(2k+1)^2} \, \cos \left( (2k+1)x \right) - \sum_{n\ge 1} \frac{(-1)^n}{n} \,\sin nx . \]
At point x=π, the function f has a jump of magnitude π≈3.1415926. Therefore, the overshoot will be around
N[Pi*(0.1789797444/2)+Pi]
3.42273
and undershoot similarly
-N[Pi*(0.1789797444/2)]
-0.281141
Now we plot partial sums to check Gibbs phenomenon around points x=-π and x=π:
F[x_]= Pi/4 -2/Pi*Sum[1/(2*k+1)^2 * Cos[(2*k+1)*x], {k,0,10}] - Sum[(-1)^n /n *Sin[n*x], {n,1,21}];
f[x_]= Piecewise[{{0,-Pi < x < 0},{x,0< x < Pi}}];
Plot[{F[x],f[x]},{x,-3.45,3.45}, PlotStyle->Thick, Ticks -> {{-3,-2,-1,1,2,3},{3.42, -0.28}}]
Plot[{F[x], f[x]}, {x, 3.0, 3.45}, PlotStyle -> Thick, Ticks -> {{Pi}, {3.42, -0.28}}]
Fourier series approximation.
Fourier series approximation near point π.
   ■
End of Example 2

Example 3: Consider the piecewise continuous function on the interval (−2,2):

\[ f(x) = \begin{cases} 1, & \ \mbox {if } -2 < x < -1 , \\ x^2 -1, & \ \mbox {if } -1< x< 2 , \end{cases} \]
We usually define a function in an open interval because its values at endpoints are dictated by its Fourier expansion. Its Fourier coefficients were found previously in section iv:
\[ f(x) = \frac{1}{4} + \sum_{n\ge 1} \frac{8(-1)^n n\pi + 4n\pi \cos \left( \frac{n\pi}{2} \right) - \left( 8 + n^2 \pi^2 \right) \sin \left( \frac{n\pi}{2} \right)}{n^3 \pi^3} \, \cos \frac{n\pi x}{2} - \sum_{n\ge 1} \frac{2 (-1)^n \left( n^2 \pi^2 -4 \right) + \left( 8 + n^2 \pi^2 \right) \cos \left( \frac{n\pi}{2} \right) + 4n\pi\,\sin \left( \frac{n\pi}{2} \right)}{n^3 \pi^3} \,\sin \frac{n\pi x}{2} \]
On the main interval (-2,2), the given function has two jumps of discontinuity:
\[ \lim_{x\to -1-0} \,f(x) = 1, \quad \lim_{x\to -1+0} \,f(x) = 0; \qquad \lim_{x\to 2-0} \,f(x) = 3, \quad \lim_{x\to 0+0} \,f(x) = 1. \]
So, at point x = -1 the function experiences a jump of 1 unit and at the point x = 2 it has a jump of 2 units. Therefore, the undershoot at point x = -1 will be closed to -w ≈ -0.089489872 (Wilbraham number) and overshoot of 1+w ≈ 1.089489872 as the number of terms in the partial Fourier sums approaches infinity. At another point of discontinuity x = 2, the given function will have the undershoot of the amount 2w, which will be about 2.82102. Correspondingly, its overshoot will approach 3.1789797. Upon making a computational experiment with plotting partial sums in neighborhoods of these two points, we see that indeed, the predicted overshoots and undershoots are very closed to the predicted ones.
Plot[F10[x], {x, -2.4, 1.4}, PlotStyle -> {Thick, Red}, Ticks -> {{-1}, {-0.089489872, 1.089489872}}]
Plot[F20[x], {x, -1.25, -0.85}, PlotStyle -> {Thick, Orange}, Ticks -> {{-1}, {-0.089489872, 1.089489872}}]
Plot[F100[x], {x, -1.05, -0.95}, PlotStyle -> {Thick, Blue}, Ticks -> {{-1}, {-0.089489872, 1.089489872}}]
   Fourier approximation with 10 terms    Fourier approximation with 20 terms    Fourier approximation with 100 terms
       
Similarly, we plot Fourier series in neighborhoods of another point x = 2.
Plot[F10[x], {x, 1.55, 2.55}, PlotRange -> {{1.55, 2.55}, {0.66, 3.2}}, PlotStyle -> {Thick, Red}, Ticks -> {{2}, {1, 2, 3, 0.8210203, 3.1789797}}]
Plot[F20[x], {x, 1.65, 2.45}, PlotRange -> {{1.65, 2.45}, {0.66, 3.2}}, PlotStyle -> {Thick, Orange}, Ticks -> {{2}, {1, 2, 3, 0.8210203, 3.1789797}}]
Plot[F100[x], {x, 1.95, 2.09}, PlotRange -> {{1.95, 2.09}, {0.66, 3.2}}, PlotStyle -> {Thick, Blue}, Ticks -> {{2}, {1, 2, 3, 0.8210203, 3.1789797}}]
   Fourier approximation with 10 terms    Fourier approximation with 20 terms    Fourier approximation with 100 terms
       
   ■
End of Example 4

Example 4: Consider the piecewise signum function on the interval (-π,π):

\[ \mbox{sign}(x) = \begin{cases} -1, & \ \mbox {if } -\pi < x < 0 , \\ \phantom{-}1, & \ \mbox {if } 0< x< \pi . \end{cases} \]
Since the signum function is odd, we expand it into sine Fourier series
\[ \mbox{sign}(x) = \sum_{n\ge 1} b_n \sin \left( n\,x \right) , \]
where
\[ b_n = \frac{2}{\pi} \int_0^{\pi} \sin \left( n\,x \right) {\text d}x = \frac{2}{n\pi} \left( 1 - (-1)^n \right) = \frac{4}{\pi} \times \begin{cases} \frac{1}{2k+1} , & \ \mbox{ if } n = 2k+1 , \\ 0, & \ \mbox{ otherwise.} \end{cases} \]
because
Integrate[Sin[n*x], {x, 0, Pi}]*2/Pi
(2 (1 - Cos[n \[Pi]]))/(n \[Pi])
First, we plot its Fourier series using standard Mathematica command.
ff0[x_] = unitstep[Mod[x, Pi, -Pi]];
f[x_] = Piecewise[{{-1, -Pi < x < 0}, {1, 0 < x < Pi}}];
ft6ff0 = FourierSinSeries[ff0[x], x, 6];
ftff0[x_, nmax_] := Sum[Sin[(2 n + 1) .9 x]/(2 n + 1), {n, 0, nmax}] 1.275
stepplot[n_] := Plot[{f[x], ftff0[x, n]}, {x, -1, 1}, Ticks -> {{1.18, 1, -1.18}}, PlotStyle -> {{Thick, Blue}, {Thick, Red}}, Exclusions -> None]
stepplot[10]
stepplot[100]
             
       Fourier approximation with n=10 terms.            Fourier approximation with n=100 terms.

Then we verify this series directly
\[ \mbox{sign}(x) = \frac{4}{\pi} \,\sum_{k\ge 0} \frac{1}{2k+1}\, \sin \left( (2k+1)\,x \right) . \]
The corresponding plot confirms:
g[x_] = 4/Pi*Sum[1/(2*k + 1)*Sin[(2*k + 1)*x], {k, 0, 30}];
Plot[g[x], {x, -2*Pi, 2*Pi}]
Since the signum function has a finite jump of 2 at the origin, the overshoot is about 1.18 (which is 2 times the Wilbraham constant) and undershoot is −1.18.    ■
End of Example 4

Lanczos σ-factors

The Lanczos σ-factorse multipliers applied directly to the Fourier coefficients. They smooth out the abrupt cutoff of high-frequency terms by scaling each coefficient using a normalized sinc function. They were introduced by Cornelius Lanczos (1893--1974) in 1956, who also invented what is now called the Fast Fourier Transform. For an N-th partial sum, the factors are defined as (independently on of the period):

\[ \sigma_k = \mbox{sinc} \left( \frac{k}{N} \right) = \dfrac{\sin \left( \frac{\pi k}{N} \right)}{\frac{\pi k}{N}} \]
The modified Fourier partial sum becomes:
\[ S_N^{\sigma} (f;x) = \sum_{k=-N+1}^{N-1} \ \sigma_k \,\hat{f} (k)\, e^{\mathbf{j}\,kx} . \]
Similar for trigometric series:
\[ S_N^{\sigma} (f;x) = \frac{1}{2}\,a_0 + \sum_{k=1}^N \sigma_k \left[ a_k\,\cos \left( k\,\frac{\pi}{\ell}\,x \right) + b_k \,\sin \left( k\,\frac{\pi}{\ell}\,x \right) \right] . \]
The mathematical magic of the σ-factor lies in how it translates to the spatial domain. Multiplying coefficients by a sinc function in the frequency domain is equivalent to taking a local moving average in the spatial domain.

Specifically, the modified sum is exactly equal to the average of the original partial sum over a small window of width 2π/N:

\[ S_N^{\sigma} (f;x) = \frac{N}{2\pi}\,\int_{x- \pi /N}^{x+ \pi /N}\ S_N (f;t)\,{\text d}t . \]
Because the window size 2π/N matches the exact period of the highest-frequency oscillation in the partial sum, the integration flattens out the rapid peaks and valleys of the Gibbs ringing.

By modifying the Fourier coefficients, the underlying convolution kernel changes. The standard Dirichlet kernel is replaced by the Lanczos Kernel (LN):

\[ L_N (x) = \sum_{k=-N+1}^{N-1} \ \mbox{sinc} \left( \frac{k}{N} \right) e^{\mathbf{j}\,kx} . \]
While it does not completely eliminate negative side lobes like the Fejér kernel (Cesàro summation) does, it reduces them significantly enough to suppress the overshoot from 9% down to less than 1.2%, while preserving a much sharper main lobe (better resolution at the jump).

Lanczos σ-factors smooth out the overshoot by acting as a low-pass window function. They implement a spatial moving average of width 2π/N, which successfully cancels out the high-frequency oscillations of the Gibbs phenomenon

Now we derive the moving average identity by starting with the spatial moving average of the standard partial sum SN(𝑓;𝑥) over an interval of width 2π/N and prove that it simplifies to the Fourier series modified by the Lanczos σ-factors.

Let A(𝑥) be the moving average of the -th partial sum SN(𝑓;𝑥) centered at 𝑥:

\[ A(x) = \frac{N}{2\pi}\,\int_{x - \pi /N}^{x + \pi /N} \ S_N (f;t)\,{\text d}t . {\text d}t . \]
Substitute the standard definition of the Fourier partial sum, \( \displaystyle \quad S_N (f;x) = \sum_{k=-N+1}^{k=N-1} \ \hat{f}(k) \,e^{\mathbf{j}\,kt/span> into the integral to obtain
\[ A(x) = \frac{N}{2\pi}\,\int_{x - \ell /N}^{x + \ell /N} \ \left( \sum_{k= -N+1}^{N-1} \ \hat{f}(k)\, e^{\mathbf{j}\,k\pi t/\ell} \right) {\text d}t . \]
Because the sum is finite, we can safely swap the order of summation and integration, pulling out the constants and coefficients:
\[ A(x) =\sum_{k= -N+1}^{N-1} \ \hat{f} (k) \cdot \left[ \frac{N}{2\pi}\,\int_{x - \pi /N}^{x + \pi /N} \ e^{\mathbf{j}\,kt} \ {\text d}t \right] \]
To evaluate the inner integral, we consider two cases for k. .

Case 1 when k = 0.

The exponential term becomes e⁰ = 1. The integral simplifies to the width of the interval:
\[ \frac{N}{2\pi}\,\int_{x - \pi /N}^{x + \pi /N} \ 1\,{\text d}t = \frac{N}{2\pi} \cdot \left( \frac{2\pi}{N} \right) = 1 . \]

Case 2 when k ≠ 0.

Find the antiderivative of exp{ⅉkx} and evaluate it at the boundaries:
\[ \int_{x - \pi /N}^{x + \pi /N} \ e^{\mathbf{j}\,kt} \ {\text d}t = \left[ \frac{e^{\mathbf{j}\,kt}}{\mathbf{j}k} \right]_{t= x- \pi /N}^{t = x + \pi /N} \ = \dfrac{e^{\mathbf{j}k \left( x + \frac{\pi}{N}\right)} - e^{\mathbf{j}k \left( x - \frac{\pi}{N}\right)}}{\mathbf{j} k} . \]
Factor out the common term exp{ⅉkx} in the numerator:
\[ \frac{1}{\mathbf{j}k} \, e^{\mathbf{j}\,kx} \,\left( e^{\mathbf{j}\,\pi k/N} - e^{-\mathbf{j}\,\pi k/N}\right) . \]
Apply Euler's identity to rewrite the exponential difference as a sine wave
\[ e^{\mathbf{j}\,kx} \,\cdot\,\frac{2\mathbf{j}\,\sin \left( \frac{\pi k}{N} \right)}{\mathbf{j}k} = \frac{2\,\sin \left( \frac{\pi k}{N} \right)}{k} \cdot e^{\mathbf{j}\,kx} . \]
Now, multiply this result by the front scaling factor N/(2π) from our moving average definition:
\[ \frac{N}{2\pi} \cdot \left[ \frac{2\,\sin \left( \frac{\pi k}{N} \right)}{k} \,e^{\mathbf{j}\,kx} \right] = \dfrac{\sin \left( \frac{\pi k}{N} \right)}{\frac{\pi k}{N}} \cdot e^{\mathbf{j}\,kx} . \]
Recognize that this fraction is exactly the normalized function used to define the Lanczos σ-factor
\[ \sigma_k = \dfrac{\sin \left( \frac{\pi k}{N} \right)}{\frac{\pi k}{N}} . \]
Substitute this simplified expression back into our summation formula from Step 2 (when k ≠ 0):
\[ A(x) = \sum_{k= -N+1}^{N-1} \ \sigma_k \,\hat{f}(k)\, e^{\mathbf{j}\,kt} . \]

The spatial moving average is mathematically identical to multiplying the Fourier coefficients by the Lanczos σ-factors:

\[ \frac{N}{2\pi} \,\int_{x - \pi /N}^{x + \pi /N} \ S_N (f;t)\,{\text d}t = S_N^{\sigma} (f;x) . \]

Example 5: We consider the odd function: \[ f(x) = \mbox{sign}\left( \sin \left( \frac{\pi x}{\ell} \right) \right) , \qquad x \in [-\ell , \ell ] . \] We use Mathematica to demonstrate specific numerical milestones to illustrate exactly how the σ-factor tames the Gibbs phenomenon for an partial sum on an interval where ℓ = 1:

  • The Discontinuity (The Jump): Located at 𝑥 = 0, where the true function abruptly jumps from −1 to +1.
  • The Classical Peak (Gibbs Overshoot): Occurs at 𝑥 ≈ 1/N = 0.02 because we use N = 50.
    • Classical Fourier Value: ≈ 1.089 (an 8,95% overshoot above the true ceiling of 1.0).
    • Lanczos -Factored Value: ≈ 1.011 (the overshoot is suppressed to just 1.1%).
  • The Ringing Valley: Occurs slightly further out at 𝑥 ≈ 2/N = 0.04.
    • Classical Fourier Value: Drops hard to ≈ 0.951 (visible dipping/ringing ripple).
    • Lanczos -Factored Value: ≈ 1.000 (perfectly flattened out).
????? to be checked

(* 1. Define Parameters *) L = 1.0; N = 50; (* 2. Define Fourier Sums and True Function *) fTrue[x_] := Sign[Sin[(Pi*x)/L]]; fClassical[x_, maxN_] := Sum[ If[OddQ[k], (4/(Pi*k))*Sin[(Pi*k*x)/L], 0], {k, 1, maxN - 1} ]; fLanczos[x_, maxN_] := Sum[ If[OddQ[k], With[{\[Sigma] = Sin[(Pi*k)/maxN]/((Pi*k)/maxN)}, \[Sigma] * (4/(Pi*k))*Sin[(Pi*k*x)/L] ], 0], {k, 1, maxN - 1} ]; (* 3. Generate Presentation Plot *) Plot[ {fTrue[x], fClassical[x, N], fLanczos[x, N]}, {x, -1.5, 1.5}, PlotStyle -> { Directive[Black, Dashed], Directive[Red, Thick], Directive[Blue, Thickness[0.007]] }, PlotLegends -> {"True Square Wave", "Classical Sum (N=50)", "Lanczos Smooth Sum (N=50)"}, Frame -> True, FrameLabel -> {"Position (x)", "f(x)"}, PlotLabel -> "Gibbs Phenomenon vs. Lanczos \[Sigma]-Factor Smoothing", GridLines -> Automatic, ImageSize -> Large ]
  1. The Core Trade-off: Notice how the blue line (Lanczos) completely deletes the ripples, but it makes the transition slope slightly wider. This explicitly demonstrates the mathematical trade-off between sharpness (frequency resolution) and smoothness (spatial control).
  2. The Amplitude Compression: The classical red spikes shoot wildly past the bounds of [−1,1]. The Lanczos series structurally stays locked within the physical reality of the original dataset.
   ■
End of Example 5

Example 6: In medical imaging, preventing ringing is a literal matter of health safety. An MRI machine does not take a standard picture. It scans the human body in the frequency domain, capturing raw mathematical data known as k-space data. To turn this raw data into a cross-section image of a brain or knee, the scanner applies an Inverse Fast Fourier Transform (IFFT).

Because scanning takes time, the machine can only capture a finite number of frequencies ( -terms) before the patient moves. At boundaries where tissue changes abruptly—such as the sharp transition between dark cerebrospinal fluid and bright brain matter—the abrupt cutoff in the raw k-space data triggers the Gibbs phenomenon.

This causes dark and bright ringing bands to ripple across the actual brain tissue in the final image. In clinical practice, these are called truncation artifacts, and they can look exactly like false lesions, small tumors, or internal bleeding, leading to a dangerous misdiagnosis.

Before passing the raw data to the IFFT algorithm, the scanner's computer multiplies the high-frequency k-space data points by the Lanczos -factors.

  • By smoothing the edge of the data collection window, the false ripples inside the brain vanish.
  • The doctor gets a clean, artifact-free image of the tissue.
To see this in action mathematically, look at a standard 1D discontinuous wave—a Square Wave of length ℓ. . The standard analytical Fourier coefficients for a square wave jumping between −1 and +1 are: \[ a_0 = 0, \quad a_k = 0, \quad b_k = \frac{4}{\pi k} \quad (\mbox{for odd }\ k). \] Below is how you would evaluate a point 𝑥 with and without the σ-factor:    ■
End of Example 6

 

  1. Hewitt, E. and Hewitt, R.E., The Gibbs--Wilbraham phenomenon: An episode in Fourier analysis, Archive for History of Exact Science, 1979, Vol. 21, No. 2, pp. 129--160.
  2. B. Kuttner, Note on the Gibbs Phenomenon, Journal of the London Mathematical Society, Volume s1-20, Issue 3, 1 June 1945, Pages 136–139, https://doi.org/10.1112/jlms/s1-20.3.136

    A correction has been published: Journal of the London Mathematical Society, Volume s1-20, Issue 4, 1 October 1945, Pages 238, https://doi.org/10.1112/jlms/s1-20.4.238-s

 

Return to Mathematica page
Return to the main page (APMA0340)
Return to the Part 1 Matrix Algebra
Return to the Part 2 Linear Systems of Ordinary Differential Equations
Return to the Part 3 Non-linear Systems of Ordinary Differential Equations
Return to the Part 4 Numerical Methods
Return to the Part 5 Fourier Series
Return to the Part 6 Partial Differential Equations
Return to the Part 7 Special Functions