Inverse central ordering for the Newton interpolation formula

An inverse central ordering of the nodes is proposed for the Newton interpolation formula. This ordering may improve the stability for certain distributions of nodes. For equidistant nodes, an upper bound of the conditioning is provided. This bound is close to the bound of the conditioning in the Lagrange interpolation formula, whose conditioning is the lowest. This ordering is related to a pivoting strategy of a matrix elimination procedure called Neville elimination. The results are illustrated with examples.


Introduction
The Lagrange interpolation operator associates to each function its Lagrange interpolating polynomial of degree less than or equal to at 1 distinct nodes 0 on the interval . Using the Lagrange interpolation formula, the Lagrange interpolation operator can be written in the form , where , 0 . In formula (2) of [5], a conditioning associated to a representation of the form 0 is introduced, cond 0 (1) where 0 are functionals belonging to the space generated by the evaluation functionals 0 and 0 is a basis of , the space of polynomials of degree not greater than . By formula (4) and Theorem 4 of [5], we have 0 cond cond (2) that is, the conditioning of the Lagrange representation coincides with the Lebesgue function and it is lower than the conditioning of any other representation, and in particular, than the Newton representation. This representation is given by In order to compare the conditionings of the Lagrange and the Newton representations, we recall the following asymptotic formula for the Lebesgue constant shown by Schönhage in [10] max 2 1 log (5) where 0.5772156649 is the Euler-Mascheroni constant. Newton's formula is sensitive to the ordering of the nodes in floating point arithmetic. An analysis of this fact can be found in [4,5]. In [5] it is shown that max cond 3 for equidistant nodes in increasing order. In [6], a central ordering of the nodes, based on the distances of the nodes to a given center was considered. For equidistant nodes and the choice of the center 2, it was proved that max cond 1 2 1 2 2 .
In this paper, we propose the inverse central ordering, also based on the distance of the nodes to a center. The bound for the conditioning of the Newton representation with this ordering is smaller than the two previous bounds. More precisely, in Theorem 3, for a sequence of equidistant nodes in 2 2 of the form 2 for even 1 2 for odd 0 with , we show that max cond 7 2 .
In Section 2, we introduce the inverse central ordering. Bounds for the conditioning of the Newton representation are provided. The bounds for the inverse central ordering are close to the maximum conditioning of the Lagrange representation, which has the best conditioning. In this sense, the inverse central ordering for equidistant nodes is near optimal. Section 3 contains numerical experiments to illustrate the properties of the conditioning of the Newton formula with nodes following an inverse central ordering. In Section 4, Neville elimination, a matrix elimination procedure alternative to Gaussian elimination, is considered. Neville elimination is especially useful when dealing with some structured matrices (see [2,7,8]) including Vandermonde matrices. We show that the inverse central ordering corresponds to partial pivoting for Neville elimination in Vandermonde matrices.

Inverse central ordering and conditioning
The inverse central order consists in arranging a set of nodes, starting with the furthest node to a given center and finishing with the closest one. .
We now consider equidistant nodes and, for the sake of simplicity, we take symmetric intervals and the center 0. Therefore, the interval has the form 2 2 , where is the degree and is the distance between neighboring nodes. This ordering may not be unique. So, we propose the following choice: in the case where two nodes lie at the same distance from the center, we take first the 1 1 .
Therefore, we have By formula (10), we have Using Pascal's identity, and formula (11) we obtain We proceed in a similar way for the divided difference functionals of odd order. By formulae (11) and (10), we have We obtain in the following result explicit formulae for the -norm of the divided difference functionals.
Theorem 1 Let 0 be the nodes given by (6) and , 0 , the divided difference functionals (8). We have Next we derive some formulae for the conditioning of the Newton formula with equidistant nodes and inverse central ordering. Hence the first formula in the statement of Lemma 1 holds. Now assume that 2 1 , for some 0 1 2 . Substituting 2 , 2 1 , we can write In order to obtain a bound for cond 0 , we will use the well-known Vandermonde identity 0 0 .
We will also use the following relation between generalized binomial coefficients: Let us obtain an upper bound for cond at the interpolation nodes.
Theorem 2 Let 0 be the nodes given by (6). Then, we have By (9), all denominators are greater than or equal to 1 in the above formula, except in the case where is even and 2. We use (9) again to show that So taking 2, we obtain Using Pascal's identity (12) and Vandermonde identity (15), we deduce that For the rest of the cases, 1 2, all the binomial coefficients in the denominators are larger than or equal to 1 by (9). Applying Pascal's identity (12) and (15) If is odd and 1 2, there is a term with denominator Then we deduce from (12) and (15)  So, the first formula follows. The second one follows from the fact that the sequence 1 1 2 , 0 , is nondecreasing by Lemma 2.

