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

Authors: Nassar H. S. Haidar
Journal: International Journal of Computational Mathematics
 (2015)

---

## Abstract

We report on the inverse problem for the truncated Fourier series representation off
x
∈
BV
-
L
,
L 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 f
x.


---

## Body

## 1. Introduction

This paper reinvestigates the Fourier series [1, 2] representation(1)
g
x
=
a
0
2
+
∑
m
=
1
∞
a
m
cos
⁡
m
π
L
x
+
b
m
sin
⁡
m
π
L
x
,of a piecewise continuous 2
L –periodic signal f
(
x
)
∈
B
V
(
-
L
,
L
), with(2)
a
0
=
1
L
∫
-
L
L
f
x
d
x
,
a
m
=
1
L
∫
-
L
L
f
x
cos
⁡
m
π
L
x
d
x
,
b
m
=
1
L
∫
-
L
L
f
x
sin
⁡
m
π
L
x
d
x
.The Gibbs phenomenon (see, e.g., [3]) is a statement of the fact that the infinite series g
(
x
) tends to overshoot the positive corner of a discontinuity of f
(
x
) 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 m
=
N. However, according to a theorem by Fejer [1], the infinite series g
(
x
) of a B
V
(
-
L
,
L
) function converges to f
(
x
) at each point x of continuity of f
(
x
).This shortcoming in the infinite Fourier series representation of piecewise continuousf
(
x
) 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) a
m and b
m happen to decay slowly with increasing m, and (iii) the global nature of the approximation of f
(
x
): 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 signalf
(
x
)
∈
B
V
(
-
L
,
L
) 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 Section2 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

Letf
(
x
)
∈
B
V
(
-
L
,
L
) be a 2
L –periodic piecewise continuously differentiable function on (
-
L
,
L
), and let x
s be a point in (
-
L
,
L
) at which f
(
x
) has a discontinuity of the first kind. Consider x
s
-
=
x
s
-
0, x
s
+
=
x
s
+
0 and f
s
-
=
f
(
x
s
-
), and f
s
+
=
f
(
x
s
+
), to assume, without loss of generality, that f
s
+
>
f
s
- and define the discontinuity jump by(3)
Δ
s
=
f
s
+
-
f
s
-
.Moreover,(4)
f
c
x
=
f
x
;
x
≠
x
s
undefined;
x
=
x
s
.Theorem 1.
Let a2
L –periodic f
(
x
)
∈
B
V
(
-
L
,
L
) be discontinuous at x
s. The Fourier series representation g
(
x
) of f
(
x
) converges over (
-
L
,
L
), as N
→
∞ so as,(5)
g
x
=
f
s
-
-
ρ
Δ
s
;
x
=
x
s
-
1
2
f
s
-
+
f
s
+
;
x
=
x
s
f
s
+
+
ρ
Δ
s
;
x
=
x
s
+
f
c
x
;
x
≠
x
s
o
r
x
s
±
,where ρ
=
u
(
0
), an uncertainty.Proof.
According to distribution theory [6] when g
(
x
) coincides a.e. with f
(
x
) over (
-
L
,
L
), it may always be differentiated; namely,(6)
g
′
x
=
f
c
′
x
+
Δ
s
δ
x
-
x
s
,with δ as Dirac’s delta. Cancellation of this differentiation is performed by sweeping the x-axis during a compensating integration from -
L to L. Indeed, integration of (6) first from -
L to x leads to(7)
∫
-
L
x
g
′
τ
-
f
c
′
τ
d
τ
=
Δ
s
∫
-
L
x
δ
τ
-
x
s
d
τ
,and since g
(
-
L
)
=
f
c
(
-
L
), then(8)
g
x
=
f
c
x
+
Δ
s
u
x
-
x
s
,with u as Heaviside’s unit step function.
By a theorem of Fejer [1], g
(
x
) should converge to f
c
(
x
) when N
→
∞, ∀
x
≠
x
s on (
-
L
,
L
). Hence the term Δ
s
u
(
x
-
x
s
) in (8) can have a compact support only on an infinitesimal interval σ
=
lim
N
→
∞
⁡
σ
N
=
0 to the right of x
s, that is, ending at x
s
+. Hence(9)
g
x
=
f
c
x
+
Δ
s
u
x
s
+
-
x
s
.Furthermore, despite the fact that u
(
x
s
+
-
x
s
)
=
1, u
(
x
s
+
-
x
s
)
≈
u
(
0
), which is an uncertainty to be denoted by ρ. Consequently, relation (9) is representable as(10)
g
x
=
f
c
x
;
x
≠
x
s
f
c
x
+
ρ
Δ
s
;
x
=
x
s
+
.
Next complete the sweeping by integrating (6) from x to L; namely,(11)
∫
x
L
g
′
τ
-
f
c
′
τ
d
τ
=
Δ
s
∫
x
L
δ
τ
-
x
s
d
τ
.Here again since g
(
L
)
=
f
c
(
L
), then the left-hand side of (11) becomes -
g
(
x
)
+
f
c
(
x
). The right-hand side requires standardization by means of the substitution τ
=
-
ν, leading to d
τ
=
-
d
ν, [
x
,
L
] becoming [
-
x
,
-
L
] and(12)
∫
x
L
δ
τ
-
x
s
d
τ
=
-
∫
-
x
-
L
δ
-
ν
-
x
s
d
ν
=
-
∫
-
x
-
L
δ
ν
+
x
s
d
ν
=
∫
-
L
-
x
δ
ν
+
x
s
d
ν
=
u
-
x
+
x
s
.It is obvious then that (11) is reducible to(13)
g
x
=
f
c
x
-
Δ
s
u
-
x
+
x
s
.
Repeated application of the same arguments employed in the derivation of (9)-(10) to (13) allows for writing it as(14)
g
x
=
f
c
x
-
Δ
s
u
x
s
-
+
x
s
,and subsequently in the form(15)
g
x
=
f
c
x
;
x
≠
x
s
f
c
x
-
ρ
Δ
s
;
x
=
x
s
-
.Taking into consideration that f
c
(
x
s
+
)
=
f
s
+ and f
c
(
x
s
-
)
=
f
s
- 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 g
(
x
) at x
s
- and x
s
+, respectively, while asserting the existence of a spiky behavior of g
(
x
) individually at x
s
+ then at x
s
-. Moreover, the uncertainty ρ had been resolved in the past (by Wilbraham; see also the Appendix) computationally to satisfy ρ
=
0.09 and long before the advance of the theory of distributions (generalized functions) by S. L. Sobolev in 1936. The number 0.09 for ρ remains until now, however, a mysterious theoretical puzzle that calls for a rigorous distribution-theoretic justification.Clearly, in the neighborhood ofx
s, g
(
x
) can only be a rudimentary approximation to f
(
x
), even when N
→
∞, with an added spike pair (16)
Q
s
x
=
+
ρ
Δ
s
;
x
=
x
s
+
-
ρ
Δ
s
;
x
=
x
s
-
0
;
x
≠
x
sto f
s
± at x
s
±. Relation (16) represents the Gibbs-Wilbraham effect which, due to the fact that x
s
+
-
x
s
-
=
0, is not practically computable. This effect degenerates, however, when N is finite, to an undershoot-overshoot pair Q
s
;
σ
N
(
x
) on an interval I
s
;
σ
N
=
[
x
s
-
σ
N
,
x
s
+
σ
N
], splitted by Δ
s at x
s. In this case, it is always possible to compute(17)
g
x
=
~
f
x
+
Q
s
;
σ
N
x
;
x
∈
I
s
;
σ
N
,
N
<
∞
~
f
x
;
x
∉
I
s
;
σ
N
.Furthermore, since lim
N
→
∞
⁡
Q
s
;
σ
N
(
x
)
=
Q
s
(
x
), then these computations with (17) over I
s
;
σ
N can reveal only a feature Q
s
;
σ
N
(
x
), associated with the Gibbs-Wilbraham effect Q
s
(
x
) but not the effect itself.


## 3. The Inverse Problem for Truncated Fourier Series

Substitution of any value ofx
∈
[
-
L
,
L
] in (1) yields a number c
=
g
(
x
), which may differ [2] from f
(
x
) by Q
s
;
σ
N
(
x
) if x is in the neighborhood of a discontinuity of f
(
x
). Obviously c satisfies, the nonlinear in x, trigonometric series equation(18)
∑
m
=
1
N
G
γ
,
x
,
m
=
∑
m
=
1
N
a
m
cos
⁡
m
π
L
x
+
b
m
sin
⁡
m
π
L
x
=
c
-
a
0
2
,which is, however, linear in the data γ
=
{
a
m
,
b
m
}
m
=
1
∞.Definition 2 (the inverse problem for the truncated Fourier series representation).
Given a numberc
=
g
(
x
) satisfying (18), what is the corresponding x?This inverse problem happens to be a discretized version of the similar deconvolution problem [4] of x from a given number q in a singular first kind integral equation of the form(19)
∫
0
∞
K
γ
,
x
,
τ
d
τ
=
q
,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 x. This situation may further prevail when the operator ∫
0
∞
K
(
γ
,
x
,
τ
)
d
τ 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 h
(
x
)
=
c
>
f
(
x
) intersects with a sharp peak (overshoot) of a truncated g
(
x
). 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 forx) is not directly evident from the structure of G
(
γ
,
x
,
m
), and despite its linearity in the data γ
=
{
a
m
,
b
m
}
m
=
1
∞.


### 3.1. The Tangent Half-Angle Stereographic Projection

The tangent half-angle substitution (the Euler-Weirstrass substitution [9]) tan
⁡
(
θ
/
2
) is widely used in integral calculus for finding antiderivatives of rational expressions of circular function pairs (
cos
⁡
θ
,
sin
⁡
θ
). Geometrically the construction goes like this: draw the unit circle, and let P be the point (
-
1,0
). A line through P (except the vertical line) passes through a point (
cos
⁡
θ
,
sin
⁡
θ
) on the circle and crosses the y-axis at some point y
=
t. Obviously t 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 P. 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 (
cos
⁡
θ
,
sin
⁡
θ
) on the unit circle to slopes t. In fact this is a stereographic projection [9] that expresses these pairs in terms of one variable: tan
⁡
(
θ
/
2
)
=
t.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 θ
m
=
m
(
π
/
L
)
x to the slope t
m
=
tan
⁡
(
θ
m
/
2
) of a line intersecting the unit circle centered around (
0,0
).Indeed asθ
m varies, the point (
cos
⁡
θ
m
,
sin
⁡
θ
m
) winds repeatedly around the unit circle centered at (
0,0
). The point(20)
1
-
t
m
2
1
+
t
m
2
,
2
t
m
1
+
t
m
2
=
cos
⁡
θ
m
,
sin
⁡
θ
mgoes only once around the circle as t
m goes from -
∞ to ∞ and never reaches the point (
-
1,0
), which is approached as a limit when t
m
→
±
∞.Clearly the stereographic projectiont
m
=
tan
⁡
(
θ
m
/
2
) is defined on the entire unit circle except at the projection point P
:
(
-
1,0
). 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(21)
d
θ
m
=
2
d
t
m
1
+
t
m
2
.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 off
(
x
).Lemma 3.
The solutionx of the inverse problem for the truncated Fourier series representation satisfies(22)
∑
m
=
1
N
H
γ
,
x
,
m
=
∑
m
=
1
N
a
m
2
+
b
m
2
-
a
m
tan
⁡
m
π
/
2
L
x
-
b
m
2
a
m
tan
2
⁡
m
π
/
2
L
x
+
1
=
c
-
a
0
2
.Proof.
By studying solvability forx of (18) Let t
m
=
tan
⁡
m
(
π
/
2
L
)
x, and then cos
⁡
m
(
π
/
L
)
x
=
(
1
-
t
m
2
)
/
(
1
+
t
m
2
), and sin
⁡
m
(
π
/
L
)
x
=
2
t
m
/
(
1
+
t
m
2
). Further consideration of these facts in (18) leads to the tangent half-angle stereographically transformed, S
[
]
=
[
]
~, form(23)
-
∑
m
=
1
N
S
G
γ
,
x
,
m
=
∑
m
=
1
N
a
m
t
m
2
-
2
b
m
t
m
-
a
m
t
m
2
+
1
=
a
0
2
-
c
.
The harmonic trinomial in the left-hand side of (23) has two harmonic zeros, ξ
m
±, for each x; namely,(24)
ξ
m
±
=
2
L
m
π
tan
-
1
⁡
t
m
±
,
t
m
±
=
b
m
a
m
±
σ
m
,
σ
m
=
a
m
2
+
b
m
2
a
m
.It is worth noting at this point that these ξ
m
± zeros should, by no means, be confused with the possible roots of the x variable. Substitution then of (24) in (23) leads to(25)
∑
m
=
1
N
a
m
t
m
-
t
m
+
t
m
-
t
m
-
t
m
2
+
1
=
a
0
2
-
c
,which is the same as(26)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
-
σ
m
t
m
-
b
m
/
a
m
+
σ
m
t
m
2
+
1
=
a
0
2
-
c
.
The quadratic degeneracy of the harmonic zerosξ
m
± of t
m can disappear only when σ
m
=
0, for all m. And this happens to be equivalent with(27)
∑
m
=
1
N
a
m
2
+
b
m
2
=
0
.Obviously, satisfaction of (27) by (26) reduces it to the equation(28)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
2
t
m
2
+
1
=
a
0
2
-
c
,which is the same as(29)
∑
m
=
1
N
a
m
tan
⁡
m
π
/
2
L
x
-
b
m
2
a
m
tan
2
⁡
m
π
/
2
L
x
+
1
=
a
0
2
-
c
.
It is interesting to remark that condition (27) is entirely different from(30)
a
0
2
+
∑
m
=
1
N
a
m
2
+
b
m
2
>
0
,which can be derived from the Parseval identity [1] (completeness relation)(31)
a
0
2
+
∑
m
=
1
∞
a
m
2
+
b
m
2
=
1
L
∫
-
L
L
f
2
x
d
x
.Hence a possible nondegenerate solution for x 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 x, 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(32)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
2
-
σ
m
2
t
m
2
+
1
=
a
0
2
-
c
.
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(33)
θ
m
=
θ
m
x
=
m
π
L
x
,in G
(
γ
,
x
,
m
) to rewrite it as(34)
G
γ
,
x
,
m
=
a
m
cos
⁡
θ
m
x
+
b
m
sin
⁡
θ
m
x
.The tangent half-angle (stereographic) transformationt
m
=
tan
⁡
(
θ
m
/
2
) or θ
m
=
2
tan
-
1
⁡
t
m is conceived as(35)
θ
m
x
=
2
t
a
n
-
1
t
m
x
=
Φ
t
m
x
,and subsequently(36)
G
=
G
θ
m
x
=
G
Φ
t
m
x
=
H
t
m
x
.It is clear then that this stereographic projection effectively leads to insertion of an additional bracket into the composition of G as a function of x. This paves the way to the definition(37)
S
G
θ
m
x
=
G
Φ
t
m
x
=
H
t
m
x
=
H
γ
,
x
,
m
.Proposition 4.
The following:(38)
S
d
d
x
∑
m
=
1
N
G
θ
m
x
=
d
d
x
S
∑
m
=
1
N
G
θ
m
x
.always holds.Proof.
SinceN is finite, then the proof of S
(
d
/
d
x
)
G
=
(
d
/
d
x
)
S
G
=
(
d
/
d
x
)
H implies the correctness of this proposition. So let us differentiate G
(
x
)
=
G
[
θ
m
(
x
)
] with respect to x by applying the chain rule; namely,(39)
d
G
d
x
=
d
G
d
Φ
d
Φ
d
t
m
d
t
m
d
x
=
d
G
d
θ
m
d
θ
m
d
t
m
d
t
m
d
x
,which is the same as(40)
d
G
d
x
=
-
a
m
sin
⁡
θ
m
+
b
m
cos
⁡
θ
m
2
1
+
t
m
2
θ
m
2
sec
2
⁡
θ
m
2
.Take then the S transformation of the previous relation to arrive, after some algebra, at(41)
S
d
d
x
G
=
-
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Invoke then(42)
H
=
H
γ
,
x
,
m
=
a
m
+
2
b
m
t
m
-
a
m
t
m
2
1
+
t
m
2
,to differentiate it with respect to x as(43)
d
d
x
H
=
d
d
x
S
G
=
-
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
,which is the required result.Let us rewrite (1) in the form(44)
φ
x
=
g
x
-
c
=
0then try to differentiate it term-by-term as(45)
φ
′
x
=
d
d
x
g
x
=
d
d
x
∑
m
=
1
∞
G
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
-
b
m
cos
⁡
m
π
L
x
+
a
m
sin
⁡
m
π
L
x
=
-
∑
m
=
1
∞
m
π
L
G
∗
γ
,
x
,
m
,where G
∗
(
γ
,
x
,
m
) is the harmonic conjugate of G
(
γ
,
x
,
m
). The trigonometric series g
∗
(
x
)
=
∑
m
=
1
∞
G
∗
(
γ
,
x
,
m
) is called the trigonometric series conjugate (see, e.g., [10]) to g
(
x
).In addition to the fact that the convergence ofg
∗
(
x
) can sometimes become worse than convergence of g
(
x
) (see, e.g., [10]) the presence of the m
(
π
/
L
) factor in (
d
/
d
x
)
g
(
x
) can aggravate the convergence problem of (45) when the coefficients a
m and b
m do not decay fast enough.Next we may stereographically transform(
∂
/
∂
x
)
g
(
x
), using the notation S
(
d
/
d
x
)
g
(
x
)
=
(
d
/
d
x
)
g
(
x
)
~, where(46)
d
d
x
g
x
~
=
-
∑
m
=
1
∞
m
π
L
G
∗
~
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Similarly, for the stereographic Fourier series representation(47)
g
x
~
=
a
0
2
+
∑
m
=
1
∞
H
γ
,
x
,
m
,in which(48)
H
γ
,
x
,
m
=
a
m
+
2
b
m
t
m
-
a
m
t
m
2
1
+
t
m
2
,we may approximate (44) by(49)
ψ
x
=
g
x
~
-
c
=
0
.Differentiate then (49) term-by-term as(50)
ψ
′
x
=
d
d
x
g
x
~
=
d
d
x
∑
m
=
1
∞
H
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
H
∗
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Note here that(51)
H
∗
γ
,
x
,
m
=
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
=
b
m
tan
⁡
m
π
/
2
L
x
-
a
m
2
-
a
m
2
+
b
m
2
b
m
tan
2
⁡
m
π
/
2
L
x
+
1is the harmonic conjugate of H
(
γ
,
x
,
m
) and that(52)
H
∗
γ
,
x
,
m
=
G
∗
~
γ
,
x
,
m
.This rather interesting fact leads to a remark that follows.Remark 5.
Away from points of discontinuity off
(
x
), the stereographically transformed differentiated Fourier series (
d
/
d
x
)
g
(
x
)
~, which defines G
∗
~
(
γ
,
x
,
m
), is the same as the derivative (
d
/
d
x
)
g
(
x
)
~, which defines H
∗
(
γ
,
x
,
m
), of the stereographic Fourier series g
(
x
)
~. Moreover, even at points of discontinuity of f
(
x
)
,
d
/
d
x
g
x
~
=
(
d
/
d
x
)
g
(
x
)
~.Theorem 6 (see [1]).
Letf
(
x
) be (i) continuous on [
-
L
,
L
] such that (ii) f
(
-
L
)
=
f
(
L
), and let (iii) f
′
(
x
) be piecewise continuous over [
-
L
,
L
]; then the corresponding Fourier series g
(
x
) of (1) or (47) is differentiable at each point where (iv) f
′
′
(
x
) 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 the2
π-periodic odd signal f
(
x
), defined over one period, by f
(
x
)
=
x. Its Fourier series is known to be(53)
g
x
=
2
∑
m
=
1
∞
-
1
m
+
1
m
sin
⁡
m
x
,that is, with a
m
=
0, ∀
m, b
m
=
2
(
(
-
1
)
m
+
1
/
m
), and m
(
π
/
L
)
=
m.
Differentiating term-by-term, we obtain(54)
g
′
x
=
2
∑
m
=
1
∞
-
1
m
+
1
cos
⁡
m
x
,which is not uniformly convergent because f
(
x
) is not only discontinuous at x
=
±
π but f
(
-
π
)
≠
f
(
π
) and f
′
′
(
±
π
) does not exist. In particular, at x
=
0, the sum of the series g
′
(
0
)
=
2
∑
m
=
1
∞
(
-
1
)
m
+
1 oscillates between 0 and 1, that is, divergent. It is also well known, nevertheless, that the series ∑
m
=
1
∞
(
-
1
)
m
+
1 is C
-
1 summable [1, 2] to 1
/
2 and that turns out to yield the right answer g
′
(
0
)
=
f
′
(
0
)
=
1. Moreover, at x
=
π, the sum of the series g
′
(
π
)
=
2
∑
m
=
1
∞
(
-
1
)
2
m
+
1
→
-
∞, that is, divergent.
Alternatively, by means of (47)(55)
g
x
~
=
∑
m
=
1
∞
b
m
2
t
m
1
+
t
m
2
=
2
∑
m
=
1
∞
-
1
m
+
1
m
sin
⁡
m
xand by means of (46)(56)
g
′
x
~
=
d
d
x
g
x
~
=
-
∑
m
=
1
∞
m
G
∗
~
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
b
m
t
m
2
-
1
1
+
t
m
2
=
2
∑
m
=
1
∞
-
1
m
+
1
1
-
tan
2
⁡
m
1
/
2
x
1
+
tan
2
⁡
m
1
/
2
x
=
2
∑
m
=
1
∞
-
1
m
+
1
cos
⁡
m
x
.Here, at x
=
π, the sum of the series g
′
(
π
)
~
=
2
∑
m
=
1
∞
(
-
1
)
2
m
+
1
→
-
∞ is also divergent.Note in general that, for even signals, that is, whenb
m
=
0, ∀
m, it follows from (22) that(57)
H
γ
,
x
,
m
b
m
=
0
=
a
m
1
-
t
m
2
1
+
t
m
2
=
a
m
cos
⁡
m
π
L
x
=
G
γ
,
x
,
m
b
m
=
0
.Moreover, for odd signals, that is, when a
m
=
0, ∀
m, the situation is not of the same clarity.Indeed, according to (22),(58)
H
γ
,
x
,
m
a
m
=
0
=
0
0
b
m
1
+
t
m
2
,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(59)
H
γ
,
x
,
m
a
m
=
0
=
b
m
2
t
m
1
+
t
m
2
=
b
m
sin
⁡
m
π
L
x
=
G
γ
,
x
,
m
a
m
=
0
.Also when b
m
=
0, ∀
m,(60)
H
∗
γ
,
x
,
m
b
m
=
0
=
0
0
a
m
1
+
t
m
2
,and this is resolvable in a similar fashion.BothG
(
γ
,
x
,
m
), of (2), and H
(
γ
,
x
,
m
), of (22), equivalent forms are nonlinear trigonometric functions of the x solution. However, while G
(
γ
,
x
,
m
) is linear in γ
=
{
a
m
,
b
m
}
m
=
1
∞, thestereographic H
(
γ
,
x
,
m
) is nonlinear in this γ. This can raise a question on a possibility for (
∂
/
∂
x
)
g
(
x
)
≠
(
∂
/
∂
x
)
g
(
x
)
~ at singular points of f
(
x
). A question that will be answered negatively in the next section. Furthermore, since the data γ
=
{
a
m
,
b
m
}
m
=
1
∞ is discrete, then the operators ∑
m
=
1
∞
G
(
γ
,
x
,
m
) and ∑
m
=
1
∞
H
(
γ
,
x
,
m
) 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 usingH
(
γ
,
x
,
m
), is a nonlinear indicator of the existence of the peak Q
s
;
σ
N
(
x
) 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 g
(
x
)
~ in integral form as(61)
g
x
~
=
a
0
2
-
∑
m
=
1
∞
∫
0
x
m
·
π
L
b
m
tan
⁡
m
π
/
2
L
ξ
-
a
m
2
-
a
m
2
+
b
m
2
b
m
tan
2
⁡
m
π
/
2
L
ξ
+
1
d
ξ
=
a
0
2
-
∑
m
=
1
∞
∫
0
x
m
·
π
L
b
m
tan
2
⁡
m
π
/
2
L
ξ
-
2
a
m
tan
⁡
m
π
/
2
L
ξ
-
b
m
tan
2
⁡
m
π
/
2
L
ξ
+
1
d
ξ
.


