Numerical Solvability and Solution of an Inverse Problem Related to the Gibbs Phenomenon

Abstract

We report on the inverse problem for the truncated Fourier series representation of in a form with a quadratic degeneracy, revealing the existence of the Gibbs-Wilbraham phenomenon. A new distribution-theoretic proof is proposed for this phenomenon. The paper studies moreover the iterative numerical solvability and solution of this inverse problem near discontinuities of .

1. Introduction

This paper reinvestigates the Fourier series [1, 2] representation of a piecewise continuous –periodic signal , with

The Gibbs phenomenon (see, e.g., [3]) is a statement of the fact that the infinite series tends to overshoot the positive corner of a discontinuity of by ~9% of the jump size and to undershoot the negative corner by the same amount. An overshoot/undershoot effect that is accompanied with spurious oscillations “ringing” near the discontinuity (as described in the appendix) when the series in (1) is truncated at . However, according to a theorem by Fejer [1], the infinite series of a function converges to at each point of continuity of .

This shortcoming in the infinite Fourier series representation of piecewise continuous was first observed by H. Wilbraham in 1848 and then analyzed in detail [4] by J. W. Gibbs in 1898. The main reasons of the Gibbs-Wilbraham effect are that (i) not all frequencies (only integer ones) are employed in (1), (ii) and happen to decay slowly with increasing , and (iii) the global nature of the approximation of : the expansion coefficients are obtained, via (2), by integration over the entire period, including the points of discontinuity.

What is unpleasant though, with all of this, is that the Gibbs-Wilbraham effect is generic and is present for any periodic signal with isolated discontinuities. The presence of this effect can in fact lead to quite negative consequences when single infinite Fourier series, multiple infinite Fourier series, or even infinite wavelet series are employed to approximate signals of various dimensions, in many fields such as radio engineering and signal transmission.

The Gibbs-Wilbraham effect can nevertheless have both positive and negative consequences in different applications. The negative consequences call for Gibbs effect reduction, and this can in principle be achieved with the use [4] of a variety of filters. This effect can also be reduced theoretically and for all purposes. In addition to classical mathematical filters, recently (in 2011) Rim and Yun [5] defined a kind of spectral series to filter off completely the Gibbs effect near a discontinuity. The construction of this series is based on the method of adding the Fourier coefficients of a Heaviside function to the given Fourier partial sums.

This paper is organized as follows. After this introduction, we propose in Section 2 a new distribution-theoretic proof for the Gibbs-Wilbraham effect. In Section 3, we advance the inverse problem of the truncated Fourier series and the new stereographic truncated Fourier series and its relation to resolution of the undershoot/overshoot pair associated to the Gibbs-Wilbraham effect. Here we illustrate how solving the inverse problem for the truncated representation (1) contains a quadratic degeneracy that indicates the existence of the Gibbs phenomenon. Section 4 deals with the numerical iterative solution of the previous inverse problem and its convergence, with a proof and demonstration that stereographic projection does not affect such a numerical solution and its posedness.

2. Distribution-Theoretic Proof for the Gibbs-Wilbraham Effect

Let be a –periodic piecewise continuously differentiable function on , and let be a point in at which has a discontinuity of the first kind. Consider , and , and , to assume, without loss of generality, that and define the discontinuity jump by Moreover,

Theorem 1. Let a –periodic be discontinuous at . The Fourier series representation of converges over , as so as, where , an uncertainty.

Proof. According to distribution theory [6] when coincides a.e. with over , it may always be differentiated; namely, with as Dirac’s delta. Cancellation of this differentiation is performed by sweeping the -axis during a compensating integration from to . Indeed, integration of (6) first from to leads to and since , then with as Heaviside’s unit step function.
By a theorem of Fejer [1], should converge to when , on . Hence the term in (8) can have a compact support only on an infinitesimal interval to the right of , that is, ending at . Hence Furthermore, despite the fact that , , which is an uncertainty to be denoted by . Consequently, relation (9) is representable as Next complete the sweeping by integrating (6) from to ; namely, Here again since , then the left-hand side of (11) becomes . The right-hand side requires standardization by means of the substitution , leading to , becoming and It is obvious then that (11) is reducible to Repeated application of the same arguments employed in the derivation of (9)-(10) to (13) allows for writing it as and subsequently in the form Taking into consideration that and in (10) and (15), respectively, then combining these equations with Dirichlet’s theorem [1] leads to the required result.

