# Series Solution of the System of Fuzzy Differential Equations

Authors: M. S. Hashemi; J. Malekinagad; H. R. Marasi
Journal: Advances in Fuzzy Systems
 (2012)

---

## Abstract

The homotopy analysis method (HAM) is proposed to obtain a semianalytical solution of the system of fuzzy differential equations (SFDE). The HAM contains the auxiliary parameterħ, which provides us with a simple way to adjust and control the convergence region of solution series. Concept of ħ-meshes and contour plots firstly are introduced in this paper which are the generations of traditional h-curves. Convergency of this method for the SFDE has been considered and some examples are given to illustrate the efficiency and power of HAM.


---

## Body

## 1. Introduction

In many cases of the modeling real world phenomena, information about the behavior of a dynamical system is uncertain. In order to obtain a more realistic model, we have to take into account these uncertainties. Since 1965, when Zadeh published his pioneering paper [1], hundreds of examples have been supplied where the nature of uncertainty in the behavior of a given system processes is fuzzy rather than stochastic nature. The concept of fuzzy derivative was first introduced by Chang and Zadeh in [2]. It was followed up by Dubois and Prade in [3], who defined and used the extension principle. Other methods have been discussed by Puri and Ralescu in [4] and Goetschel and Voxman in [5]. The initial value problem for fuzzy differential equation (FIVP) has been studied by Kaleva in [6, 7] and by Seikkala in [8].The purpose of this paper is to find the approximate solution of fuzzy differential equations system with the homotopy analysis method (HAM) introduced first by Liao in 1992 [9, 10], that is, analytic approach to get series solutions of various types of linear and nonlinear equations.Some of numerical methods have been applied to obtain the solution of fuzzy differential equations [11–15]. Sami Bataineh et al. have applied the HAM for systems of ODEs and PDEs in [16, 17]. Also, recently many types of nonlinear problems solved with HAM by others [18–22].In Section2, some basic definitions which will be used later in the paper are provided. In Sections 3 and 4, system of fuzzy differential equations and then basic ideas of HAM applied to these types of equations have been reviewed, respectively. Convergency of HAM for SFDE that shows its reliability has been considered in Section 5. The proposed method is illustrated by solving several examples in Section 6, and finally the conclusion is drawn in Section 7.


## 2. Preliminaries

In this section, the most basic notations used in this paper are introduced.Definition 1 (See [23]).
An arbitrary fuzzy numberu in parametric form is an ordered pair (u_(r),u¯(r)) of functions u_(r), u¯(r); 0≤r≤1, which satisfies the following requirements:(1)
u_(r) is a bounded left-continuous nondecreasing function over [0,1],(2)
u¯(r) is a bounded left-continuous nonincreasing function over [0,1],(3)
u_(r)≤u¯(r), 0≤r≤1.
The set of all such fuzzy numbers is represented by𝔼1.Remark 2.
For arbitraryu=(u_(r),u¯(r)),  v=(v_(r),v¯(r)), and k∈ℝ, we define addition and multiplication by k as
(1)(u+v)_(r)=u_(r)+v_(r),(u+v)¯(r)=u¯(r)+v¯(r),ku_(r)=ku_(r),ku¯(r)=ku¯(r),ifk≥0,ku_(r)=ku¯(r),ku¯(r)=ku_(r),ifk<0.Definition 3.
For arbitrary fuzzy numbersu, v∈𝔼1, we use the distance
(2)D(u,v)=sup⁡0≤r≤1max⁡{|u¯(r)-v¯(r)|,|u_(r)-v_(r)|},
and it is shown that (𝔼1,D) is a complete metric space.Definition 4.
Letf:[a,b]→𝔼1, for each partition P={t0,…,tn} of [a,b] and for arbitrary ξi∈[ti-1,ti], 1≤i≤n suppose
(3)Rp=∑i=1n‍f(ξi)(ti-ti-1),Δ:=max⁡{|ti-ti-1|,i=1,…,n}.
The definite integral of f(t) over [a,b] is
(4)∫ab‍f(t)dt=lim⁡Δ→0Rp,
provided that this limit exists in the metric D. If the fuzzy function f(t) is continuous in the metric D, its definite integral exists and also,
(5)∫ab‍f(t;r)dt_=∫ab‍f_(t;r)dt,∫ab‍f(t;r)dt¯=∫ab‍f¯(t;r)dt.Definition 5.
Letf:ℝ1→𝔼1 be a fuzzy function and let t0∈ℝ1. The derivative f′(t0) of f at the point t0 is defined by
(6)f′(t0)=lim⁡h→0f(t0+h)-f(t0)h,
provided that this limit, taken with respect to the metric D, exists.The elementsf(t0+h), f(t0) at the right-hand side of  (6) are observed as elements in the Banach space 𝔹=ℂ[0,1]×ℂ[0,1]. Thus, if f(t0+h)=(a_,a¯) and f(to)=(b_,b¯), the difference is simply f(t0+h)-f(to)=(a_-b_,a¯-b¯).Clearly[f(t0+h)-f(t0)]/h may not be a fuzzy number for all h. However, if it approaches f′(t0) (in 𝔹) and f′(t0) is also a fuzzy number (i.e., in 𝔼1), this number is the fuzzy derivative of f(t) at t0. In this case, if f=(f_,f¯), it can be easily shown that
(7)f′(t0)=(f_′(t0),f¯′(t0)),
where f_′ and f¯′ are the classic derivatives of f_ and f¯, respectively.


