# Two New Approximations for Variable-Order Fractional Derivatives

Authors: Ruilian Du; Zongqi Liang
Journal: Discrete Dynamics in Nature and Society
 (2017)

---

## Abstract

We introduced a parameterσ(t) which was related to α(t); then two numerical schemes for variable-order Caputo fractional derivatives were derived; the second-order numerical approximation to variable-order fractional derivatives α(t)∈(0,1) and 3-α(t)-order approximation for α(t)∈(1,2) are established. For the given parameter σ(t), the error estimations of formulas were proven, which were higher than some recently derived schemes. Finally, some numerical examples with exact solutions were studied to demonstrate the theoretical analysis and verify the efficiency of the proposed methods.


---

## Body

## 1. Introduction

Fractional differential equations include constant-order and variable-order equations; a great quantity of natural phenomena can be modeled by variable-order fractional differential equations; the study of such problems has attracted much attention. In recent years, profound background of physical applications for the variable-order fractional calculus has been already established; the definition of variable-order operator has been investigated in [1–4]. Numerous problems in mathematical physics and engineering have been modeled by variable-order fractional differential equations, such as successful applications in mechanics [5], in the simulation of linear and nonlinear viscoelasticity oscillators [6], and in other cases where the order of the derivative varies with time [7]. In addition, a physical experimental investigation of variable-order operators has been considered in [8]. Variable-order fractional derivatives can be used to model anomalous diffusion, as they can describe the time dependent diffusion process more specifically than fractional derivatives of constant order, just as shown in [9]. The difference schemes of fractional derivatives with constant and variable order are investigated in [10].Due to the existence of variable fractional derivative, it is usually difficult to obtain the analytical solution of such equations; therefore, it is particularly significant to give numerical solutions to these problems. Cao and Qiu [11] derived a high order numerical method for variable-order fractional ordinary differential equation by establishing a second-order numerical approximation to variable-order Riemann-Liouville fractional derivative. Fu et al. [12] adopted the method of approximate particular solutions for both constant-order and variable-order time fractional diffusion models. Several finite difference methods for variable-order fractional partial diffusion equations were proposed in [13–18].As what we can see above, there are many different definitions of the variable-order fractional derivatives. However, different from other definitions in mathematics, only initial condition is needed for the variable-order Caputo definition which can be easily used in physical field; this definition means that the memory rate of system is determined by the current time instant and changes with time. For simplicity, only variable-order Caputo definition was discussed in the whole paper.A good approximation of the variable-order Caputo derivatives was observed in [19], it was defined at the points tn+1/2 and tn; then what will happen if it is redefined in a neighborhood at the point tn? Approximation formula for the α-order (0≤α≤1) Caputo derivatives at the point t=tn-1+σ was investigated in [20], where σ=1-α/2 is a superconvergence point which ensures that the convergence order improves from the local second order to the overall 3-α order. In the current work, we applied the ideas in [20]; a parameter σ(t) was given, which was changed with α(t), and the selection of σ(t) played a vital role in the numerical results in the paper. We proposed two new approximation formulas of second-order and 3-α(t) accuracy for variable-order time fractional operator with orders 0<α(t)<1 and 1<α(t)<2, respectively. Specifically, we adopted the following definition of variable-order Caputo fractional derivatives: (1)Dtαt0Cft=1Γn-αt∫0tfnst-sαt-n+1ds,n-1<αt<n.The paper is organized as follows. In the next section, we present two new formulas, selectσ(t)=1-α(t)/2 for 0<α(t)<1 and σ(t)=3/2-α(t)/2 for 1<α(t)<2, respectively, and provide the corresponding error analysis. In Section 3, the numerical verifications are presented; as can be seen from the numerical results, minor changes to the selected parameters σ(t) will have a significant impact on the error estimates. Finally, in Section 4, the conclusion is drawn.


## 2. Approximation Formulas

In this section, we present the following lemma, which plays a vital role in the later analysis.Lemma 1 (see [19]).
Assuming that the derivatives of the functionf(x) exist to the order of n on [a,b] and to the order of n+1 on (a,b),a≤x0<x1<⋯<xn≤b,pn(x) is the nth degree interpolation polynomial of f(x) based on the points x0,x1,…,xn, then, for x∈[a,b], the following holds: (2)fkx-pnkx=fn+1ξn-k+1!x-x0kx-x1k⋯x-xn-kk,k=0,1,2,…,where ξ∈(a,b) depends on k and x, and xi<xi(k)<xi+k(i=0,1,…,n-k).
For a given functionf(x) and an integer N>0, denote tn=nτ,fn=f(tn),n=0,1,…,N, where τ=T/N is the step. Introduce the following notation: (3)tn+σ=n+στ,δtfn+1/2=1τfn+1-fn,δt2fn=δtfn+1/2-δtfn-1/2τ,and, for convenience, denote σ=σ(tn).


### 2.1. Second-Order Formula forα(t)∈(0,1)

Denoteαn-1+σ=α(tn-1+σ),σ=1-α(tn-1+σ)/2,n=0,1,…,N-1.Estimating the derivativeDtα(t)0Cf(t) with order α(t)∈(0,1) at the grid point tn-1+σ,n=0,1,…,N-1, from (1), we directly obtain (4)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∫0tn-1+σf′stn-1+σ-sαn-1+σds.Evaluate the integration on each subinterval, leading to(5)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkf′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σf′stn-1+σ-sαn-1+σds.For each interval[tk-1,tk],k=1,…,n-1, denote the second-degree interpolation polynomial in the Lagrange form as follows: (6)Lk2fs=fk-1t-tkt-tk+12τ2-fkt-tk-1t-tk+1τ2+fk+1t-tk-1t-tk2τ2, and for the last interval [tn-1,tn-1+σ], denote the first-degree interpolation polynomial in the Lagrange form as (7)Ln1fs=-fn-1t-tnτ+fnt-tn-1τ.As an approximation formula, the following result is obtained:(8)fs=Lk2fs+rks,s∈tk-1,tk,k=1,…,n-1,where (9)rks=f′′′ξk6s-tk-1s-tks-tk+1,ξk∈tk-1,tk+1,(10)fs=Ln1fs+rns,s∈tn-1,tn-1+σ,where (11)rns=f′′ξn2s-tn-1s-tn,ξn∈tn-1,tn.Substituting (8) and (10) into (5) yields (12)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkLk2fs′tn-1+σ-sαn-1+σds+∫tn-1tn-1+σLn1fs′tn-1+σ-sαn-1+σds+Rn-1+σ,and the truncation error is (13)Rn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds.Here we denote the discrete approximation formula for the variable-order derivative with orderαn-1+σ as Δtαn-1+σ0fn-1+σ; thus from (12), it is easy to get the following result: (14)Δtαn-1+σ0fn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkLk2fs′tn-1+σ-sαn-1+σds+∫tn-1tn-1+σLn1fs′tn-1+σ-sαn-1+σds=1Γ1-αn-1+σ∑k=1n-1∫tk-1tk1/2δtfk+1/2+δtfk-1/2+δt2fks-tktn-1+σ-sαn-1+σds+∫tn-1tn-1+σδtfn-1/2tn-1+σ-sαn-1+σds=1Γ1-αn-1+σ∑k=1n-1δtfk+1/2+δtfk-1/2τ1-αn-1+σ21-αn-1+σn+σ-k1-αn-1+σ-n+σ-k-11-αn-1+σ+δtfn-1/2τ1-αn-1+σ1-αn-1+σσ1-αn-1+σ+∑k=1n-1∫tk-1tkδt2fks-tktn-1+σ-sαn-1+σds=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1δtfk+1/2+δtfk-1/22n+σ-k1-αn-1+σ-n+σ-k-11-αn-1+σ+∑k=1n-1δt2fkτ12-αn-1+σn+σ-k2-αn-1+σ-n+σ-k-12-αn-1+σ-n+σ-k1-αn-1+σ+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1an-kαn-1+σδtfk+1/2+δtfk-1/2+∑k=1n-1bn-kαn-1+στδt2fk+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1an-kαn-1+σδtfk+1/2+δtfk-1/2+∑k=1n-1bn-kαn-1+σδtfk+1/2-δtfk-1/2+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=2n-1an-k+1αn-1+σ+an-kαn-1+σ+bn-k+1αn-1+σ-bn-kαn-1+σδtfk-1/2+a1αn-1+σ+a0αn-1+σ+b1αn-1+σδtfn-1/2+an-1αn-1+σ-bn-1αn-1+σδtf1/2=τ1-αn-1+σΓ2-αn-1+σ∑k=1ndn-kαn-1+σδtfk-1/2=τ-αn-1+σΓ2-αn-1+σ∑k=1ndn-kαn-1+σfk-fk-1,where (15)a0αn-1+σ=σ1-αn-1+σ,akαn-1+σ=12k+σ1-αn-1+σ-k+σ-11-αn-1+σ,1≤k≤n-1,bkαn-1+σ=12-αn-1+σk+σ2-αn-1+σ-k+σ-12-αn-1+σ-k+σ1-αn-1+σ,1≤k≤n-1.dkαn-1+σ=a1αn-1+σ+a0αn-1+σ+b1αn-1+σ,k=0,ak+1αn-1+σ+akαn-1+σ+bk+1αn-1+σ-bkαn-1+σ,1≤k≤n-2,akαn-1+σ-bkαn-1+σ,k=n-1.Next, the analysis for the approximation error of formula (12) was given in detail.Theorem 2.
Letαn-1+σ∈(0,1),f∈C3([0,tn+1]); the following holds: (16)Dtαn-1+σ0Cftn-1+σ-Δtαn-1+σ0fn-1+σ≤1Γ1-αn-1+σΘnmaxt0≤t≤tn-1⁡f′′t+σ1-αn-1+στ-αn-1+σ1-αn-1+σmaxtn-1≤t≤tn⁡f′′tτ2, where Θ1=τ-ασσ1-ασ/1-ασ,Θn=c1tn-1+σ1-αn-1+σ/1-αn-1+σ,n=2,…,N-1,c1 is a positive constant.Proof.
From (12), obviously (17)Rn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds=1Γ1-αn-1+σR1+R2.
Using Lemma1, it produces (18)rk′s=12f′′′ηks-t0ks-t1k,ηk∈tk-1,tk+1,t0k∈tk-1,tk,t1k∈tk,tk+1,1≤k≤n-1,(19)rn′s=f′′ηns-t0n,ηn∈tn-1,tn,t0n∈tn-1,tn.
Next, the error on each interval is analyzed.
Whenn=1, from (19), it follows that(20)∫0tσr1′stσ-sασds≤τmaxt0≤t≤t1⁡f′′t∫0tσ1tσ-sασds=τ-ασσ1-ασ1-ασmaxt0≤t≤t1⁡f′′tτ2.
Whenn=2,…,N-1, from (18), this leads to the following estimation: (21)R1=∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds=∑k=1n-112∫tk-1tkf′′′ηks-t0ks-t1ktn-1+σ-sαn-1+σds≤c1maxt0≤t≤tn-1⁡f′′′t∫t0tn-11tn-1+σ-sαn-1+σdsτ2=c1tn-1+σ1-αn-1+σ1-αn-1+σmaxt0≤t≤tn-1⁡f′′′tτ2,where c1 is a positive constant.
From (19), it can be obtained in the same method. (22)R2=∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds≤τmaxtn-1≤t≤tn⁡f′′t∫tn-1tn-1+σ1tn-1+σ-sαn-1+σds=σ1-αn-1+στ-αn-1+σ1-αn-1+σmaxtn-1≤t≤tn⁡f′′tτ2.
Substituting (20)–(22) into (17) leads to the theorem.