It should be pointed out that in the previous proof relations (10) and (15) say nothing specifically about the behavior of at and , respectively, while asserting the existence of a spiky behavior of individually at then at . Moreover, the uncertainty had been resolved in the past (by Wilbraham; see also the Appendix) computationally to satisfy and long before the advance of the theory of distributions (generalized functions) by S. L. Sobolev in 1936. The number for remains until now, however, a mysterious theoretical puzzle that calls for a rigorous distribution-theoretic justification.

Clearly, in the neighborhood of , can only be a rudimentary approximation to , even when , with an added spike pair to at . Relation (16) represents the Gibbs-Wilbraham effect which, due to the fact that , is not practically computable. This effect degenerates, however, when is finite, to an undershoot-overshoot pair on an interval , splitted by at . In this case, it is always possible to compute Furthermore, since , then these computations with (17) over can reveal only a feature , associated with the Gibbs-Wilbraham effect but not the effect itself.

3. The Inverse Problem for Truncated Fourier Series

Substitution of any value of in (1) yields a number , which may differ [2] from by if is in the neighborhood of a discontinuity of . Obviously satisfies, the nonlinear in , trigonometric series equation which is, however, linear in the data .

Definition 2 (the inverse problem for the truncated Fourier series representation). Given a number satisfying (18), what is the corresponding ?

This inverse problem happens to be a discretized version of the similar deconvolution problem [4] of from a given number in a singular first kind integral equation of the form in which represents the data. It is well known, additionally (see, e.g., [7, 8]) that this deconvolution problem can be computationally ill-posed if an arbitrarily small deviation of the data may cause an arbitrarily large deviation in the solution . This situation may further prevail when the operator does not continuously depend on .

Moreover, the inverse problem of the previous definition has a unique solution that satisfies (19) and is expected to accept a double root when the horizontal line intersects with a sharp peak (overshoot) of a truncated . We shall call this feature a duplicated (quadratic) degeneracy of the solution for this inverse problem.

The existence of the Fourier series truncation peak (and the corresponding two roots for ) is not directly evident from the structure of , and despite its linearity in the data .

3.1. The Tangent Half-Angle Stereographic Projection

The tangent half-angle substitution (the Euler-Weirstrass substitution [9]) is widely used in integral calculus for finding antiderivatives of rational expressions of circular function pairs . Geometrically the construction goes like this: draw the unit circle, and let be the point . A line through (except the vertical line) passes through a point on the circle and crosses the -axis at some point . Obviously determines the slope of this line. Furthermore, each of the lines (except the vertical line) intersects the unit circle in exactly two points, one of which is . This determines a function from points on the unit circle to points on a straight line. The circular functions determine therefore a map from points on the unit circle to slopes . In fact this is a stereographic projection [9] that expresses these pairs in terms of one variable: .

Numerical algorithms employing these trigonometric pairs like Fourier series (1) could in principle be influenced by such a projection. In this context, the tangent half-range transformation relates the angle to the slope of a line intersecting the unit circle centered around .

Indeed as varies, the point winds repeatedly around the unit circle centered at . The point goes only once around the circle as goes from to and never reaches the point , which is approached as a limit when .

Clearly the stereographic projection is defined on the entire unit circle except at the projection point . Where it is defined, this mapping is smooth and bijective [9]. It is conformal, meaning that it preserves angles. it is however not isometric; that is, it does not preserve distances as

The rest of this section illustrates how stereographic projection of Fourier series reveals a quadratic degeneracy that indicates the existence of the Gibbs-Wilbraham effect. It also addresses the question of solving the inverse problem for the truncated Fourier series representations. The impact of stereographic projection on the ill-posedness of the numerical solution of this inverse problem is moreover a main question that has been addressed in this paper.

Our analysis starts however by illustrating how the tangent half-angle stereographic projection can provide a new insight into the behavior of truncated Fourier series near a discontinuity of .

Lemma 3. The solution of the inverse problem for the truncated Fourier series representation satisfies

Proof. By studying solvability for of (18) Let , and then , and . Further consideration of these facts in (18) leads to the tangent half-angle stereographically transformed, , form The harmonic trinomial in the left-hand side of (23) has two harmonic zeros, , for each ; namely, It is worth noting at this point that these zeros should, by no means, be confused with the possible roots of the variable. Substitution then of (24) in (23) leads to which is the same as The quadratic degeneracy of the harmonic zeros of can disappear only when , for all . And this happens to be equivalent with Obviously, satisfaction of (27) by (26) reduces it to the equation which is the same as It is interesting to remark that condition (27) is entirely different from which can be derived from the Parseval identity [1] (completeness relation) Hence a possible nondegenerate solution for of the nonlinear equation (29) is, on one hand, fundamentally incorrect, because it involves a violation of (30) and Parseval’s identity. On the other hand, a consistent solution , in the sense of satisfaction of (31), is intrinsically degenerate because of its emergence, as a must, from the nonlinear equation (26), which rewrites as Equation (32) is easily transformable to the required result (23); here the proof completes.

3.2. The Stereographic Fourier Series

To further analyze the variability of stereographically projected Fourier series, we make the substitution in to rewrite it as

The tangent half-angle (stereographic) transformation or is conceived as and subsequently It is clear then that this stereographic projection effectively leads to insertion of an additional bracket into the composition of as a function of . This paves the way to the definition

Proposition 4. The following: always holds.

Proof. Since is finite, then the proof of implies the correctness of this proposition. So let us differentiate with respect to by applying the chain rule; namely, which is the same as Take then the transformation of the previous relation to arrive, after some algebra, at Invoke then to differentiate it with respect to as which is the required result.

Let us rewrite (1) in the form then try to differentiate it term-by-term as where is the harmonic conjugate of . The trigonometric series is called the trigonometric series conjugate (see, e.g., [10]) to .

In addition to the fact that the convergence of can sometimes become worse than convergence of (see, e.g., [10]) the presence of the factor in can aggravate the convergence problem of (45) when the coefficients and do not decay fast enough.

Next we may stereographically transform , using the notation , where Similarly, for the stereographic Fourier series representation in which we may approximate (44) by Differentiate then (49) term-by-term as Note here that is the harmonic conjugate of and that This rather interesting fact leads to a remark that follows.

Remark 5. Away from points of discontinuity of , the stereographically transformed differentiated Fourier series , which defines , is the same as the derivative , which defines , of the stereographic Fourier series . Moreover, even at points of discontinuity of .

Theorem 6 (see [1]). Let be (i) continuous on such that (ii) , and let (iii) be piecewise continuous over ; then the corresponding Fourier series of (1) or (47) is differentiable at each point where (iv) exists.

The previous remark can incidentally have useful applications when differentiating certain slowly converging Fourier series, as illustrated by the following example.

Example 7. Consider the -periodic odd signal , defined over one period, by . Its Fourier series is known to be that is, with , , , and .
Differentiating term-by-term, we obtain which is not uniformly convergent because is not only discontinuous at but and does not exist. In particular, at , the sum of the series oscillates between and , that is, divergent. It is also well known, nevertheless, that the series is summable [1, 2] to and that turns out to yield the right answer . Moreover, at , the sum of the series , that is, divergent.
Alternatively, by means of (47) and by means of (46) Here, at , the sum of the series is also divergent.

Note in general that, for even signals, that is, when , , it follows from (22) that Moreover, for odd signals, that is, when , , the situation is not of the same clarity.

Indeed, according to (22), which is an uncertainty. Luckily however this uncertainty is not an essential one and can straightforwardly be resolved by resorting to the stereographic form (23) of (22), which yields the correct result Also when , , and this is resolvable in a similar fashion.

Both , of (2), and , of (22), equivalent forms are nonlinear trigonometric functions of the solution. However, while is linear in , thestereographic is nonlinear in this . This can raise a question on a possibility for at singular points of . A question that will be answered negatively in the next section. Furthermore, since the data is discrete, then the operators and do not depend continuously on , and this is the main reason for ill-posedness of the present inverse problem.

