Successive Approximation Technique in the Study of a Nonlinear Fractional Boundary Value Problem

: We studied one essentially nonlinear two–point boundary value problem for a system of fractional differential equations. An original parametrization technique and a dichotomy-type approach led to investigation of solutions of two “model”-type fractional boundary value problems, containing some artiﬁcially introduced parameters. The approximate solutions of these problems were constructed analytically, while the numerical values of the parameters were determined as solutions of the so-called “bifurcation” equations.


Introduction
Fractional differential equations have been of high interest during recent decades. The variety of their applications in biology, physics, engineering and economics (see [1][2][3][4]) has led to development of techniques to study the qualitative behavior of solutions of these equations.
Particular attention has to be paid to the class of nonlinear fractional boundary value problems (FBVPs). Construction of their exact solutions may be impossible or one may even face computational difficulties trying to find their analytical representation. However, precise approximate methods may help to simplify and even solve this task.
In the current paper, we provide a new view on the successive approximations approach [5], recently used for the study of FBVPs for periodic, Cauchy-Nicoletti-type, and interpolation boundary conditions (see [6][7][8][9][10]). An original parametrization technique, initially suggested for the reduction of nonlinearities in boundary restrictions (see discussions [11,12]), and a dichotomy-type approach (see results in [13][14][15][16]) led to the investigation of solutions of two "model"-type FBVPs, containing some artificially introduced parameters. The approximate solutions of these problems were constructed analytically, while the numerical values of the parameters were determined as solutions of the so-called "bifurcation" equations.
The novel technique suggested in this paper has never been applied in study of FBVPs. It allowed us to improve and essentially sharpen the estimates obtained in [6][7][8][9][10]. Along with the other well-known approximate methods, dealing with the fractional differential equations and their systems (see discussions in [17][18][19][20][21]), the aforementioned approach complements the fundamental study of such essentially nonlinear problems.

Main Notations and Definitions
For a fixed n ∈ N and bounded set D ⊂ R n , the following notations apply: • For any vector x = col(x 1 , x 2 , . . . , x n ) ∈ R n and n × n real matrix A operations "| · |," "=," "≤," "≥," and "max" are understood componentwise; • I n is a unit n-dimensional matrix; • O n is a zero n-dimensional matrix; • r(A) is the maximal (in modulus) eigenvalue of matrix A. Definition 1. [2] The left and right Riemann-Liouville fractional integrals of order α ∈ R + are defined by: and: respectively, provided the right-hand sides are pointwise defined on [a, b]. Definition 2. [2] The left and right Riemann-Liouville fractional derivatives of order α ∈ R + are defined by: d n dt n t a (t − s) n−α−1 f (s)ds, t > a and: respectively, where n = [α] + 1, and [α] is the integer part of α.

Definition 4.
For any non-negative vector ρ ∈ R n under the componentwise ρ−neighborhood of a point z 0 ∈ R n , we understand: Definition 5. For a given bounded connected set D ⊂ R n , we introduce its componentwise ρ−neighborhood as follows: Definition 6. For a set D ⊂ R n , closed interval [a, b] ⊂ R, Caratheodory function f : [a, b] × D → R n and the n-dimensional square matrix K with non-negative entires, we write: if the inequality: holds for all {u, v} ⊂ D and a.e. t ∈ [a, b].

