# Chebyshev-Fourier Spectral Methods for Nonperiodic Boundary Value Problems

Authors: Bojan Orel; Andrej Perne
Journal: Journal of Applied Mathematics
 (2014)

---

## Abstract

A new class of spectral methods for solving two-point boundary value problems for linear ordinary differential equations is presented in the paper. Although these methods are based on trigonometric functions, they can be used for solving periodic as well as nonperiodic problems. Instead of using basis functions periodic on a given interval−
1,1, we use functions periodic on a wider interval. The numerical solution of the given problem is sought in terms of the half-range Chebyshev-Fourier (HCF) series, a reorganization of the classical Fourier series using half-range Chebyshev polynomials of the first and second kind which were first introduced by Huybrechs (2010) and further analyzed by Orel and Perne (2012). The numerical solution is constructed as a HCF series via differentiation and multiplication matrices. Moreover, the construction of the method, error analysis, convergence results, and some numerical examples are presented in the paper. The decay of the maximal absolute error according to the truncation number N for the new class of Chebyshev-Fourier-collocation (CFC) methods is compared to the decay of the error for the standard class of Chebyshev-collocation (CC) methods.


---

## Body

## 1. Introduction and Formulation

A standard problem in approximation theory is to compute the coefficients of a Fourier series to approximate smooth and periodic functions. This can be efficiently done by using the FFT, which is a stable and well-understood method that yields spectral convergence. Things look very different when dealing with nonperiodic or nonsmooth functions. This is due to the so called Gibbs phenomenon which causes oscillations near the points of discontinuity and/or near the boundary as well as slow decay of Fourier coefficients. There are several possibilities to overcome these difficulties and successfully deal with such functions. Some of them were analyzed by Gottlieb and Shu in [1] and by Tadmor in [2]. One possibility is to use some periodizing transformation and compute the Fourier series of the transformed function. There is a widely used transformation which yields Chebyshev polynomials of the first kind. Much about this can be found in the literature by Boyd in [3], by Fornberg and Sloan in [4], or by Trefethen in [5]. These polynomials are arranged as a Chebyshev series.Another approach was recently presented by Huybrechs in [6] where he analyzed the problem which was stated by Boyd in [7] and by Bruno et al. in [8].Problem 1.
ForT
>
1, let G
n be the space of 2
T-periodic functions of the form
(1)
g
∈
G
n
:
g
(
x
)
=
a
0
2
+
∑
k
=
1
n
‍
(
a
k
cos
⁡
π
k
x
T
+
b
k
sin
π
k
x
T
)
.
The Fourier extension of f
∈
L
2
(
-
1,1
), defined on the interval [
-
1,1
] to the interval [
-
T
,
T
], is the solution to the optimization problem
(2)
g
n
:
=
arg
⁡
min
⁡
g
∈
G
n
∥
f
-
g
∥
L
2
.In his paper, Huybrechs considered a square-integrable functionf
∈
L
2
(
-
1,1
) that is not necessarily smooth or periodic. The idea to obtain a spectrally accurate Fourier series is to extend the given function f to a function g that is periodic on a wider interval [
-
T
,
T
], T
>
1. The Fourier series of the constructed function is obviously point-wise convergent to f on the interval [
-
1,1
].Huybrechs analyzed this problem for the choiceT
=
2 and developed three numerical methods for solving it. Besides proving the existence and uniqueness of the solution, he characterized the solution with two nonclassical families of orthogonal polynomials related to Chebyshev polynomials of the first and second kind: the half-range Chebyshev polynomials of the first and second kind. The second and third methods are based on projection and collocation. The convergence is in most cases exponential rather than superalgebraic.The half-range Chebyshev polynomials of the first and second kind and the corresponding half-range Chebyshev-Fourier (HCF) series were further explored by Orel and Perne in [9]. The authors presented an efficient method for the construction of these polynomials using the modified Chebyshev algorithm for the computation of the recursion coefficients via the three-term recurrence relation formula. Moreover, the authors developed some necessary tools for the construction of spectral methods using these polynomials.In this paper, instead of approximation problems, like Problem1, which was successfully solved by Huybrechs in [6], we are interested in solving BVPs in ODEs via spectral methods, that is, solving problems as defined below.Problem 2.
Find the numerical solution of a linear two-point boundary value problem of the form(3)
L
y
(
x
)
=
f
(
x
)
,
x
∈
[
-
1,1
]
,
with boundary conditions
(4)
B
y
(
x
)
=
0
,
x
∈
{
-
1,1
}
,
where L is a linear differential operator:
(5)
L
=
α
(
x
)
d
2
d
x
2
+
β
(
x
)
d
d
x
+
γ
(
x
)
I
,
I is the identity operator, and B is a set of linear boundary differential operators.There are several different classes of numerical methods to solve boundary value problems. In the company of such methods as finite differences (FDM) and finite elements (FEM) we are interested in spectral methods (SM). It is well known that spectral methods approximate the solution in a finite dimensional subspace of a Hilbert space. The basis functions used are defined globally (on the whole interval). On the contrary, the basis functions for finite element methods are defined locally (only on a small interval).Different approaches are used if the underlying problem is periodic or nonperiodic. In the case of periodic problems, the natural basis functions are trigonometric functions. In other words, the solution is approximated with a truncated Fourier series. Equidistant mesh points are used to discretize the interval on which the problem is defined. In the case of nonperiodic problems, orthogonal polynomials are used, especially the Chebyshev polynomials of the first kind. In this case, the solution is approximated with a truncated Chebyshev series. There are two types of Chebyshev mesh points. Chebyshev points of the second kind, where the Chebyshev polynomials of the first kind reach their extreme values, are generally used in spectral methods to discretize the underlying interval. Chebyshev points of the first kind are the zeros of the Chebyshev polynomials of the first kind. Usually, Chebyshev points in nonperiodic problems outweigh equidistant points because they are denser near the boundary of the interval. Besides, such a distribution of points overwhelms problems caused by Gibbs and/or Runge phenomenon. There are different types of spectral methods depending on the method used for computing the expansion coefficients, for example, Galerkin method, Tau method, or collocation. Spectral methods based on collocation are usually called pseudospectral methods. Classical references on spectral methods include textbooks by Boyd, [3], Canuto et al. [10], Fornberg [11], Gottlieb and Orszag [12], and Trefethen [5], and more recent ones are Canuto et al. [13] and Shen et al. [14].There were several attempts to solve nonperiodic problems using a trigonometric basis. Adcock in [15, 16] solved the problem with modified Fourier series, using Galerkin method to compute expansion coefficients. Huybrechs in [6] proposed a new set of trigonometric functions, which includes sines and cosines as well as half-sines and half-cosines.Our intention is to provide a new class of spectral methods using approaches described in [6, 9], combined with collocation to compute expansion coefficients. This approach yields pseudospectral method for solving nonperiodic problems with tools used for solving periodic problems. We focus on Dirichlet boundary conditions, although it is not difficult to extend this approach to Neumann or mixed boundary conditions. The generalization to higher order linear boundary problems is also possible. The restriction to the interval [
-
1,1
] is only a matter of simplification, since it is well known how to map an arbitrary interval [
a
,
b
] to the unit interval.The paper is organized as follows. In Section2 we define the form of the exact solution to the discrete problem by introducing the set of basis functions proposed by Huybrechs in [6] and two nonclassical families of orthogonal polynomials. In Section 3 follows the developement of a spectral method for the solution of the boundary value problem, based on a pseudospectral (collocation) approach. Error analysis and convergence results based on error analysis approach in [13, 14] are addressed in Section 4. Finally, we present some numerical examples in Section 5. Section 6 concludes the paper.


