Numerical Solutions for the Three-Point Boundary Value Problem of Nonlinear Fractional Differential Equations

Abstract

We present an efficient numerical scheme for solving three-point boundary value problems of nonlinear fractional differential equation. The main idea of this method is to establish a favorable reproducing kernel space that satisfies the complex boundary conditions. Based on the properties of the new reproducing kernel space, the approximate solution is obtained by searching least value techniques. Moreover, uniformly convergence and error estimation are provided for our method. Numerical experiments are presented to illustrate the performance of the method and to confirm the theoretical results.

1. Introduction

Fractional differential equations have gained considerable importance due to their frequent applications in various fields of science and engineering including physics [1–4], bioengineering [5–7], hydrology [8–10], solid mechanics [11], chaos [12–14], control theory [15], and finance [16–18]. It has been found that fractional derivatives provide an excellent instrument for th e description of memory and hereditary properties of different substances [19]. With these features, the fractional-order models become more practical and realistic than the classical integerd-order models, in which such effects are not taken into account.

Finding exact solutions in closed forms for most differential equations of fractional order is a difficult task. As a result, a number of methods have been proposed and applied successfully to approximate fractional differential equations, such as Adomian decomposition method [20], variational iteration method [21], homotopy analysis method [22], implicit and explicit difference method [23], and collocation method [24]. Especially Momani and Odibat have applied He’s variational iteration method to fractional differential equations [25–27]. Meanwhile, various fractional order differential equation have been solved very recently including fractional advection-dispersion equations [28, 29], reaction-diffusion system with fractional derivatives [30], fractional partial differential equations fluid mechanics [31], and fractional-order two-point boundary value problem [32].

In contrast to the initial value and two-point boundary value problems, not much attention has been paid to the multipoint fractional boundary value problem (MFBVP). Ahmed and Wang [33] considered existence and uniqueness of solutions for a four-point nonlocal boundary value problem of nonlinear impulsive differential equations of fractional order. Rehman and Khan [34] also studied existence and uniqueness of solutions for a class of multipoint boundary value problems for fractional differential equations. However, the research on the numerical treatment for MFBVP has proceeded very slowly in the recent years. This motivates us to investigate computationally efficient numerical techniques for solving the MFBVP.

In the present work, we are concerned with the numerical solution of the following three-point boundary value problem of fractional differential equation [35] in a reproducing kernel space:𝐷^𝛼₀₊𝑢(𝑥)((𝑥𝑢(0)=𝑓𝑥,𝑢)),0<𝑥<1,3<𝛼≤4,=𝑢^′(0)=𝑢^′′(0)=0,𝑢^′′(1)=𝛽𝑢^′′𝜂, (1.1) where 𝐷^𝛼₀₊ is the Riemann-Liouville fractional derivative and 0<𝜂<1 satisfies 0<𝛽𝜂^𝛼−3<1, while 𝑓∶[0,1]×[0,+∞)→[0,+∞) is continuous.

Recently, the reproducing kernel space method (RKSM) has been used for obtaining approximate solutions of differential and integral equations [36–38]. However, due to the multipoint boundary value conditions in (1.1), especially for fractional differential equations, it is difficult to find the corresponding reproducing kernel space by applying traditional RKSM.

The aim of this work is to extend the RKSM to derive the numerical solutions of (1.1). One important improvement is that we successfully construct a novel reproducing kernel space so as to overcome difficulties with the nonlocal multipoint boundary value conditions. By using the new reproducing kernel functions, we present an efficient numerical algorithm to solve problem (1.1). The emphasis of the result is that uniformly convergence of the approximate solution, error estimation, and complexity analysis of our algorithm are studied.

The organization of this paper is as follows. In Section 2, we present some important definitions and preparations used in this paper. In Section 3, we construct and develop algorithms for solving nonlinear fractional differential equation with three-point boundary value conditions. In Section 4 the proposed methods are applied to several examples. Also a conclusion is given in Section 5.

2. The Construction of a New Reproducing Kernel Space

We give some basic definitions and theories which are used further in this paper.

Definition 2.1 (see [39]). The Riemann-Liouville fractional derivative of order 𝛼 is defined as 𝐷^𝛼₀₊𝑢(𝑥)=1Γ()𝑑𝑛−𝑠𝑑𝑥^𝑛𝑥₀𝑢(𝑡)(𝑥−𝑡)^{𝛼−𝑛+1}𝑑𝑡. (2.1)Γ(⋅) is the gamma function and 𝑛=[𝛼]+1, [𝛼] denotes the integerd part of number 𝛼.

Definition 2.2 (see [40]). Let 𝑊 be a Hilbert function space on a set 𝑋. 𝑊 is called a reproducing kernel space if and only if for any 𝑥∈𝑋, there exists a unique function 𝐾_𝑥(𝑦)∈𝑊, such that ⟨𝑓,𝐾_𝑥⟩=𝑢(𝑥) for any 𝑢∈𝑊. Meanwhile, 𝐾(𝑥,𝑦)≐𝐾_𝑥(𝑦) is called a reproducing kernel.

We now define a new reproducing kernel space 𝑊 which includes the nonlocal boundary value conditions and give an explicit representation formula for calculation of the reproducing kernel.

Definition 2.3. 𝑊≐𝑊[0,1]={𝑢(𝑥)∣𝑢⁽⁴⁾(𝑥) is an absolutely continuous real value function in [0,1], 𝑢⁽⁵⁾(𝑥)∈𝐿²[0,1], 𝑢(0)=𝑢^′(0)=𝑢^′′(0)=0, 𝑢^′′(1)−𝛽𝑢^′′(𝜂)=0,0<𝜂<1,0<𝛽𝜂^𝛼−3<1}. The inner product is given by⟨𝑢(𝑥),𝑣(𝑥)⟩=𝑢^′′𝜂𝑣^′′𝜂+𝑢⁽⁴⁾(0)𝑣⁽⁴⁾(0)+¹₀𝑢⁽⁵⁾(𝑥)𝑣⁽⁵⁾(𝑥)𝑑𝑥. (2.2)

Theorem 2.4. 𝑊 is a reproducing kernel space, that is, there exists a function 𝐾(𝑥,𝑦)∈𝑊, for any fixed 𝑦∈[0,1] and any 𝑢(𝑥)∈𝑊, such that 𝑢(𝑦)=⟨𝑢(𝑥),𝐾(𝑥,𝑦)⟩. Moreover, the reproducing kernel can be denoted by 𝐾=𝑘𝑥,𝑦¹𝑥,𝑦=𝑎₁+𝑎₂𝑥+𝑎₃𝑥²+⋅⋅⋅+𝑎₁₀𝑥⁹𝑘,𝑦<𝑥<𝜂,₂𝑥,𝑦=𝑏₁+𝑏₂𝑥+𝑏₃𝑥²+⋅⋅⋅+𝑏₁₀𝑥⁹𝑘,𝑦<𝜂<𝑥,₃𝑥,𝑦=𝑐₁+𝑐₂𝑥+𝑐₃𝑥²+⋅⋅⋅+𝑐₁₀𝑥⁹𝑘,𝜂<𝑦<𝑥,₄𝑥,𝑦=𝑑₁+𝑑₂𝑥+𝑑₃𝑥²+⋅⋅⋅+𝑑₁₀𝑥⁹𝑘,𝜂<𝑥<𝑦,₅𝑥,𝑦=𝑒₁+𝑒₂𝑥+𝑒₃𝑥²+⋅⋅⋅+𝑒₁₀𝑥⁹𝑘,𝑥<𝜂<𝑦,₆𝑥,𝑦=𝑓₁+𝑓₂𝑥+𝑓₃𝑥²+⋅⋅⋅+𝑓₁₀𝑥⁹,𝑥<𝑦<𝜂. (2.3)

