We analyze the continuous time evolution of a $d$-dimensional system of $N$
self propelled particles with a kinematic constraint on the velocities
inspired by the original Vicsek's one [29]. Interactions among
particles are specified by a pairwise potential in such a way that the
velocity of any given particle is updated to the weighted average velocity of
all those particles interacting with it. The weights are given in terms of the
interaction rate function. The interaction is not of mean field type and the
system is non-Hamiltonian. When the size of the system is fixed, we show the
existence of an invariant manifold in the phase space and prove its
exponential asymptotic stability. In the kinetic limit we show that the
particle density satisfies a nonlinear kinetic equation of Vlasov type, under suitable conditions on
the interaction. We study the qualitative behaviour of the solution and we
show that the Boltzmann-Vlasov entropy is strictly decreasing in time.
1.
Introduction
An important tool for describing the natural process is a system of fractional differential equations [1,2]. An ordinary system of fractional differential equations (FDEs) is described by
where $ a\in \mathbb{R}, $ $ \bf{q}:[0, a] \rightarrow \mathbb{R}^\nu $ is a $ \nu $-dimensional vector function and $ \bf{p}:[0, a]\times \mathbb{R}^\nu\rightarrow \mathbb{R}^\nu. $ Here, $ {^C_0D_t^{\alpha}} $ is the Caputo fractional derivative of order $ \alpha > 0. $
FDEs of the form (1.1) are used in an extensive amount of published papers for describing diverse evolutionary processes. Pollution model for lake [3], dynamical models of happiness [4] and Irving-Mullineux oscillaton [5] are among them.
A system of type (1.1) has at least $ n = \lceil \alpha\rceil $ (ceil of $ \alpha $) degree of freedom due to the presence of $ n $ the integer-order derivative. To have a unique solution more information of the system's states is required. In most applied problems this information on the boundary is known. These problems are referred to as boundary value problems. Initial value problems, terminal value problems (TVPs) and Sturm-Liouville problems are among boundary value problems. For initial value problems, the states of the system in the beginning time of the model are known. These states may include derivatives or more complex operators of the modeled process in the initial time. For initial value problems, the states of the system in the beginning time/point of the model are known. These states may include derivatives or more complex operators of the modeled process in the initial time. For terminal value problems, the system's states are known at the end terminal. However, for Sturm-Liouville problems states of the system are known in both initial and end boundaries. These systems have completely different dynamics and behaviors [6].
Adding boundaries information is not the only way to reduce the degree of freedom. For example, the information in the intermediate point/time can be used [7] for such a reduction. Such problems are known as intermediate value problems.
Recently, TVPs for fractional systems of the order less than one ($ 0 < \alpha < 1 $) have received some extensive attention and interest [8,9,10,11,12,13,14]. These problems can be described by
where $ D^{\alpha} $ is fractional derivative and $ \bf{q}_a $ is a given $ \nu $-dimensional vector. At the first glance, these systems seem to be well-posed with regular source functions [15,16]. Surprisingly, Cong and Tuan [17] revealed some counterexamples. The papers [14,18] pointed out that the well-posedness depends on the terminal value $ a. $ For larger $ a $ we may not obtain a unique solution.
Surprisingly, TVPs for fractional differential equations of higher-order derivatives $ \alpha\geq 1 $ are not well-studied. An order $ \alpha\in(1, 2] $ problem on an infinite interval have been studied in [19]. These problems have the form
where $ \bf{q}'(\infty)\in\mathbb{R} $ is a given number and $ {^C_{0}D_t^{\alpha}} $ stands for Caputo derivative.
This paper's important aim is to study the existence results and regularity of the higher-order terminal value problem for systems of FDEs. The other important aim of this paper is to introduce an analyzed high-order numerical method for solving such problems. Thus, the piecewise polynomial collocation method (PPCM) as a numerical solver is introduced and analyzed in detail.
The PPCMs have some superior advantages in solving differential or integral equations. For example, if we use polynomials collocation methods (known as spectral methods) we may encounter Runge-phenomena. The piecewise characteristic of the PPCMs avoids such probable divergence. In comparison with Runge-Kutta methods, the coefficient and parameters of the given methods are obtained constructively. It is remarkable to mention that for ordinary differential equations PPCMs, for some piecewise polynomial spaces are equivalent to Runge-Kutta methods. Conclusively, having more parameters (collocation parameters, degree of polynomial space, step size, meshing method) for controlling error, complexity and order of convergence makes the PPCMs competitively superior category of methods.
This paper is organized as follows: In Sections 2 and 3 we obtain the existence and uniqueness of a mild solution for higher-order FDEs. The regularity for the linear case is studied in Section 4. A numerical method is introduced in Section 5 and analysand in Section 6. Supportive numerical experiments are given in Section 7.
2.
Higher order terminal value problem
Consider the system (1.1) with terminal values
where $ n-1 < \alpha < n $ and $ n\in\mathbb{N}. $ The $ \nu $-dimensional vector function $ \bf{q} = [q_1, \ldots, q_\nu] $ is an unknown vector and $ \bf{q}^{(i)}_{a}\in\mathbb{R}^{\nu} $ ($ i = 0, \ldots, n-1 $) are known vectors.
Further, we impose the Lipschitz condition on components of the function $ \bf{p}:[0, a]\times \mathbb{R}^{\nu}\rightarrow \mathbb{R}^{\nu}. $ Thus, we suppose
where $ L_i $ ($ i = 0, \ldots, \nu $) are constants and are not depend on $ t. $ Let $ \bf{z} = \bf{q}^{(n-1)}, $ ($ \bf{z}' = \bf{q}^{(n)} $). Taking repeatedly integrals and using terminal values, we obtain
By definition of Caputo derivative, we have
Acting Riemann-Liouville integral $ {_0J^{\alpha-n+1}_t} $ on both sides of Eq (1.1), we obtain
Finally, putting $ t = a $ and using terminal conditions, we obtain
Remark 1. The solution of the coupled system of (2.2) and (2.4) can be regarded as a mild solution of the system (2.1).
3.
Existence of a unique solution
Since we use vector functions, we need a norm combined with vector norm and function norm. This combinations brings us to face complexity in calculating induced norm. In this respect, the max norm defined by
seems to be more easier. We establish well-posedness of the inverse problems (1.1) and (2.1). We define the operators $ \mathcal{P} $ and $ \mathcal{Q}:(\mathcal{C}[0, a])^{2\nu} \rightarrow (\mathcal{C}[0, a])^\nu $ by the right hand sides of systems (2.2) and (2.4), i. e.,
and
Let $ \mathcal{P}_i $ and $ \mathcal{Q}_i:(\mathcal{C}[0, a])^{2\nu} \rightarrow \mathcal{C}[0, a] $ be the $ i $th component of the operator $ \mathcal{P} $ and $ \mathcal{Q}, $ respectively. Setting $ \bf{w} = [\bf{q}, \bf{z}] $ and defining the operator $ \mathcal{T}:\mathbb{R}^{2\times \nu}\rightarrow \mathbb{R}^{2\times \nu} $ by
the Eqs (2.2) and (2.4) can be compactly written as
It is straightforward to show that for $ \bf{w}_1 $ and $ \bf{w}_2\in (\mathcal{C}[0, a])^{2\nu}. $ We can write
and
where $ L_M = \max_{i = 1, \ldots, \nu} L_i. $ Now, the Eqs (3.4) and (3.5) give
Conclusively, we can estate the following theorem.
Theorem 3.1. Let each component of the continuous function $ \bf{p}:[0, a]\times \mathbb{R}^{\nu}\rightarrow \mathbb{R}^{\nu} $ ($ \nu\in \mathbb{N} $) satisfy the Lipschitz condition with the Lipschitz constant $ L_i $ and let
Then, the problem (1.1) with terminal conditions (2.1) has a unique mild solution on $ (\mathcal{C}[0, a])^{2\nu} $ if
Proof. Suppose $ \bf{w}\in (\mathcal{C}[0, a])^{2\nu}. $ Since all integral operators involving in the definition of $ \mathcal{T} $ transform a continuous functions into a continuous functions, therefore, $ \mathcal{T}(\bf{w})\in (\mathcal{C}[0, a])^{2\nu}. $ The operator $ \mathcal{T} $ is contractor by (3.6) and the space $ (\mathcal{C}[0, a])^{2\nu} $ is a Banach space by the norm $ \|.\|_{\infty}. $ Thus, by contraction mapping theorem $ \bf{w} = \mathcal{T}(\bf{w}) $ has a unique solution which completes the proof.
The next corollary is an important result of Theorem 3.1 which shows the regularity of the solution on more smooth spaces.
Corollary 1. Suppose the hypotheses of Theorem 3.1 are satisfied. Then, the mild solution $ \bf{q} $ obtained from solving systems (2.2) and (2.4) satisfies the system (1.1) with terminal conditions (2.1). Moreover,
Proof. Taking derivative from Eq (2.2) and using Theorem 3.1 results $ \bf{q}^{(n-1)} = \bf{z}\in (\mathcal{C}[0, a])^{\nu}. $ Thus, $ \bf{q}\in (\mathcal{C}^{n-1}[0, a])^{\nu}. $ Similarly, putting $ t = a $ in Eq (2.2) and its derivatives, proves that $ \bf{q} $ satisfies terminal conditions. Replacing $ \bf{z} = \bf{q}^{(n-1)} $ in (2.2) we obtain
where $ {_a \mathcal{J}_t^{\alpha-n+1}} $ is the Riemann-Liouville fractional integral. Noting that $ \alpha-n+1 > 0 $, we can take a fractional derivative of order $ \alpha-n+1 $ from (3.8) to obtain
which shows that the component $ \bf{q} $ of the mild solutions satisfies (1.1).
4.
Regularity for linear FDEs
The regularity of the solution is important for analyzing numerical methods, especially finite difference methods and methods based on projection [20]. The regularity speaks about the order of differentiability of solutions. Usually, the regularity is investigated in the space $ \mathcal{C}^m[a, b], $ (see Corollary 1). But, some good functions such as $ \sqrt{t}\in \mathcal{C}[a, b] $ do not belong to $ \mathcal{C}^m[a, b]. $ The space $ \mathcal{C}^m(a, b]\cap \mathcal{C}[a, b] $ contains such functions. This space is not Banach or complete. Therefore, we need to introduce another Banach space such that
In this paper, we study regularity on the following weighted space $ {\mathcal{C}}^{m, \alpha}(0, a] $ introduced in [13].
Definition 4.1. [13] Let $ 0\leq \alpha < 1 $ and $ m\in \mathbb{N}. $ Then $ q\in {\mathcal{C}}^{m, \alpha}(0, a] $ if there exist functions $ q_i\in {\mathcal{C}}[0, a] $ for $ i = 0, \cdots, m+1 $ such that $ q = t^{\alpha} q_0+q_1 $ and
Theorem 4.2. The space $ {\mathcal{C}}^{m, \alpha}(0, a] $ with the norm
is a Banach space.
Proof. Let $ {q_{n}}\in {\mathcal{C}}^{m, \alpha}(0, a], $ $ n\in\mathcal{N} $ be a Cauchy sequence with norm $ \|.\|_{\alpha, m}. $ Thus, by the definition of $ {\mathcal{C}}^{m, \alpha}(0, a] $ there exists functions $ q_{n, i} $ in $ \mathcal{C}[0, a] $ for $ i = 0, \cdots, m+1 $ such that $ q_{n}(t) = t^{\alpha} q_{n, 0}(t)+q_{n, 1}(t) $ and
From the definition of the norm $ \|.\|_{\alpha, m}, $ the sequences $ q_{n, i} $ for a fixed $ i $ are Cauchy sequence in the Banach space $ \mathcal{C}[0, a] $ and thus have a unique limit say $ q_{i}. $ Let $ q $ defined by $ q(t) = t^{\alpha} q_0(t)+q_1(t) $ on $ [0, a]. $ Following items show that $ q\in {\mathcal{C}}^{m, \alpha}(0, a]. $
● $ \lim_{n\rightarrow \infty} q_n(t) = \lim_{n\rightarrow \infty} t^{\alpha} q_{n, 0}(t)+q_{n, 1}(t) = t^{\alpha} q_0(t)+q_1(t) = q $ for all $ t\in [a, b]; $
● $ q\in \mathcal{C}[0, a], $ since $ q_0\in \mathcal{C}[0, a] $ and $ q_1\in \mathcal{C}[0, a]; $
● $ \lim_{n\rightarrow \infty} D^{i}q_n(t) = \lim_{n\rightarrow \infty} t^{\alpha-i}q_{n, i+1}(t) = t^{\alpha-i}q_{i+1}(t). $ For each $ t\in(a, b] $ we can find $ \epsilon > 0 $ such that $ D_{\epsilon} = [t-\epsilon, t]\subset (a, b]. $ It is trivial that $ D^{i}q_n(t) $ converges uniformly to $ t^{\alpha-i}q_{i+1}(t) $ on $ D_{\epsilon} $. According to Theorem 7.17 of [21] $ D^{i}q(t) = \lim_{n\rightarrow \infty} D^{i}q_n(t) $ on $ D_{\epsilon}. $ Thus $ D^{i}q(t) = t^{\alpha-i}q_{i+1}(t) $ for all $ t\in D_\epsilon $ and hence for all $ t\in (a, b]. $
Let $ \nu\in \mathbb{N}. $ For a $ \nu $ dimensional vector functions, $ \bf{p} = [f_1, \ldots, f_{\nu}] $ in $ ({\mathcal{C}}^{m, \alpha}(0, a]\subset ({\mathcal{C}}^{m}(0, a])^{r}, $ we use the max norm defined by
The system (1.1) with
is a system of linear FDEs. Here $ A $ is a given $ \nu\times \nu $ dimensional matrix function and $ \bf{b}\in (\mathcal{C}[0, a])^\nu $ is a given source function. To study linear systems we introduce the operators $ \mathcal{W}_1 $ and $ \mathcal{W}_2:(\mathcal{C}[0, a])^{2\nu} \rightarrow (\mathcal{C}[0, a])^\nu $ by
and
Let us also define the vector valued functions $ \mathcal{G}_2 $ and $ \mathcal{G}_1:[0, a] \rightarrow \mathbb{R}^\nu $ by
and
Let $ \bf{T} = [\mathcal{W}_1, \mathcal{W}_2]^T $ and $ \mathcal{G} = [\mathcal{G}_1, \mathcal{G}_2]^T $. Then, the Eqs (2.2) and (2.4) can be written as an inhomogeneous system
where $ \bf{w} = [\bf{q}, \bf{z}]^T. $
Theorem 4.3. Let $ a\in\mathbb{R} $ and $ \alpha > 0. $ Assume each component of $ A $ and $ b $ are in $ {\mathcal{C}}^{m, \alpha}(0, a]) $ and the hypotheses of Theorem 3.1 are satisfied for related $ \bf{p}:[0, a]\times \mathbb{R}^{\nu}\rightarrow \mathbb{R}^{\nu} $ ($ \nu\in \mathbb{N} $).Then, the problem (1.1) with terminal conditions (2.1) has a unique mild solution on $ ({\mathcal{C}}^{m, \alpha}(0, a])^{2\nu}. $
Proof. Since $ \bf{T} $ is a combination of weakly singular integral operators, it is a compact linear operator on $ (\mathcal{C}[0, a])^{2\nu} $ and $ ({\mathcal{C}}^{m, \alpha}(0, a])^{2\nu} $. Also, it is clear that
System (25) has a unique solution on $ (\mathcal{C}[0, a])^{2\nu} $ by Theorem 3.1. Therefore, the homogeneous system $ \bf{w}-\bf{T}(\bf{w}) = 0 $ has the trivial null space in $ (\mathcal{C}[0, a])^{2\nu} $ and thus $ ({\mathcal{C}}^{m, \alpha}(0, a])^{2\nu}. $ Conclusively, alternative Fredholm theorem asserts that the system (4.7) has a unique solution on $ ({\mathcal{C}}^{m, \alpha}(0, a])^{2\nu}. $
Remark 2. Regularity of the nonlinear case needs further investigations.
Example 1. Consider the system (1.1) with $ a = 0.5, $ $ \alpha = 2.5 $ and $ \bf{p} $ defined by (4.2)
Obviously, $ L_1 = 0.5, $ $ L_2 = 0.5. $ and thus $ L_M = 0.5. $ Therefore,
and Condition (3.7) holds. By Theorem 3.1 the terminal value problem (1.1) with arbitrary terminal value (2.1) has a unique solution on $ \mathcal{C}[0, a] $. Since all components of $ A $ and $ \bf{b} $ belong to $ {\mathcal{C}}^{m, \alpha}(0, a], $ ($ m\in\mathbb{N} $) the mild solution belongs to $ ({\mathcal{C}}^{m, \alpha}(0, a])^4, $ by Theorem 4.3. Furthermore, we obtain the regularity of the solution of the system (1.1) on the closed interval $ [0, a] $ by Corollary 1 and we have
5.
Numerical method
Considering the memory, one of the best methods for solving the coupled systems (2.2)–(2.4) is to use collocation methods on piecewise polynomial spaces. To implement such methods we partition the solution interval $ [0, a] $ into sub-intervals $ \sigma_i = [t_{i}, t_{i+1}], $ $ i = 0, \ldots, N-1, $ with length $ h_i: = t_{i+1}-t_i $ where $ 0 = t_0 < \ldots < t_{N} = a $ are nodes of a chosen mesh (uniform or graded mesh) and $ N\in \mathbb{N}. $ A graded mesh with exponent $ r\geq1 $ is described by
Let $ 0 < c_1 < \ldots < c_m\leq 1 $ ($ m\in\mathbb{N} $) be collocation parameters, $ t_{n, i} = t_n+c_ih_n $ ($ n = 0, \ldots, N-1 $) be collocation points and let $ \hat{\bf{q}}_N(t) $ and $ \hat{\bf{z}}_N(t) $ be approximate solutions. The restriction of approximate functions to the $ \sigma_k $ is fully determined by Lagrange polynomials interpolation formula
and
where $ Q_{k, j} = \hat{\bf{q}}_N(t_{n, j}), $ $ Z_{k, j} = \hat{\bf{z}}_N(t_{k, j}), $ $ t_{k, j} = t_k+hc_j $ and $ L_j $ ($ j = 1, \ldots, m $) are Lagrange polynomials of degree $ m-1. $ The integrals in the operators $ \mathcal{P} $ and $ \mathcal{Q} $ of systems (3.1) and (3.2) can be discretized as:
where
For $ t\in [t_k, t_{k+1}] $ ($ k = 0, \ldots, N-1 $) we have
where
and
Similarly, we discretize integrals of the operator $ \mathcal{Q}. $ We have
By interpolating $ \bf{p}(t_l+sh_l, \bf{q}(t_l+sh_l)) $ on $ c_i $ we obtain
and thus
where
and finally
where
The discretized mappings $ \mathcal{P}_N $ and $ \mathcal{Q}_N:(\mathcal{C}[0, a])^{2\nu} \rightarrow (PC[0, a])^\nu $ ($ PC $ stands for piecewise continuous space) related to systems (3.1) and (3.2) can be defined by
and
on $ \sigma_k. $ Setting $ \mathcal{T}_N = [\mathcal{P}_N, \mathcal{Q}_N], $ unknowns vectors $ Q_{l, j} $ and $ Z_{l, j} $ can be obtained by solving the nonlinear system
for $ k = 0, \ldots, N-1 $ and $ o = 1, \ldots, m. $ Taking into account that
the dense solution can be evaluated in all points of the desired interval by Eqs (5.1) and (5.2). We note that the six parameters $ \gamma_{n, i, l, j}, $ $ \eta_{n, l, j}, $ $ \bar{\eta}_{n, k, j}, $ $ \varrho_{n, l, j}, $ $ \zeta_{n, l, j} $ and $ \bar{\zeta}_{n, k, j} $ fully determines the method for each $ m $. The equations of these six parameters for $ m = 1 $ and $ m = 2 $ are provided in Tables 1 and 2, respectively.
6.
Numerical analysis
For simplifying our analysis we recall some notations. Let $ S_{m-1}^{-1}(I_h) $ be the space of piecewise polynomials of degree less than $ m, $ ($ m\in \mathbb{N} $) on the mesh partitioning $ I_h = \{t_i:\; i = 0, \cdots, N\}. $ The projection operator $ \Pi_{m-1, N}:({\mathcal{C}}[0, a])^{2\nu}\rightarrow \left(S_{m-1}^{-1}(I_h)\right)^{2\nu} $ is uniquely determined by interpolation on collocation points such that
Since $ \hat{\bf{w}}_N(t): = [\hat{\bf{q}}_N(t), \hat{\bf{z}}_N(t)]\in \left(S_{m-1}^{-1}(I_h)\right)^{2\nu} $ we have
Thus, Eq (5.10) can be written as
A useful theorem regarding $ \Pi_{m-1, N} $ that simplifies existence and convergence analysis is the following theorem.
Theorem 6.1. Let $ m, N\in\mathbb{N}. $ Then, $ \Pi_{m-1, N} $ is a bounded linear operator. Let
be induced norm of $ \Pi_{m-1, N} $ and
Then, $ \Lambda\leq \Lambda_1 $ and
Proof. In each sub-interval $ \sigma_n $ of $ [0, a] $ by Lagrange interpolation formula, we have
The rest of the proof is straightforward by taking the max norm.
6.1. Existence results for the proposed numerical method
One of the most fundamental questions is whether the system (6.1) has a unique solution? The answer is affirmative. For $ \hat{\bf{w}}_N(t) $ in $ \left(S_{m-1}^{-1}(I_h)\right)^{2\nu} $ the operators $ \mathcal{P}_N $ of (5.8) and $ \mathcal{P} $ of (3.1) are equivalent and we have
Thus, the inequality (3.4) holds for $ \mathcal{P}_N $ and
for all $ \bf{w}_1, \bf{w}_2\in\left(S_{m-1}^{-1}(I_h)\right)^{2\nu}. $ However, $ \mathcal{Q}_N $ is different from $ \mathcal{Q} $ in this space and we need further computing. Actually, $ \mathcal{Q}_N $ can be described by
where $ \bf{q} $ is the first $ \nu $ element of $ \bf{w}. $ More precisely
Therefore,
Theorem 6.2. Let
Then, the numerical solution of system (1.1) obtained by (6.1) exists and is unique.
Proof. The operator $ \Pi_{m-1, N} \mathcal{T}_N:\left(S_{m-1}^{-1}(I_h)\right)^{2\nu}\rightarrow \left(S_{m-1}^{-1}(I_h)\right)^{2\nu} $ is a contractor by (6.7). Thus, by using contraction mapping theorem, $ \Pi_{m-1, N} \mathcal{T}_N $ admits a unique fixed-point.
Remark 3. Obviously, if
holds, then Eq (6.8) holds too. For $ m = 1 $ we have $ \Lambda_1 = 1, $ and Eq (3) matches with Eq (3.7). However, our estimate may not be optimal for higher degrees of approximations. Also, we guess using the convergence properties of $ \Pi_{m-1, N} $ when $ N\rightarrow \infty, $ we may obtain a convergence result without dependency on $ \Lambda_1. $
6.2. Detailed implementation: Solving nonlinear system
An important question is how to solve the nonlinear Eq (6.1). There are many methods available in literature that we can employ. The Newton iteration method is one of them [14]. The advantage is that it is fast. The disadvantage is that it needs computation of Jacoby as well as a good initial value for an iteration. It increases competition cost for higher dimension (computation of $ \nu^2 $ function for each iteration) and complexity of the implementation. However, strict restriction of $ \lambda_\Lambda $ and analytic discussion of the previous section is constructive and suggests the use of an iterative method. Since the operator $ \Pi_{m-1, N} \mathcal{T}_N $ is a contractor operator, beginning with an arbitrary initial value grantees the convergence of the related iterative method. In our algorithm, we initialize the iteration by
Then, we estimate the solution by the recursive formula
The computation of iteration (6.10) continues with the smallest $ \kappa $ such that
where $ Tol $ is a desired user tolerance. We report such $ \kappa $ as a number of iteration to achieve the given tolerance.
6.3. Convergence analysis
Theorem 6.3. Suppose $ \hat{\bf{w}}_N = [\hat{\bf{w}}_N, \hat{\bf{w}}_N] $ ($ N\in\mathbb{N} $) be a discretized collocation solution of (2.2) and (2.4) described by (6.1). Let the conditions of Theorems 3.1 and 6.2 be fulfilled. Then, the coupled systems (2.2) and (2.4) have a unique solution $ \bf{w} $ on $ (\mathcal{C}[0, a])^{2\nu} $ and the $ \hat{\bf{w}}_N $ converges to the $ \bf{w}. $ Furthermore, supposing $ \bf{w}\in ({\mathcal{C}}^{m, \alpha}(0, a])^{2\nu} $ implies
Proof. Let $ \bf{w} $ be the solution of the coupled systems (2.2) and (2.4). Define the residue operator by
and let $ \bf{e} = \hat{\bf{w}}_N-\bf{w}. $ Subtracting Eq (3.3) from (6.1), we obtain
Taking maximum norm, we obtain
Applying Eq (6.8), we obtain
and finally
This is a good news since we can obtain the asymptotic behavior of $ \|R(\bf{w})\|_{\infty} $ for the known $ \bf{w}. $ Thanks to Theorem 8 of [20], one can easily prove that
which completes the proof.
7.
Examples
All numerical experiments is implemented in MATLAB software. The tolerance number for solving related nonlinear equations is $ Tol = 10^{-14}. $ This tolerance is close to the floating-point relative accuracy $ 2.2204e-16 $ for data type of double-precision (64 bit) in standard machine. The numerical examples show theoretically the obtained order of convergence is sharp and can not be improved more.
Example 2. Let $ 1 < \alpha\leq 2 $ and
on $ [0, 0.15]. $ We note that $ \bf{p}(t, .)\in (\mathcal{C}[0, 0.25])^{2} $ but $ \bf{p}(t, .)\notin (\mathcal{C}^1[0, 0.25])^{2}. $ Since $ L_M = 0.1, $ one can check that the conditions of Theorem 3.1 for $ a = 0.25 $ hold and the system (2.1) with boundary condition
and
has a unique solution on $ [0, a]. $ The solution by Corollary 1 belongs to $ (\mathcal{C}^1[0, a])^2. $ We apply a numerical method based on collocation parameter $ c = [0.5, 1]. $ Thus, the dense solution on each sub-interval $ \sigma_k $ ($ k = 0, \cdots, N-1 $) of a given partition can be described by
An estimated norm of the given method can be determined by (6.2) and we have
Taking into account (6.8), we obtain
Consequently, by using Theorem 6.2 the nonlinear Eq (6.1) has a unique solution and the iteration process (6.10) converges to that solution. Theorem 6.3 implies the numerical solution converges to the exact solution by following order ($ \mathrm{O} $)
In Tables 3 and 4, we provided the maximum error
estimated order of proposed method ($ O_N $) and the number of iteration for required tolerance. As we see $ O_N = 0.5, 1, 1.5, 2 $ for $ r = 1, 2, 3, 4 $ respectively, which is in agreement with theoretical result (see Eq (7.2)).
Remark 4. Theorem 1 provides a sharp result when we consider all the components of the state of the system 1 ($ \bf{w} = [\bf{q}, \bf{z}] $). However, it can be improved for the desired state of the solution $ \bf{q}. $ According to Corollary 1, $ \bf{q} $ has better regularity for $ n\geq 1. $ Thus, we guess that the order of convergence for $ \bf{q} $ component with uniform mesh $ r = 1, $ can be greater than or equal to $ n-1. $ Thus, for that component, a greater graded mesh exponent is not necessary. The error $ E_{\bf{q}}(r, N) = \|\hat{\bf{q}}_N-\bf{q}\|_{\infty} $ and the estimated order of convergence of Example 2 is reported in Table 5. The order $ 2 $ is achieved with $ r = 1.5 $.
To show the sharpness of theoretical results, we add numerical experiments for the case $ m = 1 $ with $ c = 0.5 $ (for other choice of $ c $ on $ [0, 1] $ we have similar results). The dense solution on each sub-interval $ \sigma_k $ can be described by
In this case $ \Lambda_1 = 1 $ and by (6.8)
Consequently, by using Theorem 6.2 the nonlinear Eq (6.1) related to $ m = 1 $ has a unique solution and the iteration process (6.10) converges to the exact solution. By Theorem 6.3 we expect
The sharpness of this result can be seen in Table 6.
Remark 5. The low order method have the advantage of supporting larger class of terminal value problem with large value of $ a $. The reason is that with $ m = 1, $ we have $ \Lambda_1 = 1 $ which is three times smaller than $ \Lambda_1 $ for the case $ m = 2. $
8.
Conclusions
TVPs for higher-order (greater than one) fractional differential equations are rarely studied in the literature. In this paper, we tried to have a comprehensive study with a simple analysis to cover all general interests. Many questions in this topic deserve to study with more scrutiny. The regularity of nonlinear cases as well as optimal bound for obtaining well-posed problems are among them. We think the terminal value problem is important in applied since it monitors the past evolution of a dynamical system. Also, it can be regarded as a control problem for having the desired future by finding suitable initial value. We think this topic will catch more interest similar to the initial value problem.
Conflict of interest
The authors declare that there are no conflicts of interest.
DownLoad:
