# The Iterative Scheme and the Convergence Analysis of Unique Solution for a Singular Fractional Differential Equation from the Eco-Economic Complex System’s Co-Evolution Process

Authors: Teng Ren; Helu Xiao; Zhongbao Zhou; Xinguang Zhang; Lining Xing; Zhongwei Wang; Yujun Cui
Journal: Complexity
 (2019)

---

## Abstract

In this paper, we focus on a class of singular fractional differential equation, which arises from many complex processes such as the phenomenon and diffusion interaction of the ecological-economic-social complex system. By means of the iterative technique, the uniqueness and nonexistence results of positive solutions are established under the condition concerning the spectral radius of the relevant linear operator. In addition, the iterative scheme that converges to the unique solution is constructed without request of any monotonicity, and the convergence analysis and error estimate of unique solution are obtained. The numerical example and simulation are also given to demonstrate the application of the main results and the effectiveness of iterative process.


---

## Body

## 1. Introduction

This paper is inspired by a singular fractional differential equation with Riemann-Stieltjes integral condition(1)-Dt11/4xt=13t-1/21+sin⁡D1/4xt,t∈0,1,Dt1/4x0=Dt5/4x0=0,Dt3/2x1=∫01Dt1/4xsdAs,which arises from some complex system of economic and engineering science, where A is a bounded variation function satisfying A(t)=0,t∈[0,1/3), A(t)=3,t∈[1/3,2/3), A(t)=2,t∈[2/3,1]. For example, usually phenomenologic viscoelastic models are based on springs and dashpots which obey Hooke’s law; however the problem is that the corresponding ordinary differential equations have a relatively restricted class of solutions, which is, in general, too limited to provide an adequate description for the complex systems discussed; see [1]. To overcome this shortcoming, by defining the stress y(t) and the strain x(t), Schiessel et al. [1] introduced a fractional order viscoelasticity Kelvin-Voigt system whose stress decays after a shear jump in an algebraic manner (2)yt=Eταdαdtαxt+Eτβdβdtβxt,where α>β>0, E is a constant and dα/dtα is Riemann-Liouville derivative with an order of α. This implies that the fractional order Kelvin-Voigt model can be generalized by the following mathematical model: (3)dαdtαxt=ft,dβdtβxt,with α>β>0. Thus the equation (1) belongs to a particular fractional order viscoelasticity Kelvin-Voigt system with the stress (1/6)t-1/2, which is important to study and understand dynamic behaviour for the corresponding viscoelasticity process, also see [2]. In addition, here we also notice that the nonlinearity of (1) exhibits a blow-up behaviour at singular time variable t=0. In particular, the singular problems [3–17] as well as impulsive phenomena [18–35] often lead to some blow-up behaviour [28, 36–40] in various complex process of economic and engineering science such as in; ecological-economic-social complex system (Eco-economic System); since the energy consumption fiercely increases resulting in a rapid decrease in total energy, blow-up phenomena will happen at certain time point [41]. In mechanics process, a blow-up behaviour also occurs near the crack tip in elastic fracture like r-1/2, where r is the distance from the crack tip [42].Motivated by the above problems, in this paper, we consider the iterative scheme and the convergence analysis of unique solution for the following singular fractional viscoelasticity complex system with Riemann-Stieltjes integral condition:(4)-Dtαxt=ft,Dtγxt,t∈0,1,Dtγx0=Dtγ+1x0=0,Dtβx1=∫01DtγxsdAs,where Dtα,Dtβ,Dtγ are the standard Riemann-Liouville derivatives with 0<γ<1<β<2<α≤3 and α-γ>2,β-γ>1, ∫01x(s)dA(s) is denoted by a Riemann-Stieltjes integral, and A is a bounded variation function with a sign-changing measure dA; the nonlinearity f(t,·) may be singular at t=0,1.Mathematical models involving fractional derivatives can describe many advection-dispersion processes [43–45], viscoelasticity characteristics [46–48], thermostat model [49] and the bioprocesses with long memory [50]. Especially when one wants to describe long-term ecological-economic-social complex system phenomena and diffusive interaction, fractional differential operator possesses a higher accuracy than the traditional integer order differential model in depicting the coevolution process of economic, social, and ecological subsystems and the transport of solute in highly heterogeneous porous media [51–53]. For example, Teng et al. [54] considered the maximum and minimum solutions for a fractional order differential system, involving a p-Laplacian operator and nonlocal boundary conditions, which arises from a complex process of ecological economy phenomena and diffusive interaction.Normally, only positive solutions are meaningful for most practical problems; thus nonlinear analysis methods, such as iterative methods [37, 55–65], variational methods [66–73], the fixed point theorems [74–86], and upper and lower solution methods [87–89], play an important role for the study of various differential equtions. Recently, by means of monotone iterative technique, Zhang et al. [46] established the existence and uniqueness of the positive solution for a fractional differential equation with derivatives(5)-Dtαxt=ft,xt,-Dtβxt,t∈0,1,Dtβx0=0,Dtγx1=∑j=1p-2ajDtγxξj,where 1<α≤2,α-β>1, 0<β≤γ<1, 0<ξ1<ξ2<⋯<ξp-2<1, aj∈[0,+∞) with c=∑j=1p-2ajξjα-γ-1<1, and Dtα,Dtβ,Dtγ are the standard Riemann-Liouville derivative. f:(0,1)×[0,+∞)×(-∞,+∞)→[0,+∞) is continuous and increasing on second and third variables. However, in order to apply the iterative technique, most of works require monotonicity conditions, see [37, 55–57, 59–61, 64]. Thus the aim of this paper is to weaken the monotonicity request; that is, we establish the iterative scheme and the convergence analysis of unique solution for singular fractional differential equation (4) without any monotonicity conditions.The present paper has some new features. Firstly, both the nonlinear term and the boundary conditions involve fractional order derivatives of unknown functions. Secondly, the uniqueness and nonexistence results are established under the condition concerning the spectral radius of the relevant linear operator; that is, we do not require any monotonicity conditions for nonlinearity. Thirdly, the iterative scheme is constructed and the convergence analysis of unique solution is carried out. Finally, the nonlocal boundary condition possesses weaker positivity sincedA can be changing-sign measure.The rest of the paper is organized as follows. In Section2, we firstly recall some definitions and basic properties on Riemann-Liouville derivative and integral and then give some properties of the Green function. In Sections 3 and 4, the uniqueness and nonexistence results are established under the condition concerning the spectral radius of the relevant linear operator. In Section 5, numerical example and simulation are given to demonstrate the main results and the effectiveness of iterative process.


## 2. Preliminaries and Lemmas

For further discussion, here we briefly recall some definitions, notations, and known results, which will be found in the recent monographs.The work space of this paper is Banach spaceE=C[0,1] with the norm z=maxt∈0,1zt. Let P={z∈E:z(t)≥0,t∈[0,1]} be a cone in E and construct a subset of P as follows: (6)Q=z∈P:thereexisttwopositivenumberskz≤1≤Kzsuch  thatkztα-γ-1≤zt≤Kztα-γ-1,t∈0,1.Definition 1.
Letα>0 with α∈R; the αth Riemann-Liouville fractional integral for a function f:[a,∞)→R is defined by (7)Iαft=1Γα∫att-sα-1fsdsprovided that the right-hand side is point-wise defined on (a,+∞).Definition 2.
The Riemann-Liouville fractional derivative with orderα>0 for a function is given by (8)Dtαft=1Γn-αddtn∫att-sn-α-1fsds, where n∈N is the unique positive integer satisfying n-1≤α<n and t>a.Lemma 3.
The Riemann-Liouville fractional derivative and integral enjoy the following properties.
(1) Iff,g:(0,+∞)→R with order α>0, then (9)Dtαft+gt=Dtαft+Dtαgt.
(2) Iff∈L1(0,1),α>β>0, then (10)IαIβft=Iα+βft,DtβIαft=Iα-βft,DtαIαft=ft.
(3) Ifβ>0,α>0, then (11)Dtβtα-1=ΓαΓα-βtα-β-1.
(4) Letα>0, then (12)IαDαft=ft+c1tα-1+c2tα-2+⋯+cntα-n,where ci∈R,i=1,2,…,n(n=[α]+1).Now letx(t)=Iγz(t),z(t)∈C[0,1] and consider the following modified equation of (4): (13)-Dtα-γzt=ft,zt,z0=z′0=0,Dtβ-γz1=∫01zsdAs,and we have the following lemma.Lemma 4.
Equation (4) is equivalent to the modified boundary value problem (13). Moreover, if z∈C([0,1],[0,+∞)) is a solution of problem (13), then the function x(t)=Iγz(t) is a positive solution of (4).Proof.
Firstly, it follows from Definitions1 and 2 and Lemma 3 that (14)Dtαxt=d3dt3I3-αxt=d3dt3I3-αIγzt=d3dt3I3-α+γzt=Dtα-γzt,Dtγxt=DtγIγzt=zt,Dtγ+1xt=Dtγ+1Iγzt=z′t. Substituting (14) into (4), then (4) reduces to the following boundary value problem: (15)-Dtα-γzt=ft,zt,z0=z′0=0,Dtβ-γz1=∫01zsdAs,that is the modified boundary value problem (13).
Conversely, using (13) again, the modified boundary value problem (13) is also transformed to the form (4). Thus the problem (13) is indeed equivalent to (4) and if z solves the modified boundary value problem (13), from the monotonicity and property of Iγ that (16)Iγz∈C0,1,0,+∞.thus x(t)=Iγz(t) also solves (4).Let(17)Jt,s=1Γα-γtα-γ-11-sα-β-1-t-sα-γ-1,0≤s≤t≤1,tα-γ-11-sα-β-1,0≤t≤s≤1.Then we have the following lemma.Lemma 5.
Givenh∈L1(0,1), then the following boundary value problem(18)-Dtα-γzt=ht,0<t<1,z0=z′0=0,Dtβ-γz1=∫01zsdAs,has the unique solution zt=∫01Ht,shsds, where(19)Ht,s=tα-γ-1Γα-γ/Γα-β-AWs+Jt,s,(20)A=∫01tα-γ-1dAt,Ws=∫01Jt,sdAt.Proof.
It follows from Lemma3(4) that (18) reduces to the following equivalent integral equation: (21)zt=-Iα-γht+c1tα-γ-1+c2tα-γ-2+c3tα-γ-3,c1,c2,c3∈R.Noticing that z(0)=z′(0)=0 and (21), we have c2=c3=0, which implies that (22)zt=-Iα-γht+c1tα-γ-1.By (22) and Lemma 3 (3), one has (23)Dtβ-γzt=-Dtβ-γIα-γht+c1Dtβ-γtα-γ-1=-Iα-βht+c1Γα-γΓα-βtα-β-1=-∫0tt-sα-β-1Γα-βhsds+c1Γα-γΓα-βtα-β-1,and then (24)Dtβ-γz1=-∫011-sα-β-1Γα-βhsds+c1Γα-γΓα-β.
On the other hand, we have(25)∫01zsdAs=∫01-Iα-γhs+c1sα-γ-1dAs=∫01-∫0ss-τα-γ-1Γα-γhτdτdAs+c1∫01sα-γ-1dAs.It follows from (18) that (26)c1=∫011-sα-β-1/Γα-βhsds-∫01∫0ss-τα-γ-1/Γα-γhτdτdAsΓα-γ/Γα-β-∫01sα-γ-1dAs. Consequently, the unique solution of problem (18) is (27)zt=-∫0tt-sα-γ-1Γα-γhsds+∫011-sα-β-1/Γα-βhsds-∫01∫0ss-τα-γ-1/Γα-γhτdτdAsΓα-γ/Γα-β-∫01sα-γ-1dAstα-γ-1=-∫0tt-sα-γ-1Γα-γhsds+tα-γ-1∫011-sα-β-1Γα-γhsds+∫011-sα-β-1/Γα-βhsds-∫01∫0ss-τα-γ-1/Γα-γhτdτdAsΓα-γ/Γα-β-∫01sα-γ-1dAs-∫011-sα-β-1Γα-γhsdstα-γ-1=∫01Jt,shsds+tα-γ-1Γα-γ/Γα-β-∫01sα-γ-1dAs×-∫01∫0ss-τα-γ-1Γα-γhτdτdAs+∫011-sα-β-1Γα-γhsds∫01sα-γ-1dAs=∫01Jt,shsds+tα-γ-1Γα-γ/Γα-β-∫01sα-γ-1dAs×∫01-∫0ss-τα-γ-1Γα-γhτdτ+sα-γ-1∫011-τα-β-1Γα-γhτdτdAs=∫01Ht,shsds.The proof is completed.Lemma 6.
The functionJ(t,s) and H(t,s) have the following properties provided that 0≤A<Γ(α-γ)/Γ(α-β) and W(s)≥0 for any s∈[0,1].
(1)J(t,s) and H(t,s) are nonnegative and continuous for (t,s)∈[0,1]×[0,1].
(2)J(t,s)≤ω2(s)tα-γ-1,fort,s∈[0,1], where (28)ω2s=1-sα-β-1Γα-γ.
(3)tα-γ-1κ(s)≤J(t,s)≤κ(s),fort,s∈[0,1], where (29)κs=1-sα-β-1-1-sα-γ-1Γα-γ.
(4)(30)tα-γ-1ω1s≤Ht,s≤ω1s,t,s∈0,1, where (31)ω1s=WsΓα-γ/Γα-β-A+κs.Proof.
It follows from (17) and (19) that the properties (1) and (2) hold. Obviously, (4) also holds if (3) is satisfied. So we only need to prove (3). We divide the proof into two cases,
Case i. If 0≤t≤s≤1, we have (32)Γα-γJt,s-tα-γ-1κs=tα-γ-11-sα-β-1-tα-γ-11-sα-β-1-1-sα-γ-1=tα-γ-11-sα-γ-1≥0,Γα-γκs-Jt,s=1-sα-β-1-1-sα-γ-1-tα-γ-11-sα-β-1≥1-sα-β-1-1-sα-γ-1-sα-γ-11-sα-β-1=1-sα-β-11-sα-γ-1-1-sβ-γ.
Since(33)ddssα-γ-1+1-sβ-γ=α-γ-1sα-γ-2+β-γ1-sβ-γ-1≥0,s∈0,1,we have (34)Γα-γκs-Jt,s≥0.
Case ii. If 0≤s≤t≤1, we have (35)Γα-γJt,s-tα-γ-1κs=tα-γ-11-sα-β-1-t-sα-γ-1-tα-γ-11-sα-β-1-1-sα-γ-1=tα-γ-11-sα-γ-1-t-sα-γ-1≥0,.On the other hand, since 0<γ+1<β, we have (36)∂∂ttα-γ-11-sα-β-1-t-sα-γ-1=α-γ-1tα-γ-21-sα-β-1-t-sα-γ-2≥α-γ-1tα-γ-21-sα-β-1-1-sα-γ-2≥α-γ-1tα-γ-21-sα-β-11-1-sβ-γ-1≥0,s∈0,1.Thus (37)Γα-γJt,s=tα-γ-11-sα-β-1-t-sα-γ-1≤1-sα-β-1-1-sα-γ-1=Γα-γκs.This completes the proof of (3).Lemma 7 (Krein-Rutmann).
LetL:E→E be a continuous linear operator, P be a total cone, and L(P)⊂P. If there exist ψ∈E∖(-P) and a positive constant c such that cL(ψ)≥ψ, then the spectral radius r(L)≠0 and has a positive eigenfunction corresponding to its first eigenvalue λ=r(L)-1.Lemma 8 (Gelfand’s formula).
For a bounded linear operatorL and the operator norm ∥·∥, the spectral radius of Tn satisfies (38)rL=limn→+∞Ln1/n.


