Effect of Asymmetric Finite Difference Formulas on the Orders of Central Difference Approximations for the Second Derivative of a Periodic Function ()
1. Introduction
Numerical differentiation is an elementary and essential part in scientific modeling and numerical analysis [1] . It is widely used for the differentiation of the functions which are difficult to differentiate analytically, and for finding the derivative of the sampled data for which the generating function is not known. This is the case when solving the ordinary or partial differential equations numerically.
Many methods have been introduced and discussed to determine the derivatives numerically [2] . These methods can be classified into two approaches [2] . The first approach aims to develop formulas for calculating the derivatives numerically. This includes the Taylor expansion based finite difference method [3] [4] [5] [6] [7] , the operator method [3] [8] , the interpolating polynomial method [4] [5] [6] [8] [9] , and the lozenge diagram method [10] . The second approach does not give an explicit formula for the derivative; it just aims to evaluate it by using the function data [1] .
The Taylor series based finite difference approximation is used to numerically evaluate the derivative of a function at a grid reference point by using data samples at the other neighboring points within the domain of this function. Depending on the pattern of the samples used in calculation, the most commonly used approximation formulas for derivatives are classified as forward, backward and central difference formulas. Forward difference approximations use the samples at a mesh point and next (forward) equally spaced points of analysis for calculating the derivative at the mesh point. In contrast, backward difference approximations use the samples at a mesh point and the previous (backward) equally spaced points, whereas central difference approximations use both forward and backward samples in calculating the derivative at the specified mesh point. Forward and backward differencing are part of a larger group of what might be called lateral (as distinct from central) differencing methods, including extrapolative differencing that might be appropriate for singularity [11] .
We have to emphasize that the number of points used to approximate the derivative defines the order of the approximation and generally, the greater is this number, the more is the accuracy of the approximation. For example, a finite difference approximation of the second derivative using
points is said to be of order
for the central approximation and
for both the forward and backward approximations, which means that central difference formula gets an extra order of accuracy for free. The central differencing is therefore the most accurate of the above three methods for a given number of equidistant data samples [12] . However, central difference formulas are useful only for interior points not for end points. For example, if we have a grid of
points
given by
,
for some fixed tabular interval h, and some integer
, the
-point central difference approximation for the second derivative of a certain function
can be used only at grid points
,
but not for
or
,
being an odd integer such that
. In order to get approximations of
at the
first points (resp. last points) while preserving the accuracy of the results, one can theoretically use a
-point forward difference (resp. backward difference) formula, which must surely affect the accuracy and precision of the results.
The main purpose of this paper is to study the dependence of the orders of central difference approximations for the second derivative of a function on the oscillatory properties of this function and the value of the sampling period h in the case where it is necessary to use forward/backward approximations for this derivative at some points belonging to a given grid of equally spaced points.
The rest of this paper is organized as follows. Section 2 contains a brief presentation of Khan and Ohba’s closed-form expressions for finite difference approximations of first and higher derivative based on Taylor’s series, which have been shown to be suitable for computing derivatives of any degree with arbitrary order of accuracy over all the known function sampling points with a minimum effort [13] . Thereafter, we use these expressions to establish forward, backward and central difference formulas of several orders for the first two derivatives. In Section 3, we show how to determine numerically and graphically the order of a finite difference formula for a derivative of a certain degree q. In Section 4, we present and comment our results concerning the orders of some finite difference formulas for the second derivative of the function
defined by
, α being a positive real parameter which allows working with various oscillation frequencies of the function under consideration. We concentrate on the behavior of the order of each formula as a function of the α parameter and the sampling period h. The conclusion is given in Section 5.
2. Finite Difference Formulas of the First and Second Derivatives Based on Taylor Series
Following the notations as found in [13] , the Taylor series defines the relation between the discrete time values of a time function
sampled at
, where
, and h is the sampling period, to the value of the function and its derivatives at origin
. It can, mathematically, be written as [13]
(1)
where
denotes the value of
at
,
denotes the value of the kth derivative of f at
and
is a term of the order of
coming from the truncation of the series after
terms. The value of the kth derivative of f at
is denoted by
.
In this Section, we are mainly interested in the forward, backward and central difference approximations of the first and second derivatives of a function. We found it clarifying to use the notations of Ishtiaq Rasool Khan and Ryoji Ohba [13] . The advantage of this notation is that it is suitable for calculating derivatives of any degree with arbitrary order of accuracy over all the known function sampling points with a minimum effort.
The following paragraphs summarize a number of expressions for finite difference approximations of first and second derivatives found in [13] .
2.1. Case of the First Derivative
Based on forward values of the function, we have [13]
(2)
where
denotes the coefficient of a term
in a forward approximation of first derivative of order n.
Ishtiaq Rasool Khan and Ryoji Ohba calculated the
coefficients for different orders of approximation, and observed that they can be expressed by the following explicit formulas [13] :
(3)
and
(4)
where
is defined as
.
Similarly based on the backward values of the function, we may write
(5)
where the superscript B refers to backward difference approximations. The
coefficients are such that
(6)
and
(7)
According to Khan I. R. and Ohba R. [13] , a central difference approximation of order 2n for the first derivative is, for its part, given by the following Equation:
(8)
where
(9)
and
(10)
The approximations of the first derivative for some values of n are listed below (See Table 1).
2.2. Case of the Second Derivative
For this case we have [13] :
Table 1. Backward, forward and central difference approximations of
for some values of n.
• Forward formulas
(11)
where
(12)
and
(13)
• Backward formulas
(14)
where
, the coefficient of
, is the same as the corresponding kth coefficient in forward difference approximation given by Equations (12) and (13).
• Central difference formulas
(15)
where
(16)
and
(17)
Table 2 presents finite difference approximations of the second derivative for some values of n.
3. Numerical and Graphical Determination of the Order of a Finite Difference Formula
It is possible to calculate numerically and graphically the order of a given finite difference formula for a derivative of a function. To achieve this, it is necessary to choose a function whose successive derivatives are exactly known in a certain interval
of the independent variable. We then consider several values of h and compute for each of them the error
defined as
(18)
where
denotes the value of the qth derivative of the function f at
,
is a positive integer such that
, (19)
which means that we discretize the interval of interest by
points, i.e.,
. If
Table 2. Forward, backward and central difference approximations of
for some values of n.
(20)
for h sufficiently small (as
), then p is the order of the finite difference formula. Here, C is a constant depending on the function, the degree of the derivative and the finite difference formula. From Equation (20), we have
. (21)
Using the logarithmic scale resulting from the change of variable
and
, we obtain the expression
, (22)
which is the equation of a straight line. The slope of this line is therefore the order of the finite difference approximation.
4. Numerical Experimentation
In order to study orders of finite difference formulas for the second derivative, we consider the periodic function
defined by
. (23)
The parameter α allows working with various oscillation frequencies of the function under consideration.
Figure 1 shows the dependence of the function
on the independent variable x and the parameter α. As expected, the number of oscillations increases with α.
In Figures 1-5, we show the errors
plotted against h on a log-log scale for several values of the α parameter and the three numerical differentiation methods mentioned above. This is a good way to plot errors when we expect
Figure 1. Plots of
for several values of the parameter
.
them to behave like some power of h, since if the error
behaves like
then
. So on a log-log scale, the error behaves linearly with a slope that is equal to p, the order of accuracy.
It must be remembered that each of the three considered difference methods using np points to approximate
requires us to use at least one of the two other methods in the approximation of this derivative for at least one point. More precisely:
• A np-point forward formula requires use of a np-point backward formula to approximate
at the
last points.
• A np-point backward formula requires use of a np-point forward formula to approximate
at the
first points.
• A np-point central difference formula requires use of the
-point forward formula (resp. backward formula) for the
first points (resp. last points).
It appears from Figures 2-5 that the shape of each plot depends on the range the sampling period h belongs to, the method used to approximate
at the majority of points (those for which the method is useful) and the value of the α parameter. Indeed, the potential values of h can obviously be placed in three different regions when a
-point central difference formula is applied to the function
with α positive integer or half-integer and
(See Figure 2, Figure 3(a) and Figure 3(b), Figure 4 and Figure 5(a) and Figure 5(b)): the small-h region where the curve is a straight line with slope 2n, the medium-h region where we have a straight line of slope
and the high-h region where the slope of the curve is not where defined. The length and endpoints of each region depend on the values of the integer n and the parameter α, which implies that sizes and endpoints of the three regions are dependent on both the oscillation frequency of the function and the number of points used in the approximation of the derivative. With the same function, but values of α which are neither integers nor half-integers, only two regions are observed: there is no medium-h region (See Figure 3(c) and Figure 3(d) and Figure 5(c) and Figure 5(d)). We have to add that all results obtained from forward differencing formulas fit with those of the corresponding backward formulas. Furthermore, no cases of medium-h region were reported for these two methods.
In an effort to better understand the origin of the medium-h region, we show the errors
(Figure 2, Figure 4, Figure 3(a) and Figure 3(b) and Figure 5(a) and Figure 5(b)) defined as
(24)
for a
-point central difference approximation of the second derivative. We have to emphasize that
is obtained from Equation (18) by omitting all points at which the central difference approximation for the second derivative is not valid. It is remarkable that any case of medium-h region is observed for
, which means that this region is simply the consequence of asymmetric finite difference approximations of the derivative in question at some points (those at which the central difference formula is not valid).
Now consider the problem of estimating the slope of the line segment connecting two any consecutive points on each of the curves of Figures 2-5 for only central differencing formulas and integer or half-integer values of the parameter α. We used the expression
to denote the value of the slope for the line segment whose endpoints are
and
:
(25)
The plots of the functions highlighting
for several values of the parameter α are shown in Figure 6 and Figure 7. We see the same three regions appearing in Figures 2-5. The numerically calculated starting values of the sampling period for the small-h region (
for
and
for
) are given in Table 3 and Table 4, the values of h being listed in descending order.
A careful analysis of our numerical data presented in Table 3 leads to the conclusion that the relationship between
and α (integer or half-integer) can be described by the following simple formula:
(26)
for the 5-point central difference approximation of
, where
. (27)
From the numerical data of Table 4, we found that the small-h region associated with the 7-point central difference approximation of
is such that
(28)
with
(29)
Table 3. Starting values of the sampling period for the small-h region in the case where n = 2 (5-point central difference formula). Fifteen values of the parameter
are considered.
Table 4. Starting values of the sampling period for the small-h region in the case where n = 3 (7-point central difference formula). Six values of the parameter
are considered.
5. Conclusions
The focus of our work was to investigate the influence of forward and backward finite difference approximations on the orders of central difference approximations for the second derivative of a periodic function. First, we have briefly presented very practical closed-form expressions for the finite difference approximations of first and second derivatives based on Taylor series and found in the literature. Next, these expressions were used to obtain forward, backward and central difference formulas of several orders for the first two derivatives of any function of one variable. Thereafter, we have described one major method to determine numerically and graphically the order of a finite difference formula for the second derivative of a certain degree q. Finally, we have applied the method in question to the second derivative of the periodic function
with
and
. The results obtained from this application led to the conclusion that any
-point central difference formula for the second derivative of the considered periodic function is of order
for medium values of the sampling period h and of order 2n for small values of this parameter in the case where it is necessary to use forward and backward approximations of the derivative of interest at some points belonging to a given grid of equally spaced points. Sizes and starting points of the small-h and the medium-h regions depend on h and α. There is no medium-h region when the parameter α is neither integer nor half-integer. Careful analysis of our numerical results led us to very simple formulas to describe relationships between the starting values of h for the small-h region and the α parameter in the cases of
and
, the values of h being listed in descending order.
One has to add that errors
defined by Equation (18) become too small to be well calculated when h is of very small magnitude and the parameter α not large (
), unless the computation is performed with sufficiently high-precision. To deal with this problem, we appealed to a new software package for arbitrary precision computation, named MPFUN2020 and developed by David H. Bailey [14] .