### 2.2. Second-Order Formula forα(t)∈(1,2)

Denotef0′=f′(t0), σ=3/2-αn-1+σ/2. Computing the variable-order derivative with order α(t)(1<α(t)<2) at tn-1+σ, (23)Dtαn-1+σ0Cftn-1+σ=1Γ2-αn-1+σ∫t0tn-1+σf′′stn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σ∫t0t1f′′stn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1f′′stn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σf′′stn-1+σ-sαn-1+σ-1ds.For the interval[t0,t1], the cubic interpolation polynomial in the Hermite form is as (24)L03fs=f0+f0′s-t0+1/τf1-f0-f0′τs-t02+12τ1/τf1-f0-f0′τ-1/τf2-f1-1/τf1-f02τs-t02s-t1,and it follows that(25)fs=L03fs+r0s;the truncation error is (26)r0s=f4ξ024s-t02s-t1s-t2,ξ0∈t0,t2.For the interval[tj,tj+1],j=1,…,n-2, the cubic interpolation polynomial at points tj-1,tj,tj+1, and tj+2 in the Lagrange form is as (27)Lj3fs=∑k=j-1j+2fk∏m=j-1,m≠kj+2s-tmtk-tm, leading to (28)fs=Lj3fs+rjs;the truncation error is (29)rjs=f4ξj24∏k=-12s-tk+j,ξj∈tj-1,tj+2;using the cubic interpolation polynomial at points tn-2,tn-1, and tn in the Lagrange form to approximate f(s) on the last interval [tn-1,tn-1+σ], we get (30)Ln2fs=fn-2t-tn-1t-tn2τ2-fn-1t-tn-2t-tnτ2+fnt-tn-2t-tn-12τ2,the following holds: (31)fs=Ln2fs+rns;the truncation error is (32)rns=f′′′ξn6s-tn-2s-tn-1s-tn,ξn∈tn-2,tn.Substituting (25), (28), and (31) into (23), there is (33)Dtαn-1+σ0Cftn-1+σ=1Γ2-αn-1+σ∫t0t1L03fs+r0s′′tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1Lj3fs+rjs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs+rns′′tn-1+σ-sαn-1+σ-1ds.Denote the approximation formula of the variable-order derivative with orderαn-1+σ as Δtαn-1+σ0fn-1+σ and obtain (34)Δtαn-1+σ0fn-1+σ=1Γ2-αn-1+σ∑j=0n-2∫tjtj+1Lj3fs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs′′tn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σ∫t0t1L03fs′′tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1Lj3fs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs′′tn-1+σ-sαn-1+σ-1ds,where (35)L03fs′′=s-t0τf2-2f1+f0τ2+t1-sτ8f1-f2-7f0-6τf0′2τ2,Lj3fs′′=tj+1-sτfj-1-2fj+fj+1τ2+s-tjτfj-2fj+1+fj+2τ2,1≤j≤N-2,Ln2fs′′=δt2fn-1.Substituting the above formula into (34), thus (36)Δtαn-1+σ0fn-1+σ=1Γ2-αn-1+σ∫t0t1s-t0/τf2-2f1+f0/τ2+t1-s/τ8f1-f2-7f0-6τf0′/2τ2tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1tj+1-s/τfj-1-2fj+fj+1/τ2+s-tj/τfj-2fj+1+fj+2/τ2tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σδt2fn-1tn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σI1+I2+I3.For the first term in the right hand side of (36), this leads to (37)I1=∫t0t1s-t0/τf2-2f1+f0/τ2+t1-s/τ8f1-f2-7f0-6τf0′/2τ2tn-1+σ-sαn-1+σ-1ds=-δt2f1τ∫t0t1tn-1+σ-s2-αn-1+σds+n-1+σδt2f1∫t0t1tn-1+σ-s1-αn-1+σds+8f1-f2-7f0-6τf0′2τ3∫t0t1tn-1+σ-s2-αn-1+σds-n+σ-28f1-f2-7f0-6τf0′2τ2∫t0t1tn-1+σ-s1-αn-1+σds=3f2-12f1+9f0+6τf0′τ-αn-1+σ22-αn-1+σ3-αn-1+σn+σ-13-αn-1+σ-n+σ-23-αn-1+σ+8f1-f2-7f0-6τf0′τ-αn-1+σ22-αn-1+σn+σ-12-αn-1+σ-2f2-4f1+2f0τ-αn-1+σ22-αn-1+σn+σ-22-αn-1+σ.For the second term in the right hand side of (36), the following holds: (38)I2=∫tjtj+1tj+1-s/τfj-1-2fj+fj+1/τ2+s-tj/τfj-2fj+1+fj+2/τ2tn-1+σ-sαn-1+σ-1ds=δt2fjτ∫tjtj+1tj+1-stn-1+σ-sαn-1+σ-1ds+δt2fj+1τ∫tjtj+1s-tjtn-1+σ-sαn-1+σ-1ds=fj+2-3fj+1+3fj-fj-1τ-αn-1+σ2-αn-1+σ3-αn-1+σn+σ-j-13-αn-1+σ-n+σ-j-23-αn-1+σ+fj+1-2fj+fj-1τ-αn-1+σ2-αn-1+σn+σ-j-12-αn-1+σ-fj+2-2fj+1+fjτ-αn-1+σ2-αn-1+σn+σ-j-22-αn-1+σ.For the third term in the right hand side of (36), this yields (39)I3=∫tn-1tn-1+σδt2fn-1tn-1+σ-sαn-1+σ-1ds=fn-2fn-1+fn-2σ2-αn-1+σ2-αn-1+στ-αn-1+σ.Substituting (37)–(39) into (36), we obtain the approximate difference scheme of the Dtαn-1+σ0Cf(tn-1+σ): (40)Δtαn-1+σ0fn-1+σ=τ-αn-1+σΓ3-αn-1+σ∑j=0n-2an-jαn-1+σfj+2-2fj+1+fj+∑j=1n-2bn-jαn-1+σfj+1-2fj+fj-1+12bnαn-1+σ8f1-f2-7f0-6τf0′+σ2-αn-1+σfn-2fn-1+fn-2,where (41)ajαn-1+σ=13-αn-1+σj+σ-13-αn-1+σ-j+σ-23-αn-1+σ-j+σ-22-αn-1+σ,1≤j≤n,(42)bjαn-1+σ=j+σ-12-αn-1+σ-13-αn-1+σj+σ-13-αn-1+σ-j+σ-23-αn-1+σ,2≤j≤n-1.Similar to the proof of the Theorem2, the truncation error of formula (40) is given by following theorem.Theorem 3.
Forαn-1+σ∈(1,2),f∈C4([0,tn+1]), the following holds: (43)Dtαn-1+σ0Cftn-1+σ-Δtαn-1+σ0fn-1+σ≤1Γ3-αn-1+σtn-1+σ2-αn-1+σmaxt0≤t≤tn-1⁡f4tτ2+c2σ2-αn-1+σmaxtn-1≤t≤tn⁡f′′′tτ3-αn-1+σ,where c2 is a positive constant.Proof.
From (33), we obtain the truncation error (44)Rn-1+σ=1Γ2-αn-1+σ∑j=0n-2∫tjtj+1rj3s′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σrn2s′′tn-1+σ-sαn-1+σ-1ds.
From Lemma1, the result below is natural. (45)r03s′′=f4η02s-t00s-t10,η0∈t0,t2,t00∈t0,t1,t10∈t0,t2,rj3s′′=f4ηj2s-t0js-t1j,ηj∈tj-1,tj+2,t0j∈tj-1,tj+1,t1j∈tj,tj+2,rn2s′′=f′′′ηns-t0n,ηn∈tn-2,tn,t0n∈tn-2,tn.
Substituting (45) into (44), it yields the estimation (46)∫t0t1f4η0/2s-t00s-t10tn-1+σ-sαn-1+σ-1ds≤maxt0≤t≤t2⁡f4t∫t0t1tn-1+σ-s1-αn-1+σdsτ2≤tn-1+σ2-αn-1+σ2-αn-1+σmaxt0≤t≤t2⁡f4tτ2,∑j=1n-2∫tjtj+1f4ηj/2s-t0js-t1jtn-1+σ-sαn-1+σ-1ds≤2maxt0≤t≤tn-1⁡f4t∫t1tn-1tn-1+σ-s1-αn-1+σdsτ2≤tn-1+σ2-αn-1+σ2-αn-1+σmaxt1≤t≤tn-1⁡f4tτ2,∫tn-1tn-1+σf′′′ηns-t0ntn-1+σ-sαn-1+σ-1ds≤c2maxtn-1≤t≤tn⁡f′′′t∫tn-1tn-1+σtn-1+σ-s1-αn-1+σdsτ≤c2σ2-αn-1+σ2-αn-1+σmaxtn-1≤t≤tn⁡f′′′tτ3-αn-1+σ,where c2 is a positive constant.
Substituting (46) into (44) leads to the theorem.Remark 4.
From the second term of the right hand of the truncation error estimate in Theorem2, we can easily obtain the fact that formula (14) does not have exact second-order accuracy; it is related to the selection of α(t) and σ. This can be seen from Table 2 in the numerical example below; its accuracy is much higher than the second order.


