# Finite difference and numerical differentiation: General formulae from deferred corrections

This paper provides a new approach to derive various arbitrary high order finite difference formulae for the numerical differentiation of analytic functions. In this approach, various first and second order formulae for the numerical approximation of analytic functions are given with error terms explicitly expanded as Taylor series of the analytic function. These lower order approximations are successively improved by one or two (two order improvement for centered formulae) to give finite difference formulae of arbitrary high order. The new approach allows to recover the standard backward, forward, and centered finite difference formulae which are given in terms of formal power series of finite difference operators. Examples of new formulae suited for deferred correction methods are given.

## Authors

• 2 publications
• 4 publications
• ### High-order accurate finite difference discretisations on fully unstructured dual quadrilateral meshes

We present a novel approach for high-order accurate numerical differenti...
04/27/2021 ∙ by Yulong Pan, et al. ∙ 0

• ### High order methods for acoustic scattering: Coupling Farfield Expansions ABC with Deferred-Correction methods

Arbitrary high order numerical methods for time-harmonic acoustic scatte...
06/16/2020 ∙ by Vianey Villamizar, et al. ∙ 0

• ### High-order Differentiable Autoencoder for Nonlinear Model Reduction

This paper provides a new avenue for exploiting deep neural networks to ...
02/19/2021 ∙ by Siyuan Shen, et al. ∙ 44

• ### Meshless Hermite-HDMR finite difference method for high-dimensional Dirichlet problems

In this paper, a meshless Hermite-HDMR finite difference method is propo...
05/12/2019 ∙ by Xiaopeng Luo, et al. ∙ 0

• ### Stochastic Algorithmic Differentiation of (Expectations of) Discontinuous Functions (Indicator Functions)

In this paper we present a method for the accurate estimation of the der...
11/14/2018 ∙ by Christian P. Fries, et al. ∙ 0

• ### Spectral solutions of PDEs on networks

To solve linear PDEs on metric graphs with standard coupling conditions ...
04/30/2021 ∙ by M. Brio, et al. ∙ 0

• ### A new approach to proper orthogonal decomposition with difference quotients

In a recent work [B. Koc et al., arXiv:2010.03750, SIAM J. Numer. Anal.,...
06/18/2021 ∙ by Sarah K. Locke, et al. ∙ 0

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

Finite differences are commonly used for discrete approximations of derivatives. Large classess of schemes for the numerical approximation of ordinary differential equations (ODEs) and partial differential equations (PDEs) are derived from finite differences. Formulae for numerical differentiations are generally obtained from a linear combination of Taylor series, which leads to solving a system of linear equations

khan1999closed ; khan2000new ; khan2003taylor ; quarteroni2010

or calculating derivatives of interpolating polynomials (for instance see

quarteroni2010 ). References hildebrand1974introduction ; chung2010computational ; dahlquist2008numerical give a number of finite difference formulae, for high order approximation of derivatives, in term of formal power series of finite difference operators.

The purpose of this paper is to provide some basic results on finite difference approximations, which results are required for the numerical analysis of higher order time-stepping schemes for ODEs and PDEs. We introduce a new approach to derive arbitrary high order finite difference formulae which avoids the need for solving a system of linear equations. We provide various formulae for the discrete approximation of any order derivative of an analytic function at a point using arbitrary points evenly spread around . These discrete approximations are of order 1 or 2 (order 2 for centred formulae), with errors explicitly expanded in terms of Taylor series with the derivatives , . Substituting successively by their finite difference approximations in the error term for the discrete approximation of , we improve successively by 1 or 2 the order of the discrete approximation of . An efficient choice of the discrete points minimizes the number of points needed for a given order of accuracy of the discrete approximation of . Our approach can be used to recover the existing finite difference formulae, but it also provides various new formulae. We give three new finite difference formulae which are useful for the construction of new high order time-stepping schemes and their efficient starting procedures via the deferred correction (DC) method. In fact, the use of standard backward and central finite differences in building high order time-stepping schemes via the DC method leads to the computation of starting values for these schemes outside the solution interval while the standard forward finite difference formula leads to unstable schemes (see, e.g., daniel1967interated ; MR2058857 ; kress2002deferred ; koyaguerebo2019arbitrary ; koyaguerebo2020unconditionally ).

The paper is organized as follows: in section 2 we recall the main finite difference operators and prove some of their main properties; section 3 presents general first and second order approximations of derivatives with error terms explicitly expressed as Taylor series; section 4 gives many results for arbitrary high order finite difference approximations, and section 5 deals with a numerical test.

## 2 Properties of finite difference operators

In this section we recall the standard finite difference operators and provided some of their useful properties.

For a given spacing and a real , we denote and , for each integer . The centered, forward and backward difference operators , and , respectively, related to , and applied to a function from into a Banach space , are defined as follows:

 Du(tn+1/2)=u(tn+1)−u(tn)k,
 D+u(tn)=u(tn+1)−u(tn)k,

and

 D−u(tn)=u(tn)−u(tn−1)k.

The average operator is denoted by :

 Eu(tn+1/2)=ˆu(tn+1)=u(tn+1)+u(tn)2.

The composites of and are defined recursively. They commute, that is

 (D+D−)u(tn)=(D−D+)u(tn)=D−D+u(tn),

and satisfy the identities

 (D+D−)mu(tn)=k−2m2m∑j=0(−1)j(2mj)u(tn+m−j), (1)
 D−(D+D−)mu(tn)=k−2m−12m+1∑j=0(−1)j(2m+1j)u(tn+m−j), (2)

and

 Dm1+Dm2−u(tn)=k−m1−m2m1+m2∑j=0(−1)j(m1+m2j)u(tn+m1−j), (3)

for each nonnegative integer , , and such that these sums exist. Formulae (1)-(3) can be proven by a straightforward induction argument.

We introduce the double index such that

 Dαmu(tn)=Dαm1+Dαm2−u(tn). (4)
###### Remark 1.

If is even, then we have

 Dαmu(tn)=(D+D−)|αm|/2u(tm′), (5)

for some integer . For example,

 D+D3−u(tn)=(D+D−)2u(tn−1),

and

 D4−u(tn)=(D+D−)2u(tn−2).
###### Theorem 1 (Finite difference approximation of a product).

Suppose that is a Banach algebra. Then, for any functions , we have

 D−(fg)(tn)=D−f(tn)g(tn)+f(tn)D−g(tn)−kD−f(tn)D−g(tn), (6)
 D+(fg)(tn)=D+f(tn)g(tn)+f(tn)D+g(tn)+kD+f(tn)D+g(tn), (7)

and

 D+D−(fg)(tn)= D+D−f(tn)g(tn)+f(tn)D+D−g(tn)+D+f(tn)D−g(tn) (8) +D−f(tn)D+g(tn)+k2D+D−f(tn)D+D−g(tn).

More generally, for each integer such that exists, we have the formula

 (D+D−)m(fg)(tn)=m∑j=0(mj)k2j∑αm+βm=(m+j,m+j)Dαmf(tn)Dβmg(tn). (9)
###### Proof.

The formulae (6)-(8) can be obtained by a straightforward calculation, so we just need to establish (9). We proceed by induction on the positive integer . From the index notation introduced in (4), we can write

 D+D−f(tn)g(tn)+f(tn)D+D−g(tn)+D+f(tn)D−g(tn)+D−f(tn)D+g(tn) =∑α1+β1=(1,1)Dα1f(tn)Dβ1g(tn),

and

 D+D−f(tn)D+D−g(tn)=Dα1f(tn)Dβ1g(tn),   with α1+β1=(2,2).

These two identities combined with (8) yield

 D+D−(fg)(tn)=1∑j=0(1j)k2j∑α1+β1=(1+j,1+j)Dα1f(tn)Dβ1g(tn),

that is formula (9) holds for . Now suppose that (9) holds until some rank . We are going to show that it remains true for . By the induction hypothesis, we can write

 (D+D−)m+1(fg)(tn)=m∑j=0(mj)k2j∑αm+βm=(m+j,m+j)D+D−[Dαmf(tn)Dβmg(tn)]. (10)

Expanding as in the formula (8), we deduce that

 ∑αm+βm=(m+j,m+j)D+D−[Dαmf(tn)Dβmg(tn)]=S(j)+k2S(j+1), (11)

where

 S(j)=∑αm+1+βm+1=(m+1+j,m+1+j)Dαm+1f(tn)Dβm+1g(tn).

We have

 m∑j=0(mj)k2j[S(j)+ k2S(j+1)]=S(0) +m∑j=1k2j[(mj−1)+(mj)]S(j)+k2m+2S(m+1),

and deduce from (10), (11) and the identity that the formula (9) holds for . Finally, we conclude by induction that this formula is true for each suitable positive integer . ∎

###### Theorem 2 (Finite difference approximation of a composite).

Consider two functions and with values into Banach spaces such that the composite is defined on and the differential is integrable. Then

 D−f(u(tn))=∫10df(u(tn−1)+τkD−u(tn))(D−u(tn))dτ (12)

and

 D+f(u(tn))=∫10df(u(t)+ΔtD+u(t)τ)(D+u(t))dτ (13)
###### Proof.

As in standard mean value theorem. ∎

## 3 First and second order discrete approximation of derivatives

In this section we provide various formulae for the finite difference approximation of arbitrarry high order derivatives of analytic functions. The approximations are of order one or two, and the error terms are explicitly expanded i terms of Taylor series. We need the following lemma which proof is an easy induction.

###### Lemma 1.

For positive integers and and for any real , we have

 (14)

In particular, for any nonnegative integer , we have

 2m∑j=0(−1)j(2mj)(m−j)2p+1=0, (15)
 2m+1∑j=0(−1)j(2m+1j)(m−j+1/2)2p=0, (16)

and

 2m∑j=0(−1)j(2mj)[(m−j+1/2)2p+1+(m−j−1/2)2p+1]=0. (17)
###### Theorem 3.

Suppose that the function is analytic. Let , , be a partition of the interval . For each positive integer , we have

 Dm+u(tn)=u(m)(tn)+∞∑i=m+1ki−mi!u(i)(tn)m∑j=0(−1)j(mj)(m−j)i, (18)
 Dm−u(tn)=u(m)(tn)+∞∑i=m+1ki−mi!u(i)(tn)m∑j=0(m−1)j(mj)(−j)i, (19)
 D−(D+D−)mu(tn)=u(2m+1)(tn) (20) +∞∑i=2m+2ki−2m−1i!u(i)(tn)2m+1∑j=0(−1)j(2m+1j)(m−j)i,
 (D+D−)mu(tn)=u(2m)(tn)+∞∑i=m+1k2i−2m(2i)!u(2i)(tn)2m∑j=0(−1)j(2mj)(m−j)2i, (21)
 D(D+D−)mu(tn+1/2)=u(2m+1)(tn+1/2) (22) +∞∑i=m+1k2i−2m(2i+1)!u(2i+1)(tn+1/2)2m+1∑j=0(−1)j(2m+1j)(m−j−1/2)2i+1,

and

 (D+D−)mEu(tn+1/2)=u(2m)(tn+1/2)+∞∑i=m+1amik2i−2m(2i)!u(2i)(tn+1/2), (23)

where

 ami=122m∑j=0(−1)j(2mj)[(m−j+1/2)2i+(m−j−1/2)2i].
###### Proof.

We only prove formula (22). The other formulae can be proven similarly. By Taylor expansion series we have

 u(tn+m−j)=u(tn+s)+∞∑i=1kii!(m−s−j)iu(i)(tn+s).

Choosing in this formula, we deduce from (2) that

 D(D+D−)mu(tn+1/2)=k−2m−12m+1∑j=0(−1)j(2m+1j)u(tn+m−j) =k−2m−1∞∑i=1kii!u(i)(tn+1/2)2m+1∑j=0(−1)j(2m+1j)(m−j−1/2)i,

and (22) follows from (14) and (16). ∎

###### Theorem 4.

Let be , , and , , be a partition . Let and be two positive integers such that . Then, for each integer such that , is bounded independently of

, and we have the estimate

 ∥∥Dm1+Dm2−u(tn)∥∥≤Cmaxtn−m2≤t≤tn+m1∥∥u(m1+m2)(t)∥∥,

where is a constant depending only on the integer .

###### Proof.

According to Remark 1, it is enough to just prove the theorem for or , for suitable positive integer (the case is trivial). As in the previous proof, Taylor expansion of order with integral remainder together with formulae (1) and (14) yields

 (D+D−)pu(tn)=2p∑j=0(−1)j(2p−1)!(2pj)(p−j)2p∫10(1−s)2p−1u(2p)(tn+(p−j)ks)ds.

It follows that

 ∥(D+D−)pu(tn)∥ ≤1(2p)!2p∑j=0(2pj)(p−j)2pmaxtn−p≤t≤tn+p∥∥u(2p)(t)∥∥.

Similar reasoning can be applied in the case of .

## 4 Arbitrary high order finite difference approximations

###### Theorem 5.

There exists a sequence of real numbers such that for any function , where is a positive integer, and a partition , , of , we have

 u′(tn+1/2)=u(tn+1)−u(tn)k−p∑i=1c2i+1k2iD(D+D−)iu(tn+1/2)+O(k2p+2), (24)

and

 u(tn+1/2)=u(tn+1)+u(tn)2−p∑i=1c2ik2i(D+D−)iEu(tn+1/2)+O(k2p+2), (25)

for . The error constants for the formulae (24) and (25) are, respectively, and . Table 1 gives the first ten coefficients .

###### Proof.

By Taylor expansion we can write

 u(tn+1)=u(tn)+ku′(tn+1/2)+p∑i=1d1,2i+1(2i+1)!k2i+1u(2i+1)(tn+1/2)+O(k2p+3) (26)

and

 u(tn+1)=−u(tn)+2u(tn+1/2)+p∑i=1d1,2i(2i)!k2iu(2i)(tn+1/2)+O(k2p+2), (27)

with , for . Therefore, substituting successively the derivatives , , … and , , … by their expansion given by the formulae (22) and (23), respectively, into (26) and (27), we deduce the identities

 u(tn+1)=u(tn)+ku′(tn+12)+d1,33!k3DD+D−u(tn+12)+...+ dq,2q+1(2q+1)!k2q+1D(D+D−)qu(tn+12)+p∑i=q+1dq+1,2i+1(2i+1)!k2i+1u(2i+1)(tn+12)+O(k2p+3)

and

 u(tn+1)=−u(tn)+2u(tn+1/2)+d1,22!k2D+D−Eu(tn+1/2)+... +dq,2q(2q)!k2q(D+D−)qEu(tn+1/2)+p∑i=q+1dq+1,2i(2i)!k2iu(2i)(tn+1/2)+O(k2p+2)

where, for , and , we have

 dq+1,2i+1=dq,2i+1−dq,2q+1(2q+1)!2q+1∑j=0(−1)j(2q+1j)(q−j−1/2)2i+1,

and

 dq+1,2i=dq,2i−dq,2q(2q)!×22q∑j=0(−1)j(2qj)[(q−j−1/2)2i+(q−j−3/2)2i].

Finally, the identities (24 ) and (25) follow by setting and , for . ∎

###### Remark 2.

The approximations (24) and (25) are, from the coefficients computed in Table 1, equivalent to the central-difference approximation of the first derivative and the centered Bessel’s formulae (see (hildebrand1974introduction, , p.142 & p.183) or chung2010computational ; dahlquist2008numerical ).

###### Remark 3.

Formula (24) gives the finite difference approximations in khan2000new , writing

 u′(tn)=u(tn+1/2)−u(tn−1/2)k−p∑i=1c2i+1k2iD(D+D−)iu(tn)+O(k2p+2), (28)

where

 p∑i=1c2i+1k2iD(D+D−)iu(tn)=k−1p∑i=1[c2i+12i+1∑j=0(−1)j(2i+1j)u(tn+i−j+1/2)].

-   For we have

 u′(tn) =u(tn+1/2)−u(tn−1/2)k−124k2D(D+D−)u(tn)+O(k4) =u(tn+1/2)−u(tn−1/2)k−u(tn+3/2)−3u(tn+1/2)+3u(tn−1/2)−u(tn−3/2)24k +O(k4).

-   For we have

 u′(tn) =u(tn+1/2)−u(tn−1/2)k−124k2D(D+D−)u(tn)+18255!k4D(D+D−)2u(tn) +O(k6),

and then

 u′(tn)=u(tn+1/2)−u(tn−1/2)k+11920k[9−125330−330125−9]UTn,5+O(k6),

where

is the transpose of the vector

 Un,5=[u(tn+5/2)u(tn+3/2)u(tn+1/2)u(tn−1/2)u(tn−3/2)u(tn−5/2)].

The following theorem gives a new form of centered finite difference formulae which is useful for efficient starting procedures of high order time-stepping schemes via deferred correction strategy koyaguerebo2019arbitrary ; koyaguerebo2020unconditionally .

###### Theorem 6 (Interior centered approximations).

Let , where is a positive integer and , , is a real interval. Given a uniform partition of , that is with , and , there exist reals such that

 u′(τp+1/2)=u(b)−u(a)b−a−1b−ap∑i=1cp2i+1k2i+1D(D+D−)iu(τp+1/2)+O(k2p+2). (29)

and

 u(τp+1/2)=u(b)+u(a)2−p∑i=1cp2ik2i(D+D−)iEu(τp+1/2)+O(k2p+2), (30)

Table 2 gives the coefficients for .