In conclusion, the quadratic degeneracy of the solution of the present inverse problem, explicit when using , is a nonlinear indicator of the existence of the peak pertaining to the Gibbs phenomenon. However, the numerical ill-posedness, associated with the relative values of a trial iterative variable in the set, happens to remarkably remain invariant under the tangent half-angle stereographic projection.

Following the same arguments of Wilbraham, we may cancel the artificially induced differentiation by integrating (50) to represent the stereographic Fourier series in integral form as

4. Numerics of Solving the Inverse Problem

The solution of the present inverse problem is in fact a problem of finding the roots of (44) when is truncated at . For that reason, we shall focus on such roots in the neighborhood of a step-like discontinuity of at , say, where , . The associated with the Gibbs phenomenon overshoot of or over and undershoot of this below . The intersection of a horizontal line , corresponding to , can be at two points and on both sides of the peak location , representing two roots of . It can also be at one point , corresponding to the maximum of the truncated ; that is, a double root of .

Obviously, , where and are possible points of inflection, defined by The existence of both and or at least one of them, namely, , is however not mandatory. Moreover, for the , , and points we expect to have

Equation (44) can be solved for by a Newton-Raphson (see, e.g., [11, 12]) iterative process which is repeated until a sufficiently accurate value is reached. This process, for each series truncation number , rewrites as or in the stereographic form which is the same as and is started off with some rather free initial value . Obviously here .

The sequence will usually converge, provided that is close enough to the unknown zero and that

4.1. Computation of the Peak Location

Being a maximization problem for , the determination of by the Newton-Raphson process [11] becomes or in stereographic form.

Straightforwardly, this process can be expressed as or in the stereographic form which is the same as It is remarkable how this algorithm is distinctively independent of .

Proposition 8. The ill-posedness of the inverse problem of the truncated Fourier series is invariant under stereographic projection.

Proof. The tangent half-angle stereographic projection, though not isometric (distance preserving), is nevertheless variability preserving. Indeed, the main ingredients of this inverse problem are, namely, , , and , with . Under this projection transforms to , and transforms to . But according to Proposition 4, . Hence the variation of with in the Fourier series inverse problem and the variation of with in the stereographic form are identical.
Equivalently, this makes the ill-posedness of this inverse problem invariant under stereographic projection.

Expected manifestations of the correctness of this proposition are as follows: (i) the results of computations by the (71) and (73) algorithms for with the same must be identical and (ii) the results of computations by the (65) and (67) algorithms for (or ) with the same and must also be identical.

4.2. Convergence Analysis

Let us focus first on the convergence of the iterative process (65) or (67). In many situations where term-by-term differentiability does not hold, the Newton-Raphson method can always be replaced by some quasi-Newton method [12, 13], where or in (65) or (67) is replaced by a suitable approximation. For example, in the chord method is replaced by for all iterations.

If , and is some point between and , it can be easily proved (see [11] or [12]) that

Taking absolute value of both sides gives implying quadratic convergence; that is, is squared (the number of correct digits roughly doubles) at each step. The same arguments apply of course to the equation.

In actual fact for a distinct the convergence of in (65)–(67) is at least quadratic in the neighborhood of this . The process may face a difficulty when or diverges. This is however a special but not a general problem, and (75) always holds if the following conditions are satisfied:(i), , where for some .(ii) is finite .(iii) is sufficiently close to .The term sufficiently close in this context means the following:(a)Taylor approximation is accurate enough such that we can ignore higher order terms.(b) for some .(c) for .Consequently, (75) can be expressed in the following way: where

The initial point has to be chosen such that conditions (i) through (iii) are satisfied, where the third condition requires that . Moreover, with increasing the distance is known to decrease in a way making the satisfaction of these conditions very difficult and actually impossible as . Here it should be underlined that the iterative method (65)–(67) may fail to converge to the root in the following situations:(i)If a point like is included in , the method will terminate due to division by zero.(ii)When the initial estimate is poor, the pertaining wrong can contribute to nonconvergence of the algorithm.

As for convergence of the Newton-Raphson process towards , let , , and is some point between and , to arrive at

The convergence of in (71)–(73) is at least quadratic in the neighborhood of this . The process may face a difficulty, however, when or diverges. Moreover (79) always holds if the following conditions are satisfied:(i), , where for some .(ii) is finite .(iii) is sufficiently close to .Sufficient closeness in this context means that(a)Taylor approximation is accurate enough such that we can ignore higher order terms;(b) for some ;(c) for .Consequently, where