## 2.1. Second-Order Formula forα(t)∈(0,1)

Denoteαn-1+σ=α(tn-1+σ),σ=1-α(tn-1+σ)/2,n=0,1,…,N-1.Estimating the derivativeDtα(t)0Cf(t) with order α(t)∈(0,1) at the grid point tn-1+σ,n=0,1,…,N-1, from (1), we directly obtain (4)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∫0tn-1+σf′stn-1+σ-sαn-1+σds.Evaluate the integration on each subinterval, leading to(5)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkf′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σf′stn-1+σ-sαn-1+σds.For each interval[tk-1,tk],k=1,…,n-1, denote the second-degree interpolation polynomial in the Lagrange form as follows: (6)Lk2fs=fk-1t-tkt-tk+12τ2-fkt-tk-1t-tk+1τ2+fk+1t-tk-1t-tk2τ2, and for the last interval [tn-1,tn-1+σ], denote the first-degree interpolation polynomial in the Lagrange form as (7)Ln1fs=-fn-1t-tnτ+fnt-tn-1τ.As an approximation formula, the following result is obtained:(8)fs=Lk2fs+rks,s∈tk-1,tk,k=1,…,n-1,where (9)rks=f′′′ξk6s-tk-1s-tks-tk+1,ξk∈tk-1,tk+1,(10)fs=Ln1fs+rns,s∈tn-1,tn-1+σ,where (11)rns=f′′ξn2s-tn-1s-tn,ξn∈tn-1,tn.Substituting (8) and (10) into (5) yields (12)Dtαn-1+σ0Cftn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkLk2fs′tn-1+σ-sαn-1+σds+∫tn-1tn-1+σLn1fs′tn-1+σ-sαn-1+σds+Rn-1+σ,and the truncation error is (13)Rn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds.Here we denote the discrete approximation formula for the variable-order derivative with orderαn-1+σ as Δtαn-1+σ0fn-1+σ; thus from (12), it is easy to get the following result: (14)Δtαn-1+σ0fn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkLk2fs′tn-1+σ-sαn-1+σds+∫tn-1tn-1+σLn1fs′tn-1+σ-sαn-1+σds=1Γ1-αn-1+σ∑k=1n-1∫tk-1tk1/2δtfk+1/2+δtfk-1/2+δt2fks-tktn-1+σ-sαn-1+σds+∫tn-1tn-1+σδtfn-1/2tn-1+σ-sαn-1+σds=1Γ1-αn-1+σ∑k=1n-1δtfk+1/2+δtfk-1/2τ1-αn-1+σ21-αn-1+σn+σ-k1-αn-1+σ-n+σ-k-11-αn-1+σ+δtfn-1/2τ1-αn-1+σ1-αn-1+σσ1-αn-1+σ+∑k=1n-1∫tk-1tkδt2fks-tktn-1+σ-sαn-1+σds=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1δtfk+1/2+δtfk-1/22n+σ-k1-αn-1+σ-n+σ-k-11-αn-1+σ+∑k=1n-1δt2fkτ12-αn-1+σn+σ-k2-αn-1+σ-n+σ-k-12-αn-1+σ-n+σ-k1-αn-1+σ+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1an-kαn-1+σδtfk+1/2+δtfk-1/2+∑k=1n-1bn-kαn-1+στδt2fk+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=1n-1an-kαn-1+σδtfk+1/2+δtfk-1/2+∑k=1n-1bn-kαn-1+σδtfk+1/2-δtfk-1/2+δtfn-1/2σ1-αn-1+σ=τ1-αn-1+σΓ2-αn-1+σ∑k=2n-1an-k+1αn-1+σ+an-kαn-1+σ+bn-k+1αn-1+σ-bn-kαn-1+σδtfk-1/2+a1αn-1+σ+a0αn-1+σ+b1αn-1+σδtfn-1/2+an-1αn-1+σ-bn-1αn-1+σδtf1/2=τ1-αn-1+σΓ2-αn-1+σ∑k=1ndn-kαn-1+σδtfk-1/2=τ-αn-1+σΓ2-αn-1+σ∑k=1ndn-kαn-1+σfk-fk-1,where (15)a0αn-1+σ=σ1-αn-1+σ,akαn-1+σ=12k+σ1-αn-1+σ-k+σ-11-αn-1+σ,1≤k≤n-1,bkαn-1+σ=12-αn-1+σk+σ2-αn-1+σ-k+σ-12-αn-1+σ-k+σ1-αn-1+σ,1≤k≤n-1.dkαn-1+σ=a1αn-1+σ+a0αn-1+σ+b1αn-1+σ,k=0,ak+1αn-1+σ+akαn-1+σ+bk+1αn-1+σ-bkαn-1+σ,1≤k≤n-2,akαn-1+σ-bkαn-1+σ,k=n-1.Next, the analysis for the approximation error of formula (12) was given in detail.Theorem 2.
Letαn-1+σ∈(0,1),f∈C3([0,tn+1]); the following holds: (16)Dtαn-1+σ0Cftn-1+σ-Δtαn-1+σ0fn-1+σ≤1Γ1-αn-1+σΘnmaxt0≤t≤tn-1⁡f′′t+σ1-αn-1+στ-αn-1+σ1-αn-1+σmaxtn-1≤t≤tn⁡f′′tτ2, where Θ1=τ-ασσ1-ασ/1-ασ,Θn=c1tn-1+σ1-αn-1+σ/1-αn-1+σ,n=2,…,N-1,c1 is a positive constant.Proof.
From (12), obviously (17)Rn-1+σ=1Γ1-αn-1+σ∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds+∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds=1Γ1-αn-1+σR1+R2.
Using Lemma1, it produces (18)rk′s=12f′′′ηks-t0ks-t1k,ηk∈tk-1,tk+1,t0k∈tk-1,tk,t1k∈tk,tk+1,1≤k≤n-1,(19)rn′s=f′′ηns-t0n,ηn∈tn-1,tn,t0n∈tn-1,tn.
Next, the error on each interval is analyzed.
Whenn=1, from (19), it follows that(20)∫0tσr1′stσ-sασds≤τmaxt0≤t≤t1⁡f′′t∫0tσ1tσ-sασds=τ-ασσ1-ασ1-ασmaxt0≤t≤t1⁡f′′tτ2.
Whenn=2,…,N-1, from (18), this leads to the following estimation: (21)R1=∑k=1n-1∫tk-1tkrk′stn-1+σ-sαn-1+σds=∑k=1n-112∫tk-1tkf′′′ηks-t0ks-t1ktn-1+σ-sαn-1+σds≤c1maxt0≤t≤tn-1⁡f′′′t∫t0tn-11tn-1+σ-sαn-1+σdsτ2=c1tn-1+σ1-αn-1+σ1-αn-1+σmaxt0≤t≤tn-1⁡f′′′tτ2,where c1 is a positive constant.
From (19), it can be obtained in the same method. (22)R2=∫tn-1tn-1+σrn′stn-1+σ-sαn-1+σds≤τmaxtn-1≤t≤tn⁡f′′t∫tn-1tn-1+σ1tn-1+σ-sαn-1+σds=σ1-αn-1+στ-αn-1+σ1-αn-1+σmaxtn-1≤t≤tn⁡f′′tτ2.
Substituting (20)–(22) into (17) leads to the theorem.


## 2.2. Second-Order Formula forα(t)∈(1,2)

Denotef0′=f′(t0), σ=3/2-αn-1+σ/2. Computing the variable-order derivative with order α(t)(1<α(t)<2) at tn-1+σ, (23)Dtαn-1+σ0Cftn-1+σ=1Γ2-αn-1+σ∫t0tn-1+σf′′stn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σ∫t0t1f′′stn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1f′′stn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σf′′stn-1+σ-sαn-1+σ-1ds.For the interval[t0,t1], the cubic interpolation polynomial in the Hermite form is as (24)L03fs=f0+f0′s-t0+1/τf1-f0-f0′τs-t02+12τ1/τf1-f0-f0′τ-1/τf2-f1-1/τf1-f02τs-t02s-t1,and it follows that(25)fs=L03fs+r0s;the truncation error is (26)r0s=f4ξ024s-t02s-t1s-t2,ξ0∈t0,t2.For the interval[tj,tj+1],j=1,…,n-2, the cubic interpolation polynomial at points tj-1,tj,tj+1, and tj+2 in the Lagrange form is as (27)Lj3fs=∑k=j-1j+2fk∏m=j-1,m≠kj+2s-tmtk-tm, leading to (28)fs=Lj3fs+rjs;the truncation error is (29)rjs=f4ξj24∏k=-12s-tk+j,ξj∈tj-1,tj+2;using the cubic interpolation polynomial at points tn-2,tn-1, and tn in the Lagrange form to approximate f(s) on the last interval [tn-1,tn-1+σ], we get (30)Ln2fs=fn-2t-tn-1t-tn2τ2-fn-1t-tn-2t-tnτ2+fnt-tn-2t-tn-12τ2,the following holds: (31)fs=Ln2fs+rns;the truncation error is (32)rns=f′′′ξn6s-tn-2s-tn-1s-tn,ξn∈tn-2,tn.Substituting (25), (28), and (31) into (23), there is (33)Dtαn-1+σ0Cftn-1+σ=1Γ2-αn-1+σ∫t0t1L03fs+r0s′′tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1Lj3fs+rjs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs+rns′′tn-1+σ-sαn-1+σ-1ds.Denote the approximation formula of the variable-order derivative with orderαn-1+σ as Δtαn-1+σ0fn-1+σ and obtain (34)Δtαn-1+σ0fn-1+σ=1Γ2-αn-1+σ∑j=0n-2∫tjtj+1Lj3fs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs′′tn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σ∫t0t1L03fs′′tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1Lj3fs′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σLn2fs′′tn-1+σ-sαn-1+σ-1ds,where (35)L03fs′′=s-t0τf2-2f1+f0τ2+t1-sτ8f1-f2-7f0-6τf0′2τ2,Lj3fs′′=tj+1-sτfj-1-2fj+fj+1τ2+s-tjτfj-2fj+1+fj+2τ2,1≤j≤N-2,Ln2fs′′=δt2fn-1.Substituting the above formula into (34), thus (36)Δtαn-1+σ0fn-1+σ=1Γ2-αn-1+σ∫t0t1s-t0/τf2-2f1+f0/τ2+t1-s/τ8f1-f2-7f0-6τf0′/2τ2tn-1+σ-sαn-1+σ-1ds+∑j=1n-2∫tjtj+1tj+1-s/τfj-1-2fj+fj+1/τ2+s-tj/τfj-2fj+1+fj+2/τ2tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σδt2fn-1tn-1+σ-sαn-1+σ-1ds=1Γ2-αn-1+σI1+I2+I3.For the first term in the right hand side of (36), this leads to (37)I1=∫t0t1s-t0/τf2-2f1+f0/τ2+t1-s/τ8f1-f2-7f0-6τf0′/2τ2tn-1+σ-sαn-1+σ-1ds=-δt2f1τ∫t0t1tn-1+σ-s2-αn-1+σds+n-1+σδt2f1∫t0t1tn-1+σ-s1-αn-1+σds+8f1-f2-7f0-6τf0′2τ3∫t0t1tn-1+σ-s2-αn-1+σds-n+σ-28f1-f2-7f0-6τf0′2τ2∫t0t1tn-1+σ-s1-αn-1+σds=3f2-12f1+9f0+6τf0′τ-αn-1+σ22-αn-1+σ3-αn-1+σn+σ-13-αn-1+σ-n+σ-23-αn-1+σ+8f1-f2-7f0-6τf0′τ-αn-1+σ22-αn-1+σn+σ-12-αn-1+σ-2f2-4f1+2f0τ-αn-1+σ22-αn-1+σn+σ-22-αn-1+σ.For the second term in the right hand side of (36), the following holds: (38)I2=∫tjtj+1tj+1-s/τfj-1-2fj+fj+1/τ2+s-tj/τfj-2fj+1+fj+2/τ2tn-1+σ-sαn-1+σ-1ds=δt2fjτ∫tjtj+1tj+1-stn-1+σ-sαn-1+σ-1ds+δt2fj+1τ∫tjtj+1s-tjtn-1+σ-sαn-1+σ-1ds=fj+2-3fj+1+3fj-fj-1τ-αn-1+σ2-αn-1+σ3-αn-1+σn+σ-j-13-αn-1+σ-n+σ-j-23-αn-1+σ+fj+1-2fj+fj-1τ-αn-1+σ2-αn-1+σn+σ-j-12-αn-1+σ-fj+2-2fj+1+fjτ-αn-1+σ2-αn-1+σn+σ-j-22-αn-1+σ.For the third term in the right hand side of (36), this yields (39)I3=∫tn-1tn-1+σδt2fn-1tn-1+σ-sαn-1+σ-1ds=fn-2fn-1+fn-2σ2-αn-1+σ2-αn-1+στ-αn-1+σ.Substituting (37)–(39) into (36), we obtain the approximate difference scheme of the Dtαn-1+σ0Cf(tn-1+σ): (40)Δtαn-1+σ0fn-1+σ=τ-αn-1+σΓ3-αn-1+σ∑j=0n-2an-jαn-1+σfj+2-2fj+1+fj+∑j=1n-2bn-jαn-1+σfj+1-2fj+fj-1+12bnαn-1+σ8f1-f2-7f0-6τf0′+σ2-αn-1+σfn-2fn-1+fn-2,where (41)ajαn-1+σ=13-αn-1+σj+σ-13-αn-1+σ-j+σ-23-αn-1+σ-j+σ-22-αn-1+σ,1≤j≤n,(42)bjαn-1+σ=j+σ-12-αn-1+σ-13-αn-1+σj+σ-13-αn-1+σ-j+σ-23-αn-1+σ,2≤j≤n-1.Similar to the proof of the Theorem2, the truncation error of formula (40) is given by following theorem.Theorem 3.
Forαn-1+σ∈(1,2),f∈C4([0,tn+1]), the following holds: (43)Dtαn-1+σ0Cftn-1+σ-Δtαn-1+σ0fn-1+σ≤1Γ3-αn-1+σtn-1+σ2-αn-1+σmaxt0≤t≤tn-1⁡f4tτ2+c2σ2-αn-1+σmaxtn-1≤t≤tn⁡f′′′tτ3-αn-1+σ,where c2 is a positive constant.Proof.
From (33), we obtain the truncation error (44)Rn-1+σ=1Γ2-αn-1+σ∑j=0n-2∫tjtj+1rj3s′′tn-1+σ-sαn-1+σ-1ds+∫tn-1tn-1+σrn2s′′tn-1+σ-sαn-1+σ-1ds.
From Lemma1, the result below is natural. (45)r03s′′=f4η02s-t00s-t10,η0∈t0,t2,t00∈t0,t1,t10∈t0,t2,rj3s′′=f4ηj2s-t0js-t1j,ηj∈tj-1,tj+2,t0j∈tj-1,tj+1,t1j∈tj,tj+2,rn2s′′=f′′′ηns-t0n,ηn∈tn-2,tn,t0n∈tn-2,tn.
Substituting (45) into (44), it yields the estimation (46)∫t0t1f4η0/2s-t00s-t10tn-1+σ-sαn-1+σ-1ds≤maxt0≤t≤t2⁡f4t∫t0t1tn-1+σ-s1-αn-1+σdsτ2≤tn-1+σ2-αn-1+σ2-αn-1+σmaxt0≤t≤t2⁡f4tτ2,∑j=1n-2∫tjtj+1f4ηj/2s-t0js-t1jtn-1+σ-sαn-1+σ-1ds≤2maxt0≤t≤tn-1⁡f4t∫t1tn-1tn-1+σ-s1-αn-1+σdsτ2≤tn-1+σ2-αn-1+σ2-αn-1+σmaxt1≤t≤tn-1⁡f4tτ2,∫tn-1tn-1+σf′′′ηns-t0ntn-1+σ-sαn-1+σ-1ds≤c2maxtn-1≤t≤tn⁡f′′′t∫tn-1tn-1+σtn-1+σ-s1-αn-1+σdsτ≤c2σ2-αn-1+σ2-αn-1+σmaxtn-1≤t≤tn⁡f′′′tτ3-αn-1+σ,where c2 is a positive constant.
Substituting (46) into (44) leads to the theorem.Remark 4.
From the second term of the right hand of the truncation error estimate in Theorem2, we can easily obtain the fact that formula (14) does not have exact second-order accuracy; it is related to the selection of α(t) and σ. This can be seen from Table 2 in the numerical example below; its accuracy is much higher than the second order.