Proof. For any 𝑢(𝑥)∈𝑊, we only need to prove that there exists a 𝐾(𝑥,𝑦)∈𝑊 for any fixed 𝑦∈[0,1] and any 𝑢(𝑥)∈𝑊, 𝐾(𝑥,𝑦) must satisfy 𝑢(𝑥),𝐾𝑥,𝑦=𝑢𝑦. (2.4) By (2.2), we get 𝑢(𝑥),𝐾𝑥,𝑦=𝑢^′′𝜂𝜕²𝐾𝑥,𝑦𝜕𝑥²_𝑥=𝜂+𝑢⁽⁴⁾(0)𝜕⁴𝐾𝑥,𝑦𝜕𝑥⁴_𝑥=0+¹₀𝑢⁽⁵⁾(𝑥)𝜕⁵𝐾𝑥,𝑦𝜕𝑥⁵𝑑𝑥. (2.5) We have the following equality using the integration by parts: ¹₀𝑢⁽⁵⁾(𝑥)𝜕⁵𝐾𝑥,𝑦𝜕𝑥⁵𝑑𝑥=⁴_𝑖=0(−1)^4−𝑖𝑢^(𝑖)(𝑥)𝜕^9−𝑖𝐾𝑥,𝑦𝜕𝑥^9−𝑖1₀−¹₀𝑢(𝑥)𝜕¹⁰𝐾𝑥,𝑦𝜕𝑥¹⁰𝑑𝑥. (2.6) Substituting (2.6) in (2.5), we get 𝑢(𝑥),𝐾𝑥,𝑦=𝑢^′′𝜂𝜕²𝐾𝑥,𝑦𝜕𝑥²_𝑥=𝜂𝜕−𝛽⁷𝐾𝑥,𝑦𝜕𝑥⁷_𝑥=1+𝑢⁽⁴⁾(0)𝜕⁴𝐾𝑥,𝑦𝜕𝑥⁴_𝑥=0−𝜕⁵𝐾𝑥,𝑦𝜕𝑥⁵_𝑥=0+𝑢⁽³⁾(0)𝜕⁶𝐾𝑥,𝑦𝜕𝑥⁶_𝑥=0+¹_𝑖=0(−1)^4−𝑖𝑢^(𝑖)(1)𝜕^9−𝑖𝐾𝑥,𝑦𝜕𝑥^9−𝑖_𝑥=1+⁴_𝑖=3(−1)^4−𝑖𝑢^(𝑖)(1)𝜕^9−𝑖𝐾𝑥,𝑦𝜕𝑥^9−𝑖_𝑥=1−¹₀𝑢(𝑥)𝜕¹⁰𝐾𝑥,𝑦𝜕𝑥¹⁰𝑑𝑥. (2.7) Therefore, by (2.4) 𝐾(𝑥,𝑦) is the solution of the following generalized differential equation: 𝜕²𝐾𝑥,𝑦𝜕𝑥²_𝑥=𝜂𝜕−𝛽⁷𝐾𝑥,𝑦𝜕𝑥⁷_𝑥=1𝜕=0,₄𝐾𝑥,𝑦𝜕𝑥⁴_𝑥=0−𝜕⁵𝐾𝑥,𝑦𝜕𝑥⁵_𝑥=0𝜕=0,₆𝐾𝑥,𝑦𝜕𝑥⁶_𝑥=0𝜕=0,_9−𝑖𝐾𝑥,𝑦𝜕𝑥^9−𝑖_𝑥=1−𝜕=0,𝑖=0,1,3,4,₁₀𝐾𝑥,𝑦𝜕𝑥¹⁰=𝛿𝑥−𝑦, (2.8) where 𝛿 denotes 𝛿 function. For 𝑥≠𝑦, it is known that 𝐾(𝑥,𝑦) is the solution of the following linear homogeneous differential equation: 𝜕¹⁰𝐾𝑥,𝑦𝜕𝑥¹⁰=0, (2.9) with the boundary value conditions: 𝜕²𝐾𝑥,𝑦𝜕𝑥²_𝑥=𝜂𝜕−𝛽⁷𝐾𝑥,𝑦𝜕𝑥⁷_𝑥=1𝜕=0,₄𝐾𝑥,𝑦𝜕𝑥⁴_𝑥=0−𝜕⁵𝐾𝑥,𝑦𝜕𝑥⁵_𝑥=0𝜕=0,₆𝐾𝑥,𝑦𝜕𝑥⁶_𝑥=0𝜕=0,_9−𝑖𝐾𝑥,𝑦𝜕𝑥^9−𝑖_𝑥=1=0,𝑖=0,1,3,4. (2.10) We find that (2.9) has characteristic equation 𝜆¹⁰=0, and the eigenvalue 𝜆=0 is a root whose multiplicity is 10. Applying the feature of functions in 𝑊, the general solution of (2.9) is given by the universal representation as (2.3), in which every function 𝑘_𝑖(𝑥,𝑦)(𝑖=1,2,3,…,6) has the situation in Figure 1. Next we will calculate 60 coefficients in (2.3).
By integrating repeatedly 𝜕¹⁰𝐾(𝑥,𝑦)/𝜕𝑡¹⁰=−𝛿(𝑥−𝑦) from 𝑦−𝜀 to 𝑦+𝜀 with respect to 𝑥, we have 𝜕⁹𝐾𝑥,𝑦𝜕𝑥^9𝑦+𝜀_{𝑥=𝑦−𝜀}𝜕=−1,^𝑖𝐾𝑥,𝑦𝜕𝑥^{𝑖𝑦+𝜀}_{𝑥=𝑦−𝜀}=0,𝑖=8,7,6,…,0, (2.11) and the following 20 equations are inferred as 𝜀→0𝜕⁹𝑘₁𝑥,𝑦𝜕𝑥⁹_𝑥=𝑦−𝜕⁹𝑘₆𝑥,𝑦𝜕𝑥⁹_𝑥=𝑦𝜕=−1,₉𝑘₃𝑥,𝑦𝜕𝑥⁹_𝑥=𝑦−𝜕⁹𝑘₄𝑥,𝑦𝜕𝑥⁹_𝑥=𝑦𝜕=−1,_𝑖𝑘₁𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝑦−𝜕^𝑖𝑘₆𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝑦𝜕=0,𝑖=8,7,6,…,0,_𝑖𝑘₃𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝑦−𝜕^𝑖𝑘₄𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝑦=0,𝑖=8,7,6,…,0. (2.12) By (2.10), one can obtain 14 equations. In view of boundary value conditions, the following 8 equations can be obtained: 𝜕^𝑖𝑘₅𝑥,𝑦𝜕𝑥^𝑖_𝑥=0(),𝜕=0𝑖=0,1,2^𝑖𝑘₆𝑥,𝑦𝜕𝑥^𝑖_𝑥=0(),𝜕=0𝑖=0,1,2₂𝑘₂𝑥,𝑦𝜕𝑥²_𝑥=1𝜕−𝛽²𝑘₂𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂𝜕=0,²𝑘₃𝑥,𝑦𝜕𝑥²_𝑥=1𝜕−𝛽²𝑘₄𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂=0. (2.13) Finally, due to the smoothing of functions in 𝑊, we get 18 equations for 𝑖=0,1,2,…,8𝜕^𝑖𝑘₁𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂−𝜕^𝑖𝑘₂𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂𝜕=0,_𝑖𝑘₄𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂−𝜕^𝑖𝑘₅𝑥,𝑦𝜕𝑥^𝑖_𝑥=𝜂=0. (2.14) Hence, the unknown coefficients of 𝐾(𝑥,𝑦) are governed by solving the 60 independent equations by (2.10)–(2.14).

Figure 1 

Plot of distribution for 𝑘𝑖(𝑥,𝑦), 𝑖=1,2,…,6.

3. Description of the Proposed Numerical Method

The method consists of two steps. In the first step, a normal orthogonal basis is established in the reproducing kernel space 𝑊, and in the second step, it is used to successively obtain the approximate solution of (1.1). Let us consider these steps in detail.

3.1. A Normal Orthogonal Basis in 𝑊

Define a bounded linear operator T∶𝑊[0,1]→𝐿¹[0,1] satisfyingT𝑢(𝑥)=𝐷^𝛼₀₊𝑢(𝑥). (3.1)

The proof of existence and uniqueness of solution for (1.1) has been studied in [35]. Therefore T is reversible. Now, (1.1) is turned into the following operator equation in 𝑊:T𝑢(𝑥)=𝑓(𝑥,𝑢(𝑥)). (3.2) Choosing a countable dense subset {𝑥_𝑖}^∞_𝑖=1 on [0,1], for the reproducing kernel 𝐾(𝑥,𝑦) of 𝑊, we define a complete system in 𝑊 as 𝜓_𝑖(𝑥)=(T𝐾(𝑥,(⋅)))𝑥_𝑖=1Γ()𝜕4−𝛼𝜕𝑥⁴⁴_𝑖^𝑥𝑖₀𝐾𝑥,𝑦𝑥_𝑖−𝑦^𝛼−3𝑑𝑦. (3.3) Then the orthogonal system of 𝑊 is derived from Gram-Schmidt orthogonalization process, namely, 𝜓_𝑖(𝑥)=^𝑖_𝑘=1𝛽_𝑖𝑘𝜓_𝑘(𝑥). (3.4) Next, the complexity estimation of the orthogonal basis is discussed. We know that the orthogonal basis is obtained by orthogonalization of complete system 𝜓_𝑖(𝑥). The algorithm with time complexity may be analyzed as following.

Step 1 (Computing ⟨𝜓_𝑖(𝑥),𝜓_𝑘(𝑥)⟩). In fact, according to the properties of reproducing kernel and bounded linear operator, we have 𝜓_𝑖(𝑥),𝜓_𝑘(𝑥)=T𝜓_𝑖𝑥_𝑘. (3.5) Thus, we only need to calculate the specific function value T𝜓_𝑖(𝑥_𝑘) instead of the usual integral. Denote the computing time of the specific function value by 𝑇. It needs (((𝑛(𝑛+1)))/2)𝑇 to compute those ⟨𝜓_𝑖(𝑥),𝜓_𝑘(𝑥)⟩𝑘=1,2,…,𝑖, 𝑖=1,2,…,𝑛.

Step 2. Orthogonalization can be obtained by the following cycle.
Firstly, let 𝛽₁₁=1/‖𝜓₁‖ and 𝐶_𝑖𝑗=𝜓_𝑖,𝜓_𝑗𝐻=𝜓_𝑖,^𝑖_𝑘=1𝛽_𝑗𝑘𝜓_𝑘=^𝑖_𝑘=1𝛽_𝑗𝑘𝜓_𝑖,𝜓_𝑘. (3.6) For 1≤𝑗≤𝑖−1, it takes (𝑖(𝑖−1))/2 times multiplication operations for computing 𝐶_𝑖𝑗. Secondly, we denote 𝛽_𝑖𝑚=−∑^𝑖−1_𝑘=𝑚𝐶_𝑖𝑘𝛽_𝑘𝑚𝜓_𝑖, (3.7) where 𝜓_𝑖=𝜓_𝑖²−^𝑖−1_𝑘=1𝐶_𝑖𝑘^21/2=𝜓_𝑖,𝜓_𝑖−^𝑖−1_𝑘=1𝐶²_𝑖𝑘^1/2. (3.8) For 1≤𝑚≤𝑖−1, ((𝑖+2)(𝑖−1))/2 times multiplication operations are demanded to compute 𝛽_𝑖𝑚. Finally, for 2≤𝑖≤𝑛, it uses all (𝑛(𝑛−1)(2𝑛+3))/4 times multiplication operations. To sum up, from Step 1 and Step 2, we get orthogonal basis complexity result of (𝑛(𝑛−1)(2𝑛+3))/4times multiplication operations plus ((𝑛(𝑛+1))/2)𝑇. Therefore, the construction of orthogonal basis costs a total of 𝒪(𝑛³) operations.

3.2. The Approximate Solution of (3.2) in 𝑊