## 3.1. The Tangent Half-Angle Stereographic Projection

The tangent half-angle substitution (the Euler-Weirstrass substitution [9]) tan
⁡
(
θ
/
2
) is widely used in integral calculus for finding antiderivatives of rational expressions of circular function pairs (
cos
⁡
θ
,
sin
⁡
θ
). Geometrically the construction goes like this: draw the unit circle, and let P be the point (
-
1,0
). A line through P (except the vertical line) passes through a point (
cos
⁡
θ
,
sin
⁡
θ
) on the circle and crosses the y-axis at some point y
=
t. Obviously t 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 P. 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 (
cos
⁡
θ
,
sin
⁡
θ
) on the unit circle to slopes t. In fact this is a stereographic projection [9] that expresses these pairs in terms of one variable: tan
⁡
(
θ
/
2
)
=
t.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 θ
m
=
m
(
π
/
L
)
x to the slope t
m
=
tan
⁡
(
θ
m
/
2
) of a line intersecting the unit circle centered around (
0,0
).Indeed asθ
m varies, the point (
cos
⁡
θ
m
,
sin
⁡
θ
m
) winds repeatedly around the unit circle centered at (
0,0
). The point(20)
1
-
t
m
2
1
+
t
m
2
,
2
t
m
1
+
t
m
2
=
cos
⁡
θ
m
,
sin
⁡
θ
mgoes only once around the circle as t
m goes from -
∞ to ∞ and never reaches the point (
-
1,0
), which is approached as a limit when t
m
→
±
∞.Clearly the stereographic projectiont
m
=
tan
⁡
(
θ
m
/
2
) is defined on the entire unit circle except at the projection point P
:
(
-
1,0
). 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(21)
d
θ
m
=
2
d
t
m
1
+
t
m
2
.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 off
(
x
).Lemma 3.
The solutionx of the inverse problem for the truncated Fourier series representation satisfies(22)
∑
m
=
1
N
H
γ
,
x
,
m
=
∑
m
=
1
N
a
m
2
+
b
m
2
-
a
m
tan
⁡
m
π
/
2
L
x
-
b
m
2
a
m
tan
2
⁡
m
π
/
2
L
x
+
1
=
c
-
a
0
2
.Proof.
By studying solvability forx of (18) Let t
m
=
tan
⁡
m
(
π
/
2
L
)
x, and then cos
⁡
m
(
π
/
L
)
x
=
(
1
-
t
m
2
)
/
(
1
+
t
m
2
), and sin
⁡
m
(
π
/
L
)
x
=
2
t
m
/
(
1
+
t
m
2
). Further consideration of these facts in (18) leads to the tangent half-angle stereographically transformed, S
[
]
=
[
]
~, form(23)
-
∑
m
=
1
N
S
G
γ
,
x
,
m
=
∑
m
=
1
N
a
m
t
m
2
-
2
b
m
t
m
-
a
m
t
m
2
+
1
=
a
0
2
-
c
.
The harmonic trinomial in the left-hand side of (23) has two harmonic zeros, ξ
m
±, for each x; namely,(24)
ξ
m
±
=
2
L
m
π
tan
-
1
⁡
t
m
±
,
t
m
±
=
b
m
a
m
±
σ
m
,
σ
m
=
a
m
2
+
b
m
2
a
m
.It is worth noting at this point that these ξ
m
± zeros should, by no means, be confused with the possible roots of the x variable. Substitution then of (24) in (23) leads to(25)
∑
m
=
1
N
a
m
t
m
-
t
m
+
t
m
-
t
m
-
t
m
2
+
1
=
a
0
2
-
c
,which is the same as(26)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
-
σ
m
t
m
-
b
m
/
a
m
+
σ
m
t
m
2
+
1
=
a
0
2
-
c
.
The quadratic degeneracy of the harmonic zerosξ
m
± of t
m can disappear only when σ
m
=
0, for all m. And this happens to be equivalent with(27)
∑
m
=
1
N
a
m
2
+
b
m
2
=
0
.Obviously, satisfaction of (27) by (26) reduces it to the equation(28)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
2
t
m
2
+
1
=
a
0
2
-
c
,which is the same as(29)
∑
m
=
1
N
a
m
tan
⁡
m
π
/
2
L
x
-
b
m
2
a
m
tan
2
⁡
m
π
/
2
L
x
+
1
=
a
0
2
-
c
.
It is interesting to remark that condition (27) is entirely different from(30)
a
0
2
+
∑
m
=
1
N
a
m
2
+
b
m
2
>
0
,which can be derived from the Parseval identity [1] (completeness relation)(31)
a
0
2
+
∑
m
=
1
∞
a
m
2
+
b
m
2
=
1
L
∫
-
L
L
f
2
x
d
x
.Hence a possible nondegenerate solution for x 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 x, 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(32)
∑
m
=
1
N
a
m
t
m
-
b
m
/
a
m
2
-
σ
m
2
t
m
2
+
1
=
a
0
2
-
c
.
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(33)
θ
m
=
θ
m
x
=
m
π
L
x
,in G
(
γ
,
x
,
m
) to rewrite it as(34)
G
γ
,
x
,
m
=
a
m
cos
⁡
θ
m
x
+
b
m
sin
⁡
θ
m
x
.The tangent half-angle (stereographic) transformationt
m
=
tan
⁡
(
θ
m
/
2
) or θ
m
=
2
tan
-
1
⁡
t
m is conceived as(35)
θ
m
x
=
2
t
a
n
-
1
t
m
x
=
Φ
t
m
x
,and subsequently(36)
G
=
G
θ
m
x
=
G
Φ
t
m
x
=
H
t
m
x
.It is clear then that this stereographic projection effectively leads to insertion of an additional bracket into the composition of G as a function of x. This paves the way to the definition(37)
S
G
θ
m
x
=
G
Φ
t
m
x
=
H
t
m
x
=
H
γ
,
x
,
m
.Proposition 4.
The following:(38)
S
d
d
x
∑
m
=
1
N
G
θ
m
x
=
d
d
x
S
∑
m
=
1
N
G
θ
m
x
.always holds.Proof.
SinceN is finite, then the proof of S
(
d
/
d
x
)
G
=
(
d
/
d
x
)
S
G
=
(
d
/
d
x
)
H implies the correctness of this proposition. So let us differentiate G
(
x
)
=
G
[
θ
m
(
x
)
] with respect to x by applying the chain rule; namely,(39)
d
G
d
x
=
d
G
d
Φ
d
Φ
d
t
m
d
t
m
d
x
=
d
G
d
θ
m
d
θ
m
d
t
m
d
t
m
d
x
,which is the same as(40)
d
G
d
x
=
-
a
m
sin
⁡
θ
m
+
b
m
cos
⁡
θ
m
2
1
+
t
m
2
θ
m
2
sec
2
⁡
θ
m
2
.Take then the S transformation of the previous relation to arrive, after some algebra, at(41)
S
d
d
x
G
=
-
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Invoke then(42)
H
=
H
γ
,
x
,
m
=
a
m
+
2
b
m
t
m
-
a
m
t
m
2
1
+
t
m
2
,to differentiate it with respect to x as(43)
d
d
x
H
=
d
d
x
S
G
=
-
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
,which is the required result.Let us rewrite (1) in the form(44)
φ
x
=
g
x
-
c
=
0then try to differentiate it term-by-term as(45)
φ
′
x
=
d
d
x
g
x
=
d
d
x
∑
m
=
1
∞
G
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
-
b
m
cos
⁡
m
π
L
x
+
a
m
sin
⁡
m
π
L
x
=
-
∑
m
=
1
∞
m
π
L
G
∗
γ
,
x
,
m
,where G
∗
(
γ
,
x
,
m
) is the harmonic conjugate of G
(
γ
,
x
,
m
). The trigonometric series g
∗
(
x
)
=
∑
m
=
1
∞
G
∗
(
γ
,
x
,
m
) is called the trigonometric series conjugate (see, e.g., [10]) to g
(
x
).In addition to the fact that the convergence ofg
∗
(
x
) can sometimes become worse than convergence of g
(
x
) (see, e.g., [10]) the presence of the m
(
π
/
L
) factor in (
d
/
d
x
)
g
(
x
) can aggravate the convergence problem of (45) when the coefficients a
m and b
m do not decay fast enough.Next we may stereographically transform(
∂
/
∂
x
)
g
(
x
), using the notation S
(
d
/
d
x
)
g
(
x
)
=
(
d
/
d
x
)
g
(
x
)
~, where(46)
d
d
x
g
x
~
=
-
∑
m
=
1
∞
m
π
L
G
∗
~
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Similarly, for the stereographic Fourier series representation(47)
g
x
~
=
a
0
2
+
∑
m
=
1
∞
H
γ
,
x
,
m
,in which(48)
H
γ
,
x
,
m
=
a
m
+
2
b
m
t
m
-
a
m
t
m
2
1
+
t
m
2
,we may approximate (44) by(49)
ψ
x
=
g
x
~
-
c
=
0
.Differentiate then (49) term-by-term as(50)
ψ
′
x
=
d
d
x
g
x
~
=
d
d
x
∑
m
=
1
∞
H
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
H
∗
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
π
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
.Note here that(51)
H
∗
γ
,
x
,
m
=
b
m
t
m
2
+
2
a
m
t
m
-
b
m
1
+
t
m
2
=
b
m
tan
⁡
m
π
/
2
L
x
-
a
m
2
-
a
m
2
+
b
m
2
b
m
tan
2
⁡
m
π
/
2
L
x
+
1is the harmonic conjugate of H
(
γ
,
x
,
m
) and that(52)
H
∗
γ
,
x
,
m
=
G
∗
~
γ
,
x
,
m
.This rather interesting fact leads to a remark that follows.Remark 5.
Away from points of discontinuity off
(
x
), the stereographically transformed differentiated Fourier series (
d
/
d
x
)
g
(
x
)
~, which defines G
∗
~
(
γ
,
x
,
m
), is the same as the derivative (
d
/
d
x
)
g
(
x
)
~, which defines H
∗
(
γ
,
x
,
m
), of the stereographic Fourier series g
(
x
)
~. Moreover, even at points of discontinuity of f
(
x
)
,
d
/
d
x
g
x
~
=
(
d
/
d
x
)
g
(
x
)
~.Theorem 6 (see [1]).
Letf
(
x
) be (i) continuous on [
-
L
,
L
] such that (ii) f
(
-
L
)
=
f
(
L
), and let (iii) f
′
(
x
) be piecewise continuous over [
-
L
,
L
]; then the corresponding Fourier series g
(
x
) of (1) or (47) is differentiable at each point where (iv) f
′
′
(
x
) 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 the2
π-periodic odd signal f
(
x
), defined over one period, by f
(
x
)
=
x. Its Fourier series is known to be(53)
g
x
=
2
∑
m
=
1
∞
-
1
m
+
1
m
sin
⁡
m
x
,that is, with a
m
=
0, ∀
m, b
m
=
2
(
(
-
1
)
m
+
1
/
m
), and m
(
π
/
L
)
=
m.
Differentiating term-by-term, we obtain(54)
g
′
x
=
2
∑
m
=
1
∞
-
1
m
+
1
cos
⁡
m
x
,which is not uniformly convergent because f
(
x
) is not only discontinuous at x
=
±
π but f
(
-
π
)
≠
f
(
π
) and f
′
′
(
±
π
) does not exist. In particular, at x
=
0, the sum of the series g
′
(
0
)
=
2
∑
m
=
1
∞
(
-
1
)
m
+
1 oscillates between 0 and 1, that is, divergent. It is also well known, nevertheless, that the series ∑
m
=
1
∞
(
-
1
)
m
+
1 is C
-
1 summable [1, 2] to 1
/
2 and that turns out to yield the right answer g
′
(
0
)
=
f
′
(
0
)
=
1. Moreover, at x
=
π, the sum of the series g
′
(
π
)
=
2
∑
m
=
1
∞
(
-
1
)
2
m
+
1
→
-
∞, that is, divergent.
Alternatively, by means of (47)(55)
g
x
~
=
∑
m
=
1
∞
b
m
2
t
m
1
+
t
m
2
=
2
∑
m
=
1
∞
-
1
m
+
1
m
sin
⁡
m
xand by means of (46)(56)
g
′
x
~
=
d
d
x
g
x
~
=
-
∑
m
=
1
∞
m
G
∗
~
γ
,
x
,
m
=
-
∑
m
=
1
∞
m
b
m
t
m
2
-
1
1
+
t
m
2
=
2
∑
m
=
1
∞
-
1
m
+
1
1
-
tan
2
⁡
m
1
/
2
x
1
+
tan
2
⁡
m
1
/
2
x
=
2
∑
m
=
1
∞
-
1
m
+
1
cos
⁡
m
x
.Here, at x
=
π, the sum of the series g
′
(
π
)
~
=
2
∑
m
=
1
∞
(
-
1
)
2
m
+
1
→
-
∞ is also divergent.Note in general that, for even signals, that is, whenb
m
=
0, ∀
m, it follows from (22) that(57)
H
γ
,
x
,
m
b
m
=
0
=
a
m
1
-
t
m
2
1
+
t
m
2
=
a
m
cos
⁡
m
π
L
x
=
G
γ
,
x
,
m
b
m
=
0
.Moreover, for odd signals, that is, when a
m
=
0, ∀
m, the situation is not of the same clarity.Indeed, according to (22),(58)
H
γ
,
x
,
m
a
m
=
0
=
0
0
b
m
1
+
t
m
2
,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(59)
H
γ
,
x
,
m
a
m
=
0
=
b
m
2
t
m
1
+
t
m
2
=
b
m
sin
⁡
m
π
L
x
=
G
γ
,
x
,
m
a
m
=
0
.Also when b
m
=
0, ∀
m,(60)
H
∗
γ
,
x
,
m
b
m
=
0
=
0
0
a
m
1
+
t
m
2
,and this is resolvable in a similar fashion.BothG
(
γ
,
x
,
m
), of (2), and H
(
γ
,
x
,
m
), of (22), equivalent forms are nonlinear trigonometric functions of the x solution. However, while G
(
γ
,
x
,
m
) is linear in γ
=
{
a
m
,
b
m
}
m
=
1
∞, thestereographic H
(
γ
,
x
,
m
) is nonlinear in this γ. This can raise a question on a possibility for (
∂
/
∂
x
)
g
(
x
)
≠
(
∂
/
∂
x
)
g
(
x
)
~ at singular points of f
(
x
). A question that will be answered negatively in the next section. Furthermore, since the data γ
=
{
a
m
,
b
m
}
m
=
1
∞ is discrete, then the operators ∑
m
=
1
∞
G
(
γ
,
x
,
m
) and ∑
m
=
1
∞
H
(
γ
,
x
,
m
) 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 usingH
(
γ
,
x
,
m
), is a nonlinear indicator of the existence of the peak Q
s
;
σ
N
(
x
) 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 g
(
x
)
~ in integral form as(61)
g
x
~
=
a
0
2
-
∑
m
=
1
∞
∫
0
x
m
·
π
L
b
m
tan
⁡
m
π
/
2
L
ξ
-
a
m
2
-
a
m
2
+
b
m
2
b
m
tan
2
⁡
m
π
/
2
L
ξ
+
1
d
ξ
=
a
0
2
-
∑
m
=
1
∞
∫
0
x
m
·
π
L
b
m
tan
2
⁡
m
π
/
2
L
ξ
-
2
a
m
tan
⁡
m
π
/
2
L
ξ
-
b
m
tan
2
⁡
m
π
/
2
L
ξ
+
1
d
ξ
.