Nonlinear FBVP and Its Decomposition
Consider a system of fractional differential equations (FDEs): for some p ∈ (0, 1], subjected to the two-point boundary constraints: where C a D p t is the generalized Caputo fractional derivative with the lower limit at a, defined by (1); f : G f → R n and g : D × D → R n are continuous vector functions, and D ⊂ R n is a closed and bounded domain.
To solve this problem, we used the so-called "freezing" (or parametrization) technique (see discussions in [11,12,16]), coupled with the modification of the numerical-analytic method [5]. The aim of such an approach consists of the introduction of an appropriate parametrization with further reduction of the original FBVP (6) and (7) with nonlinear boundary conditions to two problems, containing already linear separated boundary constraints. Numerical values of the introduced parameters are then evaluated from the corresponding determining system of algebraic equations, one of which is the relation (7).
As a first step in study of the FBVP (6) and (7), let us fix three closed sets D a , D a+b 2 , D b ⊂ R n , where we look for the initial values of solutions x(·): Without loss of generality, we chose these sets to be convex. Then, we formally substitute the boundary and intermediate values of solution x(·) by vector parameters: z := x(a), where z = col(z 1 , z 2 , . . . , z n ) ∈ D a , λ = col(λ 1 , λ 2 , . . . , λ n ) ∈ D a+b 2 , and η = col(η 1 , η 2 , . . . , η n ) ∈ D b . Using relation (9), we reduce the study of the original problem (6) and (7) to an investigation of solutions of two decomposed "model" problems with separated linear boundary conditions, dependent on parameters: and: with z ∈ D a , λ ∈ D a+b 2 , and η ∈ D b defined by (9). In addition, we connect to the parametrized FBVPs (10)-(13) the correspondent sets: and their ρ u -, ρ v -neighborhoods of the form: where z ∈ D a , λ ∈ D a+b 1], and ρ u , ρ v are defined in accordance with Definition 4. Remark 1. The parametrization (9) reduces the study of the original FBVP (6) and (7) with nonlinear boundary conditions on the full interval [a, b] to the investigation of two decomposed problems (10), (11) and (12), (13), defined on the half intervals a, a+b 2 and a+b 2 , b , respectively. This approach allows to diminish some values in the qualitative analysis of the given FBVP and to essentially improve the estimates of the iteration schemes, presented in the coming sections (see also discussions in [6][7][8][9][10]).

Remark 2.
It is also worth emphasizing that the set of solutions of the FBVP (6) and (7) coincides with the set of solutions of the modified problems (10), (11), and (12), (13) under additional conditions (9).

Auxiliary Statements
In this section, we formulate some auxiliary lemmas, needed later on. In terms of fractional integrals, they were first proven in [6] for the interval [0, T] and later generalized over the interval [a, b] (see discussion in [10]).
Then, for all t ∈ [a, b], the following estimate is true: where I := b−a 2 : and Γ(·) denotes the Gamma function.
For the detailed proofs of Lemmas 1 and 2, we refer the reader to the discussions presented in [6,10].

Successive Approximation Technique on the Half Intervals
To derive the recursive sequence of functions, approximating solutions of the auxiliary FBVP (10) and (11), we first consider a perturbed system: coupled with the parametrized boundary conditions (11): where ∆ ∈ R n is an unknown vector to be defined. Direct integration shows that the general solution of (22) can be written in the form: After substituting (23) into the boundary restrictions (11), and taking into account the parametrization (9), we can obtain: and: Relation (24) is satisfied, since it corresponds to the first of the boundary conditions in (11). On the other hand, the relation (25) will hold if the perturbation term ∆ is defined as: Hence, we can obtain the exact solution to the perturbed FBVP (22), (11) in the form: Now, assume that f ∈ Lip(K u , D u ) and the ρ u -neighborhood of domain D u satisfies an inequality: where: Based on (26), we introduce an iterative sequence of functions: , associated with the parametrized BVP (10) and (11).
Using a similar approach to (22) and ρ v satisfying an inequality: we connect to the parametrized FBVP (12), (13) the corresponding sequence of functions: Remark 3. Note that the sequences (29), (30) and (33), (34) are derived in such a way that they satisfy the parametrized boundary conditions (11) and (13) beforehand.
For the sequence of functions (30), the following convergence theorem holds.

Theorem 1.
Let, for the parametrized FBVP (10) and (11), there exist a non-negative vector ρ u satisfying inequality (27), such that f ∈ Lip(K u , D u ) on interval t ∈ a, a+b 2 . In addition, assume that for the matrix: inequality: holds.
Then, for an arbitrary pair of vector parameters (z, λ) ∈ D a × D a+b 2 : 1.
All functions of the sequence are continuous on the interval a, a+b 2 and satisfy the linear parametrized boundary conditions (11).

