=Paper=
{{Paper
|id=Vol-1839/MIT2016-p36
|storemode=property
|title= Modification of Fourier approximation for solving boundary value problems having singularities of boundary layer type
|pdfUrl=https://ceur-ws.org/Vol-1839/MIT2016-p36.pdf
|volume=Vol-1839
|authors=Boris Semisalov,Georgy Kuzmin
}}
== Modification of Fourier approximation for solving boundary value problems having singularities of boundary layer type==
https://ceur-ws.org/Vol-1839/MIT2016-p36.pdf
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Modification of Fourier Approximation for
Solving Boundary Value Problems Having
Singularities of Boundary Layer Type
Boris Semisalov and Georgy Kuzmin
Institute of Computational Technologies, SB RAS,
Academic Lavrentiev av. 6, 630090, Novosibirsk, Russia,
Novosibirsk State University
Pirogova str. 1, 630090, Novosibirsk, Russia
ViBiS@ngs.ru,kubgeom@mail.ru
Abstract. A method for approximating smooth functions has been de-
veloped using non-polynomial basis obtained by mapping of Fourier se-
ries domain to the segment [β1, 1]. High rate of convergence and stability
of the method is justified theoretically for four types of coordinate map-
pings, the dependencies of approximation error on values of derivatives
of approximated functions are obtained. Algorithms for expanding of
functions into series with coupled basis composed of Chebyshev polyno-
mials and designed non-polynomial functions are implemented. It was
shown that for functions having high order of smoothness and extremely
steep gradients in the vicinity of bounds of segment the accuracy of pro-
posed method cardinally exceeds that of Chebyshevβs approximation. For
such functions method allows to reach an acceptable accuracy using only
π = 10 basis elements (relative error does not exceed 1 per cent)
Keywords: singular perturbation, small parameter, coordinate map-
ping, boundary value problem, Fourier series, Chebyshev polynomial,
non-polynomial basis, estimate of convergence rate, collocation method
1 Introduction
By now, a huge amount of urgent scientific and technological problems reduces
to boundary value problems for differential equations having pronounced singu-
larities of boundary layer type. The most popular approaches to solving them are
based on construction of computational grids with piecewise linear/polynomial
approximation of the unknown function in each cell of grid [1]. Such approaches
provide relatively low rate of convergence and lead to essential refinement of grid
in the vicinity of boundary layer and consequently to growth of computational
costs and errors. In [2, 3] methods of coordinate transformations were developed
which allow to decrease the influence of mentioned effect due to application of
special coordinate mappings eliminating the singularity. Nevertheless, the anal-
ysis of key issue concerning the smoothness of such transformations and its
influence on the quality of approximation method is absent. Frequently, authors
406
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
restrict their self by transforming uniform or more special grid and performing
numerical experiments using finite difference or spectral methods [2]β [5].
In the present paper a step aside from the traditional grid approaches is made
and the approximations based on mapping of Fourier series domain to the seg-
ment [β1, 1] are used to approximate functions having singularities of boundary
layer type. One of such mappings assigned by function cos(π₯) transforms Fourier
basis to basis composed of Chebysev polynomials. A special feature inherent to
these bases is the absence of saturation of corresponding approximations [6]. It
means that using Fourier and Chebyshev bases allows to obtain the asymptotic
of error of best polynomial approximation while approximating functions having
any order of smoothness or singularities in complex plain. The loss of effective-
ness of the mentioned approximations while solving problems having singularities
of boundary layer type is caused by degradation of asymptotical properties of
best polynomial approximations with fast increase of gradients of approximated
function. In order to eliminate this problem, trigonometric Fourier basis can be
transformed into non-polynomial algebraic one retaining all its good properties
of high convergence rate and computational stability, but specially adapted to
approximation of smooth functions having singularities of boundary layer type.
2 Preliminary information. Problem description
In framework of approximation theory of continuous and smooth functions
π (π₯) (π₯ β π· β R) that are elements of spaces πΆ(π·), πππ (π, π·) = {π β πΆ π (π·) :
|π (π) | β€ π = π (π)} the following notions are often used (see for example [6]).
(οΈ )οΈ π1
β«οΈ π
1. Norms of function βπ β = max |π |, βπ βπ = π (π₯)ππ₯ .
π₯βπ· π·
2. Finite-dimensional approximating space πΎπ with elements used for
approximation of π (π₯) that usually are series or polynomials including π sum-
mands or monomials.
3. Operator of approximation (or simply approximation) of function
π (π₯) is a continues mapping ππ performing a projection of functional space on
approximating one (ππ : πΆ(π·) β πΎπ or ππ : πππ (π, π·) β πΎπ ).
4. Best approximation of function π β πΆ(π·) is element ππ (π, πΎπ ) β πΎπ
providing the lower bound to be reached π*π (π, πΎπ ) = = π*π (π ) = inf βπ β πβ.
πβπΎπ
5. Method without saturation (loose definition) is a method of approxi-
mation of function π (π₯) having asymptotic of error of the best polynomial ap-
proximation for any order of smoothness of π (π₯). Rigorous definition based on
the analysis of asymptotic of Alexandrovβs diameters is given in [6].
The results obtained in works by Lebegue, Faber, Jackson, Bernstain al-
low to separate three classes of smooth functions with fundamentally different
asymptotical behavior of errors of best approximations in space πΎπ of algebraic
polynomials of πth power. The similar asymptotical behavior is valid for periodic
functions and trigonometrical polynomials
407
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
I. If π (π₯) β πππ (π, π·) is π-times continuously differentiable function on
segment π· and all its derivatives up to order π are bounder by value π (π), then
sup π*π (π, πΎπ ) β€ π (π)πΆπ πβπ , (1)
π βπππ (π,π·)
where πΆπ depends on π only, [7].
II. If π (π₯) β πΆ β (π·) is infinite differentiable function with a singularity (like
pole) on complex plain, then one can find a number π (0 < π < 1) and a sequence
of polynomials ππ (π₯), such that
βπ (π₯) β ππ (π₯)β β€ πΆπ π , π₯ β π·, (2)
here π < 1 is defined by location of singularity in complex plain, πΆ is constant, [8].
III. If π (π₯) β πΈππ‘ is entire function, then
π (π)(diamπ·)π 21β2π
π*π (π, πΎπ ) β€ , (3)
π!
where π (π) = βπ (π) β = π(π!). It follows from estimates of error of polynomial
interpolation (see [9]) and CauchyβHadamard inequality.
Remark 1. For smooth functions of I and III classes the accuracy of best approxi-
mations depends on maximal values of derivatives of function on π· (values of
π (π) and π (π) in (1), (3)). For functions of class II the error can be defined
through the values of function itself beyond the segment π· on complex plain.
In order to implement the properties of best approximations in this work the
Fourier and Chebyshev series are used. Note that such approximations are equal
in a specific sense. Indeed, let π β πΆ([β1, 1]), then performing the change of
variable π₯ = cos π, π β [0, 2π] one can obtain 2π-periodic even function πΛ(π) =
βοΈ
β
π (cos π). Fourier decomposition of it is πΛ(π) = ππ cos(ππ). Hence
π=0
β
βοΈ β
βοΈ
π (π₯) = ππ cos(π arccos(π₯)) = ππ ππ (π₯). (4)
π=0 π=0
In other words Chebyshev polynomials ππ (π₯) can by obtained by mapping Fourier
series domain to the segment [β1, 1], see [10]. In [6] is proved that such approx-
imations presents the methods without saturation and therefore they are ex-
tremely efficient for solving problems with smooth solutions. However, if values
π (π), πΆ, π (π) in (1)β(3) are large, then it can be wrong.
A simple example is a boundary-value problem for differential equation of
second order with small factor π of second derivative
π2 π
π β π = ππ β²β² (π₯) β π(π₯), π (β1) = 1, π (1) = β1, (5)
ππ₯2
here π₯ β [β1, 1], π(π₯) is a given smooth function. Solution to the problem is
π (π₯) = π(π₯) + π(π₯), π(π₯) = πΆ1 ππ΄(0.5π₯+0.5) + πΆ2 πβπ΄(0.5π₯+0.5) ,
408
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
1 + πβπ΄ 1 + ππ΄
πΆ1 = , πΆ2 = .
πβπ΄ β ππ΄ ππ΄ β πβπ΄
Here
βοΈ π(π₯) is an exponential boundary value component of π (π₯), πΆ1 , πΆ2 , π΄ =
1/π are constants. Table 1 shows the values of dimension of πΎπ ensuring that
the relative errors of approximation is never higher than 1 per cent. These results
were obtained in accordance with the estimates of best polynomial approxima-
tions (1), (3), while function π(π₯) has different orders of smoothness.
Table 1. Dimension π ensuring that relative error of approximation
is not higher than 1 per cent
Order of The value of small parameter π
smoothness π = 10β4 π = 10β6 π = 10β8 π = 10β10
π=1 146 430 1 170 100 19 537 000 219 560 000
π=2 5 997 65 952 714 010 7 643 800
π=3 2 265 24 251 256 530 2 691 200
π=4 1 426 15 037 157 040 1 629 500
π=5 1 090 11 390 118 000 1 215 900
entire function 136 1 359 13 591 135 910
Thus, it can be observed that the higher order of smoothness is, the less data
is necessary for recovering solution with a desired accuracy. However, even in
the case of infinitely differentiable function a space πΎπ of dimension of many
thousands can be required to reach considerably low accuracy of 1 per cent. This
effect is shown on Fig 1 where a graph of solution to a problem (5) is given with
π(π₯) β‘ 0 (i.e. a graph of function π(π₯)) and logarithm of error of approximation
of π(π₯) in space of Chebyshev polynomials
βοΈ
πβ1
π = max |π(π₯) β ππ ππ (π₯)|.
π₯β[β1,1] π=0
a b
1 0
οΏ½(x) log10οΏ½
0.5
β5
0
β10
β0.5
β1
β1 β0.5 0 0.5
x
1
β15
0 20 40 60 80
n
100
Fig. 1. Solution to the problem (5) with π = 10β3 (solid line), π = 10β4 (dash line),
π = 10β6 (dote-and-dash line): a β graph of π(π₯); b β dependance of log10 π on π
Fig. 1 shows that the values of π, given in Table 1 for function π(π₯) are
overestimated. This effect can be explained by concentration of Chebyshev nodes
409
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
in the vicinity of boundary layers, see. [10], [11]. Considering a problem with
boundary layer of size π having solution without singularities in inner part of
the segment [β1, 1], to reach an accuracy of 1 per cent one should use a basis of
β
π β 3/ π (6)
Chebyshev polynomials for approximation of unknown function. In the consid-
ered case π β 5.6/π΄ (here condition π(π) = 0.1 is used). Appearance of expres-
β
sion π in the denominator of fraction is caused by concentration of Cheby-
shev nodes. Indeed, uniformly distributed zeroes of trigonometrical monomials
cos(ππ) are concentrated in the vicinity of points Β±1 under the map π₯ = cos(π)
(see fig. 2 a). Moreover as cos(π) βΌ 1 β π2 /2 when π β 0 and the first Chebyshev
node π₯1 = cos(π1 ) = cos π/2π, then |1 β π₯1 | βΌ π 2 /8π2 . Further, the empiri-
cal requirement that even three nodes should lay on the boundary layer gives
3π
3|1 β π₯1 | β€ π or π β β β that corresponds to (6).
2 2 π
3 Description of a method
For efficient approximation of function having boundary layer component a mod-
ification of map π₯ = cos(π) that transformed Fourier basis to Chebyshev one (4)
is proposed. According to (6) a natural requirement is to use stronger concen-
tration of zeroes of basis functions in the vicinity of segment borders (see 2 b).
Let us assume that π· = [β1, 1]. Define a function π₯ = Γ¦(π¦) : [β1, 1] β [β1, 1]
satisfying the following basic requirements:
1) function Γ¦(π¦) is bijective mapping, Γ¦(1) = 1, Γ¦(β1) = β1;
2) Γ¦(π¦) is an infinite differentiable or even entire function;
3) an inverse function π¦ = Γ¦β1 (π₯) can be expressed in analytical form or easily
computed;
4) a derivative Γ¦β² (π¦) in the vicinity of points Β±1 is close to zero.
Note that concentration of zeroes of basis functions in the vicinity of segment
borders can be obtained due to the last requirement. One of possible forms of
function Γ¦(π¦) is given on fig 2 b.
For approximation of function π (π₯), π₯ β [β1, 1] let us consider the expansion
of 2π-periodic even function π (Γ¦(cos π)) (π β R) into Fourier series:
πβ1
βοΈ
π (Γ¦(π¦)) = π (Γ¦(cos π)) β ππ cos(ππ). (7)
π=0
As a result one obtains
πβ1
βοΈ
π (π₯) β ππ (π₯) = ππ cos[π arccos(Γ¦β1 (π₯))]. (8)
π=0
410
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
a
1 1
2
3
0.5 4 y=cos οΏ½
5
0 6
οΏ½
1 2 3 4 5 6 7 8 9 10 11
7
β0.5 8
9
10
β1 11
0 0.5 1 1.5 2 2.5 3
b 11
1 10
9
8
0.5 7 x=ae(y)
y
0
1 2 3 4 5 6 7 8 9 10 11
β0.5 5
4
β1 3
2
β1 β0.5 01 0.5 1
Fig. 2. Concentration of nodes of trigonometrical monomials cos(ππ) as π = 11 (indi-
cated by arrows): a β using mapping π¦ = cos(π), b β using mapping π₯ = Γ¦(π¦)
Here we denote π¦ = cos π. Thus, the basis of approximating space πΎπ can be
specified
(οΈ as π΅π (π₯))οΈ = cos[π arccos(Γ¦β1 (π₯))], where zeroes of π΅π (π₯) are π₯ππ =
(2π + 1)π
Γ¦ cos , π, π = 0, ..., π β 1. Note, that under the properties 1, 3,
2π
Γ¦(π¦), π΅π (π₯) are bijective easily computed functions.
Lemma 1. Approximation (8) is equivalent to expansion of function π (Γ¦(π¦))
into series with basis consists of Chebyshev polynomials ππ (π¦) = cos(π arccos(π¦)),
that is why
1) for approximations (8) error estimations of best approximations (1)β(3) hold;
2) elements of matrices Bπ β πππ = π΅π (π₯ππ ), π, π = 0, ..., π β 1 do not depend
on Γ¦.
The proof of Lemma 1 is obvious taking into account (8) and that Cheby-
shev and Fourier expansions are approximations without saturation. The second
condition of Lemma 1 concerns numerical approximation of π by collocation
method with nodes π₯ππ . Lemma declares that such a method is universal, its
properties do not depend on the choice of function Γ¦(π¦).
To settle a key question on rate of growth of coefficients π (π), π (π) in
estimates (1), (3) the following property was proved.
411
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Lemma 2. For derivative of composed function the followin equality is valid
π
βοΈ [οΈ βοΈ π
βοΈ ]οΈ
(π) (π) π (πΌπ )
[π (Γ¦(π¦))] = π (Γ¦(π¦)) πΆππ (Γ¦(π¦)) , (9)
π=1 πππ =(πΌ1 ,...,πΌπ ) π=1
πΌ1 +...+πΌπ =π
where the second summation goes over all possible integer partitions of number
π
π β πππ , that consist of π positive integer numbers πΌ1 ,...,πΌπ , π = 1, ..., πΏ; πΆππ are
π π
constants. Moreover βπ as π = 1 and π = π one has: πΏ = 1, πΆ11 = 1, πΆπ1 = 1.
This Lemma corresponds to the well-known Fa di Brunoβs formula and can
be proved using mathematical induction technique.
4 Analysis of four types of function Γ¦(π¦).
Now let us propose and investigate four types of function Γ¦(π¦). To this end
the following notations are necessary. Let π·π β π·, π·π β π· be neighborhoods of
boundaries of segment π· representing boundary layers, π·0 be a certain neighbor-
hood of central point of π·. Let us denote by β Β· βπΏ , β Β· βπ· the supremum norms of
continues function on π·π βͺ π·π and on π·0 correspondingly. Further, we assume
π = β1, π = 1, βπ < π, π΄π βπ (π ) β = βπ (π +1) β, where π΄π > 1 is constant, π is or-
der of smoothness of π (π₯). Typically for problems with boundary layer one has
π΄π >> 1. Let π΄ = max π΄π , π(π΄) be a size of boundary layer (as it follows from
π <π
the example with exponential boundary layer, π depends on π΄, see comments to
π(π (π ) )
formula (6)), π(π (π ) ) be such value that lim = 0.
π΄1 ,...,π΄πβ1 ββ βπ (π ) βπΏ
I. Trigonometric function
(οΈ )οΈ
ππ¦
Γ¦(π¦) = sin . (10)
2
Theorem 1. Let π2 (π΄)π΄ β 0 as π΄ β β. For the approximation (8) with
function Γ¦(π¦) of form (10) the estimate (1) is valid, but constant
π π (π/2)
π (π) = πΆπ/2 π (π/2) βπ β + π(π (π/2) ) (if π is even),
π
π (π) = πΆ[π/2]π π(π΄)(π/2)π βπ ([π/2]) β + π(π ([π/2]) ) (if π is odd),
π π
where [π ] denotes integer part of number π , πΆπ/2 π , πΆ[π/2]π are coefficients of (9),
corresponding to integer partitions ππ/2 π = (2, ..., 2), π[π/2] π = (2, ..., 2, 1).
β β β β
π/2 [π/2]
II. Polynomial of third power
Γ¦(π¦) = ππ¦ 3 + ππ¦ 2 + ππ¦ + π. (11)
After taking into account properties 1)β4) of Γ¦(π¦)-function one obtains π = π =
0, π = 1 β π, 1 β€ π β€ 1.5. Assuming π = π, one gets
Γ¦(π¦) = (1 β π)π¦ 3 + ππ¦, (12)
412
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
where 1 β€ π β€ 1.5 is free parameter equal to the value of derivative Γ¦β² (0) and
determining the value Γ¦β² (Β±1): as π β 1.5 Γ¦β² (Β±1) β 0.
Theorem 2. For approximation (8) with function Γ¦(π¦) of form (12) as π = 1.5
and π2 (π΄)π΄ β 0 as π΄ β β the estimate (1) is valid, but constant
π π/2
π (π) = πΆπ/2 π3 βπ (π/2) β + π(π (π/2) ) (if π is even),
π
π (π) = πΆ[π/2]π π(π΄)3[π/2]+1 βπ ([π/2]) β + π(π ([π/2]) ) (if π is odd),
π π
where πΆπ/2 π , πΆ[π/2]π are the same as in theorem 1.
Theorems 1,2 can be proved using Lemmas 1,2 and taking into account that
values of first derivatives of considered Γ¦-functions in points Β±1 vanishes. The
component βπ (π/2) β appears in the obtained expressions for π (π) because the
second derivative of Γ¦-functions in points Β±1 are not equal to zero. This com-
ponent grows much slower than βπ (π) β.
Remark 2. By changing in given theorems index βπβ on index βπβ, one can
obtain similar results for case of entire π (π₯), namely the similar equalities for
π (π) from estimate (3) when Γ¦(π¦) has form (10), (12).
Let us consider another approach to construction of Γ¦(π¦)-function, it does
not require first derivative of Γ¦(π¦) to be equal to zero in points Β±1, but requires
all the derivatives to be very small in vicinity of Β±1.
III. Function
Γ¦(π¦) = arctan(ππ¦)/ΜοΈπ, ΜοΈπ = arctan π. (13)
Theorem 3. If in expression (8) function Γ¦(π¦) of form (13) has parameter
βοΈ
(π) β
π+πβ2 βπ πΏ
π << (π)
= ππ , βπ = 1, 2, ..., π, (14)
βπ βπ·
then (8) satisfies the estimate of accuracy of best approximation (3) with π (π)
equals for large π and π΄ to the following
(οΈ (π) )οΈ
βπ βπΏ βπ β² βπΏ
π (π) = max , π! + π(π!) + π(π (π) ). (15)
ΜοΈπ(ΜοΈππ)πβ1 ΜοΈπ
Values βοΈ
2 πβ1 βπ (π) βπΏ 1
π β π* = (16)
π βπ β² βπΏ π!
under condition 1.56 < π << min ππ provide maximal rate of decrease of right
π=1,...,π
part of (3) with growth of π and ΜοΈπ(ΜοΈππ)πβ1 -time profit in accuracy in comparison
with the best polynomial approximation.
IV. Function
(οΈ )οΈ
2 1 + exp(βπ)
ΜοΈ
Γ¦(π¦) = π β1 , ΜοΈ =
π . (17)
1 + exp(βππ¦) 1 β exp(βπ)
413
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Theorem 4. Let π (π₯) β πππ (π, π·), π > 1 and in expression (8) function Γ¦(π¦)
of form (17) has parameter, satisfying the inequalities
βοΈ
1 βπ (π) βπΏ
π << π (πΌπ π½π ) = ππ , βπ = 1, ..., π β 1, π << π π = π0 , (18)
πΌπ βπ (π) βπ·
where π (π₯) is the Lambert
βοΈ π -function (the inverse function to π₯ = π€ exp(π€)),
π βπ (π +1) βπΏ
πΌπ = , π½π = 4πΌπ π πβπ . In this case (8) satisfies the estimate of
πβπ βπ (π +1) βπ·
accuracy of best approximation (1) with π (π) equals for large π and π΄ to the
following
(οΈ )οΈ
(π) πβ1 β²
π (π) = max βπ βπΏ (2πΜοΈ π exp(βπ)) , βπ βπΏ ππ! ΜοΈ + π(π!) + π(π (π) ). (19)
Values βοΈ
(οΈ )οΈ
1 πβ1 βπ β² βπΏ π!
π β π* = βπ β (20)
2 βπ (π) βπΏ
under condition 1.57 < π << min ππ provide maximal rate of decrease of
π =0,1,...,πβ1
π))πβ1 -time profit in accuracy
right part of (1) with growth of π and (exp(π)/(2πΜοΈ
in comparison with the best polynomial approximation. Here it is considered the
branch of function π (π₯) such that π (π₯) β ββ as π₯ β 0.
Theorems 3,4 can be proved using Lemmas 1,2 and asymptotical analysis
of growth of πth derivative of π [Γ¦(π¦)] with grows of π. These results can be
extended to case of estimates (1), (3).
5 Computation of approximate values of functions with
boundary value components
In this section one can find an experimental confirmation of results of theorems 1β
4 considering approximation of solution to the boundary value problem (5) with
exponential boundary layer. Detailed comments on correspondence of theoretical
and experimental data given on fig. 3β5 and in Tables 2β5, are absent. Such
correspondence becomes obvious if in theorems 1β4 one supposes the values
π΄1 , ..., π΄πβ1 to be equal to the value of π΄ from formula of solution π (π₯) to
the problem (5). This identification is correct because the values of derivatives
π (π) (Β±1) grow with a rate of π΄π .
To implement the approximation according to (8) four considered types of
function Γ¦(π¦) were used, expressions for inverse functions Γ¦β1 (π₯) were obtained
and the values of zeroes of basis(οΈfunctions
)οΈ of approximations π₯ππ were specified.
ππ¦
1. Sinus (sin), Γ¦π (π¦) = sin :
2
[οΈ ]οΈ
2 arcsin π₯ π cos( (2π+1)π )
Γ¦β1
π (π₯) = , π₯π
π = sin 2π
, π = 0, ..., π β 1.
π 2
414
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
2. Polynomial (pol ), Γ¦π (π¦) = (1 β π)π¦ 3 + ππ¦:
[οΈ ]οΈ
arccos π§ β arccos π§
Γ¦β1
π (π₯) = π
β cos + 3 sin ,
3 3
βοΈ β β
π β3 3π₯ π β 1
π
= , π§= ,
3(π β 1) 2π3/2
[οΈ ]οΈ [οΈ ]οΈ
(2π + 1)π (2π + 1)π
π₯ππ = (1 β π) cos3 + π cos , π = 0, ..., π β 1.
2π 2π
These expressions were obtained using Cardanoβs method and analysis of three
branches of inverse function Γ¦β1
π (π₯).
arctan(ππ¦)
3. Function inverse to tangents (tg), Γ¦π‘ (π¦) = :
arctan(π)
[οΈ ]οΈ
(2π+1)π
arctan π cos( 2π )
tan(π₯ arctan π)
Γ¦β1
π‘ (π₯) = , π₯π
π = , π = 0, ..., π β 1.
π arctan(π)
[οΈ ]οΈ
2
4. Function including exponent (exp), ΜοΈ
Γ¦π (π¦) = π β1 :
1 + exp(βππ¦)
[οΈ ]οΈ [οΈ ]οΈ
1 2 2
Γ¦β1
π (π₯) = β ln β1 , π₯ππ = π
ΜοΈ {οΈ }οΈ β1 ,
π ΜοΈ + 1
ππ₯ exp β π cos( (2π+1)π
2π ) +1
π = 0, ..., π β 1.
5.1 Approximations of boundary layer component π(π₯)
In this section results on numerical approximation of the solution π (π₯) to the
boundary-value problem (5) are given for different values of small parameter
π = 10β3 , ..., 10β10 . Let us consider first the case π(π₯) β‘ 0, then π (π₯) = π(π₯) is
an exponential boundary value component.
For searching the coefficients ππ of the expansion (8), where π = 0, ..., π β 1
a system of π equations should be composed:
πβ1
βοΈ
ππ π΅π (π₯ππ ) = π (π₯ππ ).
π=0
Matrix of this system with elements bjk = π΅π (π₯ππ ) is πΆ = (bjk ). Denoting column
vectors
πΌ = (π0 , π1 , ..., ππ )π , π½ = (π (π₯0 ), π (π₯1 ), ..., π (π₯π ))π ,
one obtains a system of linear algebraic equations (SLAE) πΆπΌ = π½. Finally the
coefficients of expansion (8) can be expressed in form πΌ = πΆ β1 π½.
415
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Remark 3. It was established by computations that the growth of number of
basis function π does not affect the condition number of matrix πΆ. The first five
digits of it 1.4142 remain invariable while π grows from 1 to 200. This means
that using the orthogonal transformations matrix πΆ can be inverted on computer
with high precision.
The algorithm of searching coefficients ππ , π = 0, ..., π β 1 was implemented
on MATLAB programming language for the following values of small parameter
π = 10β3 , 10β4 , ..., 10β10 .
For each of considered types of Γ¦(π¦) a grid of points (π§1 , π§2 , ..., π§πΎ ) was gen-
erated on the segment [β1; 1]. It was used for searching the maximal deviation
of a given function π (π₯) form its approximation ππ (π₯) by enumeration over all
points. Due to the fact that function has boundary layer in the vicinity of points
β1 and 1 a grid for enumeration was composed of large amount of Chebyshev
nodes providing significant concentration of grid near the bounds of segment.
The number of nodes π§1 , π§2 , ..., π§πΎ was chosen so that the doubling of it does not
affect first 2β3 digits of error π = max |π (π§π ) β ππ (π§π )|.
π=1,...,πΎ
β5
For example if π = 10 and function Γ¦π‘ (π¦) is used, the experiments show
that for π = 48
as πΎ = 100 000 point the error π = 2.7737 Γ 10β13 , (21)
as πΎ = 200 000 points the error π = 2.7734 Γ 10β13 . (22)
Consequently, we conclude that number πΎ = 100 000 is sufficient for esti-
mating the error π with desired accuracy.
Now let us use a single coordinate plain to represent dependencies of log10 π
on the number of basis functions π of approximation ππ for different types
of function Γ¦(π¦) together with well known Chebyshev approximations. (see.
fig. 3, 4). Number of points for enumeration πΎ here is equal to a maximal one of
all πΎ obtained for each of considered types. Note that the values of parameters
π, π, π in these results are taken from the vicinity of their βoptimalβ values (see
Table 2). Further method for searching these values is explained.
Table 2. βOptimalβ values of parameters
Parameter Value of small parameter π
of approximation 10β3 10β4 10β5 10β6 10β7 10β8 10β9 10β10
π 1.1 1.2 1.3 1.35 1.45 1.46 1.48 1.48
π 1.2 1.6 3.5 15 35 75 300 680
π 1.3 2.5 3 4.5 5.5 6.8 8.2 9.4
On given figures one can observe that the decrease of values of small pa-
rameter π results in fast decrease of convergence rate of Chebyshev expansions.
Whereas, the methods designed using functions Γ¦π‘ (π¦), Γ¦π (π¦) do not significantly
416
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
a 0
οΏ½οΏ½
log10οΏ½ οΏ½οΏ½οΏ½οΏ½
β5
- tan
β10 - exp
- Chebyshev
- polynom
- sin
β15
0 20 40 60 n
b 0
log10οΏ½ οΏ½οΏ½οΏ½οΏ½
οΏ½οΏ½
β5
β10
0 20 40 60 80 100 n
Fig. 3. Dependencies of log10 π on number π for values π = 10β4 (b), π = 10β6 (b)
a 0
log10οΏ½ οΏ½οΏ½
οΏ½οΏ½οΏ½οΏ½
β5
- tan
- exp
- Chebyshev
β10 - polynom
- sin
n
0 20 40 60 80 100 120
b0
log10οΏ½
οΏ½οΏ½οΏ½
οΏ½οΏ½οΏ½οΏ½
β5
β10 n
0 20 40 60 80 100 120
Fig. 4. Dependencies of log10 π on number π for values π = 10β8 (a), π = 10β10 (b)
417
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
change their convergence rate while decreasing small parameter. As it was al-
ready proved, increase of π, π allows one to reduce the influence of steep gradient.
Now let us estimate βoptimalβ values of parameters π, π, π for approximations
based on Γ¦π , Γ¦π‘ , Γ¦π . To this end the dependencies of logarithm of error log10 π
on π and value of parameter of function Γ¦(π¦) were represented in Cartesian
coordinate system, see those for Γ¦π‘ on fig. 5, 6.
a οΏ½οΏ½
log10οΏ½ οΏ½οΏ½οΏ½οΏ½
0
β5
β10
0
0
20
0.5
n
b 1
40
b
log10οΏ½
οΏ½οΏ½
οΏ½οΏ½οΏ½οΏ½
0
β5
β10 0
20
0
10 40
b
20
30
60
n
Fig. 5. Dependencies of value log10 π on the number π and the value of parameter π of
method with function Γ¦π‘ as π = 10β4 (a), π = 10β6 (b), π = 10β8 (c), π = 10β10 (d)
The error was computed in a similar way as before, e.g. for Γ¦π‘ and π = 10β8
parameter π goes over all integer values from 0 to 100. Maximal convergence rate
is provided by such a value π = π* that specifies maximal slope of a curve laying
in section of the represented surface by a plain π = π* = const(or π = const, or π
= const). For Γ¦π‘ and π = 10β8 π* β [60; 80]. This segment is marked with circle
on the graph of fig. 6.
5.2 Approximation of solution of inhomogeneous problem (5)
Let us consider a function π (π₯) that is solution to (5) with π(π₯) = sin(ππ₯). Here
π(π₯) provides a perturbation in inner points of domain of problem. A question
is how such perturbations can affect convergence of proposed methods.
1 1 1 1
π (π₯) = πΆ1 ππ΄( 2 π₯+ 2 ) + πΆ2 πβπ΄( 2 π₯+ 2 ) + sin(ππ₯),
1 + πβπ΄ 1 + ππ΄ βοΈ (23)
πΆ1 = βπ΄ π΄
, πΆ 2 = π΄ βπ΄
, π΄ = 1/π.
π βπ π βπ
418
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
a οΏ½οΏ½
log10οΏ½ οΏ½οΏ½οΏ½οΏ½
0
β5
0
β10
20
0
20
40
60
40
n
60
b 80
100
b
log10 οΏ½ οΏ½οΏ½οΏ½
οΏ½οΏ½οΏ½οΏ½
0
β5
β10 0
0 20
200 40
400
600
800 1000 60 n
b 1200
Fig. 6. Dependencies of value log10 π on the number π and the value of parameter π of
method with function Γ¦π‘ as π = 10β8 (a), π = 10β10 (b)
Let π be maximal values of deviation of approximation (8) from the values
of function (23). π was computed in numerical experiments by enumeration over
points π§1 , ..., π§πΎ as it was described, for different types of Γ¦(π¦) and various values
of π and π: π = 10β3 , ..., 10β10 . The values of parameters π, π, π were taken from
Table 2 excluding cases of small π. Numerical experiments showed that the
smaller π is, the stronger appears the dependence of convergence rate on π.
While using functions Γ¦π , Γ¦π , Γ¦π approximation errors for (23) agree with
those obtained for π (π₯) = π(π₯) (see fig 3, 4). βOptimalβ values of parameters π
and π in these cases retain too. However the approximation based on Γ¦π‘ loose
its convergence rate much faster. This effect is explained by inequality for π (see
estimate (14)) that turns out to be more essential than, for example, (18). Even
for moderate values of βπ (π +1) βπ· (in our case βπ βπ (π +1) βπ· β 1) one obtains
that value of π should be considerably less than those given in Table 2. Then,
according to (15), π (π) becomes large and the convergence rate decreases.
In order to obtain efficient approximation with function Γ¦π‘ a coupled basis
π΅π (π₯) βͺ ππ (π₯) (π = 0, ..., π¦, π = 0, ..., π , π¦ + π = π β 2) should be used. It is
composed of designed functions π΅π (π₯) and Chebyshev polynomials ππ (π₯):
π¦
βοΈ π
βοΈ
π (π₯) β ππ π΅π (π₯) + ππ ππ (π₯),
π=0 π=0 (24)
]οΈ)οΈ(οΈ [οΈ
ππ (π₯) = cos(π arccos(π₯)), π΅π (π₯) = cos π arccos tan(π₯ΜοΈπ)/π .
419
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
The idea of a method consists in two steps: first function π (π₯) should be
βοΈ
π
expand in basis ππ (π₯) and then difference π (π₯) β ππ ππ (π₯) should be ap-
π=0
proximated using π΅π (π₯). Let us compose matrices
[οΈ ]οΈ
(2π + 1)π
πΆ = (πΎππ ), πΎππ = ππ (π₯π ), where π₯π = cos , π, π = 0, ..., π ;
2(π + 1)
(οΈ )οΈ
[οΈ (2π + 1)π ]οΈ ΜοΈ
B = (πππ ), πππ = π΅π (π₯π ), where π₯π = arctan π cos( ) /π,
2(π¦ + 1)
π, π = 0, ..., π¦.
Introducing the notations of vectors πΌ = (π0 , ..., ππ¦ ), πΎ = (π0 , ..., ππ ), ππΌ =
(π (π₯0 ), ..., π (π₯π¦ )), ππΎ = (π (π₯0 ), ..., π (π₯π )) one can found the coefficients of ex-
pansion in Chebysev basis πΎ = πΆ β1 ππΎ . Further the following SLAE was composed
BπΌ = ππΌ β πΆπΎ,
(οΈ [οΈ {οΈ }οΈ ]οΈ)οΈ
(2π + 1)π
πΆππ = ππ (π₯π ) = cos π arccos arctan π cos( ) /ΜοΈπ .
2(π¦ + 1)
and the values of coefficients π0 , ..., ππ¦ were obtained
[οΈ ]οΈ
πΌ = Bβ1 ππΌ β πΆπΆ β1 ππΎ . (25)
The computations by formula (25) can be modified. So, for small parameters
π = 10β4 , π = 10β6 the domain of distribution of Chebyshev nodes in (24)
was narrowed from segment [β1; 1] to [β0.85; 0.85] in case π = 10β4 and to
[β0.95; 0.95] in case π = 10β6 .
Tables 3β5 show the value of approximation error π obtained for function
(23) using Chebyshev basis and the designed ones π΅π with functions Γ¦(π¦) of
four considered types. If value of parameters of Γ¦(π¦) differs from those given
in Table 2 than it is explicitly indicated in brackets. So does a number of basis
functions π΅π (π₯) plus number of basis polynomials ππ (π₯) for approximations
obtained using Γ¦π‘ , see (24).
Table 3. Values of π obtained for function (23) with π = 10β6
π Chebyshev sin(π¦) πππ(π¦) exp(π¦) tan(π¦)
10 0.997 0.168 0.208 0.069(5.5) 0.036(60, 7 + 3)
20 0.727 0.047 0.064(1.5) 0.005 5.21 Γ 10β4 (80, 14 + 6)
30 0.359 0.007 0.025 2.4418 Γ 10 7.935 Γ 10β6 (50, 23 + 7)
β5
β4
40 0.147 7.659 Γ 10 0.001 3.679 Γ 10β7 4.283 Γ 10β8 (14, 30 + 10)
50 0.051 5.419 Γ 10β5 2.216 Γ 10β5 1.936 Γ 10β9 2.637 Γ 10β10 (14, 38 + 12)
60 0.015 8.019 Γ 10β6 5.348 Γ 10β7 3.738 Γ 10β11 1.436 Γ 10β12 (9, 47 + 13)
70 0.004 7.348 Γ 10β7 2.533 Γ 10β8 5.473 Γ 10β13 2.625 Γ 10β13 (9, 57 + 13)
β4
80 7 Γ 10 3.346 Γ 10β8 7.458 Γ 10β10 3.321 Γ 10β13 2.26 Γ 10β13 (9, 67 + 13)
90 1.2 Γ 10β4 2.941 Γ 10β9 1.791 Γ 10β11 3.375 Γ 10β13 2.597 Γ 10β13 (9, 77 + 13)
420
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Table 4. Values of π obtained for function (23) with π = 10β8
π Chebyshev sin(π¦) πππ(π¦) exp(π¦) tan(π¦)
10 1.0000 0.5024 0.5903(1.5) 0.2040(7.8) 0.0463(650, 7 + 3)
20 1.0000 0.1972 0.6521 0.0366 6.0843 Γ 10β4 (180, 14 + 6)
β4
30 0.9987 0.1119 0.2743 3.5359 Γ 10 1.2506 Γ 10β5 (140, 23 + 7)
β5
40 0.9730 0.0208 0.0873 1.0721 Γ 10 5.3592 Γ 10β8 (110, 30 + 10)
50 0.8920 0.0130 0.0186 3.8726 Γ 10β7 4.9031 Γ 10β10 (90, 39 + 11)
60 0.7705 0.0032 0.0020 8.8276 Γ 10β9 3.4193 Γ 10β12 (90, 47 + 13)
70 0.6384 0.0011 4.28 Γ 10 4.0243 Γ 10β10 1.8532 Γ 10β12 (90, 57 + 13)
β4
80 0.5144 3.7 Γ 10 1.42 Γ 10β4 9.3578 Γ 10β12 1.9960 Γ 10β12 (90, 67 + 13)
β4
90 0.4059 5.7 Γ 10β5 1.27 Γ 10β5 4.3484 Γ 10β12 2.3292 Γ 10β12 (90, 77 + 13)
Table 5. Values of π obtained for function (23) with π = 10β10
π Chebyshev sin(π¦) πππ(π¦) exp(π¦) tan(π¦)
10 1.0000 1.0005 1.0003 0.3875(10.4) 0.0529(6500, 7 + 3)
20 1.0000 0.3498 0.9986 0.0391 4.7344 Γ 10β4 (4000, 13 + 7)
30 1.0000 0.2231 0.9282 0.0027 1.0203 Γ 10β5 (1200, 22 + 8)
β4
40 1.0000 0.2041 0.7389 2.1276 Γ 10 9.2063 Γ 10β8 (1100, 31 + 9)
50 1.0000 0.1418 0.5336 1.5681 Γ 10β5 3.5014 Γ 10β10 (1100, 38 + 12)
60 1.0000 0.0700 0.3642 1.0989 Γ 10β6 2.1094 Γ 10β11 (1100, 47 + 13)
70 1.0000 0.0226 0.2371 7.3964 Γ 10β8 2.2963 Γ 10β11 (1100, 57 + 13)
80 0.9999 0.0223 0.1470 4.8155 Γ 10β9 1.8534 Γ 10β11 (1100, 67 + 13)
90 0.9994 0.0124 0.0865 3.0489 Γ 10β10 2.1706 Γ 10β11 (1100, 77 + 13)
100 0.9973 0.0041 0.0481 4.0388 Γ 10β11 2.1034 Γ 10β11 (1100, 87 + 13)
6 Conclusion
In this work a method for approximating smooth functions having bound-
ary layer components was developed, justified and implemented. It is based on
expansion of function into series with basis consisting of non-polynomial func-
tions obtained from trigonometric Fourier one by special mapping operation.
Such representation ensures the estimations of accuracy of best polynomial ap-
proximations, but essentially reduces the values of coefficients in them. That is
why convergence can be observed starting from small number of basis elements.
Moreover the proposed approximations have good properties of numerical sta-
bility inherent to Fourier expansions.
Further development and successful application of the method is connected
with analysis of different forms of Γ¦-function. Proposed method can be modified
for approximation of functions having singularity in inner point of domain, or
even for problems with unknown position of singularity. From the other side,
combination of these approximations with collocation methods will allow to de-
sign efficient algorithms for solving singularly-perturbed boundary value prob-
lems for differential equations.
421
Mathematical and Information Technologies, MIT-2016 β Mathematical modeling
Acknowledgments. This work was supported by Integrated program of basic
scientific research of SB RAS No. 24, project II.2 βDesign of computational
technologies for calculation and optimal design of hybrid composite thin-walled
structuresβ.
References
1. Kadalbajoo, M. K., Gupta, V.: A brief survey on numerical methods for solving sin-
gularly perturbed problems. Applied Mathimatics and Computation. 217(8), 3641β
3716 (2010).
2. Kreiss, H.-O., Nichols, N. K., Brown, D. L. Numerical methods for stiff two-point
boundary value problems. SIAM J. Numer. Anal. 23 (2), 325β368 (1986).
3. Liseikin, V. D., Likhanova, Yu. V., Shokin, Yu. I. Numerical grids and coordinate
transformations for the solution of singularly perturbed problems. Nauka. Novosi-
birsk (2007, in Russian).
4. Liu., W., Tang, T. Error analysis for a Galerkin-spectral method with coordinate
transformation for solving singularly perturbed problems. Appl. Numer. Math. 38,
315-345 (2001).
5. Wang, Yi., Chen, S. A rational spectral collocation method for solving a class of
parameterized singular perturbation problems. J. of Comput. and Appl. Math. 233,
2652-2660 (2010).
6. Babenko, K.I.: Fundamentals of numerical analysis. Regular and chaotic dynamics,
Moscow-Izhevsk (2002, in Russian).
7. Jackson, D.: On Approximation by Trigonometric Sums and Polynomials. Trans.
Amer. Math. Soc. 13, 491β515 (1912).
8. Dzydyk, V. K.: Introduction to the Theory of Uniform Approximation of Functions
by Means of Polynomials. Nauka, Moscow (1977, in Russian).
9. Bakhvalov, N.S., Zhidkov, N.P., Kobelβkov, G.M.: Numerical methods. Textbook.
Nauka, Moscow (1987, in Russian).
10. Boyd, J.: Chebyshev and Fourier Spectral Methods. Second Edition. University of
Michigan (2000).
11. Orszag, S. A., Israeli, M.: Numerical simulation of viscouse incompressible flows.
Ann. Rev. Fluid Mech. 6, 281β318 (1974).
422