## 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 g
(
x
) is truncated at m
=
N. For that reason, we shall focus on such roots in the neighborhood of a step-like discontinuity of f
(
x
) at x
=
x
o, say, where f
x
o
-
ε
<
f
(
x
o
+
ε
), 0
<
ε
≪
1. The associated with the Gibbs phenomenon overshoot of g
(
x
)
=
a
0
/
2
+
∑
m
=
1
N
G
(
γ
,
x
,
m
) or a
0
/
2
+
∑
m
=
1
N
H
(
γ
,
x
,
m
) over f
(
x
o
+
ε
) and undershoot of this g
(
x
) below f
(
x
o
-
ε
). The intersection of a horizontal line h
(
x
)
=
c
≥
f
(
x
o
+
ε
), corresponding to φ
(
x
)
=
0, can be at two points α
1 and α
2 on both sides of the peak location α
0, representing two d
i
s
t
i
n
c
t roots of φ
(
x
)
=
0. It can also be at one point α
0, corresponding to the maximum of the truncated g
(
x
); that is, a double root of φ
(
x
)
=
0.Obviously,α
1
<
α
1
~
<
α
0
<
α
2
~
<
α
2, where α
1
~ and α
2
~ are possible points of inflection, defined by(62)
g
′
′
α
1
~
=
φ
′
′
α
1
~
=
g
′
′
α
2
~
=
φ
′
′
α
2
~
=
0
.The existence of both α
1
~ and α
2
~ or at least one of them, namely, α
2
~, is however not mandatory. Moreover, for the α
1, α
0, and α
2 points we expect to have(63)
g
′
α
1
=
φ
′
α
1
>
0
,
φ
′
′
α
1
>
0
,
φ
′
′
′
α
1
<
0
,
g
′
α
0
=
φ
′
α
0
=
0
,
φ
′
′
α
0
<
0
,
φ
′
′
′
α
0
=
0
,
g
′
α
2
=
φ
′
α
2
<
0
,
φ
′
′
α
2
>
0
,
φ
′
′
′
α
2
>
0
.Equation (44) can be solved for x
:
φ
(
x
)
=
0 by a Newton-Raphson (see, e.g., [11, 12]) iterative process(64)
x
n
+
1
=
x
n
-
φ
x
n
φ
′
x
n
,which is repeated until a sufficiently accurate value is reached. This process, for each series truncation number N, rewrites as(65)
x
n
+
1
=
x
n
+
∑
m
=
1
N
G
γ
,
x
n
,
m
+
a
0
/
2
-
c
∑
m
=
1
N
m
π
/
L
G
∗
γ
,
x
n
,
mor in the stereographic form(66)
x
n
+
1
=
x
n
+
∑
m
=
1
N
H
γ
,
x
n
,
m
+
a
0
/
2
-
c
∑
m
=
1
N
m
π
/
L
H
∗
γ
,
x
n
,
m
,which is the same as(67)
x
n
+
1
=
x
n
+
∑
m
=
1
N
a
m
+
2
b
m
t
m
-
a
m
t
m
2
/
1
+
t
m
2
+
a
0
/
2
-
c
∑
m
=
1
N
m
π
/
L
b
m
t
m
2
+
2
a
m
t
m
-
b
m
/
1
+
t
m
2
,and is started off with some rather free initial value x
0. Obviously here t
m
=
tan
⁡
m
(
π
/
2
L
)
x.The sequence{
x
n
} will usually converge, provided that x
0 is close enough to the unknown zero α and that(68)
φ
′
x
0
=
-
∑
m
=
1
∞
m
π
L
G
∗
γ
,
x
0
,
m
≠
0
or
ψ
′
x
0
=
-
∑
m
=
1
∞
m
π
L
H
∗
γ
,
x
0
,
m
≠
0
.


### 4.1. Computation of the Peak Locationα
0

Being a maximization problem forφ
(
x
)
=
0, the determination of α
0 by the Newton-Raphson process [11] becomes(69)
x
n
+
1
=
x
n
-
φ
′
x
n
φ
′
′
x
nor(70)
x
n
+
1
=
x
n
-
ψ
′
x
n
ψ
′
′
x
n
,in stereographic form.Straightforwardly, this process can be expressed as(71)
x
n
+
1
=
x
n
-
∑
m
=
1
N
m
G
∗
γ
,
x
n
,
m
∑
m
=
1
N
m
2
π
/
L
G
γ
,
x
n
,
mor in the stereographic form(72)
x
n
+
1
=
x
n
-
-
∑
m
=
1
N
m
π
/
L
H
∗
γ
,
x
n
,
m
∑
m
=
1
N
m
2
π
2
/
L
2
a
m
t
m
2
-
2
b
m
t
m
-
a
m
/
1
+
t
m
2
,which is the same as(73)
x
n
+
1
=
x
n
+
∑
m
=
1
N
m
b
m
t
m
2
+
2
a
m
t
m
-
b
m
/
1
+
t
m
2
∑
m
=
1
N
m
2
π
/
L
a
m
t
m
2
-
2
b
m
t
m
-
a
m
/
1
+
t
m
2
.It is remarkable how this algorithm is distinctively independent of c.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,γ
=
{
a
m
,
b
m
}
m
=
1
N, c, and x
0, with (
d
/
d
x
)
G
[
θ
m
(
x
)
]
x
0. Under this projection θ
m
(
x
0
) transforms to t
m
(
x
0
)
=
tan
⁡
⁡
[
θ
m
(
x
0
)
/
2
], and G
[
θ
m
(
x
0
)
] transforms to H
[
t
m
(
x
0
)
]. But according to Proposition 4, (
d
/
d
x
)
H
[
t
m
(
x
)
]
x
0
=
(
d
/
d
x
)
G
[
θ
m
(
x
)
]
x
0. Hence the variation of G
(
γ
,
x
0
,
m
) with x
0 in the Fourier series inverse problem and the variation of H
(
γ
,
x
0
,
m
) with x
0 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 α
0 with the same x
0
must be identical and (ii) the results of computations by the (65) and (67) algorithms for α
1 (or α
2) with the same x
0 and c
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 φ
′
(
x
n
) or ψ
′
(
x
n
) in (65) or (67) is replaced by a suitable approximation. For example, in the chord method φ
′
(
x
n
) is replaced by φ
′
(
x
0
) for all iterations.Ifϵ
n
=
α
-
x
n, ϵ
n
+
1
=
α
-
x
n
+
1 and β
n is some point between x
n and α, it can be easily proved (see [11] or [12]) that(74)
ϵ
n
+
1
=
-
φ
′
′
β
n
2
φ
′
x
n
ϵ
n
2
.Taking absolute value of both sides gives(75)
ϵ
n
+
1
=
φ
′
′
β
n
2
φ
′
x
n
ϵ
n
2
,implying quadratic convergence; that is, ϵ
n is squared (the number of correct digits roughly doubles) at each step. The same arguments apply of course to the ψ
(
x
)
=
0 equation.In actual fact for a distinctα the convergence of {
x
n
} in (65)–(67) is at least quadratic in the neighborhood of this α. The process may face a difficulty when φ
′
(
x
n
)
=
-
∑
m
=
1
N
m
(
π
/
L
)
G
∗
(
γ
,
x
n
,
m
) or ψ
′
(
x
n
)
=
-
∑
m
=
1
N
m
(
π
/
L
)
H
∗
(
γ
,
x
n
,
m
) diverges. This is however a special but not a general problem, and (75) always holds if the following conditions are satisfied:(i)
φ
′
(
x
)
≠
0, ∀
x
∈
I, where I
=
[
α
-
η
,
α
+
η
] for some η
≥
(
α
-
x
0
).
(ii)
φ
′
′
(
x
) is finite ∀
x
∈
I.
(iii)
x
0 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)
(
1
/
2
)
φ
′
′
(
x
n
)
/
φ
′
(
x
n
)
<
K
φ
′
′
(
α
)
/
φ
′
(
α
) for some K
<
∞.
(c)
K
φ
′
′
(
α
)
/
φ
′
(
α
)
ϵ
n
<
1 for n
∈
N.Consequently, (75) can be expressed in the following way:(76)
ϵ
n
+
1
=
M
ϵ
n
2
,where(77)
M
=
sup
x
∈
I
⁡
1
2
φ
′
′
x
φ
′
x
.The initial pointx
0 has to be chosen such that conditions (i) through (iii) are satisfied, where the third condition requires that M
ϵ
0
<
1. Moreover, with increasing N the distance α
1
-
α
2 is known to decrease in a way making the satisfaction of these conditions very difficult and actually impossible as N
→
∞. 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α
0 is included in I, the method will terminate due to division by zero.
(ii)
When the initial estimatex
0 is poor, the pertaining wrong I can contribute to nonconvergence of the algorithm.As for convergence of the Newton-Raphson process towardsα
0, let ε
n
=
α
0
-
x
n, ε
n
+
1
=
α
0
-
x
n
+
1, and θ
n is some point between x
n and α
0, to arrive at(78)
ε
n
+
1
=
-
φ
′
′
′
θ
n
2
φ
′
′
x
n
ε
n
2
,
(79)
ε
n
+
1
=
φ
′
′
′
θ
n
2
φ
′
′
x
n
ε
n
2
.The convergence of{
x
n
} in (71)–(73) is at least quadratic in the neighborhood of this α
0. The process may face a difficulty, however, when φ
′
′
(
x
n
)
=
-
∑
m
=
1
∞
m
2
π
/
L
2
G
(
γ
,
x
n
,
m
) or ψ
′
′
(
x
n
) diverges. Moreover (79) always holds if the following conditions are satisfied:(i)
φ
′
′
(
x
)
≠
0, ∀
x
∈
I, where I
=
[
α
0
-
η
,
α
0
+
η
] for some η
≥
(
α
0
-
x
0
).
(ii)
φ
′
′
′
(
x
) is finite ∀
x
∈
I.
(iii)
x
0 is sufficiently close to α
0.Sufficient closeness in this context means that(a)
Taylor approximation is accurate enough such that we can ignore higher order terms;
(b)
(
1
/
2
)
φ
′
′
′
(
x
n
)
/
φ
′
′
(
x
n
)
<
K
φ
′
′
′
(
α
0
)
/
φ
′
′
(
α
0
) for some K
<
∞;
(c)
K
φ
′
′
′
(
α
0
)
/
φ
′
′
(
α
0
)
ε
n
<
1 for n
∈
N.Consequently,(80)
ε
n
+
1
=
M
ε
n
2
,where(81)
M
=
sup
x
∈
I
⁡
1
2
φ
′
′
′
x
φ
′
′
x
.