## 2. Approximation with Half-Range Chebyshev Polynomials

Let us first review the exact solution of the discrete Problem1 using the half-range Chebyshev polynomials. The idea is to extend a nonperiodic function on the interval [
-
1,1
] to an interval [
-
T
,
T
] on which the given function is periodic and to use the set of trigonometric functions, proposed by Huybrechs in [6]. He analyzed the extension for T
=
2 and proposed to use the set of basis functions
(6)
D
n
:
=
C
n
∪
S
n
,
where
(7)
C
n
=
{
1
2
}
∪
{
cos
⁡
π
k
x
2
}
k
=
1
n
,
S
n
=
{
sin
π
k
x
2
}
k
=
1
n
.
Note that the set C
n consists of even and the set S
n of odd functions. The function space is the span of D
n.Huybrechs showed thatD
∞ is not a basis of L
2
(
-
1,1
) but a tight frame with frame bound 2, that is, the frame obeys a generalized Parseval’s identity, and that the set D
∞ consists of all eigenfunctions of the Laplace operator on [
-
1,1
] subject to either homogeneous Dirichlet or Neumann boundary conditions. That is due to the fact that the functions in D
∞ are linearly dependent. On the contrary, the set D
n is a basis for a finite dimensional subspace of L
2
(
-
1,1
), for any finite n. All the functions in D
n are linearly independent.It makes perfect sense to look for an orthonormal basis on the interval[
-
1,1
]. Since the even functions in C
n and the odd functions in S
n are mutually orthogonal, the orthonormalization problem divides into two problems. Let us define two spaces, denoted by
(8)
C
n
:
=
span
⁡
C
n
,
S
n
:
=
span
⁡
S
n
.
The following two theorems are stated and proved by Huybrechs in [6].Theorem 3.
LetT
k
h
(
y
) be the unique normalized sequence of orthogonal polynomials satisfying
(9)
4
π
∫
0
1
‍
T
k
h
(
y
)
y
l
1
-
y
2
d
y
=
δ
k
,
l
,
l
=
0
,
…
,
k
-
1
.
Then the set {
T
k
h
(
cos
⁡
⁡
(
π
x
/
2
)
)
}
k
=
0
n is an orthonormal basis for C
n on [
-
1,1
].Theorem 4.
LetU
k
h
(
y
) be the unique normalized sequence of orthogonal polynomials satisfying
(10)
4
π
∫
0
1
‍
U
k
h
(
y
)
y
l
1
-
y
2
d
y
=
δ
k
,
l
,
l
=
0
,
…
,
k
-
1
.
Then the set {
U
k
h
(
cos
⁡
⁡
(
π
x
/
2
)
)
sin
⁡
(
π
x
/
2
)
}
k
=
0
n
-
1 is an orthonormal basis for S
n on [
-
1,1
].The polynomialsT
n
h
(
x
) and U
n
h
(
x
) are calledhalf-range Chebyshev polynomials of the first and second kind, respectively. They have the same weight functions as the Chebyshev polynomials of the first and second kind but are orthogonal on the interval [
0,1
] rather than on the interval [
-
1,1
]. The orthogonal polynomials are guaranteed to exist, because the weight functions are positive and integrable. The construction and additional properties of these polynomials were studied by Orel and Perne in [9]. An arbitrary function f
∈
L
2
(
-
1,1
) can be then expanded as ahalf-range Chebyshev-Fourier series:
(11)
f
(
x
)
=
∑
k
=
0
∞
‍
a
k
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
∞
‍
b
k
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
.Huybrechs in [6] proved the existence and uniqueness of the exact solution to Problem 1 as a truncated half-range Chebyshev-Fourier (HCF) series.Theorem 5.
For a givenf
∈
L
2
(
-
1,1
), the solution to Problem 1 is
(12)
g
n
(
x
)
=
∑
k
=
0
n
‍
a
k
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
n
-
1
‍
b
k
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
,
where
(13)
a
k
=
∫
-
1
1
‍
f
(
x
)
T
k
h
(
cos
⁡
π
x
2
)
d
x
,
b
k
=
∫
-
1
1
‍
f
(
x
)
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
d
x
.Convergence of the HCF series (11) was extensively studied by Huybrechs in [6], where the following theorem and its corollary were stated and proved.Theorem 6.
Letf
~
(
y
)
=
f
(
x
), where x
=
(
2
/
π
)
cos
⁡
-
1
y, be analytic in the region bounded by an ellipse with major semiaxis length R and with foci 0 and 1. The corresponding domain of analyticity of f is denoted by D
(
R
). If f is analytic in the domain D
(
R
), with R
>
1
/
2, then the solution g
n to the Problem 1 satisfies
(14)
∥
f
-
g
n
∥
~
ρ
-
n
,
with
(15)
ρ
=
min
⁡
(
3
+
2
2
,
2
R
+
4
R
2
-
1
)
,
unless f is analytic and periodic on [
-
2,2
].Corollary 7.
Under the conditions of Theorem6, the coefficients a
k and b
k of g
n in the form of (12) and (13) satisfy a
k
,
b
k
~
ρ
-
n.


