A method for approximating smooth functions has been developed using non-polynomial basis obtained by mapping of Fourier series domain to the segment [−1, 1]. High rate of convergence and stability of the method is justified theoretically for four types of coordinate mappings, 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 polynomials 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 proposed 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)
By now, a huge amount of urgent scientific and technological problems reduces to boundary value problems for differential equations having pronounced singularities 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 analysis of key issue concerning the smoothness of such transformations and its influence on the quality of approximation method is absent. Frequently, authors 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 segment [−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 effectiveness 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. Norms of function ‖ ‖ = max | |, ‖ ‖ = ∈ ︂( ︀∫ ( ) ︂) 1 .
2. Finite-dimensional approximating space with elements used for approximation of ( ) that usually are series or polynomials including summands 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 providing the lower bound to be reached *(, ) = = *( ) = ∈ ( ) is element (, ) → ).
∈ inf ‖ ) ∈ − ‖.
5. Method without saturation (loose definition) is a method of approximation of function ( ) having asymptotic of error of the best polynomial approximation 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 allow 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
) is -times continuously differentiable function on segment and all its derivatives up to order are bounder by value ( ), then sup ∈ (, ) *(, ) ≤ ( ) − , 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 here <
1 is defined by location of singularity in complex plain, is constant, [8]. III. If ( ) ∈
is entire function, then ∞ =0 ︁∑ ∞ =0 ‖ ( ) − ( )‖ ≤ , ∈ , *(, ) ≤ ( )(diam ) 21−2 ! , interpolation (see [9]) and Cauchy–Hadamard inequality. where ( ) = ‖ ( )
‖ = ( !). It follows from estimates of error of polynomial Remark 1. For smooth functions of I and III classes the accuracy of best approximations depends on maximal values of derivatives of function on (values of ( ) and
( ) in (
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( )) = ( ).
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
approximations presents the methods without saturation and therefore they are
extremely efficient for solving problems with smooth solutions. However, if values
( ), , ( ) in (
A simple example is a boundary-value problem for differential equation of second order with small factor of second derivative
2
2 − = ′′( ) − ( ), (−1) = 1, (
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
approximations (
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 (
∈[−1,1] =0 a 0
x 0.5 1 0 −5 −10 −150 b
log10 20 40 60 80
−0.5 1 0 0.5 −0.5 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/√ 3|1 − 1| ≤ or ≈ 2√2√
3 3
Description of a method Chebyshev polynomials for approximation of unknown function. In the considsion √ ered case ≈ 5.6/
in the denominator of fraction is caused by concentration of Chebyshev 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 empirical requirement that even three nodes should lay on the boundary layer gives (here condition ( ) = 0.1 is used). Appearance of expres
that corresponds to (
Let us assume that
= [−1, 1]. Define a function = ae( ) : [−1, 1] → [−1, 1]
1) function ae( ) is bijective mapping, ae(
For approximation of function ( ), ∈ [−1, 1] let us consider the expansion As a result one obtains (ae( )) = (ae(cos )) ≈
cos( ). ︁∑ −1 =0 ( ) ≈ ( ) =
cos[ arccos(ae−1( ))].
︁∑
−1
=0
(
, ,
= 0, ..., − 1. Note, that under the properties 1, 3,
Lemma 1. Approximation (
= ( ), , = 0, ..., − 1 do not depend
The proof of Lemma 1 is obvious taking into account (
To settle a key question on rate of growth of coefficients ( ), ( ) in
estimates (
( )(ae( )) ︁∑ constants. Moreover ∀ as = 1 and = one has: where the second summation goes over all possible integer partitions of number be proved using mathematical induction technique. 4
Analysis of four types of function ae( ).
< Now let us propose and investigate four types of function ae( ). To this end the following notations are necessary. Let ∈ , ∈ boundaries of segment
representing boundary layers, 0 be a certain neighborhood 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
order of smoothness of ( ). Typically for problems with boundary layer one has
>> 1. Let
= max , ( ) be a size of boundary layer (as it follows from
be neighborhoods of
the example with exponential boundary layer, depends on , see comments to
formula (
(
. ae( ) = 3 + 2 +
+ .
ae( ) = (1 − ) 3 + ,
Theorem 1. Let 2( )
function ae( ) of form (
‖ + ( (/ 2)) (if is even), ( ) = [/ 2] ( )(/ 2) ‖ ([/ 2])‖ + ( ([/ 2])) (if is odd), where [ ] denotes integer part of number , / 2 , [/ 2] are coeficients of corresponding to integer partitions / 2 = (2, ..., 2), [/ 2] = (2, ..., 2, 1). ⏟
⏞ / 2 ⏟ [/ ⏞2]
After taking into account properties 1)–4) of ae( )-function one obtains = =
0, = 1 − , 1 ≤ ≤ 1.5. Assuming = , one gets
determining the value ae′(±1): as → 1.5 ae′(±1) → 0.
where 1 ≤ ≤ 1.5 is free parameter equal to the value of derivative ae′(0) and
Theorem 2. For approximation (
Theorems 1,2 can be proved using Lemmas 1,2 and taking into account that component ‖ (/ 2)
‖ appears in the obtained expressions for
values of first derivatives of considered ae-functions in points ±1 vanishes. The
( ) because the
second derivative of ae-functions in points ±1 are not equal to zero. This
comRemark 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 (
Let us consider another approach to construction of ae( )-function, it does all the derivatives to be very small in vicinity of ±1. not require first derivative of ae( ) to be equal to zero in points ±1, but requires
Theorem 3. If in expression (
︃√ + −2 ‖
( )‖ ‖ ( )‖
= , ∀ = 1, 2, ..., , ae( ) = ︀̃ ︂( 2
︂) 1 + exp(− ) − 1 , = ︀̃ 1 + exp(− ) 1 − exp(− ) .
Values
then (
︂) ‖ , ‖ ′‖ ! + ( !) + ( ( )).
︀̃ (̃︀ ) −1 ≈ * = ︃√ ︀̃ 2 −1 ‖ ( )
‖ 1 ‖ ′‖ ! under condition 1.56 < << =1,..., provide maximal rate of decrease of right
min
with the best polynomial approximation.
part of (
Theorem 4. Let ( ) ∈ (, of form (17) has parameter, satisfying the inequalities
), > 1 and in expression (
(18) ︃√ ‖ ( )‖
‖ ( )‖ ︃√ ‖ ‖ ( +1)‖ where ( ) is the Lambert -function (the inverse function to = exp( )), , = 4 −
( +1)‖ . In this case (
− ︂) ︂( (20) ( ) = max ‖ =0,1,..., −1 provide maximal rate of decrease of
min
right part of (
Theorems 3,4 can be proved using Lemmas 1,2 and asymptotical analysis
of growth of th derivative of [ae( )] with grows of . These results can be
extended to case of estimates (
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 (
To implement the approximation according to (
: ︂(
2 ae−1( ) = 2 arcsin , = sin ︂[ cos( (2 +1) ) ]︂
2 2 2. Polynomial (pol ),
ae ( ) = (1 − ) 3 + : ae−1( ) = =
︂√ = (1 − ) cos3︂[ (2 + 1) ]︂ 3( − 1) − cos arccos 3 + √ 3 sin arccos ]︂ 3
, , = −3√3 √ − 1
, branches of inverse function ae−1( ).
3. Function inverse to tangents (tg ),
ae ( ) = ︀̃ 4. Function including exponent (exp), ae ( ) = arctan cos( (2 +1) )︂] ︂[ arctan( ) ae−1( ) = − = 0, ..., − 1.
a system of equations should be composed: an exponential boundary value component.
In this section results on numerical approximation of the solution ( ) to the
boundary-value problem (
( ) = ( ).
vectors
Matrix of this system with elements bjk = ( ) is = (bjk). Denoting column
= ( 0, 1, ..., ) , = ( ( 0), ( 1), ..., ( )) ,
one obtains a system of linear algebraic equations (SLAE)
coefficients of expansion (
The algorithm of searching coefficients , = 0, ..., − 1 was implemented on MATLAB programming language for the following values of small parameter 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, ..., affect first 2–3 digits of error = was chosen so that the doubling of it does not =1,..., | ( ) − ( )|.
max
For example if = 10−5 and function ae ( ) is used, the experiments show that for = 48 as as = 100 000 point = 200 000 points the error = 2.7737 × 10−13, the error = 2.7734 × 10−13.
(21) (22) Consequently, we conclude that number = 100 000 is sufficient for estimating 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 ae( ) 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
On given figures one can observe that the decrease of values of small parameter results in fast decrease of convergence rate of Chebyshev expansions. Whereas, the methods designed using functions ae ( ), ae ( ) do not significantly 0 20 40 60 80 100 n change their convergence rate while decreasing small parameter. As it was already proved, increase of , allows one to reduce the influence of steep gradient.
Now let us estimate ”optimal” values of parameters , , for approximations based on ae , ae , ae . To this end the dependencies of logarithm of error log10 on and value of parameter of function ae( ) were represented in Cartesian coordinate system, see those for ae on fig. 5, 6.
a b 0 −5 −10 0 0 −5
The error was computed in a similar way as before, e.g. for ae 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 ae and = 10−8 * ∈ [60; 80]. This segment is marked with circle on the graph of fig. 6.
Approximation of solution of inhomogeneous problem (
0 b 0 −5 −10 log10 20 log10 40 b 60 80 100 60
20 0 200 400 600 800 1000 1200 b 60 40 n method with function ae as = 10−8 (a), = 10−10 (b)
Let be maximal values of deviation of approximation (
as it was described, for different types of ae( ) and various values of and : = 10−3, ..., 10−10. The values of parameters , , were taken from smaller is, the stronger appears the dependence of convergence rate on .
While using functions ae , ae , ae
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 ae 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 ae 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( ̃︀ )/ ]︀) . The idea of a method consists in two steps: first function ( ) should be expand in basis ( ) and then difference ( ) − ∑︀ proximated using ( ). Let us compose matrices = ( ),
= ( ), where = cos B = ( ), = ( ), where = arctan [︀ cos(
︂( , = 0, ..., . =0 pansion in Chebysev basis = −1 . Further the following SLAE was composed B = − , = ( ) = cos arccos arctan cos( ︂( ︂[ ︂{ and the values of coefficients 0, ..., were obtained (2 + 1) }︂ 2( + 1) ) / ︀̃ ︂]) . = B−1[︀ − −1 ︀] .
(25) [−0.95; 0.95] in case = 10−6.
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 (23) using Chebyshev basis and the designed ones with functions ae( ) of four considered types. If value of parameters of ae( ) 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 ae , see (24).
( ) 100 0.9973 6
In this work a method for approximating smooth functions having boundary layer components was developed, justified and implemented. It is based on expansion of function into series with basis consisting of non-polynomial functions obtained from trigonometric Fourier one by special mapping operation. Such representation ensures the estimations of accuracy of best polynomial approximations, 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 stability inherent to Fourier expansions.
Further development and successful application of the method is connected with analysis of different forms of ae-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 design efficient algorithms for solving singularly-perturbed boundary value problems for differential equations. 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”.