Theorem 3.1. If 𝑢(𝑥) is the solution of (3.2), then 𝑢(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝑓𝑥_𝑘,𝑢𝑥_𝑘𝜓_𝑖(𝑥). (3.9)

Proof. According to the orthogonal basis {𝜓_𝑖(𝑥)}^∞_𝑖=1 of 𝑊, we have 𝑢(𝑥)=^∞_𝑖=1𝑢(⋅),𝜓_𝑖(𝑥)𝜓_𝑖(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝑢(⋅),𝜓_𝑘(𝑥)𝜓_𝑖(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝑢(⋅),T𝐾𝑥_𝑘,⋅𝜓_𝑖(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘T𝑢(⋅),𝐾𝑥_𝑘,⋅𝜓_𝑖(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘(T𝑢)𝑥_𝑘𝜓_𝑖(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝑓𝑥_𝑘,𝑢𝑥_𝑘𝜓_𝑖(𝑥). (3.10) Denoting 𝛼_𝑘=𝑓(𝑥_𝑘,𝑢(𝑥_𝑘)), (3.9) can be rewritten as𝑢(𝑥)=^∞𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝛼_𝑘𝜓_𝑖(𝑥). (3.11) To obtain 𝛼_𝑘, we truncate the series of the left-hand side of (3.11) 𝑢_𝑛(𝑥)=^𝑛𝑖_𝑖=1𝑘=1𝛽_𝑖𝑘𝛼_𝑘𝜓_𝑖(𝑥). (3.12) Then we get 𝛼_𝑘 based on the minimum point of function 𝐽𝛼₁,𝛼₂,…,𝛼_𝑛=^𝑛_𝑘=1𝑓𝑥_𝑘,𝑢_𝑛𝑥_𝑘−𝛼_𝑘². (3.13) Consequently, the approximate solution of (3.2) can be obtained.

Lemma 3.2. The approximate solution 𝑢_𝑛(𝑥) and its derivatives 𝑢_𝑛^(𝑖)(𝑥), (𝑖=1,2,3,4) uniformly converge to exact solution 𝑢(𝑥) and its derivatives, respectively.

Proof. For any 𝑥∈[0,1], according to the boundedness of ‖𝜕^𝑖_𝑥𝐾(𝑥,𝑦)‖(𝑖=0,1,2,3,4) and reproducing property of 𝐾(𝑥,𝑦), the following conclusion is obtained as 𝑛→∞𝑢_𝑛^(𝑖)(𝑥)−𝑢^(𝑖)(𝑥)=𝑢_𝑛(𝑥)−𝑢(𝑥)^(𝑖)=𝜕^𝑖_𝑥𝑢_𝑛𝑦𝑦=𝑢−𝑢,𝐾𝑥,𝑦_𝑛𝑦−𝑢𝑦,𝜕^𝑖_𝑥𝐾≤𝑢𝑥,𝑦_𝑛−𝑢𝜕^𝑖_𝑥𝐾𝑥,𝑦≤𝑀𝑢_𝑛−𝑢→0. (3.14) Finally, we give the process to obtain 𝛼_𝑘:(1^∘) pick any initial value 𝛼⁰_𝑘;(2^∘) substitute 𝛼⁰_𝑘 into (3.12) and get 𝑢_𝑛(𝑥);(3^∘) substitute 𝛼⁰_𝑘, 𝑢_𝑛(𝑥) into (3.13) and compute the value 𝐽(𝛼⁰₁,𝛼⁰₂,…,𝛼⁰_𝑛);(4^∘) if 𝐽(𝛼⁰₁,𝛼⁰₂,…,𝛼⁰_𝑛)<10⁻²⁰ then the computations terminate, otherwise substituting 𝑢_𝑛(𝑥) into (3.13) yields a new minimum point 𝛼¹_𝑘;(5^∘) calculate 𝐽(𝛼¹₁,𝛼¹₂,…,𝛼¹_𝑛);(6^∘) if 𝐽(𝛼¹₁,𝛼¹₂,…,𝛼¹_𝑛)<𝐽(𝛼⁰₁,𝛼⁰₂,…,𝛼⁰_𝑛) then replace 𝛼⁰_𝑘 with 𝛼¹_𝑘 and return to 2°, otherwise give up 𝛼¹_𝑘 and return to 1° and pick another initial value.

Theorem 3.3. |𝑢_𝑛(𝑥)−𝑢(𝑥)|=𝑜(1/𝑛).

Proof. Firstly, because of denseness, for any 𝑥∈[0,1] and 𝑛∈𝑁, we take 𝑥_𝑖∈{𝑥₁,𝑥₂,…}, 𝑖≤𝑛 such that |𝑥−𝑥_𝑖|<1/𝑛. Then due to the reproducing property and the property of projector 𝑃_𝑛, it follows that T𝑢_𝑛𝑥_𝑖=𝑢_𝑛(⋅),T𝐾𝑥_𝑖,⋅=𝑃_𝑛𝑢(⋅),𝜓_𝑖(⋅)=𝑢(⋅),𝑃_𝑛𝜓_𝑖(⋅=𝑢(⋅)),𝜓_𝑖(⋅)=T𝑢(⋅),𝐾𝑥_𝑖,⋅=T𝑢𝑥_𝑖. (3.15) This implies that T𝑢_𝑛(𝑥)−T𝑢(𝑥)=T𝑢_𝑛(𝑥)−T𝑢_𝑛𝑥_𝑖−T𝑢(𝑥)−T𝑢𝑥_𝑖=𝑢_𝑛(⋅),T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖,⋅−𝑢(⋅),T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖=𝑢,⋅_𝑛(⋅)−𝑢(⋅),T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖,⋅. (3.16) By the mean value theorem, we have T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖,⋅=𝑥−𝑥_𝑖𝜕.𝜕𝑦T𝐾𝑦,⋅ (3.17) Finally, the following conclusion follows from above: 𝑢_𝑛(𝑥)−𝑢(𝑥)=𝑢_𝑛−𝑢,T⁻¹T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖≤T,⋅₋₁𝑢_𝑛−𝑢T𝐾(𝑥,⋅)−T𝐾𝑥_𝑖≤T,⋅₋₁𝑢_𝑛−𝑢𝑥−𝑥_𝑖𝜕.𝜕𝑦T𝐾𝑦,⋅ (3.18) Thus, according to ‖𝑢_𝑛−𝑢‖→0, |𝑥−𝑥_𝑖|<1/𝑛 and the boundedness of ‖(𝜕/𝜕𝑦)T𝐾(𝑦,⋅)‖, we get |𝑢_𝑛(𝑥)−𝑢(𝑥)|=𝑜(1/𝑛).

4. Numerical Experiments

In this section, we give some computational results of three numerical experiments with methods based on preceding sections, to support our theoretical discussion.

Example 4.1. Consider the following three-point boundary value problem of nonlinear fractional differential equation [35]: 𝐷₀^7/2₊𝑢(𝑥)((𝑥𝑢(0)=𝑓𝑥,𝑢)),0<𝑥<1,=𝑢^′(0)=𝑢^′′(0)=0,𝑢^′′(1)=𝑢^′′14, (4.1) with an exact solution 𝑢(𝑥)=(𝑥³/100)(425+170𝑥−128𝑥³). Here 𝑓(𝑥,𝑢(𝑥))=ln(2+((17𝑥³)/4) + ((17𝑥⁴)/10)+((32𝑥⁶)/25))+(1275+4080𝑥−24576𝑥³)/50√𝑥Γ(0.5)−ln(2+𝑢(𝑥)).
For this example, we solved the three-point boundary value problem, by applying the technique described in preceding section as following.
Step 1. By using the representation formula for calculation of reproducing kernel in Section 2, the concrete expression of 𝐾(𝑥,𝑦) for the three-point boundary value conditions in (4.1) is given as 𝑘₁=𝑦𝑥,𝑦³×5095276462080−21𝑥³(−3141613209+64𝑥(16138815+𝑥(3227763+4𝑥(−219394+3007𝑥))))+20160𝑥³−1075921+4𝑥108017−336𝑥+64𝑥³𝑦+4032𝑥³−1075921+4𝑥108017−336𝑥+64𝑥³𝑦²−768𝑥²−658182+𝑥21049+64𝑥−105−21𝑥+4𝑥³𝑦⁴−126370944𝑥𝑦⁵[phone]𝑦⁶;𝑘₂=𝑦𝑥,𝑦³×(((((()𝑥(((((()𝑥𝑦((((((()𝑥𝑦815244233932801051889+64𝑦+4−21𝑥66115+4𝑥−198345+𝑥−3140290909+16𝑥64224685+652592𝑥21+−7+𝑥)))))−6720𝑥7+4𝑥−21+𝑥3227903+4𝑥−324191+64𝑥21+−7+𝑥)))))−336−1+4𝑥7+4𝑥−21+𝑥3227903+4𝑥−324191+64𝑥21+−7+𝑥))))))²((((((()𝑥𝑦+64−1+4𝑥7+4𝑥1974525+𝑥−63007+16𝑥1225+16𝑥21+−7+𝑥))))))⁴−505483𝑥𝑦⁵+56164𝑦⁶;𝑘₃=1𝑥,𝑦×𝑥(((((()𝑥(((((()𝑥𝑦((((((()𝑥𝑦3913172322877440−1+28𝑦+47+4𝑥−21+𝑥595175+16𝑥1225+16𝑥21+−7+𝑥)))))−28𝑥7+4𝑥−21+𝑥595175+16𝑥1225+16𝑥21+−7+𝑥)))))+84−1+4𝑥7+4𝑥−21+𝑥59575+16𝑥1225+16𝑥21+−7+𝑥))))))²((((((()𝑥𝑦−28−85025+4𝑥595175+4𝑥−1785525+𝑥−28250714681+16𝑥5780465+5873648𝑥21+−7+𝑥))))))³((((((()𝑥𝑦−15680−5+4𝑥35+4𝑥−105+𝑥16515619+256𝑥−26062+5𝑥21+−7+𝑥))))))⁴−21504−1+4𝑥7+4𝑥−21+𝑥25697𝑦⁵+516677𝑦⁶+16𝑥(1225+(−7+𝑥)𝑥)𝑦⁷−242632𝑥𝑦⁸+269591𝑦⁹;𝑘₄𝑥,𝑦=𝑘₃𝑦,𝑥;𝑘₅𝑥,𝑦=𝑘₂𝑦,𝑥;𝑘₆𝑥,𝑦=𝑘₁𝑦,𝑥. (4.2)Step 2. According to the numerical algorithm in Section 3, we get 𝜓_𝑖(𝑥)=1Γ()𝜕0.5𝜕𝑥⁴⁴_𝑖^𝑥𝑖₀𝐾𝑥,𝑦√𝑥_𝑖−𝑦𝑑𝑦. (4.3) By (3.4) we get the orthogonalization coefficients 𝛽_𝑖𝑘 and 𝜓_𝑖(𝑥). Then the approximate solution 𝑢_𝑛(𝑥) can be obtained iteratively by (3.12).
The obtained numerical results are displayed in Table 1. Furthermore, the graph of 𝑢(𝑥) and 𝑢_𝑛(𝑥) for 𝑛=128 is plotted in Figure 2. It can be shown that the numerical solutions agree with exact solution by means of the proposed method.
Due to the multipoint boundary value conditions, to our knowledge, there is no the same example as (1.1) in the literature about numerical method. For the purpose of comparison, we compare the approximate solution of our method, together with the approximate solution by Adomian Decomposition (AD) method given in [41].

Table 1 

Numerical results for 𝑛=128 for Example 4.1.

Figure 2 

Comparison of the approximate solution 𝑢𝑛(𝑥) and 𝑢(𝑥) at 𝑛=128 for Example 4.1.

Example 4.2. In this example, we consider the following nonlinear fourth-order fractional integrodifferential equation of the form: 𝐷^𝛼₀₊𝑦(𝑥)=1+^𝑥₀𝑒^−𝑡𝑦²(𝑡)𝑦(0)𝑑𝑡,0<𝑥<1,3<𝛼≤4,=𝑦^′′(0)=1,𝑦(1)=𝑦^′′(1)=𝑒. (4.4) In this case, a similar numerical method can be proposed as (1.1). In brief, the conditions 𝑢(0)=𝑢^′(0)=𝑢^′′(0)=0,𝑢^′′(1)=𝛽𝑢^′′(𝜂) in the reproducing kernel space 𝑊 are replaced by 𝑦(0)=𝑦^′′(0)=1,𝑦(1)=𝑦^′′(1)=𝑒. There are some slight changes in the process of derivation. The obtained numerical results for 𝛼=3.25 are displayed in Table 2.

Table 2 

The approximate solution in our method and in [41] for Example 4.2.

Example 4.3. Consider the nonlinear fifth-order fractional differential equation with three-point boundary value conditions as following: 𝐷₀^9/2₊𝑢(𝑥)=sin19𝑥⁵−(760/93)𝑥⁷+1520√𝑥93−448𝑥²()(𝑥)𝑢(0)31Γ0.5−sin𝑢,0<𝑥<1,=𝑢^′(0)=𝑢^′′(0)=𝑢^′′′(0)=0,𝑢^′′(1)=𝑢^′′12. (4.5) According to the numerical scheme above, one can obtain the approximation 𝑢_𝑛(𝑥) of 𝑢(𝑥) for 𝑛=128. The numerical results are displayed in Table 3 and Figure 3, which show that the numerical solutions agree with exact solution.

Table 3 

Numerical results for 𝑛=128 for Example 4.3.

Figure 3 

Comparison of the approximate solution 𝑢𝑛(𝑥) and 𝑢(𝑥) at 𝑛=128 for Example 4.3.

5. Conclusion

In this paper, the RKSM is applied to derive approximate analytical solution of nonlinear fractional-order differential equations with three-point boundary value conditions. We have constructed a novel reproducing kernel space and give the way to express the reproducing kernel function, while traditional RKSM still can not be mentioned. The explicit series solution is obtained using the orthogonal basis established in the new reproducing kernel space. The numerical results given in the previous section demonstrate the better accuracy of our algorithms. Moreover, the numerical algorithms introduced in this paper can be well suited for handling general linear and nonlinear fractional-order differential equations with multipoint boundary conditions. We note that the corresponding analytical and numerical solutions are obtained using Mathematica 7.0.

Acknowledments

This paper was supported by Youth Foundation of Heilongjiang Province under grant QC2010036 and also supported by Fundamental Research Funds for the Central Universities under Grant no.HIT.NSRIF.2009050 and Academic Foundation for Youth of Harbin Normal University 11KXQ-04 and 10KXQ-05.