## 3. Construction of Chebyshev-Fourier-Collocation Methods

In this section we construct and analyze a new class of spectral methods for the solution of Problem2, which we will then call Chebyshev-Fourier-collocation (CFC) methods. The numerical solution is sought as a half-range Chebyshev-Fourier series defined in (12). The series is expanded in terms of trigonometric functions and rearranged in terms of half-range Chebyshev polynomials defined in Theorems 3 and 4. Collocation is used for the computation of expansion coefficients defined in (13). Let us turn back to Problem 2 and rewrite it as a two-point boundary value problem (3) with Dirichlet boundary conditions (4):
(16)
α
(
x
)
d
2
y
d
x
2
+
β
(
x
)
d
y
d
x
+
γ
(
x
)
y
=
f
(
x
)
,
x
∈
[
-
1,1
]
,
y
(
-
1
)
=
A
,
y
(
1
)
=
B
.
The numerical solution of the above problem is sought as a truncated HCF series, introduced in Theorem 5:
(17)
y
(
x
)
≈
y
N
(
x
)
=
∑
k
=
0
N
‍
a
k
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
N
-
1
‍
b
k
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
.
In other words, we seek for coefficients a
k and b
k, so that the numerical solution (17) is as good as possible. Here, N is the truncation number and for a given value of N we use 2
N
+
1 orthogonal polynomials: N
+
1 half-range Chebyshev polynomials of the first and N half-range Chebyshev polynomials of the second kind.We commence, by dividing the interval[
-
1,1
] with 2
N
+
1 collocation points, the Chebyshev points of the second kind into 2
N subintervals:
(18)
-
1
=
x
0
<
x
1
<
⋯
<
x
2
N
=
1
,
where
(19)
x
i
=
-
cos
⁡
(
π
i
2
N
)
,
i
=
0,1
,
…
,
2
N
.Besides, we compute the first derivative of the truncated series (17), where a
k
′ and b
k
′ denote the coefficients of the first derivative of the truncated HCF series:
(20)
d
y
N
d
x
(
x
)
=
∑
k
=
0
N
‍
a
k
′
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
N
-
1
‍
b
k
′
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
.We obtain a similar form (20) for the first derivative as for the approximation of the solution (17). The coefficients a
k
′ and b
k
′ are linearly dependent on the coefficients a
k and b
k and are computed via the differentiation matrix D:
(21)
u
′
=
D
u
,
where u
=
(
a
0
,
…
,
a
N
,
b
0
,
…
,
b
N
-
1
)
T and u
′
=
(
a
0
′
,
…
,
a
N
′
,
b
0
′
,
…
,
b
N
-
1
′
)
T. Since the coefficients a
j
′ depend only on b
k and b
j
′ only on a
k, the differentiation matrix D
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) is block-antidiagonal:
(22)
D
=
π
2
[
0
H
1
H
2
0
]
,
where H
1
∈
R
(
N
+
1
)
×
N and H
2
∈
R
N
×
(
N
+
1
). This matrix was constructed and studied in detail by Orel and Perne in [9]. Furthermore, we proceed by computing the second derivative to obtain a similar representation:
(23)
d
2
y
N
d
x
2
(
x
)
=
∑
k
=
0
N
‍
a
k
′′
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
N
-
1
‍
b
k
′′
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
,
where the coefficients of the second derivative of the truncated series are denoted by a
k
′′ and b
k
′′. Again, we compute these coefficients using differentiation matrix D:
(24)
u
′′
=
D
2
u
,
where u
′′
=
(
a
0
′′
,
…
,
a
N
′′
,
b
0
′′
,
…
,
b
N
-
1
′′
)
T.Solving BVPs via collocation assumes that the numerical solution of the BVP exactly solves the BVP in the interior collocation points (19) x
i,   i
=
1,2
,
…
,
2
N
-
1. After inserting the truncated series (17), (20), and (23) into the differential equation (3) we obtain a system of linear equations:
(25)
α
(
x
i
)
d
2
y
N
d
x
2
(
x
i
)
+
β
(
x
i
)
d
y
N
d
x
(
x
i
)
+
γ
(
x
i
)
y
N
(
x
i
)
=
f
(
x
i
)
,
i
=
1,2
,
…
,
2
N
-
1
.
Besides, we have two additional equations, originating in boundary value conditions:
(26)
∑
k
=
0
N
‍
a
k
T
k
h
(
0
)
-
∑
k
=
0
N
-
1
‍
b
k
U
k
h
(
0
)
=
A
,
∑
k
=
0
N
‍
a
k
T
k
h
(
0
)
+
∑
k
=
0
N
-
1
‍
b
k
U
k
h
(
0
)
=
B
.Let us denote byC
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) the collocation matrix. The entries of this matrix are the values of the basis functions computed in the collocation points (19). If the set of all basis functions is denoted by {
ϕ
j
}
j
=
0
2
N, then the entries of the matrix C
=
[
c
i
j
] are
(27)
c
i
j
=
ϕ
j
(
x
i
)
.
In our case the set of basis function is composed by two sets, (half-)sines and (half-)cosines, rearranged as half-range Chebyshev polynomials of the first and second kind: T
k
h
(
cos
⁡
⁡
(
π
x
/
2
)
) and U
k
h
(
cos
⁡
⁡
(
π
x
/
2
)
)
sin
⁡
(
π
x
/
2
), respectively.We have already denoted byu the set of sought expansion coefficients a
k and b
k, arranged as a vector. Besides, we denote by   v
=
(
A
,
f
(
x
1
)
,
…
,
f
(
x
2
N
-
1
)
,
B
)
T the vector of values of the right-hand side function f at interior Chebyshev collocation points x
1
,
x
2
,
…
,
x
2
N
-
1. The first and last elements of vector v are Dirichlet boundary conditions at the endpoints of the interval [
-
1,1
].Since the coefficients of the differential equation are not necessarily constant, but functions of the independent variablex, it is in general necessary to expand these coefficient functions into truncated HCF series; for example, for η
∈
{
α
,
β
,
γ
}, we have
(28)
η
(
x
)
≈
η
N
(
x
)
=
∑
k
=
0
N
‍
μ
k
T
k
h
(
cos
⁡
π
x
2
)
+
∑
k
=
0
N
-
1
‍
ν
k
U
k
h
(
cos
⁡
π
x
2
)
sin
π
x
2
.
Now we need to multiply the truncated HCF series α
N, β
N, and γ
N with the truncated HCF series y
N (17) of the numerical solution and its derivatives y
N
′ (20) and y
N
′′ (23). In order to perform these operations, we introduce the multiplication matrices F
α, F
β, and F
γ, which were constructed and studied in detail by Orel and Perne in [9]. As before, for η
∈
{
α
,
β
,
γ
}, the multiplication matrix F
η
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) is a transformation matrix between the coefficients a
k and b
k of the truncated HCF series (17) and the coefficients a
~
j and b
~
j of the truncated multiplication of the truncated HCF series (28) and (17) using the coefficients μ
k and ν
k of the truncated HCF series (28):
(29)
u
~
=
F
η
u
,
where u
~
=
(
a
~
0
,
…
,
a
~
N
,
b
~
0
,
…
,
b
~
N
-
1
)
T denotes the vector of the expansion coefficients of the truncated HCF series of the multiplication η
N
(
x
)
·
y
N
(
x
). The multiplication matrix F
η
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) is a block matrix:
(30)
F
η
=
[
G
1
G
2
G
3
G
4
]
,
where G
1
∈
R
(
N
+
1
)
×
(
N
+
1
), G
2
∈
R
(
N
+
1
)
×
N, G
3
∈
R
N
×
(
N
+
1
), and G
4
∈
R
N
×
N. In the case of a differential equation with constant coefficients, the matrix F
η is scalar; otherwise the matrix is dense.Let us now denote byL
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) the differential operator matrix
(31)
L
=
F
α
D
2
+
F
β
D
+
F
γ
,
where D is the differentiation matrix (21) and F
α, F
β, and F
γ are multiplication matrices (29). We are now ready to rewrite the linear system of (25) and (26) into a matrix form
(32)
U
u
=
v
,
where the matrix U
∈
R
(
2
N
+
1
)
×
(
2
N
+
1
) is obtained by replacing the first and the last row of the matrix C
L with the first and the last row of the collocation matrix C to satisfy the boundary conditions (26). The solution of the system (32) gives the spectral coefficients of the truncated HCF series for the numerical solution of the two-point boundary value problem (16). Multiplication with collocation matrix C yields the solution values at collocation (Chebyshev) points (19).