Remark 2
Recall that the Lagrange representation is optimal in the sense that it provides the least conditioning (see (2)). Corollary 6 of [5] shows that the ratio between the conditioning of the Newton representation and the conditioning of the Lagrange representation is attained at a node. Using this result and Theorem 2, we obtain the following inequality for equidistant nodes arranged according to the inverse central ordering (6): We observe that the maximum bound for the condition number at the nodes is attained for the nodes closest to the origin.
In the following theorem we bound cond .
Theorem 3 Let 0 be nodes given by (6). Then Proof First assume that 2 2 , 0 2 , and let us define Observe that, by (9), the binomial coefficients 1 2 are always larger than or equal to 1, except if 2 and is even. Using (16), the corresponding term can be bounded from above by By (9), all binomial coefficients in the denominators are larger than or equal to 1 and we derive the following inequality The terms corresponding to binomial coefficients in the denominators smaller than 1 in the first sum can be bounded again by 2 . For the second sum, we first consider the case where is odd and 1 2, 1 2. By (16) and (18) Since, by (9), all binomial coefficients in the denominators are larger than or equal to 1, we deduce that Observe that the bound obtained in Theorem 3 is close to the formula (5) for the Lebesgue constant. Therefore, the inverse central ordering is a convenient strategy to order the nodes in the Newton formula, in addition to easy implementation.

Numerical experiments
We have computed the conditioning of the Newton formula for the inverse central ordering. Figures 1 and 2 show a comparison between the Lebesgue function and cond at equidistant nodes in the inverse central ordering. We see that cond is close to the Lebesgue function at the extremities of the intervals but the difference between both functions grows in a neighborhood of the origin. We can also see that the maximum of the conditioning is attained in a neighborhood of the center for 10 and that the largest value is attained in a neighborhood of the extremities of the interval for 19. Let us compare the condition of the Newton formula for equidistant nodes for the increasing ordering and an inverse central ordering. Figure 3 (left) shows the condition for nodes in increasing order for degree 70. Figure 3 (right) shows the condition for nodes following an inverse central ordering for the same degree. We see that the condition is much lower for the inverse central ordering. In each case, we have also computed the numerical error for interpolating the smooth function 1 sin 2 , where 2 24 is the unit roundoff in single precision. The numerical error has been computed subtracting the value of and the computed value of by Newton's formula in single precision arithmetic. We see in both cases that the error is lower than the bound and even imitates its shape. We also see that the highest numerical error is found for the increasing order in a neighborhood of the right extremity of the interval.
We have also checked other distributions of nodes. Generally the conditioning of the Newton formula near the extremities of the interval is lower with the inverse central ordering that with other orderings. In Fig. 4, we have compared the conditioning of the Newton formula for Chebyshev nodes using different orderings. The worst condition corresponds to the nodes in increasing order. The inverse central ordering leads to a relatively worse conditioning in the center of the interval as compared with the Lebesgue function. In contrast, the Leja order remains closer to the Lebesgue function.

Inverse central ordering and Neville elimination
In this section, we interpret the inverse central ordering in terms of a matrix elimination with a pivoting strategy. Specifically, we show a connection between the inverse central ordering for equidistant nodes and Neville elimination with partial pivoting. Neville elimination (NE) is an elimination algorithm alternative to Gaussian elimination where one subtracts to each row a multiple of the previous one (see [7]). Let is the -multiplier of the NE of , for 0 . The NE with partial pivoting is a strategy (see [1]) in which the pivots corresponding to the -th column of form a nonincreasing sequence: In the following we show that, for certain ordering of the nodes, NE with partial pivoting does not require row interchanges of the associated Vandermonde matrix. The Leja order satisfies a similar property for the Gaussian elimination with partial pivoting. In fact, Gaussian elimination with partial pivoting leads to an ordering of rows corresponding essentially to the Leja order (see Section 22.3.3 of [9]).
Taking into account that the nodes form an increasing sequence, as well as (25)  Taking into account that the distances form a nonincreasing sequence, we have . Therefore, by (25), we deduce that 1.
If the distances between adjacent nodes are greater in the center than close to the extremities, as in the case of Chebyshev nodes, then the hypotheses of Theorem 5 do not hold. This fact suggests that the inverse central ordering is not a good ordering for dealing with Chebyshev nodes, in agreement with the observations on the conditioning illustrated by Fig. 4. If the distances between adjacent nodes decrease as we approach the center of the interval, then the hypotheses of Theorem 5 hold, which implies more stable computations.
In the following result we consider nodes following the inverse central ordering (6). We show that NE with partial pivoting does not require row interchanges. Furthermore, we prove that the double sequence 0 is bimonotonically nonincreasing or, equivalently, Clearly, the 1 form a nondecreasing sequence in each of the indices and and, so, 0 is a bimonotonically nonincreasing sequence.
Funding Open Access funding provided thanks to the CRUE-CSIC agreement with Springer Nature. This work has been partially supported by the PGC2018-096321-B-I00 Spanish Research Grant, by Gobierno de Aragón E41 17R and Feder 2014-2020 "Construyendo Europa desde Aragón".

Conflict of interest The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.