## 3. Existence Results

To obtain the existence result for (4), we use the following assumptions.(H0)
A is a function of bounded variation such that W(s)≥0 for s∈[0,1] and 0≤A<Γ(α-γ)/Γ(α-β).(H1)
f∈C((0,1)×[0,+∞),[0,+∞)) and there exists a function p∈L1[0,1] such that (39)ft,u1-ft,u2≤ptu1-u2,t∈0,1,ui∈0,+∞,i=1,2.(H2)
​(40)0<∫01ωspsds<+∞,where(41)ωs=ω1s+ω2s.(H3)
f(t,0)≢0,t∈(0,1).Remark 9.
Clearly, if (H2) holds, then we have(42)0<∫01ω1spsds<+∞,0<∫01ω2spsds<+∞.Now let us define a nonlinear operatorT:P→E and a linear operator L:P→E as follows: (43)Tzt=∫01Ht,sfs,zsds,t∈0,1,(44)Lzt=∫01Ht,spszsds,t∈0,1.Clearly, if z solves the operator equation z=Tz, then z is a solution of the boundary value problem (13).Lemma 10.
Assume that(H0) and (H2) hold. Then the linear operator L:P→Q is a completely continuous operator, and the spectral radius r(L)≠0; moreover L has a positive eigenfunction ψ corresponding to its first eigenvalue λ1=(r(L))-1.Proof.
For anyz∈P, it follows from Lemma 6 that(45)Lzt=∫01Ht,spszsds≤tα-γ-1∫01ω2spszsds,(46)Lzt≥tα-γ-1∫01ω1spszsds, which implies that(47)kztα-γ-1≤Lzt≤Kztα-γ-1, where (48)kz=min∫01ω1spszsds,1,Kz=max∫01ω2spszsds,1.Thus we have that L:P→Q. By the uniform continuity of H(t,s),(t,s)∈[0,1]×[0,1], we know that the linear operator L is a completely continuous operator.
Next we show thatL has the first eigenvalue λ1>0. In fact, it follows from Lemma 6 and (H2) that there exists a constant t0∈(0,1) such that H(t0,t0)p(t0)>0, which implies that there exists [a,b]⊂(0,1) such that t0∈(a,b) and H(t,s)p(s)>0 for all t,s∈[a,b]. Take z∈Q such that z(t0)>0 and z(t)=0 for all t∉[a,b]. So for any t∈[a,b], (49)Lzt=∫01Ht,spszsds≥∫abHt,spszsds>0.Thus there exists c>0 such that c(Lz)(t)≥z(t) for t∈[0,1]. By using Lemma 7, the spectral radius r(L) of linear operator L is not zero; moreover L has a positive eigenfunction ψ corresponding to its first eigenvalue λ1=(r(L))-1 such that (50)Lψ=rLψ.The proof is completed.Lemma 11.
If(H0)-(H2) hold, then T:P→Q is completely continuous.Proof.
Similar to Lemma10, we have T(P)⊂Q. By means of the Arzela-Ascoli Theorem, T:P→Q is completely continuous.Theorem 12.
Suppose that(H0)-(H3) hold. If the spectral radius of the linear operator r(L)∈(0,1), then (4) has a unique positive solution x∗, and there exist two constants μ2>μ1>0 such that (51)μ1tα-1≤x∗t≤μ2tα-1. Moreover, for any initial z0∈P∖θ, construct successively a sequence (52)zm=∫01Ht,sfs,zm-1sds,m=1,2,…, and then the iterative sequence zm(t) converges uniformly to Dtγx∗(t) on [0,1] as n→∞, i.e., (53)um-Dtγx∗→0,as m→∞. Furthermore, there exists an error estimation (54)um-Dtγx∗≤Kz1-z0ψ1-rLkψrmL with the rate of convergence(55)um-Dtγx∗=OrmL,where ψ is the positive eigenfunction of the linear operator L.Proof.
Firstly, it follows from Lemma11 that T:P→Q is completely continuous. Since f(t,0)≢0, we know that T does not have zero fixed point. Thus we only need show that T has a unique fixed point in P.
Step 1. We shall prove that T has fixed points in P.
In fact, for anyz∈P∖{θ}, it follows from Lemma 10 that there exist two positive numbers Kz≥1≥rz>0 such that(56)kztα-γ-1≤Lzt≤Kztα-γ-1. On the other hand, by Lemma 10, L has a positive eigenfunction ψ, i.e, (57)Lψ=rLψ.It follows from (56) that(58)r-1Lkψtα-γ-1≤ψt≤r-1LKψtα-γ-1.
Now letz0∈P∖{θ} be given; we construct an iterative sequence (59)zn=Tzn-1,n=1,2,3,⋯. Without loss of generality, suppose that |z1-z0|≠θ (otherwise, the proof is completed), and then it follows from (56) and (58) that(60)Lz1t-z0t≤rLKz1-z0kψψt. So we have(61)z2t-z1t=∫01Ht,sfs,z1s-fs,z0sds≤∫01Ht,spsz1s-z0sds≤Lz1-z0t≤rLKz1-z0kψψt.(62)z3t-z2t=∫01Ht,sfs,z2s-fs,z1sds≤∫01Ht,spsz2s-z1sds≤Lz2-z1t≤rLKz1-z0kψLψt=r2LKz1-z0kψψt.⋮(63)zn+1t-znt=∫01Ht,sfs,zns-fs,zn-1sds≤∫01Ht,spszns-zn-1sds≤Lzn-zn-1t≤rn-1LKz1-z0kψLψt=rnLKz1-z0kψψt.
Consequently, for anyn,m∈N, it follows from (61)–(63) that(64)zn+mt-znt=zn+mt-zn+m-1t+zn+m-1t-zn+m-2t+⋯+zn+1t-znt≤rn+m-1L+rn+mL⋯+rnLKz1-z0kψψt=rnL1-rmLKz1-z01-rLkψψt≤rnLKz1-z01-rLkψψt.Thus we have(65)zn+m-zn≤rnLKz1-z01-rLkψψ. Noticing that 0<r(L)<1, for any m∈N, one gets(66)zn+m-zn→0,n→+∞, which implies that {zn} is a Cauchy sequence. Consequently, there exists z∗∈Q such that zn converges to z∗. Thus z∗∈Q is a fixed point of T. It follows from Lemma 4 that x∗(t)=Iγz∗(t) is a positive solution of (4). Moreover, since z∗∈Q, there exist two positive constants Kz∗>kz∗>0 such that(67)kz∗tα-γ-1≤z∗t≤Kz∗tα-γ-1. Thus (68)Γα-γkz∗Γαtα-1≤x∗t=Iγz∗t=1Γγ∫0tt-sγ-1z∗sds≤Γα-γKz∗Γαtα-1.That is, there exists a constant μ2>μ1>0 such that (69)μ1tα-1≤x∗t≤μ2tα-1.
Step 2. Next we shall show that the fixed point of T is unique. In fact, for any positive fixed point z¯≠z∗∈Q of T, similar to (56) and (60), there exists a constant Kz∗-z¯>0 such that(70)Lz1t-z0t≤rLKz∗-z¯kψψt. Thus we have(71)z∗t-z¯t=Tz∗t-Tz¯t=∫01Ht,sfs,z∗s-fs,z¯sds≤∫01Ht,spsz∗s-z¯sds≤Lz∗t-z¯t≤rLKz∗-z¯kψψt.By induction, we have(72)Tnz∗t-Tnz¯t≤rnLKz∗-z¯kψψt. Thus it follows from (71) and r(L)<1 that(73)z∗-z¯=Tnz∗-Tnz¯≤rnLKz∗-z¯kψψ→0,n→+∞,which implies that z¯=z∗, a contradiction. So the positive fixed point of T is unique.
Step 3. In the following, we consider the convergence analysis of solution. Similar toStep 1, for any initial z0∈P∖θ, construct successively a sequence (74)zm=∫01Ht,sfs,zm-1sds,m=1,2,…, and then the iterative sequence zm(t) converges uniformly to the unique fixed point z∗ of T satisfying x∗(t)=Iγz∗(t); i.e., zm(t) converges uniformly to z∗(t)=Dtγx∗(t) on [0,1] as n→∞; that is, (75)zm-Dtγx∗→0,as m→∞. Furthermore, we have error estimation (76)zm-Dtγx∗≤Kz1-z0ψ1-rLkψrmL with the rate of convergence(77)um-Dtγx∗=OrmL.Remark 13.
In Theorem12, we not only establish the condition of existence of unique positive solution for (4), but also construct an iterative sequence which converges uniformly to γ order derivative of the unique solution of (4). In particular, the error estimation of between exact solution and approximate solution and the corresponding rate of convergence are also obtained.Remark 14.
Here we also briefly state how the design parameters affect the control performance and how to choose these parameters.
( 1 ) To design parameters effect of the condition (H1) for system, we can choose nonlinearity f as linear functions or sine (cosine) function which shall very easily satisfy the control condition (H1).
( 2 ) For control condition (H2), we can choose a function p∈L1(0,1) such that ω(s)p(s) is a power function satisfying the power exponent larger than -1, and then the control condition (H2) will naturally hold.
( 3 ) For the selection of the control condition r(L)∈(0,1), according to the Gelfand’s formula, the spectral radius of Ln satisfies (78)rL=limn→+∞Ln1/n,rL≤Ln1/n. In particular, taking n=1, we have (79)rL≤L=sup0≤t≤1∫01Ht,spsds≤∫01ωspsds. Thus (H2) can be replaced: (80)0<∫01ωspsds<1. That is, we can choose the coefficient of function p such that r(L)∈(0,1). Consequently in Theorem 12, the restrictions r(L)∈(0,1) can be omitted.Remark 15.
Noticing thattα-γ-1∈P∖{θ}, the iterative process can be started from the initial value z0=tα-γ-1, which will simplify the whole iterative process.


## 4. Nonexistence Results

In this section, we focus on the nonexistence results of positive solution of (4).Theorem 16.
Assume that(H0)-(H2) hold. Then (4) has no positive solution provided that f(t,0)≡0,t∈(0,1) and the spectral radius r(L)∈(0,1).Proof.
It is sufficient to prove that the operatorT has no positive fixed point in P. In fact, if not, there exists z∈P such that Tz=z. It follows from Lemmas 10 and 11 that z∈Q and there exists a positive eigenfunction ψ∈Q such that Lψ=r(L)ψ. Notice that z,ψ∈Q, and there exists a constant b>0 such that (81)z≤bψ.
On the other hand, it follows from(H1) and f(t,0)≡0,t∈(0,1) that f(t,z)≤p(t)z,t∈(0,1),z∈[0,+∞). Thus we have z=Tz≤Lz. By induction, we can get (82)z≤Lnz,n=1,2,3,…, which implies that (83)z≤Lnz=Lnbψ≤brnLψ→0,asn→+∞, since the spectral radius r(L)∈(0,1). Therefore, z=0, which contradict with z∈Q. The proof is completed.Theorem 17.
Assume that(H0) and (H2) hold and there exists p∈L[0,1] such that (84)ft,u≥ptu,t∈0,1,u∈0,+∞. Then (4) has no positive solution provided that the spectral radius r(L)>1.Proof.
The proof is similar to Theorem16, so here we omit the proof.


## 5. Numerical Result and Simulation

In this section, we recall the example (1) in introduction and consider the existence of positive solutions for (1).Conclusion.Equation (1) has a unique positive solution x∗, and there exist two constants μ2>μ1>0 such that (85)μ1t5/2≤x∗t≤μ2t5/2.Moreover, for any initial z0∈P∖θ, construct successively a sequence (86)zm=∫01t3/2Ws1.4334+Jt,s1+sin⁡zm-1s3sds,m=1,2,…,and then the iterative sequence zm(t) converges uniformly to Dt1/4x∗(t) on [0,1] as n→∞, i.e., (87)zm-Dt1/4x∗→0,as m→∞. Furthermore, there exists a constant such that error estimation satisfies (88)zm-Dt1/4x∗≤a×0.7832mwith the rate of convergence(89)zm-Dt1/4x∗=O0.7832m.Proof.
Letα=11/4,β=3/2,γ=1/4, then we have 0<γ<1<β<2<α≤3 with α>γ+2,β>γ+1 and(90)Jt,s=1Γ5/2t3/21-s1/4-t-s3/2≕J1t,s,0≤s≤t≤1,t3/21-s1/4≕J2t,s,0≤t≤s≤1. Thus, we have (91)Ws=∫01Jt,sdAt=1Γ5/23J113,s-J123,s,0≤s<13,3J213,s-J123,s,13≤s<23,3J212,s-J223,s,23≤s≤1.=1Γ5/233-2691-s1/4+23-s3/2-313-s3/2,0≤s<13,33-2691-s1/4+23-s3/2,13≤s<23,33-2691-s1/4,23≤s≤1.which implies that W(s)≥0 for s∈[0,1].
On the other hand, we have(92)A=∫01t3/2dAt=2-32∫01Att1/2dt=0.0332<Γ5/2Γ5/4=1.4666. Therefore, the condition (H0) is satisfied.
Let(93)ft,u=1+sin⁡u3t,pt=13t,and then f∈C((0,1)×[0,+∞),[0,+∞)) with f(t,0)≢0, and (94)ft,u1-ft,u2≤ptu1-u2,t∈0,1,ui∈0,+∞,i=1,2.Thus the condition (H1) is also satisfied.
Now we check the condition(H2). In fact, since (95)ω1s≤34.3002Γ5/2+1-s1/4-1-s3/2Γ5/2,ω2s=1-s1/4Γ5/2, we have (96)0<∫01ωspsds≤∫0134.3002Γ5/2+1-s1/4-1-s3/2Γ5/2+1-s1/4Γ5/213sds=0.7832<1.
Thus according to Theorem12, (1) has a unique positive solution x∗, and there exist two constants μ2>μ1>0 such that (97)μ1t5/2≤x∗t≤μ2t5/2. Moreover, for any initial z0∈P∖θ, construct successively a sequence (98)zm=∫01t3/2Ws1.4334+Jt,s1+sin⁡zm-1s3sds,m=1,2,…,and the iterative sequence zm(t) converges uniformly to Dt1/4x∗(t) on [0,1] as n→∞, i.e., (99)zm-Dt1/4x∗→0,as m→∞. Furthermore, there exists a constant a>0 such that error estimation satisfies (100)zm-Dt1/4x∗≤a×0.7832m with the rate of convergence(101)zm-Dt1/4x∗=O0.7832m.
In the following, we are going to show the effectiveness of the proposed iterative method. Since(102)
z
m
=
∫
0
1
t
3
/
2
W
s
1.4334
+
J
t
,
s
1
+
sin
⁡
z
m
-
1
s
3
s
d
s
=
∫
0
t
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
-
3
1
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
×
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
t
1
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
-
3
1
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
1
/
3
2
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
1
+
sin
⁡
z
m
-
1
s
3
s
d
s
+
∫
2
/
3
1
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
t
3
/
2
1
-
s
1
/
4
1
+
s
i
n
z
m
-
1
s
3
s
d
s
,
t
∈
0
,
1
3
,
∫
0
1
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
-
3
1
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
×
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
1
/
3
t
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
t
2
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
2
/
3
1
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
t
3
/
2
1
-
s
1
/
4
1
+
s
i
n
z
m
-
1
s
3
s
d
s
,
t
∈
1
3
,
2
3
,
∫
0
1
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
-
3
1
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
×
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
1
/
3
2
/
3
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
2
3
-
s
3
/
2
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
1
+
sin
⁡
z
m
-
1
s
3
s
d
s
+
∫
2
/
3
t
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
t
3
/
2
1
-
s
1
/
4
-
t
-
s
3
/
2
1
+
s
i
n
z
m
-
1
s
3
s
d
s
+
∫
t
1
t
3
/
2
1.4334
3
3
-
2
6
9
1
-
s
1
/
4
+
t
3
/
2
1
-
s
1
/
4
1
+
s
i
n
z
m
-
1
s
3
s
d
s
,
t
∈
2
3
,
1
,let z0=-t3/2, and then the simulations for the iterative process (see Figure 1 and Table 1) show that iterative process is very effective and the iterative convergence speed (see Table 2) is robust.Table 1
The numerical approximation of the solution for the Eq. (1).
m t 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 m=1 0.0282 0.0714 0.1163 0.1134 0.1336 0.1442 0.1440 0.1329 0.1119 m=2 0.0299 0.0840 0.1552 0.1683 0.2314 0.2987 0.3677 0.4366 0.5033 m=3 0.0300 0.0850 0.1606 0.1765 0.2510 0.3378 0.4366 0.5469 0.6669 m=4 0.0300 0.0849 0.1598 0.1763 0.2519 0.3420 0.4484 0.5725 0.7144Table 2
Maximum errors between iterative values and exact solution.
m=2 m=3 m=4 m=5 m=6 m=7 Error 0.3680 0.0299 4.7955 e - 007 1.5217 e - 032 7.4495 e - 192 0Figure 1
The approximation of the solution for (1).Remark 18.
The maximum errors between iterative values and exact solution are exposed in Table2, which validates the convergence of the maximum errors to zero at a fast speed. Thus by using the suggested iterative process, the simulation results are robust against nonlinearities and singularity at time variable and also exhibit a reasonable performance for exploring the uniqueness of solution for the singular fractional order viscoelasticity complex system.


## 6. Conclusion

Mathematical models involving fractional derivatives can describe many chaotic systems, advection-dispersion process, viscoelasticity characteristics, and the bioprocesses with long memory. In this paper, we offer some control conditions to govern a complex system arising from ecological economy phenomena and diffusive interaction. These control conditions make us easily judge the existence of unique solution of system and further obtain the approximate solution of system according to the required precision of physical problem. Numerical results show that the iterative convergence speed is very fast, which also implies that the suggested method is more accurate to depict the coevolution process of economic, social, and ecological subsystems and the transport of solute in highly heterogeneous porous media compared to what is mentioned before. The extension of the infinite time controller design based on noncompactness measure for complex viscoelasticity systems with multiple stresses will be presented in future work.



---

*Source: 9278056-2019-09-02.xml*