## 3. Numerical Verification

In this section, the validity and numerical accuracy of the new presented formula (14) and (40) are demonstrated, respectively. Meanwhile, the corresponding computational results with the formula of [19] are given for contrast.Take a positive integerN, let T=tn-1+σ,τ=T/N, and denote (47)Fn-1+σ=Dtαn-1+σ0Cftn-1+σ,fn-1+σ=Δtαn-1+σ0fn-1+σ,0≤n≤N,ENτ=Fn-1+σ-fn-1+σ,Rate=log2⁡ENτENτ/2.Example 1 (accuracy of formula (14) and [19]).
Takeft=t5-t3,0<t≤T. Compute the Caputo fractional derivative of f(t) at T=1,2 numerically.
The exact solution is given by(48)Dtαn-1+σ0Cftn-1+σ=Γ6Γ6-αn-1+σt5-αn-1+σ-Γ4Γ4-αn-1+σt3-αn-1+σ.
From the results presented in Table1, taking σ=1-α1/2=1/2, we found that although both formula (14) and [19] have second-order approximation accuracy, the numerical accuracy by the formula (14) is significantly higher than that by the formula of [19]; furthermore, the computational errors are also obviously smaller than the formula of [19]. Thus, the new formula (14) is valid for solving complex problems that require high accuracy, and better computational results can be obtained from the formula.
In Table2, taking σ=1-α2/2=0.9323, it shows that the convergence order is higher as compared with the theoretical result, and the error is smaller than the theoretical result, but why this happens is not clear yet; it requires further investigation.
Moreover, on the basis of Table1, if a slight disturbance to σ is added to Tables 3 and 4, the experimental results will be greatly different. It can be easily found that there is a very close relationship between the convergence order and the selection of σ. The experimental results show that the convergence order is the ideal result only when σ=1-α1/2. Tables 2, 5, and 6 also illustrate this fact.Table 1
Maximum errors and convergence rates for Example1 at T=1 and σ=1-α1/2.
α
(
t
)
τ
E
N
(
τ
) (14)
Rate (14)
E
N
(
τ
) [19]
Rate [19]
t
2
+
1
2
1/40
1.400
e
-
3
2.0753
1.726
e
-
3
1.9354
1/80
3.560
e
-
4
1.9755
4.512
e
-
4
1.9421
1/120
1.576
e
-
5
2.0121
2.053
e
-
4
1.9438
1/160
8.844
e
-
5
2.0091
1.173
e
-
4
1.9468
1/200
5.653
e
-
5
—
7.600
e
-
5
—Table 2
Maximum errors and convergence rates for Example1 at T=2 and σ=1-α2/2.
α
(
t
)
τ
E
N
(
τ
) (14)
Rate (14)
E
N
(
τ
)
Rate [19]
e
-
t
1/40
5.830
e
-
4
2.5859
4.834
e
-
2
1.9724
1/80
9.138
e
-
5
2.6734
1.232
e
-
2
1.9843
1/120
3.021
e
-
5
2.7145
5.509
e
-
3
1.9889
1/160
1.362
e
-
5
2.7462
3.109
e
-
3
1.9914
1/200
7.285
e
-
6
—
1.993
e
-
3
—Table 3
Maximum errors and convergence rates for Example1 at T=1 and σ1=σ-0.1.
α
(
t
)
τ
E
N
(
τ
)
Rate
t
2
+
1
2
1/40
3.720
e
-
2
1.0847
1/80
1.800
e
-
2
1.0473
1/120
1.190
e
-
2
1.0300
1/160
8.900
e
-
3
1.0161
1/200
7.100
e
-
3
1.0101Table 4
Maximum errors and convergence rates for Example1 at T=1 and σ2=σ+0.1.
α
(
t
)
τ
E
N
(
τ
)
Rate
t
2
+
1
2
1/40
3.380
e
-
2
0.9412
1/80
1.720
e
-
2
0.9746
1/120
1.150
e
-
2
0.9874
1/160
8.700
e
-
3
0.9833
1/200
7.000
e
-
3
0.9792Table 5
Maximum errors and convergence rates for Example1 at T=2 and σ3=σ+0.1.
α
(
t
)
τ
E
N
(
τ
)
Rate
e
-
t
1/40
3.070
e
-
2
1.8040
1/80
8.600
e
-
3
1.8358
1/120
4.100
e
-
3
1.8323
1/160
2.400
e
-
3
1.8413
1/200
1.600
e
-
3
1.8329Table 6
Maximum errors and convergence rates for Example1 at T=2 and σ4=σ-0.23.
α
(
t
)
τ
E
N
(
τ
)
Rate
e
-
t
1/40
5.390
e
-
2
1.8808
1/80
1.470
e
-
2
1.8745
1/120
6.900
e
-
3
1.8688
1/160
4.000
e
-
3
1.8777
1/200
2.600
e
-
3
1.8995Example 2 (accuracy of formula (40)).
Takeft=t5,f′(0)=0,0<t≤T. Compute the Caputo fractional derivative of f(t) at T=1 numerically.
The exact solution is given by(49)Dtαn-1+σ0Cftn-1+σ=Γ6Γ6-αn-1+σt5-αn-1+σ.
Table7 lists the maximum errors and convergence rates for the variable-order sin⁡(t)+1. In particular, when T=1, α(tn-1+σ)=sin⁡(1)+1=1.8415, and σ=0.5793, the numerical results demonstrate that formula (40) leads to 3-α(t) accuracy.
Similar to Example1, we also give Tables 8 and 9 as a comparison; the results show that σ=3/2-α1/2 is the best choice.Table 7
Maximum errors and convergence rates for Example2 at T=1 and σ=3/2-α1/2.
α
(
t
)
τ
E
N
(
τ
)
Rate
sin
⁡
(
t
)
+
1
1/40
4.023
e
-
1
1.1310
1/80
1.821
e
-
1
1.1435
1/160
8.210
e
-
2
1.1493
1/320
3.690
e
-
2
1.1538
1/640
1.660
e
-
2
1.1524
1/1280
7.400
e
-
3
1.1656Table 8
Maximum errors and convergence rates for Example2 at T=1 and σ5=σ+0.2.
α
(
t
)
τ
E
N
(
τ
)
Rate
sin
⁡
(
t
)
+
1
1/40
5.644
e
-
1
1.1219
1/80
2.562
e
-
1
1.1394
1/160
1.156
e
-
1
1.1481
1/320
5.200
e
-
2
1.1526
1/640
2.33
e
-
2
1.1582
1/1280
1.050
e
-
2
1.1499Table 9
Maximum errors and convergence rates for Example2 at T=1 and σ6=σ-0.3.
α
(
t
)
τ
E
N
(
τ
)
Rate
sin
⁡
(
t
)
+
1
1/40
1.705
e
-
1
1.1216
1/80
7.760
e
-
2
1.1356
1/160
3.510
e
-
2
1.1446
1/320
1.580
e
-
2
1.1515
1/640
7.100
e
-
3
1.1540
1/1280
3.200
e
-
3
1.1497Remark 3.
Due to the fact that the selection ofσ is very important, it is worth noting that different selected σ leads to the different numerical accuracy. From considerable trials, we observed that only if we take σ=1-α(t)/2 for α(t)∈(0,1) and σ=3/2-α(t)/2 for α(t)∈(1,2), will numerical results be closer to theoretical analysis. In addition, carefully observing Table 9, although the convergence order is not higher than in Table 7, the error is lower than in Table 7; further investigation is needed on this phenomenon.


## 4. Conclusion

In this paper, two approximations to variable-order Caputo fractional derivatives were developed, and the analysis for the truncation errors of new formulas was made. In addition, numerical examples support theoretical results. In further work, we also want to give a strict proof of the selection ofσ and apply these new methods to solve more complex cases within a reasonable accuracy.



---

*Source: 1586249-2017-07-31.xml*