2.
The sequence of functions (30) converges uniformly as m → ∞ to its limit function: for all t ∈ a, a+b 2 . 3.
The limit function (37) satisfies boundary conditions: and is a unique solution of integral equation: in domain D u . In other words, it is a solution of the corresponding Cauchy problem for a perturbed system of FDEs: where ∆ : D a × D a+b 2 → R n is a mapping given by formula: 4.
The following error estimation holds: where Q u and M u are defined by (35) and (28), respectively.
Proof. Assertion 1 follows from the theorem's assumptions and by direct substitution of (30) into the parametrized boundary conditions (11). Now, we prove that, for all m ∈ N functions u m of the sequence (30) remain in their domain of definition D u , and that (30) is the Cauchy sequence in the Banach space C a, a+b 2 , R n . For this purpose, let us estimate the differences: where functions u 0 (·, z, λ) and u m (·, z, λ) are given by formulas (29) and (30). For every m ∈ N, we get: where α 1 t, a, b−a 2 and M u are defined by (19) and (28), respectively. To prove the error estimate (43), we analyse the difference: where u m (·, z, λ) are functions of the sequence (30).
For m = 0 from the inequality (44), we conclude that: Using the method of mathematical induction, we derive, that for the general case of the iteration step m: In view of inequality (45), we obtain an estimate: for α m+k (t), Q u , and M u defined by (20), (35), and (28), respectively. Under condition (36), the maximal eigenvalue r(Q u ) of matrix Q u does not exceed 1. This means that: Passing in (46) to the limit for j → ∞, we can get the estimate (43). Moreover, according to the Cauchy criteria, the sequence of functions {u m (·, z, λ)}, defined by the iterative formula (30), is uniformly convergent in the domain G u to its limit function u ∞ (·, z, λ).
Since all functions of sequence (43) satisfy the two-point parametrized boundary conditions (11), the limit function u ∞ (·, z, λ) also satisfies them.
Analogically to Theorem 1 in [6], by letting m in relation (30) tend to ∞, it is easy to show that the limit function (37)  Similarly to the Theorem 1 conditions, one can prove convergence of the sequence of functions v m (·, λ, η), i.e., the theorem holds. (12) and (13), there exist a non-negative vector ρ v satisfying an inequality (32), such that f ∈ Lip(K v , D v ), ∀t ∈ a+b 2 , b . In addition assume, that matrix:
All functions of the sequence (34) are continuous on the interval a+b 2 , b and satisfy the separated boundary conditions (13).

2.
The sequence of functions (34) converges uniformly as m → ∞ to its limit function: for all t ∈ a+b 2 , b .

The limit function (49) satisfies the boundary conditions
and is a unique solution of an integral equation: in domain D v . In other words, it is a solution of the corresponding Cauchy problem for a perturbed system of FDEs: where Θ : D a+b 2 × D b → R n is a mapping, given by formula:

4.
The following error estimation holds: where Q v and M b are defined by (47) and (31), respectively.
Proof. The proof is similar to Theorem 1.

Remark 4. Theorems 1 and 2 guarantee that under assumed conditions functions:
are well-defined for all pairs of artificially introduced parameters (z, λ) ∈ ×D a × D a+b 2 and (λ, η) ∈ D a+b 2 × D b . Then, by putting: we can obtain a well-defined continuous function x ∞ (·, z, λ, η), which, at point t = a+b 2 , attains the value:
Proof. The proof can be carried out using a similar approach described in Theorem 2 of [6].
Let us now state the main results of the paper.

Theorem 4.
Assume that the conditions of Theorems 1 and 2 hold. Then: (6) and (7), if and only if the triplet (z, λ, η) satisfies the system of determining equations:

is a continuous solution of the original nonlinear FBVPs
where ∆ and Θ are the mappings defined by formulas (42) and (54), respectively, and Ξ : D u × D v → R n is given by: z, λ, η)).