4.3. Applications

The standard iterative algorithms, (71) or (73), for finding the location of overshoots in of the Gibbs effect employ the fact that the first derivative of a peak, in or , respectively, has a downward-going-zero-crossing at the peak maximum. But the presence of many neighboring smaller peaks in the ringing close to the overshoot will cause many undesirable zero-crossings to be invoked if the initial trial root is not correctly chosen. This fact makes the numerical problem of peak finding in the Gibbs effect of pathological difficulty. A difficulty that aggravates with the increase in the number of truncated terms in the employed Fourier series. Indeed, with the increase of the number of undesirable zero-crossings increases in the neighborhood of , as illustrated in Section 1. This can increasingly lead to oscillatory and even unstable values in sequence , the more the distance contains unwanted zero-crossings.

These problems have been resolved in the two examples that follow, which simultaneously illustrate the applicability of Proposition 8.

Example 9. Consider the -periodic odd signal , of Example A.1, defined over one period, by which has, for each truncation number , the truncated Fourier series approximations This signal has a jump discontinuity at of magnitude 2. We shall focus our computations on the Gibbs phenomenon overshoot near this discontinuity.
Here is a printout list of results for the parameters of the posing inverse problem, via computations by both the direct and stereographic algorithms.
It is remarkable how results by both algorithms are practically identical and share the same following features that are observed by results of the direct algorithm, discussed below.

Case I (). The computations of for the previous example, using (71), with lead to ; this varies quickly with (only beyond the 4th decimal) to stabilize towards . Here using (1) gives , which is quite close to the theoretically expected value, by (A.7), of ~1.18.
Then using (65) for and leads to ; this varies quickly with (only beyond the -d decimal) to stabilize quickly towards . Also using (65) for and leads to ; this varies quickly with (beyond the 2nd decimal) to stabilize quickly towards .
For both and , using (1) yields , as expected.

Case II (). Computations of , using (71), with lead to ; this varies quickly with (only beyond the decimal) to stabilize towards . Here using (1) gives . This result appears to be slightly excessive, perhaps due to increased rounding off errors with increasing to 100.
Then using (65) for and leads after a negative value to ; this varies with (beyond the 2nd decimal) to stabilize quickly towards . The used value of is apparently not appropriate, and this illustrates the generic ill-posedness of this inverse problem.
The pertaining computations of , with , lead to ; this varies with (beyond the 3-d decimal) to stabilize quickly towards .
For both and , using (1), yields , as expected.

Case III (). The computations of , using (71), with lead to ; this varies insignificantly with (only beyond the 4th decimal) to stabilize towards . Here using (1) gives , which does not significantly differ from when .
Then using (65) for and leads after a negative value to ; and this varies (beyond the 3-d decimal) to stabilize quickly towards .
The pertaining computations of , with , lead to ; this varies (beyond the decimal) to stabilize towards .
For both and , using (1) yields .
During the previous computations, summarized in Table 1, it has been observed that varying in the present Newton-Raphson process, for all , can lead to strong oscillations and even completely unstable values for , , and of the Fourier series inverse problem, which though ill-posed is practically regularizible. Apart from this drawback, the (71) with (65) iterative algorithms (i) converge quickly to the theoretically expected values and are (ii) of demonstrated stability of the iterations, for properly chosen (which act as effective regularizing parameters). Moreover, computations by the stereographic algorithm, invoking (73) and (67), share the same magnitude and convergence rates with results by the direct algorithm. Finally, it is demonstrated that the present signal exhibits an overshoot, when , of ~9% of a jump discontinuity, of magnitude 2, that is, located at a distance of ~0.0105.
To illustrate robustness of both numerical algorithms and their possible effective regularizibility we shall study also the numerical inverse problem for the Fourier series representation of an arbitrary (nonsymmetric) periodic signal.

Table 1 

Results of computations for Example 9.