### 4.3. Applications

The standard iterative algorithms, (71) or (73), for finding the location of overshoots in Q
s
;
σ
N
(
x
) of the Gibbs effect employ the fact that the first derivative of a peak, in g
(
x
) or g
(
x
)
~, 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 x
0 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 N of truncated terms in the employed Fourier series. Indeed, with the increase of N the number of undesirable zero-crossings increases in the neighborhood of α
0, as illustrated in Section 1. This can increasingly lead to oscillatory and even unstable values in sequence {
x
n
}, the more the distance x
0
-
α
0 contains unwanted zero-crossings.These problems have been resolved in the two examples that follow, which simultaneously illustrate the applicability of Proposition8.Example 9.
Consider the2
π-periodic odd signal f
(
x
), of Example A.1, defined over one period, by(82)
f
x
=
-
1
,
-
π
<
x
<
0
1
,
0
≤
x
<
π
,which has, for each truncation number N, the truncated Fourier series approximations(83)
g
x
≃
4
π
∑
m
=
1
N
1
2
m
-
1
s
i
n
2
m
-
1
x
,
g
x
~
≃
8
π
∑
m
=
1
N
b
m
t
m
1
+
t
m
2
=
8
π
∑
m
=
1
N
1
2
m
-
1
t
m
1
+
t
m
2
=
8
π
∑
m
=
1
N
1
2
m
-
1
tan
⁡
2
m
-
1
/
2
x
1
+
tan
⁡
2
m
-
1
/
2
x
2
.
This signal has a jump discontinuity atx
=
0 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 (N
=
40).
The computations ofα
0 for the previous example, using (71), with x
0
=
0.02 lead to x
1
=
0.039
109; this varies quickly with n (only beyond the 4th decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.039
299
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.039
299
)
=
1.179
030, which is quite close to the theoretically expected value, by (A.7), of ~1.18.
Then using (65) for c
=
1.09 and x
0
=
0.000381 leads to x
1
=
0.021
697; this varies quickly with n (only beyond the 3-d decimal) to stabilize quickly towards x
5
=
x
6
=
x
7
=
⋯
=
0.028
368
=
α
1. Also using (65) for c
=
1.09 and x
0
=
0.04 leads to x
1
=
0.067
699; this varies quickly with n (beyond the 2nd decimal) to stabilize quickly towards x
5
=
x
6
=
x
7
=
⋯
=
0.052
873
=
α
2.
For bothα
1 and α
2, using (1) yields g
(
α
1
)
=
g
(
0.028368
)
=
g
(
α
2
)
=
g
(
0.052873
)
=
1.09
=
c, as expected.Case II (N
=
100).
Computations ofα
0, using (71), with x
0
=
0.02 lead to x
1
=
0.014
221; this varies quickly with n (only beyond the 4
t
h decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.015
708
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.015708
)
=
1.178
990.This result appears to be slightly excessive, perhaps due to increased rounding off errors with increasing
N
to 100.
Then using (65) for c
=
1.09 and x
0
=
0.015
237 leads after anegative value to x
2
=
0.014
462; this varies with n (beyond the 2nd decimal) to stabilize quickly towards x
7
=
x
8
=
x
9
=
⋯
=
0.011
3348
=
α
1.The used value of
x
0
is apparently not appropriate,and this illustrates the generic ill-posedness of this inverse problem.
The pertaining computations ofα
2, with x
0
=
0.017000, lead to x
1
=
0.020
434; this varies with n (beyond the 3-d decimal) to stabilize quickly towards x
4
=
x
5
=
x
6
=
⋯
=
0.021
149
=
α
2.
For bothα
1 and α
2, using (1), yields g
(
α
1
)
=
g
(
0.0113348
)
=
g
(
α
2
)
=
g
(
0.021149
)
=
1.09
=
c, as expected.Case III (N
=
150).
The computations ofα
0, using (71), with x
0
=
0.005 lead to x
1
=
0.010
595; this varies insignificantly with n (only beyond the 4th decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.010
472
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.010472
)
=
1.178
980,which does not significantly differ from
g
(
α
0
)
when
N
=
100.
Then using (65) for c
=
1.09 and x
0
=
0.010
158 leads after anegative value to x
2
=
0.009
643; and this varies (beyond the 3-d decimal) to stabilize quickly towards x
7
=
x
8
=
x
9
=
⋯
=
0.007
565
=
α
1.
The pertaining computations ofα
2, with x
0
=
0.010
800, lead to x
1
=
0.025
998; this varies (beyond the 2
n
d decimal) to stabilize towards x
14
=
x
15
=
x
16
=
x
17
=
x
18
=
⋯
=
0.014
099
=
α
2.
For bothα
1 and α
2, using (1) yields g
(
α
1
)
=
g
(
007
565
)
=
g
(
α
2
)
=
g
(
0.014
099
)
=
1.09
=
c.
During the previous computations, summarized in Table1, it has been observed that varying x
0 in the present Newton-Raphson process, for all N, can lead to strong oscillations and even completely unstable values for α
0, α
1, and α
2 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 x
0 (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 N
=
150, 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 Example9.
N
α
1
α
0
α
2
40
0.028 368
0.039 299
0.052 873
100
0.011 335
0.015 705
0.021 149
150
0.007 565
0.010 472
0.014 099Example 10.
Consider the2
π-periodic signal f
(
x
), defined over one period, by(84)
f
x
=
-
x
π
-
1
,
-
π
<
x
<
0
1
,
0
≤
x
<
π
.It has, for each truncation number N, the truncated Fourier series approximations(85)
g
x
≃
1
4
+
∑
m
=
1
N
-
1
m
-
1
π
2
m
2
cos
⁡
m
x
+
2
-
-
1
m
π
m
sin
⁡
m
x
.
g
x
~
≃
1
4
+
∑
m
=
1
N
-
1
m
-
1
/
π
2
m
2
+
2
2
-
-
1
m
/
π
m
tan
⁡
2
m
-
1
/
2
x
-
-
1
m
-
1
/
π
2
m
2
tan
⁡
2
m
-
1
/
2
x
2
1
+
tan
⁡
2
m
-
1
/
2
x
2
.
This signal has, like the signal of Example9, a jump discontinuity at x
=
0 of magnitude 2. For the sake of comparison, we shall also focus our computations on the Gibbs phenomenon overshoot part of Q
s
;
σ
N
(
x
) near this discontinuity.
It should be noted here that the results by both algorithms, summarized in Table2, 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 Example10.
N
α
1
α
0
α
2
40
0.0564765
0.0781391
0.105155
100
0.0226546
0.0313516
0.0422024
150
0.0151124
0.0209153
0.028156Case I (N
=
40).
The computations ofα
0 for this example, using (71) and (73), with x
0
=
0.03 lead to an identical sequence {
x
n
} that converges by both algorithms after the x
4 iteration to 0.078
139 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
76.Case II (N
=
100).
The computations ofα
0, with x
0
=
0.01, lead to an identical sequence {
x
n
} that converges by both algorithms after the x
5 iteration to 0.031
351 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
87.Case III (N
=
150).
The computations ofα
0, with x
0
=
0.007, lead also to an identical sequence {
x
n
} that converges by both algorithms after the x
5 iteration to 0.020
915 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
91.Therefore the present signal exhibits an overshoot whenN
=
150, 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 samex
=
0 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 α
0’s are however different. Indeed the nonsymmetric signal has a significantly larger α
0 (0.020
9) than the α
0 (0.0105) of the symmetric signal.A further analysis has also been performed of Example10 with N
=
100,200,300, to compute the corresponding α
0 with α
1, α
2 that correspond simultaneously to c
=
1.08,1.10,1.12 as illustrated in Figure 1. Results of these computations are summarized in Table 3 and then used for the partial Q
s
;
σ
N
(
x
) plots, of Figure 2, that interpolate α
1, α
2 (for c
=
1.08,1.10,1.12) with α
0, for each N
=
100,200,300. In this table x
∗ stands for the first convergence step of the {
x
n
} sequence of Newton Raphson iterations.Table 3
Summary of computations of all overshoot parameters in Example10.
α
0 results
N
x
0
x
∗
α
0
g
(
α
0
)
100
0.030
x
3
0.031 3516
1.178 87
200
0.015
x
3
0.015 6919
1.178 92
300
0.008
x
3
0.010 4648
1.178 94
α
1 results
N
x
0
x
∗
α
1
c
0.025
x
4
0.022 2027
1.08
100
0.025
x
3
0.023 1344
1.10
0.025
x
3
0.024 2027
1.12
0.010
x
3
0.011 1117
1.08
200
0.010
x
3
0.011 5779
1.10
0.010
x
4
0.012 1123
1.12
0.008
x
3
0.007 410 08
1.08
300
0.008
x
3
0.007 720 95
1.10
0.008
x
2
0.008 077 32
1.12
α
2 results
N
x
0
x
∗
α
2
c
100
0.35
x
4
0.042 9439
1.08
0.35
x
4
0.041 4445
1.10
0.35
x
4
0.039 8513
1.12
200
0.018
x
3
0.021 4960
1.08
0.018
x
3
0.020 7456
1.10
0.018
x
3
0.019 9483
1.12
300
0.012
x
3
0.014 3360
1.08
0.012
x
3
0.013 8356
1.10
0.012
x
4
0.013 3039
1.12Figure 1
Plots ofg
(
x
) near a discontinuity of f
(
x
) in Example 10.Figure 2
Plots of theg
(
x
) curves near a discontinuity of f
(
x
) in Example 10 interpolating α
1, α
2 (for c
=
1.08,1.10,1.12) with α
0, all for N
=
100,200,300.These results turn out to be extremely stable against variations of the trial initial rootx
0. Indeed, in the computationally most difficult situation, namely, when N
=
300, a change of x
0 from 0.008 to 0.020 leads to a variation of x
∗ from x
2 to just x
6 with the same asymptotic α
0
=
0.0104648 and α
1
=
0.00807732 (for c
=
1.12).


## 4.1. Computation of the Peak Locationα
0

Being a maximization problem forφ
(
x
)
=
0, the determination of α
0 by the Newton-Raphson process [11] becomes(69)
x
n
+
1
=
x
n
-
φ
′
x
n
φ
′
′
x
nor(70)
x
n
+
1
=
x
n
-
ψ
′
x
n
ψ
′
′
x
n
,in stereographic form.Straightforwardly, this process can be expressed as(71)
x
n
+
1
=
x
n
-
∑
m
=
1
N
m
G
∗
γ
,
x
n
,
m
∑
m
=
1
N
m
2
π
/
L
G
γ
,
x
n
,
mor in the stereographic form(72)
x
n
+
1
=
x
n
-
-
∑
m
=
1
N
m
π
/
L
H
∗
γ
,
x
n
,
m
∑
m
=
1
N
m
2
π
2
/
L
2
a
m
t
m
2
-
2
b
m
t
m
-
a
m
/
1
+
t
m
2
,which is the same as(73)
x
n
+
1
=
x
n
+
∑
m
=
1
N
m
b
m
t
m
2
+
2
a
m
t
m
-
b
m
/
1
+
t
m
2
∑
m
=
1
N
m
2
π
/
L
a
m
t
m
2
-
2
b
m
t
m
-
a
m
/
1
+
t
m
2
.It is remarkable how this algorithm is distinctively independent of c.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,γ
=
{
a
m
,
b
m
}
m
=
1
N, c, and x
0, with (
d
/
d
x
)
G
[
θ
m
(
x
)
]
x
0. Under this projection θ
m
(
x
0
) transforms to t
m
(
x
0
)
=
tan
⁡
⁡
[
θ
m
(
x
0
)
/
2
], and G
[
θ
m
(
x
0
)
] transforms to H
[
t
m
(
x
0
)
]. But according to Proposition 4, (
d
/
d
x
)
H
[
t
m
(
x
)
]
x
0
=
(
d
/
d
x
)
G
[
θ
m
(
x
)
]
x
0. Hence the variation of G
(
γ
,
x
0
,
m
) with x
0 in the Fourier series inverse problem and the variation of H
(
γ
,
x
0
,
m
) with x
0 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 α
0 with the same x
0
must be identical and (ii) the results of computations by the (65) and (67) algorithms for α
1 (or α
2) with the same x
0 and c
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 φ
′
(
x
n
) or ψ
′
(
x
n
) in (65) or (67) is replaced by a suitable approximation. For example, in the chord method φ
′
(
x
n
) is replaced by φ
′
(
x
0
) for all iterations.Ifϵ
n
=
α
-
x
n, ϵ
n
+
1
=
α
-
x
n
+
1 and β
n is some point between x
n and α, it can be easily proved (see [11] or [12]) that(74)
ϵ
n
+
1
=
-
φ
′
′
β
n
2
φ
′
x
n
ϵ
n
2
.Taking absolute value of both sides gives(75)
ϵ
n
+
1
=
φ
′
′
β
n
2
φ
′
x
n
ϵ
n
2
,implying quadratic convergence; that is, ϵ
n is squared (the number of correct digits roughly doubles) at each step. The same arguments apply of course to the ψ
(
x
)
=
0 equation.In actual fact for a distinctα the convergence of {
x
n
} in (65)–(67) is at least quadratic in the neighborhood of this α. The process may face a difficulty when φ
′
(
x
n
)
=
-
∑
m
=
1
N
m
(
π
/
L
)
G
∗
(
γ
,
x
n
,
m
) or ψ
′
(
x
n
)
=
-
∑
m
=
1
N
m
(
π
/
L
)
H
∗
(
γ
,
x
n
,
m
) diverges. This is however a special but not a general problem, and (75) always holds if the following conditions are satisfied:(i)
φ
′
(
x
)
≠
0, ∀
x
∈
I, where I
=
[
α
-
η
,
α
+
η
] for some η
≥
(
α
-
x
0
).
(ii)
φ
′
′
(
x
) is finite ∀
x
∈
I.
(iii)
x
0 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)
(
1
/
2
)
φ
′
′
(
x
n
)
/
φ
′
(
x
n
)
<
K
φ
′
′
(
α
)
/
φ
′
(
α
) for some K
<
∞.
(c)
K
φ
′
′
(
α
)
/
φ
′
(
α
)
ϵ
n
<
1 for n
∈
N.Consequently, (75) can be expressed in the following way:(76)
ϵ
n
+
1
=
M
ϵ
n
2
,where(77)
M
=
sup
x
∈
I
⁡
1
2
φ
′
′
x
φ
′
x
.The initial pointx
0 has to be chosen such that conditions (i) through (iii) are satisfied, where the third condition requires that M
ϵ
0
<
1. Moreover, with increasing N the distance α
1
-
α
2 is known to decrease in a way making the satisfaction of these conditions very difficult and actually impossible as N
→
∞. 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α
0 is included in I, the method will terminate due to division by zero.
(ii)
When the initial estimatex
0 is poor, the pertaining wrong I can contribute to nonconvergence of the algorithm.As for convergence of the Newton-Raphson process towardsα
0, let ε
n
=
α
0
-
x
n, ε
n
+
1
=
α
0
-
x
n
+
1, and θ
n is some point between x
n and α
0, to arrive at(78)
ε
n
+
1
=
-
φ
′
′
′
θ
n
2
φ
′
′
x
n
ε
n
2
,
(79)
ε
n
+
1
=
φ
′
′
′
θ
n
2
φ
′
′
x
n
ε
n
2
.The convergence of{
x
n
} in (71)–(73) is at least quadratic in the neighborhood of this α
0. The process may face a difficulty, however, when φ
′
′
(
x
n
)
=
-
∑
m
=
1
∞
m
2
π
/
L
2
G
(
γ
,
x
n
,
m
) or ψ
′
′
(
x
n
) diverges. Moreover (79) always holds if the following conditions are satisfied:(i)
φ
′
′
(
x
)
≠
0, ∀
x
∈
I, where I
=
[
α
0
-
η
,
α
0
+
η
] for some η
≥
(
α
0
-
x
0
).
(ii)
φ
′
′
′
(
x
) is finite ∀
x
∈
I.
(iii)
x
0 is sufficiently close to α
0.Sufficient closeness in this context means that(a)
Taylor approximation is accurate enough such that we can ignore higher order terms;
(b)
(
1
/
2
)
φ
′
′
′
(
x
n
)
/
φ
′
′
(
x
n
)
<
K
φ
′
′
′
(
α
0
)
/
φ
′
′
(
α
0
) for some K
<
∞;
(c)
K
φ
′
′
′
(
α
0
)
/
φ
′
′
(
α
0
)
ε
n
<
1 for n
∈
N.Consequently,(80)
ε
n
+
1
=
M
ε
n
2
,where(81)
M
=
sup
x
∈
I
⁡
1
2
φ
′
′
′
x
φ
′
′
x
.


## 4.3. Applications

The standard iterative algorithms, (71) or (73), for finding the location of overshoots in Q
s
;
σ
N
(
x
) of the Gibbs effect employ the fact that the first derivative of a peak, in g
(
x
) or g
(
x
)
~, 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 x
0 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 N of truncated terms in the employed Fourier series. Indeed, with the increase of N the number of undesirable zero-crossings increases in the neighborhood of α
0, as illustrated in Section 1. This can increasingly lead to oscillatory and even unstable values in sequence {
x
n
}, the more the distance x
0
-
α
0 contains unwanted zero-crossings.These problems have been resolved in the two examples that follow, which simultaneously illustrate the applicability of Proposition8.Example 9.
Consider the2
π-periodic odd signal f
(
x
), of Example A.1, defined over one period, by(82)
f
x
=
-
1
,
-
π
<
x
<
0
1
,
0
≤
x
<
π
,which has, for each truncation number N, the truncated Fourier series approximations(83)
g
x
≃
4
π
∑
m
=
1
N
1
2
m
-
1
s
i
n
2
m
-
1
x
,
g
x
~
≃
8
π
∑
m
=
1
N
b
m
t
m
1
+
t
m
2
=
8
π
∑
m
=
1
N
1
2
m
-
1
t
m
1
+
t
m
2
=
8
π
∑
m
=
1
N
1
2
m
-
1
tan
⁡
2
m
-
1
/
2
x
1
+
tan
⁡
2
m
-
1
/
2
x
2
.
This signal has a jump discontinuity atx
=
0 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 (N
=
40).
The computations ofα
0 for the previous example, using (71), with x
0
=
0.02 lead to x
1
=
0.039
109; this varies quickly with n (only beyond the 4th decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.039
299
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.039
299
)
=
1.179
030, which is quite close to the theoretically expected value, by (A.7), of ~1.18.
Then using (65) for c
=
1.09 and x
0
=
0.000381 leads to x
1
=
0.021
697; this varies quickly with n (only beyond the 3-d decimal) to stabilize quickly towards x
5
=
x
6
=
x
7
=
⋯
=
0.028
368
=
α
1. Also using (65) for c
=
1.09 and x
0
=
0.04 leads to x
1
=
0.067
699; this varies quickly with n (beyond the 2nd decimal) to stabilize quickly towards x
5
=
x
6
=
x
7
=
⋯
=
0.052
873
=
α
2.
For bothα
1 and α
2, using (1) yields g
(
α
1
)
=
g
(
0.028368
)
=
g
(
α
2
)
=
g
(
0.052873
)
=
1.09
=
c, as expected.Case II (N
=
100).
Computations ofα
0, using (71), with x
0
=
0.02 lead to x
1
=
0.014
221; this varies quickly with n (only beyond the 4
t
h decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.015
708
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.015708
)
=
1.178
990.This result appears to be slightly excessive, perhaps due to increased rounding off errors with increasing
N
to 100.
Then using (65) for c
=
1.09 and x
0
=
0.015
237 leads after anegative value to x
2
=
0.014
462; this varies with n (beyond the 2nd decimal) to stabilize quickly towards x
7
=
x
8
=
x
9
=
⋯
=
0.011
3348
=
α
1.The used value of
x
0
is apparently not appropriate,and this illustrates the generic ill-posedness of this inverse problem.
The pertaining computations ofα
2, with x
0
=
0.017000, lead to x
1
=
0.020
434; this varies with n (beyond the 3-d decimal) to stabilize quickly towards x
4
=
x
5
=
x
6
=
⋯
=
0.021
149
=
α
2.
For bothα
1 and α
2, using (1), yields g
(
α
1
)
=
g
(
0.0113348
)
=
g
(
α
2
)
=
g
(
0.021149
)
=
1.09
=
c, as expected.Case III (N
=
150).
The computations ofα
0, using (71), with x
0
=
0.005 lead to x
1
=
0.010
595; this varies insignificantly with n (only beyond the 4th decimal) to stabilize towards x
3
=
x
4
=
x
5
=
⋯
=
0.010
472
=
α
0. Here using (1) gives g
(
α
0
)
=
g
(
0.010472
)
=
1.178
980,which does not significantly differ from
g
(
α
0
)
when
N
=
100.
Then using (65) for c
=
1.09 and x
0
=
0.010
158 leads after anegative value to x
2
=
0.009
643; and this varies (beyond the 3-d decimal) to stabilize quickly towards x
7
=
x
8
=
x
9
=
⋯
=
0.007
565
=
α
1.
The pertaining computations ofα
2, with x
0
=
0.010
800, lead to x
1
=
0.025
998; this varies (beyond the 2
n
d decimal) to stabilize towards x
14
=
x
15
=
x
16
=
x
17
=
x
18
=
⋯
=
0.014
099
=
α
2.
For bothα
1 and α
2, using (1) yields g
(
α
1
)
=
g
(
007
565
)
=
g
(
α
2
)
=
g
(
0.014
099
)
=
1.09
=
c.
During the previous computations, summarized in Table1, it has been observed that varying x
0 in the present Newton-Raphson process, for all N, can lead to strong oscillations and even completely unstable values for α
0, α
1, and α
2 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 x
0 (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 N
=
150, 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 Example9.
N
α
1
α
0
α
2
40
0.028 368
0.039 299
0.052 873
100
0.011 335
0.015 705
0.021 149
150
0.007 565
0.010 472
0.014 099Example 10.
Consider the2
π-periodic signal f
(
x
), defined over one period, by(84)
f
x
=
-
x
π
-
1
,
-
π
<
x
<
0
1
,
0
≤
x
<
π
.It has, for each truncation number N, the truncated Fourier series approximations(85)
g
x
≃
1
4
+
∑
m
=
1
N
-
1
m
-
1
π
2
m
2
cos
⁡
m
x
+
2
-
-
1
m
π
m
sin
⁡
m
x
.
g
x
~
≃
1
4
+
∑
m
=
1
N
-
1
m
-
1
/
π
2
m
2
+
2
2
-
-
1
m
/
π
m
tan
⁡
2
m
-
1
/
2
x
-
-
1
m
-
1
/
π
2
m
2
tan
⁡
2
m
-
1
/
2
x
2
1
+
tan
⁡
2
m
-
1
/
2
x
2
.
This signal has, like the signal of Example9, a jump discontinuity at x
=
0 of magnitude 2. For the sake of comparison, we shall also focus our computations on the Gibbs phenomenon overshoot part of Q
s
;
σ
N
(
x
) near this discontinuity.
It should be noted here that the results by both algorithms, summarized in Table2, 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 Example10.
N
α
1
α
0
α
2
40
0.0564765
0.0781391
0.105155
100
0.0226546
0.0313516
0.0422024
150
0.0151124
0.0209153
0.028156Case I (N
=
40).
The computations ofα
0 for this example, using (71) and (73), with x
0
=
0.03 lead to an identical sequence {
x
n
} that converges by both algorithms after the x
4 iteration to 0.078
139 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
76.Case II (N
=
100).
The computations ofα
0, with x
0
=
0.01, lead to an identical sequence {
x
n
} that converges by both algorithms after the x
5 iteration to 0.031
351 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
87.Case III (N
=
150).
The computations ofα
0, with x
0
=
0.007, lead also to an identical sequence {
x
n
} that converges by both algorithms after the x
5 iteration to 0.020
915 with g
(
α
0
)
=
g
(
α
0
)
~
=
1.178
91.Therefore the present signal exhibits an overshoot whenN
=
150, 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 samex
=
0 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 α
0’s are however different. Indeed the nonsymmetric signal has a significantly larger α
0 (0.020
9) than the α
0 (0.0105) of the symmetric signal.A further analysis has also been performed of Example10 with N
=
100,200,300, to compute the corresponding α
0 with α
1, α
2 that correspond simultaneously to c
=
1.08,1.10,1.12 as illustrated in Figure 1. Results of these computations are summarized in Table 3 and then used for the partial Q
s
;
σ
N
(
x
) plots, of Figure 2, that interpolate α
1, α
2 (for c
=
1.08,1.10,1.12) with α
0, for each N
=
100,200,300. In this table x
∗ stands for the first convergence step of the {
x
n
} sequence of Newton Raphson iterations.Table 3
Summary of computations of all overshoot parameters in Example10.
α
0 results
N
x
0
x
∗
α
0
g
(
α
0
)
100
0.030
x
3
0.031 3516
1.178 87
200
0.015
x
3
0.015 6919
1.178 92
300
0.008
x
3
0.010 4648
1.178 94
α
1 results
N
x
0
x
∗
α
1
c
0.025
x
4
0.022 2027
1.08
100
0.025
x
3
0.023 1344
1.10
0.025
x
3
0.024 2027
1.12
0.010
x
3
0.011 1117
1.08
200
0.010
x
3
0.011 5779
1.10
0.010
x
4
0.012 1123
1.12
0.008
x
3
0.007 410 08
1.08
300
0.008
x
3
0.007 720 95
1.10
0.008
x
2
0.008 077 32
1.12
α
2 results
N
x
0
x
∗
α
2
c
100
0.35
x
4
0.042 9439
1.08
0.35
x
4
0.041 4445
1.10
0.35
x
4
0.039 8513
1.12
200
0.018
x
3
0.021 4960
1.08
0.018
x
3
0.020 7456
1.10
0.018
x
3
0.019 9483
1.12
300
0.012
x
3
0.014 3360
1.08
0.012
x
3
0.013 8356
1.10
0.012
x
4
0.013 3039
1.12Figure 1
Plots ofg
(
x
) near a discontinuity of f
(
x
) in Example 10.Figure 2
Plots of theg
(
x
) curves near a discontinuity of f
(
x
) in Example 10 interpolating α
1, α
2 (for c
=
1.08,1.10,1.12) with α
0, all for N
=
100,200,300.These results turn out to be extremely stable against variations of the trial initial rootx
0. Indeed, in the computationally most difficult situation, namely, when N
=
300, a change of x
0 from 0.008 to 0.020 leads to a variation of x
∗ from x
2 to just x
6 with the same asymptotic α
0
=
0.0104648 and α
1
=
0.00807732 (for c
=
1.12).


## 5. Conclusions

This paper demonstrates that the quadratic degeneracy of the solution of the present inverse problem, explicit when usingH
(
γ
,
x
,
m
), 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∑
m
=
1
∞
m
(
π
/
L
)
G
∗
(
γ
,
x
,
m
)
=
∑
m
=
1
∞
m
(
π
/
L
)
H
∗
(
γ
,
x
,
m
) converge continuously in x and are nonzero at α
1 (or α
2
), there exists a neighborhood (defined by I) of α such that, for all starting values x
0 in the neighborhood, {
x
n
} will converge to α
1 (or α
2).
(ii)
When∑
m
=
1
∞
m
(
π
/
L
)
G
∗
(
γ
,
x
,
m
)
=
∑
m
=
1
∞
m
(
π
/
L
)
H
∗
(
γ
,
x
,
m
) converge continuously in x and are nonzero at α
1 and α
2 (but not at α
0) and if φ
(
x
) has a second derivative, which might be zero, at α
1 and α
2 (and at α
0) then the convergence of {
x
n
} is quadratic or faster. If φ
′
′
(
x
)
≠
0 at α
1 and α
2, then the convergence of {
x
n
} is merely quadratic.
(iii)
The convergence of the Newton-Raphson iterative process toα
0 can also be quadratic, but restrictions on the derivatives of φ
(
x
) have to be observed up to the third order. Otherwise, the convergence of {
x
n
} to α
0 can slow down to become only linear [12] with a rate l
o
g
10
2.
(iv)
The reported computations clearly illustrate thatα
0
→
0 and α
1
-
α
2
→
0 when N
→
∞; that is, the infinite Fourier series will sum to f
(
x
) in the neighborhood of a discontinuity, except at the discontinuity itself.



---

*Source: 212860-2015-06-14.xml*