2.
For every function X(·) of the FBVP (6), (7) with values X(a), X a+b Proof. We refer to the proofs of Theorem 3 (see discussion in [6]) and Theorem 3 (see [12]) and note that the equations (40), (52), (58), (67), and (68) lead straightforward to the continuity of the function x ∞ (·, z, λ, η) at the point t = a+b 2 . Moreover, according to the definition (57) of the aforementioned function, its continuity at all other points of the interval [a, b] holds as well.

Some Remarks
Using our conclusions about function x ∞ (·, z, λ, η), which is given by (57), it is natural that its m-th approximation will be defined as: where the sequences of function u m (·, z, λ), v m (·, λ, η) have the form (30) and (34) accordingly. Moreover, it is more convenient to consider an approximate determining system: Proof. Since the functions u m (·, z, λ) and v m (·, λ, η), defined by the successive approximations (30) and (34), satisfy the consistency condition: Under the relation (70), this proves the continuity of function x m (·, z, λ, η) at the point t = a+b 2 . The fact that this function is also continuous at all the other points follows straightforward from its definition.

Example
In this section, we demonstrate the applicability and improvement of the numericalanalytic technique, previously presented in this paper.
Direct computations show that the vector function f (t, x(t)) in FDS (79) satisfies the Lipschitz condition (5) with the matrix: In addition, the maximum eigenvalue of the matrix: Validity of the last relation means that we cannot apply the classical version of the numerical-analytic method (see discussions in [6,11]), since one of the sufficient conditions in the convergence theorem fails.
However, the dichotomy-type approach, described in Section 3, allows to decompose the original problems (79) and (80) in such a way that the value r(Q) in the inequality (82) will be reduced.
Following our steps in Section 3 of the paper, we first introduced the following parameters: where z = col(z 1 , z 2 ), λ = col(λ 1 , λ 2 ), η = col(η 1 , η 2 ). Then, in alignment with (10)-(13), we decompose the original BVPs (79) and (80) onto two "model"-type parametrized FBVPs: and:   Let us now choose the sets D a and D a+b Thus, we can determine the convex sets D a, a+b 2 and D a+b 2 ,b , originally given by formulae (14) and (15), as: Analogically, we can conclude that the vector function f (t, v(t)) of FDS (85) satisfies the Lipschitz condition (5) with the matrix: and the maximum eigenvalue of the matrix: We want to recall that initially, the maximal eigenvalue of the matrix Q, correspondent to the original problems (79) and (80), failed to be less than 1. This issue was resolved in FBVPs (84)-(87).
Vector functions f (·, u(·)) and f (·, v(·)) in FDSs (84) and (86) are such that: Moreover, we define the ρ u -and ρ v -neighborhoods of the sets D u and D v as: Thus, all conditions of the convergence Theorems 1 and 2 hold, and we are able to proceed with construction of the approximate solution of the parametrized FBVPs (84)-(85) and (86)-(87).
The graphs of the first approximations (93) and (94) and the exact solution (81) are given in Figures 1 and 2.

Discussion
We want to highlight that the obtained results open new possibilities for future developments in the field of differential equations of an arbitrary order and their applications. In particular, one may study differential systems of a mixed order, subjected to multipoint or integral boundary constraints, restrictions containing values of the fractional derivative of solution, etc. This would essentially complement the already existing results in the study of nonlinear BVPs for ordinary differential equations of the real order.

Conclusions
This paper disclosed the recent results in the study of a system of nonlinear FDEs of the real order, subjected to essentially nonlinear two-point boundary constraints. For the analytical representation of a solution, we suggested a modified successive approximation technique (see earlier results in [6][7][8][9][10]), based on the so-called dychotomy-type approach ( [13][14][15][16]). The modification aimed to reduce the a priori error of the method for its more efficient application to the nonlinear problems discussed herein.

Acknowledgments:
The author is grateful for the reviewers' comments, which helped to improve the paper.

Conflicts of Interest:
The author declares no conflicts of interest.

Abbreviations
The following abbreviations are used in this manuscript: FBVP Fractional boundary value problem FDE Fractional differential equation BVP Boundary value problem FIVP Fractional initial value problem