Example 10. Consider the -periodic signal , defined over one period, by It has, for each truncation number , the truncated Fourier series approximations This signal has, like the signal of Example 9, a jump discontinuity at of magnitude 2. For the sake of comparison, we shall also focus our computations on the Gibbs phenomenon overshoot part of near this discontinuity.
It should be noted here that the results by both algorithms, summarized in Table 2, are also identical, a feature indicating that the degree of ill-posedness of the posing inverse problem is not affected, as predicted by Proposition 8, by the tangent half-angle stereographic projection.

Table 2 

Results of computations for Example 10.

Case I (). The computations of for this example, using (71) and (73), with lead to an identical sequence that converges by both algorithms after the iteration to with .

Case II (). The computations of , with , lead to an identical sequence that converges by both algorithms after the iteration to with .

Case III (). The computations of , with , lead also to an identical sequence that converges by both algorithms after the iteration to with .

Therefore the present signal exhibits an overshoot when , of ~9% of a jump discontinuity of magnitude 2, that is, located at a distance of ~0.020 9.

Comparison of the computed Gibbs overshoot parameters of the two considered signals with jump discontinuities at the same of the same magnitude illustrates that the relative size of the overshoot depends essentially exclusively on the magnitude of the stimulating jump discontinuity. For instance, based on the present numerical results, one can suggest here that nonsymmetric signals have equal overshoots with symmetric signals. Their respective ’s are however different. Indeed the nonsymmetric signal has a significantly larger () than the () of the symmetric signal.

A further analysis has also been performed of Example 10 with , to compute the corresponding with , that correspond simultaneously to as illustrated in Figure 1. Results of these computations are summarized in Table 3 and then used for the partial plots, of Figure 2, that interpolate , (for ) with , for each . In this table stands for the first convergence step of the sequence of Newton Raphson iterations.

Table 3 

Summary of computations of all overshoot parameters in Example 10.

Figure 1 

Plots of near a discontinuity of in Example 10.

Figure 2 

Plots of the curves near a discontinuity of in Example 10 interpolating , (for ) with , all for .

These results turn out to be extremely stable against variations of the trial initial root . Indeed, in the computationally most difficult situation, namely, when , a change of from 0.008 to 0.020 leads to a variation of from to just with the same asymptotic and (for ).

5. Conclusions

This paper demonstrates that the quadratic degeneracy of the solution of the present inverse problem, explicit when using , is a nonlinear indicator of the existence of Gibbs’ phenomenon. It provides also a distribution-theoretic proof for the existence of this phenomenon. The reported analysis of numerical solvability of this problem illustrates the following:(i)When converge continuously in and are nonzero at (or , there exists a neighborhood (defined by ) of such that, for all starting values in the neighborhood, will converge to (or ).(ii)When converge continuously in and are nonzero at and (but not at ) and if has a second derivative, which might be zero, at and (and at ) then the convergence of is quadratic or faster. If at and , then the convergence of is merely quadratic.(iii)The convergence of the Newton-Raphson iterative process to can also be quadratic, but restrictions on the derivatives of have to be observed up to the third order. Otherwise, the convergence of to can slow down to become only linear [12] with a rate .(iv)The reported computations clearly illustrate that and when ; that is, the infinite Fourier series will sum to in the neighborhood of a discontinuity, except at the discontinuity itself.

Appendix

To illustrate the Gibbs-Wilbraham effect we shall use here the same arguments employed by Wilbraham in 1848 in analyzing a similar square wave.

Example A.1. Let be a square wave periodic signal, which is defined over one period by Its Fourier series representation is This defines the truncated series which may be differentiated as Since then It is possible to cancel the previous artificially induced differentiation by integrating the last result as This integral starts at zero when in agreement with (A.2) and increases until . Then defines a turning point . The integral starts to decrease, at which point the numerator goes negative. For large the denominator remains positive.
The form (A.7), distinctively from (A.3), reveals that has turning points at the zeros of . These occur when ; , that is, at . Clearly then are turning points for . In particular defines the overshoot peak location and is its magnitude. To evaluate numerically, let us change the variable of integration; namely, , or , to rewrite the previous integral as . Clearly when , , and correct to two decimal points.
As for decreases as we move away from the discontinuity, but it should be remarked however that the location of the overshoot moves towards the discontinuity according to

Conflict of Interests

The author declares that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

The numerical computations of Examples 9 and 10 have been performed, using Mathematica 8, by Mounir Abou Yassine, of the Department of CCE of AUL.