## 4. Error Analysis

In this section we focus on error analysis and the rate of convergence of the new class of Chebyshev-Fourier-collocation spectral methods introduced in Section3. In order to prove convergence of the method and to estimate the error of the numerical solution, we follow error estimation techniques presented by Canuto et al. in [13] and by Shen et al. in [14]. Since the HCF series is a generalized trigonometric series reorganized in terms of half-range Chebyshev polynomials of the first and second kind, it is convenient to follow the steps for error estimating of Fourier spectral methods.Let us restrict to the problem of solving linear BVPs of second order with homogeneous Dirichlet boundary conditions on the interval[
-
1,1
]:
(33)
L
u
=
f
,
(34)
u
(
-
1
)
=
u
(
1
)
=
0
,
where L is a linear differential operator defined in (5). Let us assume that the coefficient function α
(
x
)
≡
-
1 on the interval [
-
1,1
] and let us further assume that the coefficient function β is differentiable and both functions β and γ are bounded and strictly positive on the interval [
-
1,1
]. Finally, let us assume that the condition γ
(
x
)
-
β
′
(
x
)
/
2
>
0 holds for every x
∈
[
-
1,1
].The required Hilbert space isL
2
(
-
1,1
), the space of all square integrable functions on the interval [
-
1,1
]. In this space the operator L is unbounded. We denote by
(35)
(
u
,
v
)
=
∫
-
1
1
‍
u
(
x
)
v
(
x
)
d
x
the appropriate inner product in L
2
(
-
1,1
) and by ∥
u
∥
=
(
u
,
u
)
1
/
2 the associated norm. Moreover, let T
N
⊂
L
2
(
-
1,1
) denote the space of all trigonometric polynomials of degree ≤
N that satisfy the boundary conditions (34) and let us further denote by (
u
,
v
)
N the appropriate discrete inner product with the associated norm ∥
u
∥
N
=
(
u
,
u
)
N
1
/
2. The collocation solution u
N
∈
T
N of (33) and (34) satisfies the equations
(36)
L
N
u
N
(
x
k
)
=
f
(
x
k
)
,
u
N
(
x
0
)
=
u
N
(
x
2
N
)
=
0
.
Here, the nodes x
k, k
=
1,2
,
…
,
2
N
-
1, are the interior collocation points defined in (19) and the operator L
N is an approximation to the operator L, obtained by replacing exact derivatives by interpolation derivatives (21) and (24) and expanding coefficient functions in HCF series.Equations (33) and (34) can be equivalently written in a weak form as a bilinear form:
(37)
(
L
u
,
v
)
=
(
f
,
v
)
,
∀
v
∈
L
2
(
-
1,1
)
,
where u satisfies the boundary conditions (34). The collocation method (36) can be then rewritten as
(38)
(
L
N
u
N
,
v
)
N
=
(
f
,
v
)
N
,
v
∈
T
N
,
where u
N
∈
T
N. Analysis of convergence properties requires the existence of a dense Hilbert subspace of L
2
(
-
1,1
). An appropriate choice in our analysis is the Sobolev space
(39)
H
1
(
-
1,1
)
=
{
v
∈
L
2
(
-
1,1
)
;
d
v
d
x
∈
L
2
(
-
1,1
)
}
,
where the derivative d
v
/
d
x in the sense of distributions belongs to L
2
(
-
1,1
). In the following we denote the first derivative by v
′
≡
d
v
/
d
x and the second derivative by v
′′
≡
d
2
v
/
d
x
2. This Sobolev space is equipped with the Sobolev norm
(40)
∥
v
∥
H
1
=
(
∥
v
∥
L
2
2
+
∥
v
′
∥
L
2
2
)
1
/
2
,
such that ∥
u
∥
L
2
≤
∥
u
∥
H
1 for all u
∈
H
1
(
-
1,1
). Note that for every N
>
0, the space T
N is contained in H
1
(
-
1,1
). Moreover, our analysis requires that the operator L, or more exactly the bilinear form (
L
u
,
v
), satisfies the coercivity condition
(41)
∃
α
*
>
0
:
(
L
u
,
u
)
≥
α
*
∥
u
∥
H
1
2
,
u
∈
T
N
,
and the continuity condition
(42)
∃
A
>
0
:
|
(
L
u
,
v
)
|
≤
A
∥
u
∥
H
1
∥
v
∥
H
1
,
u
,
v
∈
T
N
.
Since
(43)
(
L
u
,
u
)
=
∫
-
1
1
‍
(
-
u
′′
+
β
(
x
)
u
′
+
γ
(
x
)
u
)
u
d
x
=
∫
-
1
1
‍
(
u
′
)
2
d
x
+
∫
-
1
1
‍
(
γ
(
x
)
-
β
′
(
x
)
2
)
u
2
d
x
≥
∫
-
1
1
‍
(
u
′
)
2
d
x
+
λ
∫
-
1
1
‍
u
2
d
x
≥
α
*
∥
u
∥
H
1
2
,
where we use integration by parts in the second row, the coercivity condition (41) for the bilinear form (37) is satisfied with
(44)
α
*
=
min
⁡
{
1
,
λ
}
,
where
λ
=
inf
⁡
x
∈
[
-
1,1
]
(
γ
(
x
)
-
β
′
(
x
)
2
)
>
0
.
Similarly, since
(45)
|
(
L
u
,
v
)
|
=
|
∫
-
1
1
‍
(
-
u
′′
+
β
(
x
)
u
′
+
γ
(
x
)
u
)
v
d
x
|
≤
max
⁡
x
∈
[
-
1,1
]
{
1
,
β
(
x
)
,
γ
(
x
)
}
(
|
∫
-
1
1
‍
u
′′
v
d
x
|
m
m
m
m
m
m
m
m
m
m
m
m
m
m
+
|
∫
-
1
1
‍
u
′
v
d
x
|
+
|
∫
-
1
1
‍
u
v
d
x
|
)
=
max
⁡
x
∈
[
-
1,1
]
{
1
,
β
(
x
)
,
γ
(
x
)
}
(
|
∫
-
1
1
‍
u
′
v
′
d
x
|
m
m
m
m
m
m
m
m
m
m
m
m
m
m
+
|
∫
-
1
1
‍
u
v
′
d
x
|
+
|
∫
-
1
1
‍
u
v
d
x
|
)
≤
max
⁡
x
∈
[
-
1,1
]
{
1
,
β
(
x
)
,
γ
(
x
)
}
(
∥
u
′
∥
L
2
∥
v
′
∥
L
2
m
m
m
m
m
m
m
m
m
m
m
m
m
m
+
∥
u
∥
L
2
∥
v
′
∥
L
2
+
∥
u
∥
L
2
∥
v
∥
L
2
)
≤
A
∥
u
∥
H
1
∥
v
∥
H
1
,
where we use integration by parts in the second row, the Cauchy-Schwartz inequality in the third row, and inequalities ∥
u
∥
L
2
≤
∥
u
∥
H
1 and ∥
u
′
∥
L
2
≤
∥
u
∥
H
1 in the fourth row, the continuity condition (42) is satisfied with
(46)
A
=
3
max
⁡
x
∈
[
-
1,1
]
⁡
{
1
,
β
(
x
)
,
γ
(
x
)
}
.
Both constants α
* and A are independent of N. Let us note that the coercivity and continuity conditions are sufficient but not necessary as is shown in the forthcoming examples.Let us further denotee
=
u
N
-
R
N
u, where u
∈
L
2
(
-
1,1
) is the exact, u
N
∈
T
N is the numerical solution of (33) and (34), and R
N is a projection operator from L
2
(
-
1,1
) to T
N. Under the assumptions of the first Strang lemma (see Theorem 1.3, page 14, [14]), that is, the coercivity (41) and the continuity (42) condition are satisfied with constants α
* and A defined in (44) and (46), the problem (38) admits a unique numerical solution u
N
∈
T
N, satisfying
(47)
∥
u
N
∥
L
2
≤
1
α
*
sup
⁡
0
≠
v
∈
T
N
|
(
f
,
v
)
N
|
∥
v
∥
L
2
.
Moreover, the first Strang lemma states that if the coercivity and continuity conditions are fulfilled, the error estimate of the numerical solution u
N reads as follows:
(48)
∥
u
-
u
N
∥
H
1
≤
∥
u
-
R
N
u
∥
H
1
+
∥
e
∥
H
1
≤
(
1
+
A
α
*
)
∥
u
-
R
N
u
∥
H
1
+
1
α
*
|
(
Q
N
f
,
e
)
N
-
(
f
,
e
)
|
∥
e
∥
H
1
+
1
α
*
|
(
L
R
N
u
,
e
)
-
(
Q
N
L
N
R
N
u
,
e
)
N
|
∥
e
∥
H
1
.Here,Q
N is a projection operator from L
2
(
-
1,1
) to T
N according to the discrete inner product and Q
N
v is then a trigonometric polynomial of degree N, matching v at the interior collocation points (19) and vanishing at the boundary points. The method is convergent, if all three parts of the inequality (48) converge to 0 with N
→
∞. Since (see (5.5
.
15), page 294, [13])
(49)
∥
u
-
R
N
u
∥
H
1
≤
C
1
N
1
-
m
∥
u
(
m
)
∥
L
2
,
where m is the order of smoothness of the right-hand side function f and the solution u, and
(50)
|
(
Q
N
f
,
e
)
N
-
(
f
,
e
)
|
∥
e
∥
H
1
≤
∥
Q
N
f
-
f
∥
H
1
∥
e
∥
H
1
∥
e
∥
H
1
≤
D
2
N
1
-
m
∥
f
(
m
)
∥
L
2
the above is true for the first two parts of inequality (48). For the third part of that inequality we use the fact that ∥
v
∥
H
1
≤
∥
v
∥
H
2 for every v
∈
H
2
(
-
1,1
), where
(51)
H
2
(
-
1,1
)
=
{
v
∈
L
2
(
-
1,1
)
;
v
′
,
v
′′
∈
L
2
(
-
1,1
)
}
,
is the Sobolev space equipped with the norm
(52)
∥
v
∥
H
2
=
(
∥
v
∥
L
2
2
+
∥
v
′
∥
L
2
2
+
∥
v
′′
∥
L
2
2
)
1
/
2
.
Moreover, we observe that the differential operator L is bounded in this space; see Leoni [17] or Ziemer [18]. Finally, we obtain the estimate
(53)
|
(
L
R
N
u
,
e
)
-
(
Q
N
L
N
R
N
u
,
e
)
N
|
∥
e
∥
H
1
≤
∥
L
R
N
u
-
Q
N
L
N
R
N
u
∥
H
1
∥
e
∥
H
1
∥
e
∥
H
1
≤
∥
L
R
N
u
-
L
u
∥
H
1
+
∥
L
u
-
Q
N
L
N
u
∥
H
1
+
∥
Q
N
L
N
u
+
Q
N
L
N
R
N
u
∥
H
1
≤
∥
L
∥
H
2
∥
R
N
u
-
u
∥
H
2
+
∥
f
-
Q
N
f
∥
H
1
+
∥
Q
N
L
N
∥
H
2
∥
u
-
R
N
u
∥
H
2
≤
C
2
N
2
-
m
∥
u
(
m
)
∥
L
2
+
D
2
N
1
-
m
∥
f
(
m
)
∥
L
2
.
All constants C
1, C
2, and D are independent of N and m. We have proved the following theorem which states the error estimate and the rate of convergence for solutions and right-hand side functions that are continuously differentiable to a certain order m.Theorem 8.
Letu
∈
L
2
(
-
1,1
) be the exact solution of the problem (33) with boundary conditions (34) and let u
N
∈
T
N be the numerical solution obtained by the class of Chebyshev-Fourier-collocation (CFC) methods, constructed in Section 3. Let one assume that the solution and the right-hand side function are m-times continuously differentiable. Then the estimated error of the approximation of the solution for the constructed class of CFC methods is
(54)
∥
u
-
u
N
∥
H
1
≤
(
1
+
A
α
*
)
C
1
N
1
-
m
∥
u
(
m
)
∥
L
2
+
C
2
N
2
-
m
∥
u
(
m
)
∥
L
2
+
D
N
1
-
m
∥
f
(
m
)
∥
L
2
,
where the constants α
* and A are defined in (44) and (46).If the solutionu and the right-hand side function f are smooth or analytic functions in some domain, containing the interval [
-
1,1
], then under the conditions of Theorem 6, the Chebyshev-Fourier-collocation method has a spectral rate of convergence.Theorem 9.
Letu
∈
L
2
(
-
1,1
) be the exact solution of the problem (33) with boundary conditions (34) and let u
N
∈
T
N be the numerical solution obtained by the class of Chebyshev-Fourier-collocation (CFC) methods, constructed in Section 3. Let one assume that the solution and the right-hand side function are analytic functions in the domain D
(
R
), defined in Theorem 6. Then the estimated error of the approximation of the solution for the constructed class of CFC methods is
(55)
∥
u
-
u
N
∥
~
ρ
-
N
,
where ρ is defined in (15).


## 5. Numerical Examples

In the following examples we compare numerical solutions obtained with two different classes of collocation spectral methods. One is a standard Chebyshev-collocation approach, where we approximate the solution using Chebyshev series; the other one is the new class of methods constructed in Section3, where approximation of the solution with half-range Chebyshev-Fourier series is used. Much about classical Chebyshev polynomials can be found, for example, in [19]. Both methods were implemented in MATLAB by the authors.We observe that the performance of the Chebyshev-Fourier-collocation (CFC) method is comparable with the standard Chebyshev-collocation (CC) method. However, since the absence of a fast technique for the computation of the expansion coefficients, comparable with the FFT, the performance in terms of computational costs is worse for the new class of methods. Throughout this section we use the abbreviationsy
′
≡
d
y
/
d
x and y
′′
≡
d
2
y
/
d
x
2.Example 1.
As a first example we consider a second order linear differential equation with nonconstant coefficients:(56)
y
′′
+
x
y
′
=
(
2
+
x
2
)
cos
⁡
x
,
y
(
-
1
)
=
sin
(
1
)
,
y
(
1
)
=
sin
(
1
)
.
The exact solution is
(57)
y
(
x
)
=
x
sin
x
.In Figure1(a), the comparison of the decay of the maximal absolute error of the two numerical solutions, obtained with CC and CFC methods with respect to the truncation number N of terms in the Chebyshev and Chebyshev-Fourier expansions, is depicted. For both methods we have to compute N
+
1 coefficients. The exact solution of this problem is a smooth and analytic function and both methods reach machine accuracy. However, as seen from the figure, the CC numerical solution converges more rapidly, since the maximal absolute error reaches machine accuracy at N
=
14, instead of N
=
34 for the CFC numerical solution. Note the exponential decay of the maximal absolute error for both methods. In this case we obtain spectral accuracy as predicted by Theorem 9 in (55).Plots of the maximal absolute error of the numerical solution with respect to the truncation numberN. Error of the Chebyshev-Fourier-collocation (CFC) method is shown with a plus (+) mark; error of the Chebyshev-collocation (CC) method is shown with a circle (∘) mark. (a) Error for Example 1. (b) Error for Example 2. (c) Error for Example 3. (d) Error for Example 4.
(a)
Maximal error dependent onN
(b)
Maximal error dependent onN
(c)
Maximal error dependent onN
(d)
Maximal error dependent onNExample 2.
As a second example we consider the Airy differential equation(58)
y
′′
-
(
x
-
1000
)
y
=
0
,
y
(
-
1
)
=
1
,
y
(
1
)
=
1
.
The exact solution is a highly oscillatory function
(59)
y
(
x
)
=
(
A
i
(
x
-
1000
)
(
B
i
(
-
1001
)
-
B
i
(
-
999
)
)
m
m
-
B
i
(
x
-
1000
)
(
A
i
(
-
999
)
-
A
i
(
-
1001
)
)
)
×
(
A
i
(
-
999
)
B
i
(
-
1001
)
-
A
i
(
-
1001
)
B
i
(
-
999
)
)
-
1
,
where A
i and B
i are Airy functions; see, for example, [19].Figure1(b) shows the comparison of the decay of the maximal absolute error of the two numerical solutions obtained with CC and CFC methods with respect to the truncation number N. The exact solution is a smooth but oscillatory function. For this reason the maximal absolute error begins to decay not until N is big enough to overcome problems with the resolution. However, both methods again yield exponential decay of the maximal absolute error with respect to N, which gives spectral accuracy. We note that the convergence of the Chebyshev-Fourier-collocation method is faster in comparison with the Chebyshev-collocation method.Example 3.
As a third example we consider a second order linear differential equation with nonconstant and nonsmooth coefficients:(60)
y
′′
-
|
x
|
y
′
+
2
y
=
30
x
3
|
x
|
+
2
x
5
|
x
|
-
6
x
6
,
y
(
-
1
)
=
-
1
,
y
(
1
)
=
1
.
The exact solution is
(61)
y
(
x
)
=
x
5
|
x
|
.Figure1(c) shows the comparison of the decay of the maximal absolute error for the two numerical solutions with CC and CFC methods with respect to the truncation number N. In this example, the exact solution is nonsmooth, actually it is only six-times continuously differentiable. For this reason the decay of the maximal absolute error is slow for both methods; however, it is faster for the CFC method. The error decays as O
(
1
/
N
5
). This is in accordance with the result (54) of Theorem 8. We note that in this case the Chebyshev-Fourier-collocation method outweighs the Chebyshev-collocation method.Example 4.
As a last example we again consider a second order linear differential equation with nonconstant and nonsmooth coefficients:(62)
y
′′
+
2
|
x
|
y
′
+
3
y
=
56
x
5
|
x
|
+
3
x
7
|
x
|
+
16
x
8
,
y
(
-
1
)
=
-
1
,
y
(
1
)
=
1
.
The exact solution is
(63)
y
(
x
)
=
x
7
|
x
|
.Figure1(d) shows the comparison of the decay of the maximal absolute error for the two numerical solutions with CC and CFC methods with respect to the truncation number N. In this example, the exact solution is again nonsmooth; actually it is only eight times continuously differentiable. For this reason the decay of the maximal absolute error is slow for both methods; however, it is faster than in the previous example; actually the error decays as O
(
1
/
N
7
). This is again in accordance with the result (54) of Theorem 8. We note that in this case the Chebyshev-Fourier-collocation method again outweighs the Chebyshev-collocation method.


## 6. Conclusions

Based on the papers by Huybrechs [6] and by Orel and Perne [9], we first discussed two nonclassical families of orthogonal polynomials: the half-range Chebyshev polynomials of the first and second kind and the associated half-range Chebyshev-Fourier series. We have briefly recalled the solution of the approximation problem stated in Problem 1, using the HCF series, and analyzed by Huybrechs in [6]. The main focus of the paper is on analyzing Problem 2 via collocation spectral methods using truncated HCF series. These series are generalized Fourier series, where the basic trigonometric functions are rearranged in terms of the half-range Chebyshev polynomials of the first and second kind. The usage of such reorganization which yields orthogonality of the expansion basis functions allows solving nonperiodic problems with tools, otherwise reserved for periodic problems.In the paper we have constructed a new class of spectral methods, based on collocation, the Chebyshev-Fourier-collocation methods. The idea is similar to the widely used Chebyshev spectral methods, where instead of using classical Chebyshev series for the approximation of the solution, we expand the numerical solution as a half-range Chebyshev-Fourier series. We deal with linear two-point boundary value problems of the second order. Generalization to higher order and/or different intervals is also possible. Error analysis and convergence theory, presented in Section4, show that for problems that are smooth enough, that is, analytic functions, we obtain spectral accuracy, that is, the maximal absolute error depending on the truncation number N decays exponentially. Otherwise, we obtain algebraic convergence. These results are comparable with the ones for the standard Chebyshev-collocation methods. Examples shown at the end of the paper demonstrate exponential convergence for smooth and analytic functions and comparability with standard classes of spectral methods. Furthermore, for some problems the CFC method outweighs the CC method.As far as we are interested in convergence theory, things are more or less comparable with standard classes of spectral methods. This is regrettably not the case if we are concerned with the computational costs. For computing spectral coefficients for the HCF series we do not yet have such a marvelous tool as it is the FFT for Fourier or Chebyshev series. This is an open problem that needs further investigation. Moreover, as a part of future work with HCF series, we are interested in evolutive time-dependent partial differential equations, for example, heat or wave equations. More of this theme will be discussed in a subsequent paper.



---

*Source: 572694-2014-06-01.xml*