## 3. System of Fuzzy Differential Equations

In this section, we will review system of fuzzy differential equations of the forms(8)∑i=0n‍Ai(t)Y(i)(t)=F(t),Y(i)(0)=Gi,i=0,…,n-1,
where t is a scaler and A0(t),A1(t),…,An-1(t) are s×s matrixes and every component of them is a real function of t and An(t)=Is denotes the s×s identity matrix. Y, Gi, and F are fuzzy s-dimensional vectors. The rth component of Y∈Es will be denoted by yr, so that we may write
(9)Y=[y1(t,r),y2(t,r),…,ys(t,r)]T,yi(t,r)=((yl)i(t,r),(yu)i(t,r)),i=1,…,s,F=[f1(t,r),f2(t,r),…,fs(t,r)]T,fi(t,r)=((fl)i(t,r),(fu)i(t,r)),i=1,…,s,Gi=[gi1(r),gi2(r),…,gis(r)]T,
where
(10)gij(r)=((gl)ij(r),(gu)ij(r)),i=1,…,n-1,j=1,…,s.The superscriptT denotes transpose. The component that is in rth row, kth column, 1≤r,  k≤s of matrix Ap, 0≤p≤n will be denoted by ar,kp. Then, (8) can be replaced by the following equivalent system in parametric form:
(11)∑i=0n∑k=1sad,ki(t)yk(i)(t,r)_=fd(t,r)_,d=1,…,s,∑i=0n∑k=1sad,ki(t)yk(i)(t,r)¯=fd(t,r)¯,d=1,…,s,yj(i)_(0,r)=gij_(r),i=0,…,n-1,j=1,…,s,yj(i)¯(0,r)=gij¯(r),i=0,…,n-1,j=1,…,s,
where r∈[0,1].We define functionsbd,ki(t) and cd,ki(t);  1≤d, k≤s,  0≤i≤n, in the following form:
(12)bd,ki(t)={0,ad,ki(t)<0,ad,ki(t),ad,ki(t)≥0,cd,ki(t)={ad,ki(t),ad,ki(t)<0,0,ad,ki(t)≥0.
Now (11) becomes as follows:
(13)∑i=0n∑k=1s[bd,ki(t)(yl)k(i)(t,r)+cd,ki(t)(yu)k(i)(t,r)]=(fl)d(t,r),d=1,…,s∑i=0n∑k=1s[cd,ki(t)(yl)k(i)(t,r)+bd,ki(t)(yu)k(i)(t,r)]=(fu)d(t,r),d=1,…,s(yl)j(i)(0,r)=(gl)ij(r),i=0,…,n-1,j=1,…,s,(yu)j(i)(0,r)=(gu)ij(r),i=0,…,n-1,j=1,…,s.


## 4. Basic Ideas of HAM

We consider the following differential equations:(14)∑i=0n‍Ai(t)Z(i)(t,r)-F(t,r)=0,
where
(15)Z(t,r)=[z1(t,r),z2(t,r),…,zs(t,r)]T,F(t,r)=[f1(t,r),f2(t,r),…,fs(t,r)]T,
and A0(t),A1(t),…,An-1(t) are s×s matrixes and every component of them is a real function of t and An(t)=Is, where Is denotes the s×s identity matrix.The above system of equations can be written in following form:(16)Ni[z1(t,r),z2(t,r),…,zs(t,r)]=0,i=1,2,…,s,
where 𝒩i are nonlinear operators that represent the whole equations, t and r denote the independent variables, and zi(t,r),  i=1,2,…,s are unknown functions, respectively. By means of generalizing the traditional homotopy method constructed the so-called zero-order deformation equations
(17)(1-q)Li[ϕi(t,r;q)-zi,0(t,r)]=qℏNi[ϕ1(t,r;q),…,ϕn(t,r;q)],i=1,…,n,
where q∈[0,1] is an embedding parameter, ℏ is nonzero auxiliary parameter, ℒi,  i=1,2,…,s are auxiliary linear operators, zi,0(t,r) are initial guesses of zi(t,r), and ϕi(t,r;q) are unknown functions. It is important to note that one has great freedom to choose auxiliary objects such as ℏ and ℒi in HAM. Obviously, when q=0 and q=1, both
(18)ϕi(t,r;0)=zi,0(t,r),ϕi(t,r;1)=zi(t,r)
hold. Thus, as q increases from 0 to 1, the solutions ϕi(t,r;q) varies from the initial guesses zi,0(t,r) to the solutions zi(t,r). Expanding ϕi(t,r;q) in Taylor’s series with respect to q, one has
(19)ϕi(t,r;q)=zi,0(t,r)+∑m=1+∞‍zi,m(t,r)qm,
where
(20)zi,m=1m!∂mϕi(t,r;q)∂qm∣q=0.
If auxiliary linear operator, initial guesses, and auxiliary parameter ℏ are properly chosen, then the series equations (19) converges at q=1 and
(21)ϕi(t,r;1)=zi,0(t,r)+∑m=1+∞‍zi,m(t,r),i=1,…,s,
which must be solution of the original nonlinear equations.According to (20), the governing equations can be deduced from the zero-order deformation equations (17). Define the vectors
(22)z→i,m={zi,0(t,r),zi,1(t,r),…,zi,m(t,r)},i=1,…,s.
Differentiating (17) m time with respect to the embedding parameter q and then setting q=0 and finally dividing them by m!, we have the so-called mth order deformation equations
(23)Li[zi,m(t,r)-χmzi,m-1(t,r)]=ℏRi,m(z→1,m-1,…,z→s,m-1)i=1,…,s,
where
(24)Ri,m(z→1,m-1,…,z→s,m-1)=1(m-1)!∂m-1Ni[ϕ1(t,r;q),…,ϕs(t,r;q)]∂qm-1∣q=0,i=1,…,s,χm={0,m≤1,1,m>1.Now, to simplify and solve (11), we need to know the sign of components of Ai(t);  0≤i≤n and use of Remark 2. Also we should be noted that in case that some components of A0(t),A1(t),…,An(t) for some t are negative, we must divide the defined interval for t into small intervals, so that sign of each component of A0(t),A1(t),…,An(t) in each small interval be unchanged. Now for each small interval we have a separate system 2s×s(n+1).


## 5. Convergency of HAM for System of Fuzzy Differential Equations

In this section, we prove a theorem which shows the convergency of approximate HAM solution applied for (8).First we defineℛk,m[(yl→)1,m-1,…,(yl→)s,m-1, (yu→)1,m-1,…,(yu→)s,m-1];  1≤k≤s, in the following form:(25)Rk,m=1(m-1)!∂m-1Nk[(yl)1,m-1,…,(yl)s,m-1,(yu)1,m-1,…,(yu)s,m-1]∂qm-1∣q=0,ifkth equation of (13) is inclusive bd,ki(t)(yl)k(i)(t,r);  1≤i≤n, 1≤d≤s, and otherwise(26)Rk,m=1(m-1)!∂m-1Nk+s[(yl)1,m-1,…,(yl)s,m-1,(yu)1,m-1,…,(yu)s,m-1]∂qm-1∣q=0.Theorem 6.
Let the series∑m=1∞(y)k,m(t,r);1≤k≤s be uniformly converge to (y)k(t,r); 1≤k≤s. Then ((yl)k(t,r),(yu)k(t,r)),1≤k≤s is exact solution of (13).Proof.
Since∑m=1∞‍(y)k,m(t,r),1≤k≤s is convergent, we must have
(27)lim⁡m→∞(y)k,m(t,r)=0,1≤k≤s,0≤t≤T,0≤r≤1.
On the other hand, since ℒi are linear operators, thus
(28)∑m=1n‍Lk[(y)k,m(t,r)-χm(y)k,m-1(t,r)]=∑m=1n‍[Lk(y)k,m(t,r)-χmLk[(y)k,m-1(t,r)]]=Lk[(y)k,1(t,r)]+(Lk[(y)k,2(t,r)]-Lk[(y)k,1(t,r)])+⋯+(Lk[(y)k,n(t,r)]-Lk[(y)k,n-1(t,r)])=Lk[(y)k,n(t,r)],k=1,2,…,s,
then from (27) we can write
(29)∑m=1∞‍Lk[(y)k,m(t,r)-χm(y)k,m-1(t,r)]=lim⁡n→∞Lk[(y)k,n(t,r)]=Lk[lim⁡n→∞(y)k,n(t,r)]=0,k=1,2,…,s,
hence,
(30)ℏ∑m=1∞‍Rk,m[(yl→)1,m-1,…,(yl→)s,m-1,(yu→)1,m-1,…,(yu→)s,m-1]=∑m=1∞‍Lk[(y)k,m(t,r)-χm(y)k,m-1(t,r)]=0,k=1,2,…,s.
Since ℏ≠0, thus from the above equations for k=1,2,…,s, we can write
(31)∑m=1∞‍Rk,m[(yl→)1,m-1,…,(yl→)s,m-1,(yu→)1,m-1,…,(yu→)s,m-1]=0.
From uniform convergency we have
(32)∑m=1∞‍∂i∂ti(y)k,m-1(t,r)=∂i∂ti(y)k(t,r),i=0,1,…,n,k=1,2,…,s,
thus
(33)0=∑m=1∞‍[∑i=0n‍Ai(t)Y(i)(t)-(1-χm)(F)d(t,r)]=∑i=0n‍Ai(t)∑m=1∞‍Y(i)(t)-∑m=1∞‍(1-χm)(F)d(t,r)=∑i=0n‍Ai(t)∑m=1∞‍∂i∂ti(Y)k,m-1(t,r)-(F)d(t,r)=∑i=0n‍Ai(t)∂i∂ti(Y)k(t,r)-(F)d(t,r),d=1,…,s.
Therefore, {yk(t,r)=((yl)k(t,r),(yu)k(t,r));k=1,2,…s} is the exact solution of (11) and Y=[y1(t,r),y2(t,r),…,ys(t,r)]T  is the exact solution of (8) and proof is completed.


## 6. Illustrative Examples

Example 7.
Consider the following second-order fuzzy linear differential equation:(34)Y′′+Y=-t,t∈[0,1]Y(0)=(0.1r-0.1,0.1-0.1r),Y′(0)=(0.088+0.1r,0.288-0.1r).
The exact solution is as follows:
(35)Y_(t,r)=(0.1r-0.1)cos⁡(t)+(1.088+0.1r)sin⁡(t)-t,Y¯(t,r)=(0.1-0.1r)cos⁡(t)+(1.288-0.1r)sin⁡(t)-t.
According to (11), we may replace (34) by the following equivalent system:
(36)Y_′′+Y_=-t,t∈[0,1],Y¯′′+Y¯=t,t∈[0,1],Y_(0)=0.1r-0.1,Y_′(0)=0.088+0.1r,Y¯(0)=0.1-0.1r,Y¯′(0)=0.288-0.1r.
We first construct the zero-order deformation equations
(37)(1-q)Li[ϕi(t,r;q)-zi,0(t,r)]=qℏNi[ϕi(t,r;q)],i=1,2,
subject to the initial conditions
(38)Y_0=(0.1r-0.1)+(0.088+0.1r)t,Y0¯=(0.1-0.1r)+(0.288-0.1r)t,
and the linear operator
(39)Li[ϕi(t,r;q)]=σ2ϕi(t,r;q)σt2,i=1,2,
with the property
(40)Li[ci0+tci1]=0,i=1,2,
where ci0, ci1(i=1,2) are integral constants. Also from (36), we can define
(41)N1[ϕ1(t,r;q)]=σ2ϕ1(t,r;q)σt2+ϕ1(t,r;q)+t,N2[ϕ2(t,r;q)]=σ2ϕ2(t,r;q)σt2+ϕ2(t,r;q)+t.
Obviously, when q=0 and q=1,
(42)ϕ1(t,r;0)=z1,0(t,r)=Y_0,ϕ1(t,r;1)=Y_,ϕ2(t,r;0)=z2,0(t,r)=Y¯0,ϕ2(t,r;1)=Y¯.
Therefore, when the embedding parameterq increases from 0 to 1, the homotopy solutions ϕi(t,r,q) vary from zi,0(t,r) to the solutions zi(t,r) for i=1,2. Now, by expanding ϕi(t,r;q) in Taylor’s series with respect to q, we have
(43)ϕi(t,r;q)=zi,0(t,r)+∑m=1+∞‍zi,m(t,r)qm,i=1,2,
where
(44)zi,m(t,r)=1m!∂mϕi(t,r;q)∂qm|q=0,i=1,2.
Assuming that auxiliary parameter ℏ, the initial guesses and the auxiliary linear operator are properly chosen, then the above series is convergent at q=1, and
(45)Y_(t,r)=z1,0(t,r)+∑m=1∞‍z1,m(t,r),Y¯(t,r)=z2,0(t,r)+∑m=1∞‍z2,m(t,r).
The mth order deformation equations are
(46)Li[zi,m(t,r)-χmzi,m-1(t,r)]=ℏRi,m(z→i,m-1),i=1,2,
with the initial conditions
(47)zi,m(0,r)=0,i=1,2,
where
(48)Ri,m(z→i,m-1)=∂2zi,m-1∂t2+zi,m-1+(1-χm)t,i=1,2.
Therefore, we recursively obtain
(49)z1,1(t,r)=-ℏt220+ℏrt220+68ℏt3375+ℏrt360,z1,2(t,r)=-ℏt220--ℏ2t220+ℏrt220+ℏ2rt220+68ℏt3375+68ℏ2t3375+ℏrt360+ℏ2rt360-ℏ2t4240+ℏ2rt4240+17ℏ2t51875+ℏ2rt51200,z2,1(t,r)=ℏt220-ℏrt220+161ℏt3750-ℏrt360,z2,2(t,r)=ℏt220+-ℏ2t220-ℏrt220-ℏ2rt220+161ℏt3750+161ℏ2t3750-ℏrt360-ℏ2rt360+ℏ2t4240-ℏ2rt4240+161ℏ2t515000-ℏ2rt51200,⋮.
Then the solutions obtained by HAM are as follows:
(50)Y_(t,r)=z1,0(t,r)+z1,1(t,r)+z1,2(t,r)+z1,3(t,r)+⋯,Y¯(t,r)=z2,0(t,r)+z2,1(t,r)+z2,2(t,r)+z2,3(t,r)+⋯.
Figures1 and 3 show the ℏ-mesh of Y_ and Y¯ to get a proper interval for convergency. ℏ-mesh is a generalization of traditional Contour plots, and for connection between ℏ-meshes and Contour plots we plot Figures 2 and 4. Also to find the best quantity of ℏ that lies in the convergency interval, we use the residual error of norm 2 as follows:
(51)Res[Y_n]=[∫∫01‍(Y_n(t,r)-Y_(t,r))2dtdr]1/2,
which is a function with respect to ℏ. Now, by minimizing the Res[Y_10] we obtain the best choice for auxiliary parameter to approximate of Y_10(t,r) as follows:
(52)ℏ=-1.02698,
and in this case absolute error for the 10th order approximation by HAM for Y_n(t,r) is plotted in Figure 5.
By minimizing the residual error defined by(53)Res[Y¯n]=[∫∫01‍(Y¯n(t,r)-Y¯(t,r))2dtdr]1/2,
it is clear that the best choice for auxiliary parameter to approximate Y¯10(t,r) is
(54)ℏ=-1.02683
which in this case absolute error for the 10th order approximation by HAM for Y¯n(t,r) is plotted in Figure 6.Figure 1
ℏ-mesh of (∂3/∂t3)Y_10(0,r).Figure 2
Contour plot of(∂3/∂t3)Y_10(0,r).Figure 3
ℏ-mesh of (∂3/∂t3)Y¯10(0,r).Figure 4
Contour plot of(∂3/∂t3)Y¯10(0,r).Figure 5
Absolute error ofY_10(t,r) for ℏ=-1.02698.Figure 6
Absolute error ofY¯10(t,r) for ℏ=-1.02683.Example 8.
Consider the system of fuzzy differential equation(55)Y′′+e2tG′=((4r2+5r-1)e2t,(-4r3-2rt+2t+12)e2t),Y′+G′′=((2r2+2r)e2t,(-2r3+6)e2t-2r+2),Y(0)=(r2+r,3-r3),Y′(0)=(2r2+2r,-2r3+6),G(0)=(0,0),G′(0)=(r-1,0),t∈[0,1].
The exact solution is as follows:(56)Y(t,r)=((r2+r)e2t,(3-r3)e2t),G(t,r)=((r-1)t,(1-r)t2).
We may replace (55) by the equivalent system
(57)Y_′′+e2tG_′=(4r2+5r-1)e2t,Y¯′′+e2tG¯′=(-4r3-2rt+2t+12)e2t,Y_′+G_′′=(2r2+2r)e2t,Y¯′+G¯′′=(-2r3+6)e2t-2r+2,Y_(0,r)=r2+r,Y¯(0,r)=3-r3,Y_′(0,r)=2r2+2r,Y¯′(0,r)=-2r3+6,G_(0,r)=0,G¯(0,r)=0,G_′(0,r)=r-1,G¯′(0,r)=0.
We select the following initial guesses:
(58)Y_0=r2+r+(2r2+2r)t,Y¯0=3-r3+(-2r3+6)t,G_0=(r-1)t,G¯0=2-r,
and the auxiliary linear operators
(59)Li[ϕi(t,r;q)]=∂2ϕi(t,r;q)∂t2,i=1,2,3,4,
with the property
(60)Li[ci0+tci1]=0,i=1,2,3,4,
where ci0 and ci1(i=1,2,3,4) are integral constants. Also the nonlinear operators are
(61)N1[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂2ϕ1(t,r;q)∂t2+e2t∂ϕ3(t,r;q)∂t-(4r2+5r-1)e2t,N2[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂2ϕ2(t,r;q)∂t2-∂ϕ4(t,r;q)∂t-(-4r3-2rt+2t+12)e2t,N3[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂2ϕ3(t,r;q)∂t2+∂ϕ1(t,r;q)∂t-(2r2+2r)e2t,N4[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂2ϕ4(t,r;q)∂t2+∂ϕ2(t,r;q)∂t-((-2r3+6)e2t-2r+2).
Then mth order deformation equations are
(62)Li[zi,m(t,r)-χmzi,m-1(t,r)]=ℏRi,m(z→1,m-1,…,z→4,m-1),i=1,2,3,4,
with the initial conditions
(63)zi,m(0,r)=0,i=1,2,3,4,R1,m(z→1,m-1,…,z→4,m-1)=∂2z1,m-1∂t2+e2t∂z3,m-1∂t-(1-χm)(4r2+5r-1)e2t,R2,m(z→1,m-1,…,z→4,m-1)=∂2z2,m-1∂t2-∂z4,m-1∂t-(1-χm)(-4r3-2rt+2t+12)e2t,R3,m(z→1,m-1,…,z→4,m-1)=∂2z3,m-1∂t2+∂z1,m-1∂t-(1-χm)(2r2+2r)e2t,R4,m(z→1,m-1,…,z→4,m-1)=∂2z4,m-1∂t2+∂z2,m-1∂t-(1-χm)((-2r3+6)e2t-2r+2).
Therefore, we recursively obtain
(64)z1,1(t,r)=ℏr-e2tℏr+ℏr2-e2tℏr2+2ℏrt+2ℏr2t,z1,2(t,r)=ℏr-e2tℏr+21ℏ2r16-54e2tℏ2r-116e4tℏ2r+ℏr2-e2tℏr2+21ℏ2r216-54e2tℏ2r2-116e4tℏ2r2+2ℏrt+94ℏ2rt+12e2tℏ2rt+2ℏr2t+94ℏ2r2t+12e2tℏ2r2t,z2,1(t,r)=5ℏ2-52e2tℏ+ℏr2-12e2tℏr-ℏr3+e2tℏr3+11ℏt2-12e2tℏt+ℏrt2+12e2tℏrt-2ℏr3t,z2,2(t,r)=5ℏ2-52e2tℏ+47ℏ216-114e2tℏ2-316e4tℏ2+ℏr2-12e2tℏr+ℏ2r-e2tℏ2r-ℏr3+e2tℏr3-21ℏ2r316+54e2tℏ2r3+116e4tℏ2r3+11ℏt2-12e2tℏt+23ℏ2t4+12e2tℏ2t+ℏrt2+12e2tℏrt+ℏ2rt+e2tℏ2rt-2ℏr3t-94ℏ2r3t-12e2tℏ2r3t,z3,1(t,r)=ℏr2-12e2tℏr+ℏr22-12e2tℏr2+ℏrt+ℏr2t+ℏrt2+ℏr2t2,z3,2(t,r)=ℏr2-12e2tℏr+ℏ2r-e2tℏ2r+ℏr22-12e2tℏr2+ℏ2r2-e2tℏ2r2+ℏrt+2ℏ2rt+ℏr2t+2ℏ2r2t+ℏrt2+2ℏ2rt2+ℏr2t2+2ℏ2r2t2,z4,1(t,r)=3ℏ2-32e2tℏ-ℏr32+12e2tℏr3+3ℏt-ℏr3t+2ℏt2+ℏrt2-ℏr3t2,z4,2(t,r)=3ℏ2-32e2tℏ+21ℏ28-218e2tℏ2+3ℏ2r8-38e2tℏ2r-ℏr32+12e2tℏr3-ℏ2r3+e2tℏ2r3+3ℏt+11ℏ2t2-14e2tℏ2t+12ℏ2rt+14e2tℏ2rt-ℏr3t-2ℏ2r3t+2ℏt2+19ℏ2t24+ℏrt2+54ℏ2rt2-ℏr3t2-2ℏ2r3t2,⋮.
Then the solutions obtained by HAM are as follows:
(65)Y_(t,r)=z1,0(t,r)+z1,1(t,r)+z1,2(t,r)+z1,3(t,r)+⋯,Y¯(t,r)=z2,0(t,r)+z2,1(t,r)+z2,2(t,r)+z2,3(t,r)+⋯,G_(t,r)=z3,0(t,r)+z3,1(t,r)+z3,2(t,r)+z3,3(t,r)+⋯,G¯(t,r)=z4,0(t,r)+z4,1(t,r)+z4,2(t,r)+z4,3(t,r)+⋯,
and Figures 7, 8, 9, and 10 show the ℏ-meshes of (∂3/∂t3)Y_10(0,r),  (∂3/∂t3)Y¯10(0,r),  (∂/∂t)G_10(0.1,r), and (∂/∂t)G¯10(0.1,r) respectively. Also, Figures 11, 12, 13, and 14 exhibit their related Contour plots.Figure 7
ℏ-mesh of (∂3/∂t3)Y_10(0,r).Figure 8
ℏ-mesh of (∂3/∂t3)Y¯10(0,r).Figure 9
ℏ-mesh of (∂/∂t)G_10(0.1,r).Figure 10
ℏ-mesh of (∂/∂t)G¯10(0.1,r).Figure 11
Contour plot of(∂3/∂t3)Y_10(0,r).Figure 12
Contour plot of(∂3/∂t3)Y¯10(0,r).Figure 13
Contour plot of(∂/∂t)G_10(0.1,r).Figure 14
Contour plot of(∂/∂t)G¯10(0.1,r).By minimizing the residual error, we find that the best choices for auxiliary parameter to approximate ofY_15(t,r), Y¯15(t,r), G_15(t,r), and G¯15(t,r) are
(66)ℏ=-0.972656,ℏ=-0.97265,ℏ=-0.972969,ℏ=-0.972964,
respectively. Absolute errors for the 15th order approximation by HAM for Y_n(t,r) and Y¯n(t,r), G_n(t,r), G¯n(t,r) are plotted in Figures 15, 16, 17, and 18, respectively.Figure 15
Absolute error ofY_15(t,r) for ℏ=-0.972656.Figure 16
Absolute error ofY¯15(t,r) for ℏ=-0.97265.Figure 17
Absolute error ofG_15(t,r) for ℏ=-0.972969.Figure 18
Absolute error ofG¯15(t,r) for ℏ=-0.972964.Example 9.
Let us consider the first-order system of fuzzy differential equation(67)Y′-G=((-2tr+2t+2r+3)e2t,(2tr-2t-2r+7)e2t),Y+G′=((-rt+t+3r+2)e2t,(rt-t-3r+8)e2t),Y(0)=(2+r,4-r),G(0)=(r,2-r),t∈[0,1].
The exact solution is as follows:
(68)Y(t,r)=((2+r)e2t+(1-r)te2t,(4-r)e2t+(r-1)te2t),G(t,r)=(re2t,(2-r)e2t).
We may replace (67) by the equivalent system
(69)Y_′-G¯=(-2tr+2t+2r+3)e2t,Y¯′-G_=(2tr-2t-2r+7)e2t,Y_+G_′=(-rt+t+3r+2)e2t,Y¯+G¯′=((rt-t-3r+8)e2t),Y_(0,r)=2+r,Y¯(0,r)=4-r,G_(0,r)=r,G¯(0,r)=2-r.
We choose the initial guesses of system (69) as follows:
(70)Y_0=2+r,Y¯0=4-r,G_0=r,G¯0=2-r,
and the auxiliary linear operators
(71)Li[ϕi(t,r;q)]=∂ϕi(t,r;q)∂t,i=1,2,3,4,
with the properties
(72)Li[ci]=0,i=1,2,3,4,
where ci(i=1,2,3,4) are integral constants, and the nonlinear operators,
(73)N1[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂ϕ1(t,r;q)∂t-ϕ4(t,r;q)-(-2tr+2t+2r+3)e2t,N2[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂ϕ2(t,r;q)∂t-ϕ3(t,r;q)-(2tr-2t-2r+7)e2t,N3[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂ϕ3(t,r;q)∂t+ϕ1(t,r;q)-(-rt+t+3r+2)e2t,N4[ϕ1(t,r;q),…,ϕ4(t,r;q)]=∂ϕ4(t,r;q)∂t+ϕ2(t,r;q)-(rt-t-3r+8)e2t.
The mth order deformation equations are
(74)Li[zi,m(t,r)-χmzi,m-1(t,r)]=ℏRi,m(z→1,m-1,…,z→4,m-1),i=1,2,3,4,
with the initial conditions
(75)zi,m(0,r)=0,i=1,2,3,4,
where
(76)R1,m(z→1,m-1,…,z→4,m-1)=∂z1,m-1∂t-z4,m-1-(1-χm)(-2tr+2t+2r+3)e2t,R2,m(z→1,m-1,…,z→4,m-1)=∂z2,m-1∂t-z3,m-1-(1-χm)(2tr-2t-2r+7)e2t,R3,m(z→1,m-1,…,z→4,m-1)=∂z3,m-1∂t+z1,m-1-(1-χm)(-rt+t+3r+2)e2t,R4,m(z→1,m-1,…,z→4,m-1)=∂z4,m-1∂t+z2,m-1-(1-χm)(rt-t-3r+8)e2t.
Therefore, we recursively obtain
(77)z1,1(t,r)=ℏ-e2tℏ+3ℏr2-3e2tℏr2-2ℏt-e2tℏt+ℏrt+e2tℏrt,z1,2(t,r)=ℏ-e2tℏ-5ℏ24+5e2tℏ24+3ℏr2-3e2tℏr2+5ℏ2r2-5e2tℏ2r2-2ℏt-e2tℏt-25ℏ2t4-5e2tℏ2t4+ℏrt+e2tℏrt+11ℏ2rt4+5e2tℏ2rt4-2ℏ2t2+ℏ2rt22,z2,1(t,r)=4ℏ-4e2tℏ-3ℏr2+3e2tℏr2+e2tℏt-ℏrt-e2tℏrt,z2,2(t,r)=4ℏ-4e2tℏ+15ℏ24-15e2tℏ24-3ℏr2+3e2tℏr2-5ℏ2r2+5e2tℏ2r2+e2tℏt-3ℏ2t4+5e2tℏ2t4-ℏrt-e2tℏrt-11ℏ2rt4-5e2tℏ2rt4-ℏ2t2-ℏ2rt22,z3,1(t,r)=3ℏ4-3e2tℏ4+7ℏr4-7e2tℏr4+2ℏt-e2tℏt2+ℏrt+e2tℏrt2,z3,2(t,r)=3ℏ4-3e2tℏ4+ℏ2-e2tℏ2+7ℏr4-7e2tℏr4+11ℏ2r4-11e2tℏ2r4+2ℏt-e2tℏt2+3ℏ2t-e2tℏ2t+ℏrt+e2tℏrt2+5ℏ2rt2+e2tℏ2rt-ℏ2t2+ℏ2rt22,z4,1(t,r)=17ℏ4-17e2tℏ4-7ℏr4+7e2tℏr4+4ht+e2tℏt2-ℏrt-e2tℏrt2,z4,2(t,r)=17ℏ4-17e2tℏ4+13ℏ22-13e2tℏ22-7ℏr4+7e2tℏr4-11ℏ2r4+11e2tℏ2r4+4ℏt+e2tℏt2+8ℏ2t+e2tℏ2t-ℏrt-e2tℏrt2-5ℏ2rt2-e2tℏ2rt-ℏ2rt22,⋮.
Then obtained solutions are as follows:
(78)Y_(t,r)=z1,0(t,r)+z1,1(t,r)+z1,2(t,r)+z1,3(t,r)+⋯,Y¯(t,r)=z2,0(t,r)+z2,1(t,r)+z2,2(t,r)+z2,3(t,r)+⋯,G_(t,r)=z3,0(t,r)+z3,1(t,r)+z3,2(t,r)+z3,3(t,r)+⋯,G¯(t,r)=z4,0(t,r)+z4,1(t,r)+z4,2(t,r)+z4,3(t,r)+⋯.
Figures 19, 20, 21, and 22 show the ℏ-meshes of (∂3/∂t3)Y_10(0,r), (∂3/∂t3)Y¯10(0,r), (∂3/∂t3)G_10(0,r) and (∂3/∂t3)G¯10(0,r), respectively. Also, Figures 23, 24, 25, and 26 show the related Contour plots.Figure 19
ℏ-mesh of (∂3/∂t3)Y_10(0,r).Figure 20
ℏ-mesh of (∂3/∂t3)Y¯10(0,r).Figure 21
ℏ-mesh of (∂3/∂t3)G_10(0,r).Figure 22
ℏ-mesh of (∂3/∂t3)G¯10(0,r).Figure 23
Contour plot of(∂3/∂t3)Y_10(0,r).Figure 24
Contour plot of(∂3/∂t3)Y¯10(0,r).Figure 25
Contour plot of(∂3/∂t3)G_10(0,r).Figure 26
Contour plot of(∂3/∂t3)G¯10(0,r).Best choices for auxiliary parameterℏ to approximate Y_15(t,r),Y¯15(t,r),G_15(t,r), and G¯15(t,r) from minimizing the residual errors are
(79)ℏ=-0.977581,ℏ=-0.968521,ℏ=-1.01756,ℏ=-0.978537,
respectively. Finally, absolute errors for the 15th order approximation by HAM for Y_n(t,r),Y¯n(t,r),G_n(t,r) and G¯n(t,r), are plotted in Figures 27, 28, 29, and 30.Figure 27
Absolute error ofY_15(t,r) for ℏ=-0.977581.Figure 28
Absolute error ofY¯15(t,r) for ℏ=-0.968521.Figure 29
Absolute error ofG_15(t,r) for ℏ=-1.01756.Figure 30
Absolute error ofG¯15(t,r) for ℏ=-0.978537.


## 7. Conclusion

In this paper, homotopy analysis method has been implemented to derive approximate analytical solutions for the system of fuzzy differential equations. Obtained results show that we can control the convergence region of HAM series solution by the auxiliary parameterℏ. The concept of traditional h-curves has been generalized to ℏ-meshes and then contour plot firstly has been introduced in HAM. Convergence theorem and given illustrative examples show the efficiency and accuracy of the HAM.



---

*Source: 407647-2012-08-28.xml*