Article Search
닫기

## Original Article

Split Viewer

International Journal of Fuzzy Logic and Intelligent Systems 2022; 22(1): 23-47

Published online March 25, 2022

https://doi.org/10.5391/IJFIS.2022.22.1.23

© The Korean Institute of Intelligent Systems

## A New Approach to Solving Fuzzy Quadratic Riccati Differential Equations

1Department of Mathematics, Faculty of Science, Al-Balqa Applied University, Salt, Jordan
2Department of Basic Engineering Sciences, College of Engineering, Imam Abdulrahman Bin Faisal University, Dammam, Saudi Arabia
3Department of Applied Science, Ajloun College, Al-Balqa Applied University, Ajloun, Jordan
4Nonlinear Dynamics Research Centre (NDRC), Ajman University, Ajman, UAE

Correspondence to :

Received: June 4, 2021; Revised: September 5, 2021; Accepted: November 1, 2021

This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (http://creativecommons.org/licenses/by-nc/3.0/) which permits unrestricted noncommercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

In this paper, we present in detail the power series solutions to fuzzy quadratic Riccati differential equations (QRDEs) along with a suitable fuzzy constant through an interactive derivative, more specifically, the Hukuhara-strongly generalized differentiability (H-SGD) based on our new technique. This technique is called the Laplace residual power series (LRPS) method, and it mainly depends on a new expansion and the combination of the Laplace transform technique with the residual power series method. To validate the accuracy of our proposed algorithm, numerous examples were examined numerically and graphically, and we compared the results of the optimal homotopy asymptotic (OHA), multiagent neural network (MNN), and fourth-order Runge-Kutta (RK-4) methods with the LRPS method at γ = 1.

Keywords: Fuzzy-valued function, Strongly generalized differentiability, Quadratic Riccati differential equation, Laplace residual power series method, Laplace and inverse transforms

Fuzzy set theory is an important topic in the study and modeling of many real-life problems related to numerous physical phenomena and dynamical processes, such as astronomy, quantum mechanics, chromodynamics, quantum optics, electronic mechanisms, electronic-controlled systems, the modeling of anomalous hydraulic diffusion systems, and population dynamics models [114]. Moreover, several problems in applied mathematics that are associated with biology, medicine, physics, and technology can be modeled using fuzzy differential equations [16,1517].

Fuzzy set theory has been explored since the 1920s; however, the term fuzzy derivative was introduced by Chang and Zadeh [15]. Subsequently, Dubois and Prade [16] proposed the fuzzy differential calculus concept. Many researchers have recently presented other results and notations for fuzzy mapping [1821]. Moreover, Bede and his colleagues [22,23] defined strongly generalized differentiability (SGD) as a fuzzy-valued function (FVF). In general, it is not straightforward to obtain the exact solution for these types of problems because of the difficulties involved; therefore, reliable numerical techniques are needed, including a reliable numerical algorithm for handling fuzzy integral equations of the second kind in the Hilbert spaces [4], a residual power series (RPS) for solving uncertain Riccati differential equations

[24], and a fuzzy Picard method for solving fuzzy quadratic Riccati and fuzzy Painlevl equations [25]. Other examples include a novel extended approach for obtaining numerical solutions to fuzzy conformable fractional differential equations [22], a reproducing kernel algorithm for solving second-order two-point fuzzy boundary value problems and second-order fuzzy Volterra integro-differential equations [3,18,26], a Runge-Kutta method applied to the SI epidemiological model [27], a homotopy analysis method (HAM) for solving fuzzy initial value problems (FIVPs) [28], and bidimensional IVP with interactive fuzzy numbers [14]. Additional solutions are a logistic model of population growth through an interactive derivative [29], an HIV viral dynamics model for individuals under antiretroviral treatment [17], a numerical approach in the Hilbert space to solve a fuzzy Atangana-Baleanu fractional hybrid system [30], a novel technique for solving fuzzy fractional delay differential equations [31], and a homotopy perturbation Sumudu transform method to find the analytical fuzzy solution of nonlinear fuzzy integro-differential equations [32].

The existence, uniqueness, and controllability of fuzzy solutions for such nonlinear fuzzy (fractional) differential equations have been studied and discussed in [3339]. Kwun et al. [3336] and Lee et al. [37] respectively studied the existence and uniqueness of fuzzy (periodic) solutions for the nonlinear fuzzy integro-differential equations on ENn, semi linear fuzzy integro-differential equation with an initial nonlocal condition, fuzzy differential equations in an n-dimension fuzzy vector space (EN)n, fuzzy differential equations with three types of forcing terms using a Hukuhara derivative, and semilinear fuzzy integro-differential equations using an integral contractor. For the controllability, Yoon et al. [38] and Lee et al. [39] studied the controllability of fuzzy solutions for semilinear fuzzy integro-differential equations with non-local conditions and a forcing term with memory, and abstract fuzzy differential equations in the credibility space from the perspective of the Liu process, respectively.

Quadratic Riccati differential equations (QRDEs) form a precise nonlinear model for defining a certain type of physical scheme with applications in the diffusion process, optimal control, optical networks, and artificial intelligence, and usually have the following general form [4043] of

u(t)=A(t)u2(t)+B(t)u(t)+C(t),t0,

subject to

u(0)=λ.

Many numerical and analytical techniques have been employed to obtain QRDEs solutions, including the HAM [40], homotopy perpetuation method (HPM) [41], differential transformation method [42], and Adomian decomposition method [43].

Several authors [24,26,28,4045] have reformulated the QRDE as a fuzzy differential equation along with suitable fuzzy constants under the SGD concept, the fuzzy form of which is as follows [24]:

u^(t)=Au^2(t)+Bu^(t)+C,0tT,

subject to the fuzzy initial condition (FIC):

u^(0)=λ^,

where A, B, and C are constants of t, û(t): [0, T] → ℝ is an unknown fuzzy function of the crisp variable t, and λ̂ ∈ ℝ (the set of all fuzzy numbers). Several techniques have been applied to solve (fractional) differential equations [713,4648], and other techniques have been used to find the solution to an FIVP in Eq. (3), including the HPM [44], Euler method [45], fuzzy Picard technique [25], and RPS method [24].

The main objective of the research described in this study is to extend the application of the LRPS method and provide a power series solution to the FIVP in Eq. (3), as well as study the behavior of the fuzzy solution obtained.

The LRPS method is a new method that combines the Laplace transform technique with the RPS method, which is used to solve the integral and differential equations [713]. The Laplace transform method involves transforming the equation from the space within it into a new space, presenting its solution in the new space, and then returning the solution to the original space and thus obtaining the solution to the equation. However, it is sometimes difficult to offer a solution to an equation in a new space. The RPS method assumes a power series solution to the equation and provides an easy technique for determining the series coefficients. The LRPS method uses the two methods together to provide a series solution to the equation in the new space and uses the principle of the RPS method to determine the series coefficients. However, through a new approach, it is easier and faster than the user in the original technique. In fact, our proposed method is an analytical and efficient method that can be easily applied to solve different types of fuzzy problems.

Two fuzzy numbers are considered non-interactive when the joint possibility distribution involving both is given by the minimum t-norm. Thus, any operation between two fuzzy numbers, aside from the one obtained by the minimum t-norm, must be considered to be interactive. Following this reasoning, when we use the Hukuhara difference (H-difference) for two fuzzy numbers, we assume that they are interactive. Thus, all derivatives obtained by a fuzzy arithmetic operation that differs from the standard operation must be considered interactive. In this study, we find analytical and numerical solutions for FIVPs along with a suitable fuzzy constant through an interactive derivative, more specifically, the Hukuhara-strongly generalized differentiability (H-SGD) based on our new LRPS technique. The solution to the FIVP obtained here using H-SGD is more suitable than the classical derivative because the latter uses an interactive derivative obtained by the H-difference, whereas the standard sum is used to compute from the non-interactive operation, that is, the sum operation is based on the minimum t-norm. This is because the arithmetical operations compatible with the H-difference are not known. Thus, it makes no sense to combine interactive and non-interactive H-SGD. Furthermore, with the H-SGD theory, it is possible to extend and generalize our analysis and results in this study to all fuzzy interactive arithmetic operations, and not only to the H-difference operation in future studies. In addition, in the present paper, we introduce a new expansion to Taylor’s formula, which is required for creating the LRPS method and formulating the FIVP. Moreover, several examples are solved and examined numerically and graphically using this new technique, and the results of the optimal homotopy asymptotic (OHA), multiagent neural network (MNN), and fourth-order Runge-Kutta (RK-4) methods are compared with those of the LRPS method at γ = 1. For the calculations and results, Maple 2018 was used to conduct the computations. The results show that the proposed algorithm is efficient, accurate, and robust.

This paper consists of five sections. In Section 2, we study some basic concepts of fuzzy theory, establish a new expansion to Taylor’s formula, and formulate the FIVP. In Section 3, the LRPS algorithm is created and successfully tested to obtain analytical and numerical solutions of the nonlinear FIVP. In Section 4, we test the accuracy of our proposed algorithm (LRPS) by solving three important examples numerically and graphically, and comparing the results of the OHA, MNN, and RK-4 methods with the LRPS method at γ = 1. Finally, the conclusions are presented in the last section..

### 2. Basic Concepts, A New Formula and Drafting

This section reviews some theories and concepts necessary for our study and is divided into three subsections: the first relates to fuzzy calculus, the second relates to a new expansion needed to construct the LRPS solution, and we formulate the FIVP in third.

### 2.1 Basic Concepts and Notations on Fuzzy Theory

In this subsection, we study the most important definitions and notations related to fuzzy analysis and fuzzy calculus [3,4,12, 15,16,1823,4863] used throughout this paper. The provided ℝ denotes the set of all the fuzzy numbers.

Definition 2.1 [19,48,59,60]

Letting v: ℝ → [0, 1], v is then called a fuzzy number if the following conditions hold:

(i) v is normal (i.e., ∃ t0 ∈ ℝ such that v(t0) = 1),

(ii) v is convex (i.e., ∀ t, s ∈ ℝ, and 0 ≤ λ ≤ 1, and it holds that

v(λt+(1-λ)s)min{v(t),v(s)}),

(iii) v is upper semi-continuous (i.e., {t ∈ ℝ: v(t) ≥ r} is closed for all r ∈ ℝ), and

(iv) v is the bounded support (i.e., [v]0=supp(v)¯ is a compact subset of ℝ).

Definition 2.2 [14,17,19,20,53]

Let v ∈ ℝ and γ ∈ [0, 1], the γ-cut (or level) of v is then the crisp set

[v]γ={t:v(t)γ}.

Therefore, if v ∈ ℝℱ, then γ-levels are closed intervals in ℝ and are denoted by

[v]γ=[v1(γ),v2(γ)]=[v1γ,v2γ];

v1γ=v1(γ)=min{r:r[v]γ},v2γ=v2(γ)=max{r:r[v]γ}.

In addition, if v ∈ ℝ, the length of the γ-level set of v is defined as Len[v]γ = v2γv1γ for all γ ∈ [0, 1]. Note that if γ = 0, then Len[v]0 = Diam(v), where Diam(v) is the diameter of a number fuzzy v [14,17]. However, there are several common forms of a fuzzy number v, one of which is the triangular form defined by an ordered triple v = (μ1, μ2, μ3) ∈ ℝ3 with μ1< μ2< μ3 and whose γ-level is

[v]γ=[μ1+(μ2-μ1)γ,μ3-(μ3-μ2)γ],γ[0,1].

Note that a fuzzy set v of ℝ is said to be a fuzzy number if [v]γ is a non-empty, bounded, and closed interval of ℝ for every γ ∈ [0, 1]. The fuzzy numbers satisfy the following property: vω if and only if [v]γ ⊆ [ω]γ for all γ ∈ [0, 1] [14,27,29,60].

Some conditions must be satisfied by two functions: v1γ, v2γ: [0, 1] → ℝ, such that [v1γ, v2γ] can be a parametric form of a fuzzy number v for all γ ∈ [0, 1]. These conditions are presented in the following theorem.

Theorem 2.1 [57]

Suppose that v1γ, v2γ: [0, 1] → ℝ satisfy the following three conditions:

(i) v1γ is a bounded monotonic non-decreasing left-continuous function on γ ∈ (0, 1],

(ii) v2γ is a bounded monotonic non-increasing left-continuous function on γ ∈ (0, 1],

(iii) v1γ and v2γ are right-continuous at γ = 0, and

(iv) v11v21.

Subsequently, v: ℝ → [0, 1] defined by v(t) = sup{γ: v1γtv2γ} is a fuzzy number with a γ-level set: [v]γ = [v1γ, v2γ]. Conversely, if v is a fuzzy number with [v]γ = [v1γ, v2γ], then the functions v1γ, v2γ: [0, 1] → ℝ satisfy the conditions (i)–(iv).

Definition 2.3 [1016,1821,53]

Let v = [v1γ, v2γ], ω = [ω1γ, ω2γ] ∈ ℝ, μ ∈ ℝ \ {0}, and γ ∈ [0, 1]. If v and ω are non-interactive (i.e., the joint possibility distribution is given by the minimum t-norm), then the following operations for all γ ∈ [0, 1] in the levels are given by the following:

(i) standard sum,

[v+ω]γ=[v1γ+ω1γ,v2γ+ω2γ];

(ii) standard differences,

[v-ω]γ=[v1γ-ω2γ,v2γ-ω1γ];

(iii) scalar multiplication,

[μv]γ=μ[v]γ=[min{μv1γ,μv2γ},max{μv1γ,μv2γ}];

(iv) multiplication,

[v.ω]γ=[min{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ},   max{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ}],

(v) quotient,

[min{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ},max{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ}],and

(vi) equality, the two fuzzy numbers v and ω are equal if [v]γ = [ω]γ, that is, v1γ = ω1γ and v2γ = ω2γ.

Definition 2.4 [1021,53]

Letting v = [v1γ, v2γ] and ω = [ω1γ, ω2γ] be fuzzy numbers, and γ ∈ [0, 1], then z is called the H-difference of v and ω and is denoted by vω if there exists an element z ∈ ℝ such that v = ω + z for all γ ∈ [0, 1], and is given in levels by the following:

[vω]γ=[v]γ-[ω]γ=[v1γ-ω1γ,v2γ-ω2γ].

This definition was modified by Bede and Stefanini [23] to the generalized H-difference (gH-difference), as follows:

vgHω=zv=ω+zor ω=v+(-1)z,

with γ-levels:

[vgHω]γ=[min{v1γ-ω1γ,v2γ-ω2γ},max{v1γ-ω1γ,v2γ-ω2γ}].

Very recently, Wasques et al. [17,27,29,59] generalized the Hukuhara and H-differences and used interactivities for all arithmetic operations and not only for the difference operation. We summarize this work as follows.

In general, a fuzzy subset v of ℝn is described by its membership function ρv: ℝn → [0, 1], where ρv(t) means the degree to which t belongs to v. The r-levels of fuzzy subset v are classical subsets defined as follows:

[v]γ={tn:ρv(t)γ},γ(0,1],

and

[v]0={tn:ρv(t)>0}.

Note that the fuzzy subset v of ℝ is a fuzzy number if its γ-levels are closed and nonempty intervals of ℝ and the support of v, supp(v) = {t ∈ ℝ: ρv(t) > 0}, is limited. The family of fuzzy subsets of ℝn with nonempty compact and convex γ-levels is denoted by n, whereas the family of fuzzy numbers is denoted by ℝ.

Let v and ω ∈ ℝ and ℐ ∈ F(ℝ2). Then, the fuzzy relation ℐ is a joint possibility distribution of v and ω if

maxtρ(s,t)=ρv(t)and maxsρ(s,t)=ρω(s),

for all t and s ∈ ℝ.

maxtρ(s,t)=ρv(t)and maxsρ(s,t)=ρω(s),

for all t and s ∈ ℝ.

Here, v and ω are the marginal possibility distributions of ℐ. The fuzzy numbers v and ω are said to be non-interactive if and only if their joint possibility distribution ℐ is given by

ρ(s,t)=min{ρv(t),ρω(s)},

for all t and s ∈ ℝ.

Otherwise, the fuzzy numbers are said to be interactive.

Definition 2.5 [17]

The Pompeiu-Hausdorff distance d:n×n+{0} is defined by

d(v,ω)=supγ[0,1]dH([v]γ,[ω]γ),

where dH is the Pompeiu-Hausdorff distance for the compact subsets of ℝn. If v and ω are fuzzy numbers, that is, v and ω ∈ ℝ, then Eq. (20) becomes

d(v,ω)=supγ[0,1]max{v1γ-ω1γ,v2γ-ω2γ},

and gives the Hausdorff distance between [v]γ and [ω]γ.

It is provided that (ℝ, d) is a complete metric space, which processes the following properties [53]:

(i) d(v + z, ω + z) = d(v, ω) for all v, z, ω ∈ ℝ,

(ii) d(μv, μω) = |μ|d(v, ω) for all v, ω ∈ ℝ and μ ∈ ℝ \ {0}, and

(iii) d(v + z, ω + e) ≤ d(v, ω) + d(z, e) for all v, z, ω, e ∈ ℝ.

Note that the metric space d in ℝ has many nice and interesting properties for an addition operation and scalar multiplication, some of which can be found in [19,5356,61].

Definition 2.6 [58,59]

Let [a, b] ⊆ ℝ be a vector space and ℝ be a set of fuzzy numbers. Then, the function û: [a, b] → ℝ is called an FVF on [a, b]. Corresponding to such a function û and γ ∈ [0, 1], we denote two real-valued functions û1γ(t) and û2γ(t) for all t ∈ [a, b] as

[u^(t)]γ=[u^1γ(t),u^2γ(t)],

which is called the γ-level representation of an FVF.

Definition 2.7 [20,53]

The complete metric structure on ℝ is given by the Hausdorff distance mapping d: × → ℝ+ ∪ {0} such that

d(v,ω)=supt[0,1]max{v1(t)-ω1(t),v2(t)-ω2(t)}.
Definition 2.8 [53,55]

Let û1 and û2: [a, b] → ℝ be two FVFs. Then, the uniform distance between û1 and û2 is defined as

d*(u^1(t),u^2(t))=supt[a,b]d(u^1(t),u^2(t)).
Definition 2.9 [53,54,59]

Let û: [a, b] → ℝ be an FVF and t0 ∈ [a, b]. If ∀ ɛ > 0, ∃ δ > 0 such that d*(û (t), L) < ɛ, whenever |tt0| < δ, then we state that L ∈ ℝ is limit of û at t0, which is denoted by: limtt0û(t) = L. In addition, the FVF û is said to be continuous for t0 ∈ [a, b] if limtt0û(t) = û(t0). Indeed, û is continuous on [a, b] if it is continuous ∀ t ∈ [a, b]. In addition, if û: [a, b] → ℝ is a fuzzy continuous function, then we have [5356] the following:

d(u^(t),0)=supγ[0,1]max{u^1γ(t),u^2γ(t)},t[a,b].
Definition 2.10 [20]

Let û: [a, b] → ℝ be an FVF. Then, for a fixed point t0 ∈ [a, b], û is called an SGD at t0 if there exists an element û′(t0) ∈ ℝ such that either of the following hold:

(i) the H-differences û (t0 +ɛ) ⊝ û (t0) and û (t0) ⊝ û (t0ɛ) exist, and

u^(t0)=D(1)u^(t0)=limɛ0+u^(t0+ɛ)u^(t0)ɛ=limɛ0+u^(t0)u^(t0-ɛ)ɛ,

or all ɛ > 0 are sufficiently near 0 and the limits are in a metric d.

(ii) The H-differences û (t0) ⊝ û (t0+ɛ) and û (t0ɛ) ⊝ û (t0) exist, and

u^(t0)=D(2)u^(t0)=limɛ0+u^(t0)u^(t0+ɛ)-ɛ=limɛ0+u^(t0-ɛ)u^(t0)-ɛ,

for all ɛ > 0 sufficiently near 0 and the limits are in a metric d.

Remark 2.1 [20]

(i) If û is differentiable for all t0 ∈ (a, b), then û is an SGD on (a, b).

(ii) The limits in Eqs. (26) and (27)) are considered in (ℝ, d), and at the edge-point of (a, b), we consider one direction derivative and recalling that as in Eq. (15).

Theorem 2.2 [22,53]

Let û: [a, b] → ℝ be an FVF, and let [û(t)]γ = [û1γ (t), û2γ (t)], ∀t ∈ [0, 1], and D(1)û(t) or D(2)û(t) exists. Then,

(i) if û is (1)-differentiable, then û1γ (t) and û2γ (t) are differentiable functions and

[u^(t)]γ=[u^1γ(t),u^2γ(t)],

and (ii) if û is (2)-differentiable, then û1γ (t) and û2γ (t) are differentiable functions and

[u^(t)]γ=[u^2γ(t),u^1γ(t)].
Theorem 2.3 [20,22]

Let û: [a, b] → ℝ be an FVF such that û is (1)-differentiable or û is (1)-differentiable at t0 ∈ (a, b). Subsequently, û is continuous at t0.

Definition 2.11 [17,62]

A strongly measurable and limited integrable fuzzy function is known as an integrable function. The fuzzy integral of an FVF û: [a, b] → ℝ with [û(t)]γ = [û1γ (t), û2γ (t)] is defined as follows:

[Ju^(t)]γ=[abu^(t)dt]γ=ab[u^(t)]γdt=ab[u^1γ(t),u^2γ(t)]dt,

for all γ ∈ [0, 1]. Provided Eq. (30), a fuzzy number is defined, and the integral exists in ℝ.

Similarly, the fuzzy integral of an FVF û: [a, b] → ℝ with [û(t)]γ = [û1γ (t), û2γ (t)] is defined with reference point t0 as

[Jt0u^(t)]γ=[t0tu^(z)dz]γ=t0t[u^(z)]γdz=t0t[u^1γ(z),u^2γ(z)]dz,

for all γ ∈ [0, 1]. Provided Eq. (31) define a fuzzy number and the integral exists in ℝ.

Theorem 2.4 [20]

Let û: [a, b] → ℝ be a FVF. Thus, the following are working:

(i) If û(t) is (1)-differentiable, then

D(1)(Jt0u^)(t)=u^(t),

and

Jt0(D(1)u^(t))=u^(t)u^(t0).

(ii) If f is û(t) is (2)-differentiable, then.

D(2)(Jt0u^)(t)-u^(t)=0,

and

u^(t0)=u^(t)-Jt0(D(2)u^(t)).

### 2.2 A New Form to Taylor’s Formula

In this subsection, we present a new form of Taylor’s formula, which we will mainly use in our proposed method for solving FIVPs.

Theorem 2.5

Let u(t) be a piecewise continuous function on [0,∞) and be of exponential order. Suppose that the function U(s) = ℒ[u(t)] (Laplace transform of u(t)) can be expanded as follows:

U(s)=n=0cnsn+1,s>0.

Then cn = u(n)(0).

Proof

We assume that U(s) is a function that can be represented by the expansion in Eq. (34). First, note that if we multiply Eq. (34) by s and taking the limit for both sides as s → ∞, each term in the expansion after the first vanishes; thus, we obtain c0 = lims→∞ sU(s) = u(0), and Eq. (44) becomes

U(s)=u(0)s+n=1cnsn+1,s>0.

For the other aspect as well, multiplying Eq. (35) by s2 leads to the following expansion of c1:

c1=s2U(s)-su(0)-c2s-c3s2-c4s3-,s>0.

It is clear that

c1=lims(s2U(s)-su(0))=lims(s(sU(s)-u(0))=limss([u(t)](s))=u(0).

If we multiply Eq. (35) by s3 and taking the limit as s → ∞, we obtain

c2=lims(s3U(s)-s2u(0)-su(0))=limss(s2U(s)-su(0)-u(0))=limss([u(t)](s))=u(0).

Now, we can see the patterns and find the general formula for cn. However, if we continue to multiply Eq. (35) by sn+1 and compute the limit of the expansion of cn as s → ∞, we can obtain cn = u(n)(0), n = 0, 1, 2, · · · . This completes the proof of Theorem 2.5.

### Remark 2.2

The inverse Laplace transform of the expansion in Theorem 2.5 has the following form:

u(t)=n=0u(n)(0)n!tn,t0,

which is the classical Taylor’s formula.

The following theorem explains and determines the conditions for convergence of the new form of Taylor’s formula, which is presented in Theorem 2.5, as mentioned above.

### Theorem 2.6

Let U(s) = ℒ[u(t)](s). If |sℒ[u(n+1)(t)]|M, for 0 < sd, where M and d are positive numbers, then the remaining Rn(s) of the new form of Taylor’s formula in Eq. (34) satisfies the following inequality:

Rn(s)Msn+2,0<sd.
Proof

First, suppose that ℒ[u(j)(t)](s) is defined as 0 < sd for j = 0, 1, 2, · · ·, n + 1. As given, assume also that

s[u(n+1)(t)]M,0<sd.

From the definition of the reminder Rn(s)=U(s)-j=0nu(j)(0)sj+1, one can obtain the following:

sn+2Rn(s)=sn+2U(s)-j=0nsn+1-j(u(j))(0)=s(sn+1U(s)-j=0nsn-j(u(j))(0))=s[u(n+1)(t)].

It follows from Eqs. (41) and (42) that |sn+2Rn(s)|M. Hence,

-Msn+2Rn(s)M,0<sd.

Thus, by reformulating Eq. (43), we can find the inequality Rn(s)Msn+2, which completes the proof.

### 2.3 Drafting of FIVP

In this subsection, we state the conditions and assumptions required to formulate the FIVP, as in Eq. (3). First, we point out that the FIVP is expressed in Eq. (3) yields a unique solution [25,45,63]. Second, we assume the γ-level depiction of û′(t), û(t), û2(t) and û(0) as [u^1γ(t),u^2γ(t)], [û1γ (t), û2γ (t)], [u^1γ2(t),u^2γ2(t)] and [û0,1γ (t), û0,2γ (t)], respectively. Thus, the FIVP, as in Eq. (3) can be reformulated as

[u^(t)]γ=A[u^2(t)]γ+B[u^(t)]γ+C,t>0,

with the FIC

[u^(0)]γ=[u^0]γ.

### Definition 2.12 [1]

Let û: [a, b] → ℝ, where D(1)û(t) or D(2)û(t) exist. If û and D(1)û(t) satisfy FIVP in Eq. (3), we state that û is (1)-solution for the FIVP in Eq. (3). Otherwise, if û and D(2)û(t) satisfy FIVP, as in Eq. (3), we state that û is (2)-solution for the FIVP, as in Eq. (3).

Finally, to construct the LRPS solution for the target equation in a γ-level depiction, we alter the FIVP, as in Eqs. (44) and (45), into crisp ordinary differential equations (ODEs) systems and assume the following two cases to create the fuzzy solution û(t) for this FIVP regarding to the type of differentiability, where û(t) is either (1)-differentiable or (2)-differentiable.

Case 1

If û(t) is (1)-differentiable, then the IVP in Eq. (40)

u1γ(t)=Au^1γ2(t)+Bu^1γ(t)+C,u2γ(t)=Au^2γ2(t)+Bu^2γ(t)+C,

with FICs,

u^1γ(0)=u^0,1γ,u^2γ(0)=u^0.2γ.

We then have the following procedures:

(i) Construct the LRPS solution for the system in Eqs. (46) and (47).

(ii) Ensure that the solutions [û1γ (t), û2γ (t)] and [u^1γ(t),u^2γ(t)] are valid γ-level sets, ∀ γ ∈ [0, 1].

(iii) Then, [û1γ (t), û2γ (t)] are the (1)-solutions û in a γ-level depiction.

Case 2

If û(t) is (2)-differentiable, then the FIVP, as in Eqs. (44) and (45), can be transformed into the following ODEs system:

u^1γ(t)=Au^2γ2(t)+Bu^2γ(t)+C,u^2γ(t)=Au^1γ2(t)+Bu^1γ(t)+C,

with FICs,

u^1γ(0)=u^0,1γ,u^2γ(0)=u^0.2γ.

In addition, we have the following procedures:

(i) Construct the LRPS solution for systems (48) and (49).

(ii) Ensure that the solutions [û1γ (t), û2γ (t)] and [u^1γ(t),u^2γ(t)] are valid γ-level sets, ∀ γ ∈ [0, 1].

(iii) Then, [û2γ (t), û1γ (t)] are the (2)-solutions û in a γ-level depiction.

In this section, we present the (1)-solution for the FIVP, as in Eqs. (44) and (45), using the LRPS method. In addition, the same process can be pursued each time û(t) is (2)-differentiable to create the (2)-solution for the FIVP, as in Eqs. (44) and (45). To achieve this, we assumed that û(t) is (1)-differentiable. First, the technique involves applying a Laplace transformation to both sides of the system in Eq. (46); hence,

[u^1γ(t)]=A[u^1γ2(t)]+B[u^1γ(t)]+[C],[u^2γ(t)=A[u^2γ2(t)+B[u^2γ(t)]+[C].

Using the properties of the Laplace transform and the initial conditions in Eq. (47), we can rewrite the system in Eq. (50) into the following algebraic system:

U^1γ(s)=u^0,1γs+As[(-1[U^1γ(s)])2]+BsU^1γ(s)+Cs2,U^2γ(s)=u^0,2γs+As[(-1[U^2γ(s)])2]+BsU^2γ(s)+Cs2,

where Û1γ (s) = ℒ[û1γ (t)], and Û2γ (s) = ℒ[û2γ (t)].

Next, the LRPS technique represents the solution as an infinite series, and the series solution to Eq. (51) is given by the following expansions:

U^1γ(s)=n=0ann!sn+1,U^2γ(s)=n=0bnn!sn+1.

Because

limssU^1γ(s)=u^1γ(0)=u^0,1γ,

and

limssU^2γ(s)=u^2γ(0)=u^0,2γ,

this leads to a0 = û0,1γ and b0 = û0,2γ.

Therefore, the expansions in Eq. (52) can be written as

U^1γ(s)=u^0,1γs+n=1ann!sn+1,U^2γ(s)=u^0,2γs+n=1bnn!sn+1.

Consequently, the kth-truncated series of Û1γ (s) and Û2γ (s) can be given by

U^k,1γ(s)=u^0,1γs+n=1kann!sn+1,U^k,2γ(s)=u^0,2γs+n=1kbnn!sn+1.

Prior to applying the LRPS technique to determine the values of the coefficients an and bn in the expansions in Eq. (54), we must define the Laplace residual function for Eq. (51) as

LRes1γ(s)=U^1γ(s)-u^0,1γs-As[(-1[U^1γ(s)])2]-BsU^1γ(s)-Cs2s,LRes2γ(s)=U^2γ(s)-u^0,2γs-As[(-1[U^2γ(s)])2]-BsU^2γ(s)-Cs2,

and the kth-Laplace residual function as follows:

LResk,1γ(s)=U^k,1γ(s)-u^0,1γs-As[(-1[U^k,1γ(s)])2]-BsU^k,1γ(s)-Cs2,LResk,2γ(s)=U^k,2γ(s)-u^0,2γs-As[(-1[U^k,2γ(s)])2]-BsU^k,2γ(s)-Cs2.

It is clear that limjLResj,iγ(s)=LRes,iγ(s)=LResiγ(s) for each s > 0 and i = 1, 2.

Because limssLResiγ(s)=0, it is clear that

limssLResk,iγ(s)=0,

for i = 1, 2 and j = 1, 2, 3, · · ·, k, and thus lims(sj+1LRes,iγ(s))=lims(sj+1LResj,iγ(s))=0,         i=1,2,j=1,2,3,,k.

For the other aspects as well, to determine the value of the first unknown coefficients a1and b1 in Eq. (54), we should substitute Û1,1γ (s) = û0,1γ/s + a1/s2 and Û1,2γ (s) = û0,2γ/s + b1/s2 into the kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), at k = 1 in Eq. (56) to obtain the following:

LRes1,1γ(s)=   U^1,1γ(s)-u^0,1γs-As[(-1[U^1,1γ(s)])2]-BsU^1,1γ(s)-Cs2=   u^0,1γs+a1s2-u^0,1γs-As[(-1[u^0,1γs+a1s2])2]-Bs(u^0,1γs+a1s2)-Cs2=   a1s2-As(u^0,1γ2s+2u^0,1γa1s2+2a12s3)-Bs(u^0,1γs+a1s2)-Cs2,LRes1,2γ(s)=   U^1,2γ(s)-u^0,2γs-As[(-1[U^1,2γ(s)])2]-BsU^1,2γ(s)-Cs2=   u^0,2γs+b1s2-u^0,2γs-As[(-1[u^0,2γs+b1s2])2]-Bs(u^0,2γs+b1s2)-Cs2=   b1s2-As(u^0,2γ2s+2u^0,2γa1s2+2b12s3)-Bs(u^0,2γs+b1s2)-Cs2.

Now, multiplying both sides of Eq. (56) when k = 1 by s2 gives

s2LRes1,1γ(s)=a1-A(u^0,1γ2+2u^0,1γa1s+2a12s2)-B(u^0,1γ+a1s)-C,s2LRes1,2γ(s)=b1-A(u^0,2γ2+2u^0,2γa1s+2b12s2)-B(u^0,2γ+b1s)-C.

Computing the limits of both sides of Eq. (59) as s→∞gives the following expressions:

limss2LRes1,1γ(s)=a1-Au^0,1γ2-Bu^0,1γ-C,limss2LRes1,2γ(s)=b1-Au^0,2γ2-Bu^0,2γ-C.

Considering j = 1 in Eq. (57) using the expressions in Eq. (60) and solving the resulting algebraic system for a1 and b1, it is easy to obtain a1=Au^0,1γ2+Bu^0,1γ+C and b1=Au^0,2γ2+Bu^0,2γ+C. Hence, we obtain the following:

U^1,1γ(s)=u^0,1γs+(Au^0,1γ2+Bu^0,1γ+C)s2,U^1,2γ(s)=u^0,2γs+(Au^0,2γ2+Bu^0,2γ+C)s2.

Similarly, to determine the values of the second unknown coefficients a2 and b2, we substitute U^2,1γ(s)=a0/s+(Aa02+Ba0+C)/s2+2a2/s3 and U^2,2γ(s)=b0/s+(Ab02+Bb0+C)/s2+2b2/s3 into LRes2,1γ (s) and LRes2,2γ (s) in Eq. (56) to obtain the following discretized form:

LRes2,1γ(s)=   a0s+Aa02+Ba0+Cs2+2a2s3-a0s-As[(-1[a0s+Aa02+Ba0+Cs2+2a2s3])2]-Bs(a0s+Aa02+Ba0+Cs2+2a2s3)-Cs2,LRes2,2γ(s)=   b0s+Ab02+Bb0+Cs2+2b2s3-b0s-As[(-1[b0s+Ab02+Bb0+Cs2+2b2s3])2]-Bs(b0s+Ab02+Bb0+Cs2+2b2s3)-Cs2.

By multiplying both sides of Eq. (56) for k = 2 by s3, we obtain

s3LRes2,1γ(s)=   2a2-A(2(Aa03+Ba02+Ca0)+4a0a2s+2(Aa02+Ba0+C)2s+12(Aa02a2+Ba0a2+Ca2)s2+24a22s3)-B(Aa02+Ba0+C+2a2s),s3LRes2,2γ(s)=   2b2-A(2(Ab03+Bb02+Cb0)+4b0b2s+2(Ab02+Bb0+C)2s+12(Ab02b2+Bb0b2+Cb2)s2+24b22s3)-B(Ab02+Bb0+C+2b2s).

Using the fact that limssj+1LResj,iγ(s)=0 in Eq. (57) for j = 2, i = 1, 2, we obtain the following:

limss3LRes2,1γ(s)=2a2-A(2(Aa03+Ba02+Ca0))   -B(Aa02+Ba0+C)=2a2-2Aa0(Aa02+Ba0+C)   -B(Aa02+Ba0+C)=0,limss3LRes2,2γ(s)=2b2-A(2(Ab03+Bb02+Cb0))   -B(Ab02+Bb0+C)=2b2-2Ab0(Ab02+Bb0+C)   -B(Ab02+Bb0+C)=0.

By solving the resultant algebraic system in Eq. (64) for a2 and b2, we obtain a2=Aa0a1+12Ba1 and b2=Ab0b1+12Bb1. Thus,

U^1,1γ(s)=a0s+(Aa02+Ba0+C)s2+2(Aa0a1+12Ba1)s3,U^1,2γ(s)=b0s+(Ab02+Bb0+C)s2+2(Ab0b1+12Bb1)s3.

For k = 3, if we substitute Û3,1γ (s) and Û3,2γ (s) into the Laplace residual functions LRes3,1γ (s) and LRes3,2γ (s) of Eq. (65) and the utilize the facts limss4LRes3,1γ(s)=0 and lims→∞ s4LRes3,2γ (s) = 0, then the third coefficients a3 and b3 are given by the following:

a3=13A(2a0a2+a12)+13Ba2,b3=13A(2b0b2+b12)+13Bb2.

Therefore, the third LRPS estimation can also be given.

By continuing the same process until random coefficients order k = n, and using the fact that

limssn+1LResn,1γ(s)=limssn+1LResn,2γ(s)=0,

subsequently an and bn (unknown coefficients) can be achieved. Nevertheless, more precise solutions can be obtained using more iterations. Likewise, if û(t) is (2)-differentiable, the (2)-solution for the FIVP, as in Eqs. (44) and (45), can be obtained.

Finally, the exact or approximate series solution obtained must be transformed into the original space by using the inverse Laplace transform to obtain the series solution to the original FIVP, as in Eqs. (44) and (45).

To validate the accuracy of our proposed algorithm, numerous examples are examined numerically in this section. For calculations and results, Maple 2018 is used to perform the computations.

### Example 4.1

Consider the following FIVP:

u^(t)=2u^(t)-u^2(t)+1,t>0,

with FIC,

[u^(0)]γ=[γ-1,1-γ],γ[0,1].

For γ = 1, the exact solution of the FIVP, as in Eqs. (67) and (68), can be located as follows:

u(t)=1+2tanh (2t+12log(2-12+1)).

Using Definition 2.10, the FIVP, as in Eqs. (67) and (68), can be represented by the set of ODEs corresponding to their parametric forms as

u^1γ(t)=2u^1γ(t)-u^1γ2(t)+1,u^2γ(t)=2u^2γ(t)-u^2γ2(t)+1,

subject to

u^1γ(0)=γ-1,u^2γ(0)=1-γ.

Applying the LRPS method as below, taking the Laplace transform of both sides of Eq. (70) gives

s[u^1γ(t)]-u^1γ(0)=2[u^1γ(t)]-[u^1γ2(t)]+[1],s[u^2γ(t)]-u^2γ(0)=2[u^2γ(t)]-[u^2γ2(t)]+[1],

and Eq. (71) then leads to

U^1γ(s)=γ-1s+2sU^1γ(s)-1s[(-1[U^1γ(s)])2]+1s2,U^2γ(s)=1-γs+2sU^2γ(s)-1s[(-1[U^2γ(s)])2]+1s2,

where Û1γ (s) = ℒ[û1γ (t)] and Û2γ (s) = ℒ[û2γ (t)].

Now, we assume that an infinite series solution of the algebraic system in Eq. (73) has the following form:

U^1γ(s)=n=0ann!sn+1,U^2γ(s)=n=0bnn!sn+1

In addition, the kth-truncated series of Eq. (74) is expressed as follows:

U^k,1γ(s)=u^0,1γs+n=1kann!sn+1,U^k,2γ(s)=u^0,2γs+n=1kbnn!sn+1.

Depending on the initial data, limssU^1γ(s)=u^1γ(0)=γ-1 and limssU^1γ(s)=u^2γ(0)=1-γ. Therefore, the series solution to Eq. (73) becomes

U^k,1γ(s)=γ-1s+n=1kann!sn+1,U^k,2γ(s)=1-γs+n=1kbnn!sn+1.

The kth Laplace residual functions, LResk,1γ (s) and LResk,2γ (s) of (73), are as follows:

LResk,1γ(s)=U^k,1γ(s)-γ-1s-2sU^k,1γ(s)+1s[(-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-1-γs-2sU^k,2γ(s)+1s[(-1[U^k,2γ(s)])2]-1s2.

To determine the first unknown coefficients, a1 and b1, in the expansion (76), we substitute the first-truncated series Û1,1γ (s) = (γ – 1)/s + a1/s2 and Û1,2γ (s) = (1 – γ)/s + b1/s2 into the first Laplace residual functions, LRes1,1γ (s) and LRes1,2γ (s) of Eq. (77), and obtain

LRes1,1γ(s)=   a1s2-2s(γ-1s+a1s2)+1s((γ-1)2s+2(γ-1)a1s2+2a12s3)-1s2,LRes1,2γ(s)=   b1s2-2s(1-γs+b1s2)+1s((1-γ)2s+2(1-γ)b1s2+2b12s3)-1s2.

Now, by remultiplying both sides of Eq. (77) by sk+1, for k = 1, we obtain

s2LRes1,1γ(s)=   a1-2(γ-1+a1s)+((γ-1)2+2(γ-1)a1s+2a12s2)-1,s2LRes1,2γ(s)=   b1-2(1-γ+b1s)+((1-γ)2+2(1-γ)a1s+2b12s2)-1.

Using the fact that limssk+1LResk(s)=0 for k = 1 in Eq. (79) yields a1 = −(γ – 1)2 + 2(γ – 1) + 1 and b1 = −(1 – γ)2 +2(1–γ)+1. Hence, the first-approximate LRPS solution to Eq. (73) can be expressed as

U^1,1γ(s)=(γ-1)s+(-2+4γ-γ2)s2,U^1,2γ(s)=(1-γ)s+(2-γ2)s2.

Similarly, to determine the values of the second unknown coefficients a2 and b2, we substitute the second-truncated series Û2,1γ (s) = (γ – 1)/s + (−2 + 4γγ2)/s2 + 2a2/s3, and Û2,2γ (s) = (1–γ)/s+(2–γ2)/s2 +2b2/s3 into the second-Laplace residual functions LRes2,1γ (s) and LRes2,2γ (s) of Eq. (77) to obtain the following discretized form:

LRes2,1γ(s)=   (-2+2γ-γ2)s2+2a2s3-2s(γ-1s+(-2+4γ-γ2)s2+2a2s3)+1s(24a22s5+(4γ-4s3+-12γ2+48γ-24s4)a2+γ2-2γ+1s+-2γ3+10γ2-12γ+4s2+2γ4-16γ3+40γ2-32γ+8s3)-1s2,LRes2,2γ(s)=   (2-γ2)s2+2b2s3-2s(1-γs+(2-γ2)s2+2b2s3)+1s(24b22s5+(-4γ+4s3+-12γ2+24s4)b2+γ2-2γ+1s+2γ3-2γ2-4γ+4s2+2γ4-8γ2-32γ+8s3)-1s2.

Now, re-multiplying s3 by LRes2,1γ (s) and LRes2,2γ (s) in Eq. (81) yields.

s3LRes2,1γ(s)=   24a22s3+(2+4γ-8s+-12γ2+48γ-24s2)a2-2γ3+12γ2-20γ+8+2γ4-16γ3+40γ2-32γ+8s,s3LRes2,2γ(s)=   24b22s3+(2-4γs+-12γ2+24s2)b2+2γ3-4γ+2γ4-8γ2+8s.

According to the fact: limssk+1LResk(s)=0 for k = 2, then we have the following algebraic equations:

limss3LRes2,1γ(s)=2a2-2γ3+12γ2-20γ+8=0,limss3LRes2,2γ(s)=2b2+2γ3-4γ=0.

Solving Eq. (83) for a2 and b2, the following discretized forms are obtained:

a2=γ3-6γ2+10γ-4,b2=-γ3+2γ.

Hence, the second-approximate LRPS solution to Eq. (73) is expressed as

U^2,1γ(s)=(γ-1)s+(-2+4γ-γ2)s2+2(γ3-6γ2+10γ-4)s3,U^2,2γ(s)=(1-γ)s+(2-γ2)s2+2(-γ3+2γ)s3.

The kth-LRPS approximation can be obtained sequentially by using the same procedure. For further documentation, we will report the first two coefficients of ak and bk and compute the new component forms for future use in the succeeding approximation to obtain the following:

a1=-2+4γ-γ2,a2=γ3-6γ2+10γ-4,a3=13(-20+64γ-64γ2+24γ3-3γ4),a4=13(-32+128γ-180γ2+110γ3-30γ4+3γ5),a5=115(-256+1232γ-2228γ2+1920γ3-840γ4+180γ5-15γ6),a6=145(-1232+6920γ-15288γ2+17108γ3-10500γ4+3570γ5-630γ6+45γ7),a7=1315(-13840+88832γ-231872γ2+319872γ3-255024γ4+120960γ5-33600γ6+5040γ7-315γ8),

and

b1=2-γ2,b2=-γ3+2γ,b3=13(-4+8γ2-3γ4),b4=13(-8γ+10γ3-3γ5),b5=115(16-68γ2+60γ4-15γ6),b6=145(136γ-308γ3+210γ5-45γ7),b7=1315(-272+1984γ2-3024γ4+1680γ6-315γ8).

The LRPS solution to Eq. (73) in an infinite series is as follows:

U^1γ(s)=(γ-1)s+(-2+4γ-γ2)s2+2(γ3-6γ2+10γ-4)s3+,U^2γ(s)=(1-γ)s+(2-γ2)s2+2(-γ3+2γ)s3+.

By applying the inverse Laplace transform to Eq. (87), we determine the LRPS solution of Eqs. (70) and (71) as follows:

u^1γ(t)=(γ-1)+(-2+4γ-γ2)t+(γ3-6γ2+10γ-4)t2+,u^2γ(t)=(1-γ)+(2-γ2)t+(-γ3+2γ)t2+.

As a special case, when γ = 1, the LRPS solutions to Eqs. (87) and (88) have a general pattern form, which coincides with the exact solution in terms of the infinite series

u(t)=t+t2+13t3-12t4-115t5-745t6+53315t7+71315t8+.

Table 1 shows the numerical results of the exact and approximate solutions using the LRPS method for Example 4.1 at γ = 1. Moreover, the relative and absolute errors at some chosen grid points ti in [0,1] with a step-size of 0.1 are also given in Table 1 for n = 51 and γ = 1. The results indicate that the LRPS approximate solutions are almost as good as the exact solutions. However, a more accurate solution can be obtained using more iterations.

Table 2 shows a numerical comparison between the tenth-LRPS solution and the OHA method [2], the MNN method [3], and the RK-4 method [27]. We concluded that the results obtained using the LRPS method were similar to those obtained using other methods.

Figure 1 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8 with diverse values of (t ∈ {0.0, 0.25, 0.5, 0.75}). Figure 2 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

### Example 4.2

u^(t)=2u^(t)-u^2(t)+1,t>0,

with FIC,

[u^(0)]γ=[0.1+0.1γ,0.3-0.1γ],γ[0,1].

For γ = 1, the initial condition in Eq. (91) becomes u(0) = 0.2. Therefore, the exact solution to Eq. (90) can be expressed as

u(t)=1+2tanh (2t+12log(52-452+4)).

According to Definition 2.10, the FIVPs, as in Eqs. (90) and (91), can be condensed to the set of ODEs corresponding to their parametric forms as follows:

u^1γ(t)=2u^1γ(t)-u^1γ2(t)+1,u^2γ(t)=2u^2γ(t)-u^2γ2(t)+1.

subject to

u^1γ(0)=0.1+0.1γ,u^2γ(0)=0.3-0.1γ.

According to the procedure of the LRPS algorithm presented in Section 3, the Laplace transform of Eq. (90) can be expressed as

U^1γ(s)=0.1+0.1γs+2sU^1γ(s)-1s[(-1[U^1γ(s)])2]+1s2,U^2γ(s)=0.3-0.1γs+2sU^2γ(s)-1s[(-1[U^2γ(s)])2]+1s2,

where Û1γ (s) = ℒ[û1γ (t)], and Û2γ (s) = ℒ[û2γ (t)].

Depending on the initial conditions in Eq. (94), the proposed LRPS solutions Û1γ (t) and Û2γ (t) of the algebraic system in Eq. (95) can be given by

U^k,1γ(s)=0.1+0.1γs+n=1kann!sn+1,U^k,2γ(s)=0.3-0.1γs+n=1kbnn!sn+1,

and the kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), as

LResk,1γ(s)=U^k,1γ(s)-0.1+0.1γs-2sU^k,1γ(s)+1s[(-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-0.3-0.1γs-2sU^k,2γ(s)+1s[(-1[U^k,2γ(s)])2]-1s2.

By utilizing lims→∞(sk+1LResk,1γ (s))=0 and lims→∞(sk+1LResk,2γ (s)) = 0, for k = 1, 2, · · ·, the first few terms ak and bk are

a1=1102(119+18γ-γ2),a2=1103(1071+43γ-27γ2+γ3),a3=13×104(5117-5652γ-658γ2+108γ3-3γ4),a4=13×105(-168147-64585γ+5130γ2+1430γ3-135γ4+3γ5),a5=13×106(-1537123+11682γ+151955γ2+540γ3-2445γ4+162γ5-3γ6),a6=19×107(695079+18187783γ+2825739γ2-719285γ3-40635γ4+11109γ5-567γ6+9γ7),

and

b1=1102(151-14γ-γ2),b2=1103(1057+53γ-21γ2-γ3),b3=13×104(-8003+7084γ-82γ2-84γ3-3γ4),b4=13×105(-267421+30985γ+10710γ2-470γ3-105γ4-3γ5),b5=13×106(-935747-560126γ+108755γ2+13020γ3-1005γ4-126γ5-3γ6),b6=19×107(42289513-20342887γ-1706523γ2+685685γ3+38955γ4-5061γ5-441γ6-9γ7).

We can express the LRPS solution of Eq. (95) in an infinite series as follows:

U^1γ(s)=(0.1+0.1γ)s+1102(119+18γ-γ2)s2+2103(1071+43γ-27γ2+γ3)s3+,U^2γ(s)=(0.3-0.1γ)s+1102(151-4γ-γ2)s2+2103(1057+53γ-21γ2-γ3)s3+.

Applying the inverse Laplace transform to Eq. (99), the LRPS solutions of Eqs. (93) and (94) are then given by

u^1γ(t)=(0.1+0.1γ)+1102(119+18γ-γ2)t+1103(1071+43γ-27γ2+γ3)t2+,u^2γ(t)=(0.3-0.1γ)+1102(151-4γ-γ2)t+1103(1057+53γ-27γ2-γ3)t2+.

As a special case, when γ = 1, the LRPS solution of Eqs. (90) and (91) have a general pattern form, which coincides with the exact solution in terms of infinite series:

u(t)=   15+3425t+136125t2-681875t3-70729365t4-2148846875t5+163744703125t6+1145963224609375t7+16217728123046875t8+.

Tables 3 and 4 present the numerical results of the fuzzy solution upper and lower bounds of the system in Eqs. (93) and (94) at the selected grid points with γ = 1 and n = 51. We observed that for γ = 1, the LRPS upper- and lower-bound solutions are the same and coincide with the crisp solution.

Figure 3 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8with diverse values of (t ∈ {0.25, 0.5, 0.75, 1}), where the center (crisp) solution is represented by the midline for γ = 1. Figure 4 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

### Example 4.3

Consider the following FIVP:

u^(t)=-u^2(t)+1,t>0,

with FIC,

[u^(0)]γ=[γ-1,1-γ],γ[0,1].

For γ = 1 and the crisp initial condition u(0) = 0, the exact solution of (102) is located as

u(t)=(e2t-1)(e2t+1)-1.

The FIVP, as in Eqs. (102), and (103), can be reduced to the set of ODEs corresponding to their parametric forms as

u^1γ(t)=-u^1γ2(t)+1,u^2γ(t)=-u^2γ2(t)+1,

subject to

u^1γ(0)=γ-1,u^2γ(0)=1-γ.

Similar to the steps used in previous examples, the Laplace transform of Eq. (102) can be expressed as follows:

U^1γ(s)=γ-1s-1s[(-1[U^1γ(s)])2]+1s2,U^2γ(s)=1-γs-1s[(-1[U^2γ(s)])2]+1s2.

According to the initial data in Eq. (106), the LRPS solutions, Û1γ (t) and Û2γ (t) of the system in Eq. (107) can be represented as

U^k,1γ(s)=γ-1s+n=1kann!sn+1,U^k,2γ(s)=1-γs+n=1kbnn!sn+1.

The kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), can be defined as follows:

LResk,1γ(s)=U^k,1γ(s)-γ-1s+1s[(-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-1-γs+1s[(-1[U^k,2γ(s)])2]-1s2.

By utilizing the following facts: lims(sk+1LResk,1γ(s))=0 and lims→∞(sk+1LResk,2γ (s)) = 0, for k = 1, 2, · · ·, the first few coefficients ak and bk are

a1=(2γ-γ2),a2=(-2γ+3γ2-γ3),a3=13(4γ-14γ2+12γ3-3γ4),a4=13(-2γ+15γ2-25γ3+15γ4+3γ5),a5=115(4γ-62γ2+180γ3-195γ4+90γ5-15γ6),a6=145(-4γ+126γ2-602γ3+1050γ4-840γ5+315γ6-45γ7),

and

b1=(2γ-γ2),b2=(2γ-3γ2+γ3),b3=13(4γ-14γ2+12γ3-3γ4),b4=13(2γ-15γ2+25γ3-15γ4+3γ5),b5=115(4γ-62γ2+180γ3-195γ4+90γ5-15γ6),b6=145(4γ-126γ2+602γ3-1050γ4+840γ5-315γ6+45γ7).

We can express the LRPS solution of Eq. (107) as in an infinite series by

U^1γ(s)=(γ-1)s+(2γ-γ2)s2+2(-2γ+3γ2-γ3)s3+,U^2γ(s)=(1-γ)s+(2γ-γ2)s2+2(2γ-3γ2+γ3)s3+.

Applying the inverse Laplace transform to Eq. (111), we then obtain the LRPS solution of Eqs. (105) and (106) as

u^1γ(t)=(γ-1)+(2γ-γ2)t+(-2γ+3γ2-γ3)t2+,u^2γ(t)=(1-γ)+(2γ-γ2)t+(2γ-3γ2+γ3)t2+.

As a special case, when γ = 1, the LRPS solutions to Eqs. (102) and (103) have a general pattern form, which coincides with the exact solution in terms of infinite series:

u(t)=t-13t3+215t5-17315t7+622835t9+.

Table 5 presents the error analysis of the proposed algorithm for the system in Eqs. (105) and (106) at some nodes t in [0, 1] with a step size of 0.16 for γ = 1 and n = 51. Table 6 presents a numerical comparison between the tenth-LRPS solution and RK-4, OHAM, and MNN given in the literature. We observed that the numerical results of LRPS were as good as those obtained using other techniques.

Figure 5 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8 with diverse values of t (t ∈ {0.25, 0.5, 0.75, 1}), where the center (crisp) solution is represented by the midline for γ = 1. Figure 6 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

Many numerical and analytical methods have been applied to solve fuzzy differential equations, some of which have advantages over others. Some of them are accurate and effective; however, they require mathematical operations that can be difficult, long, or sometimes failed. Others are simple and fast; however, they may not provide precise solutions. Our new method, called LRPS, is characterized by its accuracy, speed, and simplicity in finding exact and accurate approximate solutions of linear and nonlinear fuzzy differential equations. In this study, we succeeded in providing accurate approximate and sometimes exact solutions to FIVPs using the newly proposed technique. The proposed expansion of our study enabled us to obtain a serial solution for equations in the Laplace transform space. The LRPS method provides an easy and fast technique for determining the proposed series coefficients as a solution to the equation. Unlike the traditional RPS method, which specifies series coefficients, the derivative must be calculated each time, whereas LRPS only requires the concept of the limit at infinity in determining series coefficients, which distinguishes this method. It is worth noting that the LRPS method can be applied to solve other types of fuzzy ordinary and partial differential equations for an integer or fractional order, such as fuzzy Poisson, KdV, Schrödinger, and apple equations. Furthermore, based on the H-SGD theory, it is possible to extend and generalize our analysis and results to all fuzzy interactive arithmetic operations mentioned in Section 2 in future studies.

Fig. 1.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.1. Solution plots for (a) t = 0, (b) t = 0.25, (c) t = 0.5, and (d) t = 0.75.

Fig. 2.

The 3-dim plot for Example 4.1. Blue and yellow are the lower and upper bounds, respectively, of the tenth-LRPS fuzzy solution.

Fig. 3.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.2. Solution plots for (a) t = 0.25, (b) t = 0.5, (c) t = 0.75, and (d) t = 1.

Fig. 4.

3-dim plot for Example 4.2. Blue and yellow are the lower and upper bounds, respectively, of the tenth-LRPS fuzzy solution.

Fig. 5.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.3. Solution plots for (a) t = 0, (b) t = 0.25, (c) t = 0.5, and (d) t = 0.75.

Fig. 6.

3-dim plot for Example 4.3. Blue and yellow are lower and upper bounds, respectively, of tenth-LRPS fuzzy solution.

Table. 1.

Table 1. Numerical results of the LRPS solutions for Example 4.1 at γ = 1.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.1102951969160.1102951969171.52656 × 10−161.38406 × 10−15
0.20.2419767996210.2419767996215.55112 × 10−172.29407 × 10−16
0.30.3951048486600.3951048486602.77556 × 10−167.02486 × 10−16
0.40.5678121662930.5678121662933.33067 × 10−165.86579 × 10−16
0.50.7560143934310.7560143934313.33067 × 10−164.40556 × 10−11
0.60.9535662164720.9535662164721.11022 × 10−161.16429 × 10−16
0.71.1529489669791.1529489669805.32907 × 10−144.62212 × 10−14
0.81.3463636553681.3463636554225.39830 × 10−114.00954 × 10−11
0.91.5269113132801.5269113368462.35658 × 10−81.54336 × 10−8)

Table. 2.

Table 2. Approximate solutions of Example 4.1 by various methods at γ = 1.

tiExact10th LRPSOHAMNNRK-4
0.10.1102950.1102950.1103280.1102950.100000
0.20.2419770.2419770.2422730.2419760.219000
0.30.3951050.3951050.3961750.3950890.358004
0.40.5678120.5678120.5702310.5676600.516788
0.50.7560140.7560140.759550.7551340.693439
0.60.9535660.9535660.9550490.9499640.884041
0.71.1529491.1529491.1424441.1414231.082696
0.81.3463641.3463581.3005691.3157231.282012
0.91.5269111.5268141.4004441.4565451.474059

Table. 3.

Table 3. Numerical results of lower bound LRPS-solution of Example 4.2.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.346763995050.3467639950591.11022 × 10−163.20167 × 10−16
0.20.513897275840.5138972758411.11022 × 10−162.16040 × 10−16
0.30.697990872510.6979908725142.22045 × 10−163.18120 × 10−16
0.40.893472160240.8934721602443.33067 × 10−163.72778 × 10−16
0.51.093134993061.0931349930602.22045 × 10−162.03126 × 10−16
0.61.289137171941.2891371719472.22045 × 10−161.72243 × 10−16
0.71.474194376961.4741943769631.05871 × 10−127.18161 × 10−13
0.81.642601362551.6426016335771.02182 × 10−96.22073 × 10−10
0.91.790797918331.7907983475064.29173 × 10−72.39655 × 10−7

Table. 4.

Table 4. Numerical results of upper bound LRPS-solution of Example 4.2.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.346763995050.3467639950591.66533 × 10−164.80250 × 10−16
0.20.513897275840.5138972758411.11022 × 10−162.16040 × 10−16
0.30.697990872510.6979908725142.22045 × 10−163.18120 × 10−16
0.40.893472160240.8934721602443.33067 × 10−163.72778 × 10−16
0.51.093134993061.0931349930602.22045 × 10−162.03126 × 10−16
0.61.289137171941.2891371719472.22045 × 10−161.72243 × 10−16
0.71.474194376961.4741943769631.05871 × 10−127.18161 × 10−13
0.81.642601362551.6426016335771.02182 × 10−96.22073 × 10−10
0.91.790797918331.7907983475064.29173 × 10−72.39655 × 10−7

Table. 5.

Table 5. Numerical results of the LRPS solutions of Example 4.3 at γ = 1.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.160.1586485042970.1586485042975.55112 × 10−173.49900 × 10−16
0.3203095069212130.3095069212135.55112 × 10−171.79354 × 10−16
0.480.4462436102490.4462436102495.55112 × 10−171.24397 × 10−16
0.640.5648995528460.5648995528460.000000.00000
0.800.6640367702680.6640367702683.33067 × 10−165.01579 × 10−16
0.960.7442768673620.7442768673584.29645 × 10−125.77265 × 10−12

Table. 6.

Table 6. Numerical comparison of Example 4.3 at γ = 1.

tiExact10th-LRPSOHAMMNNRK-4
0.10.0996680.0996680.0996680.0996680.100000
0.20.1973750.1973750.1973760.1973750.199000
0.30.2913130.2913130.2913150.2913130.295040
0.40.3799490.3799490.3799490.3799490.386335
0.50.4621170.4621210.4620920.4621210.471410
0.60.5370500.5370780.5369100.5370770.549187
0.70.6043680.6045140.6038150.6045130.619026
0.80.6640370.6646410.6622450.6646400.680707
0.90.7162980.7183920.7112870.7183900.734371

1. Oberguggenberger, M, and Pittschmann, S (1999). Differential equations with fuzzy parameters. Mathematical and Computer Modelling of Dynamical Systems. 5, 181-202. https://doi.org/10.1076/mcmd.5.3.181.3683
2. Gumah, G, Naser, MFM, Al-Smadi, M, Al-Omari, SKQ, and Baleanu, D (2020). Numerical solutions of hybrid fuzzy differential equations in a Hilbert space. Applied Numerical Mathematics. 151, 402-412. https://doi.org/10.1016/j.apnum.2020.01.008
3. Alaroud, M, Al-Smadi, M, Rozita Ahmad, R, and Salma Din, UK (2019). An analytical numerical method for solving fuzzy fractional Volterra integro-differential equations. Symmetry. 11. article no. 205
4. Al-Smadi, M (2019). Reliable numerical algorithm for handling fuzzy integral equations of second kind in Hilbert spaces. Filomat. 33, 583-597. https://doi.org/10.2298/FIL1902583A
5. Hasan, S, El-Ajou, A, Hadid, S, Al-Smadi, M, and Momani, S (2020). Atangana-Baleanu fractional framework of reproducing kernel technique in solving fractional population dynamics system. Chaos, Solitons & Fractals. 133. article no. 109624
6. Ahmadian, A, Salahshour, S, and Chan, CS (2017). Fractional differential systems: a fuzzy solution based on operational matrix of shifted Chebyshev polynomials and its applications. IEEE Transactions on Fuzzy Systems. 25, 218-236. https://doi.org/10.1109/TFUZZ.2016.2554156
7. El-Ajou, A, Oqielat, MAN, Ogilat, O, Al-Smadi, M, Momani, S, and Alsaedi, A (2019). Mathematical model for simulating the movement of water droplet on artificial leaf surface. Frontiers in Physics. 7. article no. 132
8. El-Ajou, A, Oqielat, MN, Al-Zhour, Z, and Momani, S (2019). Analytical numerical solutions of the fractional multi-pantograph system: two attractive methods and comparisons. Results in Physics. 14. article no. 102500
9. El-Ajou, A, Al-Zhour, Z, Oqielat, MA, Momani, S, and Hayat, T (2019). Series solutions of nonlinear conformable fractional KdV-Burgers equation with some applications. The European Physical Journal Plus. 134. article no. 402
10. El-Ajou, A, Oqielat, MAN, Al-Zhour, Z, Kumar, S, and Momani, S (2019). Solitary solutions for time-fractional nonlinear dispersive PDEs in the sense of conformable fractional derivative. Chaos: An Interdisciplinary Journal of Nonlinear Science. 29. article no. 093102
11. Oqielat, MAN, El-Ajou, A, Al-Zhour, Z, Alkhasawneh, R, and Alrabaiah, H (2020). Series solutions for nonlinear time-fractional Schrödinger equations: comparisons between conformable and Caputo derivatives. Alexandria Engineering Journal. 59, 2101-2114. https://doi.org/10.1016/j.aej.2020.01.023
12. El-Ajou, A (2020). Taylor’s expansion for fractional matrix functions: theory and applications. Journal of Mathematics and Computer Science. 21, 1-17. https://doi.org/10.22436/jmcs.021.01.01
13. Shqair, M, El-Ajou, A, and Nairat, M (2019). Analytical solution for multi-energy groups of neutron diffusion equations by a residual power series method. Mathematics. 7. article no. 633
14. Santo Pedro, F, de Barros, LC, and Esmi, E (2019). Population growth model via interactive fuzzy differential equation. Information Sciences. 481, 160-173. https://doi.org/10.1016/j.ins.2018.12.076
15. Chang, SSL, and Zadeh, LA (1972). On fuzzy mapping and control. IEEE Transactions on Systems, Man, and Cybernetics. 2, 30-34. https://doi.org/10.1109/TSMC.1972.5408553
16. Dubois, D, and Prade, H (1982). Towards fuzzy differential calculus. Part 1: integration of fuzzy mappings. Fuzzy Sets and Systems. 8, 1-17. https://doi.org/10.1016/0165-0114(82)90025-2
17. Wasques, V, Laiate, B, Santo Pedro, F, Esmi, E, and Barros, LCD 2020. Interactive fuzzy fractional differential equation: application on HIV dynamics. Information Processing and Management of Uncertainty in Knowledge-Based Systems, Array, pp.198-211. https://doi.org/10.1007/978-3-030-50153-2_15
18. Gumah, GN, Naser, MF, Al-Smadi, M, and Al-Omari, SK (2018). Application of reproducing kernel Hilbert space method for solving second-order fuzzy Volterra integro-differential equations. Advances in Difference Equations. 2018. article no. 475
19. Goetschel, R, and Voxman, W (1986). Elementary fuzzy calculus. Fuzzy Sets and Systems. 18, 31-43. https://doi.org/10.1016/0165-0114(86)90026-6
20. Arqub, OA, and Al-Smadi, M (2020). Fuzzy conformable fractional differential equations: novel extended approach and new numerical solutions. Soft Computing. 24, 12501-12522. https://doi.org/10.1007/s00500-020-04687-0
21. Seikkala, S (1987). On the fuzzy initial value problem. Fuzzy Sets and Systems. 24, 319-330. https://doi.org/10.1016/0165-0114(87)90030-3
22. Bede, B, and Gal, SG (2004). Almost periodic fuzzy-number-valued functions. Fuzzy Sets and Systems. 147, 385-403. https://doi.org/10.1016/j.fss.2003.08.004
23. Bede, B, and Stefanini, L (2013). Generalized differentiability of fuzzy-valued functions. Fuzzy Sets and Systems. 230, 119-141. https://doi.org/10.1016/j.fss.2012.10.003
24. Alshammari, M, Al-Smadi, M, Alshammari, S, Arqub, OA, Hashim, I, and Alias, M (2020). An attractive analytic-numeric approach for the solutions of uncertain Riccati differential equations using residual power series. Applied Mathematics and Information Sciences. 14, 177-190. https://doi.org/10.18576/amis/140202
25. Behzadi, SS, Vahdani, B, Allahviranloo, T, and Abbasbandy, S (2016). Application of fuzzy Picard method for solving fuzzy quadratic Riccati and fuzzy Painlevé I equations. Applied Mathematical Modelling. 40, 8125-8137. https://doi.org/10.1016/j.apm.2016.05.003
26. Arqub, OA, Al-Smadi, M, Momani, S, and Hayat, T (2017). Application of reproducing kernel algorithm for solving second-order, two-point fuzzy boundary value problems. Soft Computing. 21, 7191-7206. https://doi.org/10.1007/s00500-016-2262-3
27. Wasques, VF, Esmi, E, Barros, LC, and Bede, B 2019. Comparison between numerical solutions of fuzzy initial-value problems via interactive and standard arithmetics. Fuzzy Techniques: Theory and Applications, Array, pp.704-715. https://doi.org/10.1007/978-3-030-21920-8_62
28. Arqub, OA, El-Ajou, A, Momani, S, and Shawagfeh, N (2013). Analytical solutions of fuzzy initial value problems by HAM. Applied Mathematics & Information Sciences. 7. article no. 1903
29. Wasques, VF, Esmi, E, Barros, LC, and Sussner, P (2018). Numerical solutions for bidimensional initial value problem with interactive fuzzy numbers. Fuzzy Information Processing. Cham, Switzerland: Springer, pp. 84-95 https://doi.org/10.1007/978-3-319-95312-0_8
30. Hasan, S, Al-Smadi, M, El-Ajou, A, Momani, S, Hadid, S, and Al-Zhour, Z (2021). Numerical approach in the Hilbert space to solve a fuzzy Atangana-Baleanu fractional hybrid system. Chaos, Solitons & Fractals. 143. article no. 110506
31. Padmapriya, V, Kaliyappan, M, and Manivannan, A (2020). Numerical solutions of fuzzy fractional delay differential equations. International Journal of Fuzzy Logic and intelligent Systems. 20, 247-254. https://doi.org/10.5391/IJFIS.2020.20.3.247
32. Ahmad, J, Iqbal, A, and Ul Hassan, QM (2021). Study of nonlinear fuzzy integro-differential equations using mathematical methods and applications. International Journal of Fuzzy Logic and Intelligent Systems. 21, 76-85. https://doi.org/10.5391/IJFIS.2021.21.1.76
33. Kwun, YC, Han, CW, Kim, SY, and Park, JS (2004). Existence and uniqueness of fuzzy solutions for the nonlinear fuzzy integro-differential equation on E. International Journal of Fuzzy Logic And Intelligent Systems. 4, 40-44. https://doi.org/10.5391/IJFIS.2004.4.1.040
34. Kwun, YC, Park, JS, Kim, SY, and Park, JH (2006). Existence and uniqueness of solutions for the semilinear fuzzy integrodifferential equations with nonlocal conditions and forcing term with memory. International Journal of Fuzzy Logic and Intelligent Systems. 6, 288-292. https://doi.org/10.5391/IJFIS.2006.6.4.288
35. Kwun, YC, Kim, WH, Nakagiri, SI, and Park, JH (2009). Existence and uniqueness of solutions for the fuzzy differential equations in n-dimension fuzzy vector space. International Journal of Fuzzy Logic and Intelligent Systems. 9, 16-19. https://doi.org/10.5391/IJFIS.2009.9.1.016
36. Kwun, YC, Kim, JS, Hwang, JS, and Park, JH (2010). Existence of periodic solutions for fuzzy differential equations. International Journal of Fuzzy Logic and Intelligent Systems. 10, 184-193. https://doi.org/10.5391/IJFIS.2010.10.3.184
37. Lee, BY, Kwun, YC, Ahn, YC, and Park, JH (2009). The existence and uniqueness of fuzzy solutions for semilinear fuzzy integrodifferential equations using integral contractor. International Journal of Fuzzy Logic and Intelligent Systems. 9, 339-342.
38. Yoon, JH, Kwun, YC, Park, JS, and Park, JH (2007). Controllability for the semilinear fuzzy integrodifferential equations with nonlocal conditions and forcing term with memory. International Journal of Fuzzy Logic and Intelligent Systems. 7, 34-40. https://doi.org/10.5391/IJFIS.2007.7.1.034
39. Lee, BY, Youm, HE, and Kim, JS (2014). Exact controllability for fuzzy differential equations in credibility space. International Journal of Fuzzy Logic and Intelligent Systems. 14, 145-153. https://doi.org/10.5391/IJFIS.2014.14.2.145
40. El-Ajour, A, and Alawneh, A (2008). Modified Homotopy Analysis Method: Application to Linear and Nonlinear Ordinary Differential Equations of Fractional Order. Amman, Jordan: University of Jordan
41. Odibat, Z, and Momani, S (2008). Modified homotopy perturbation method: application to quadratic Riccati differential equation of fractional order. Chaos, Solitons & Fractals. 36, 167-174. https://doi.org/10.1016/j.chaos.2006.06.041
42. Biazar, J, and Eslami, M (2010). Differential transform method for quadratic Riccati differential equation. International Journal of Nonlinear Science. 9, 444-447.
43. Abbasbandy, S (2004). Solving Riccati differential equation using Adomian decomposition method. Appl Math Computers. 157, 54-63.
44. Tapaswini, S, and Chakraverty, S (2013). Approximate solution of fuzzy quadratic Riccati differential equation. Coupled Systems of Mechanics. 23, 255-269. https://doi.org/10.12989/csm.2013.2.3.255
45. Hooshangian, L (2016). Numerical solution for fuzzy Riccati equation under generalized derivation. International Journal of Applications of Fuzzy Sets and Artificial Intelligence. 6, 133-143.
46. Gnitchogna, R, and Atangana, A (2018). New two step Laplace Adam-Bashforth method for integer a noninteger order partial differential equations. Numerical Methods for Partial Differential Equations. 34, 1739-1758. https://doi.org/10.1002/num.22216
47. Atangana, A, and Doungmo Goufo, EF (2014). Extension of matched asymptotic method to fractional boundary layers problems. Mathematical Problems in Engineering. 2014. article no. 107535
48. Bede, B, and Gal, SG (2005). Generalizations of the differentiability of fuzzy-number-valued functions with applications to fuzzy differential equations. Fuzzy Sets and Systems. 151, 581-599. https://doi.org/10.1016/j.fss.2004.08.001
49. Diamond, P, and Kloeden, PE (1994). Metric Spaces of Fuzzy Sets: Theory and Applications. Singapore: World Scientific Publishing
50. Dubois, D, and Prade, H (1982). Towards fuzzy differential calculus. Part 3: Differentiation. Fuzzy Sets and Systems. 8, 225-233. https://doi.org/10.1016/s0165-0114(82)80001-8
51. Zimmermann, HJ (1992). Fuzzy Set Theory and its Applications. Boston, MA: Kluwer Academic Publisher
52. Ghaemi, F, Yunus, R, Ahmadian, A, Salahshour, S, Suleiman, M, and Saleh, SF (2013). Application of fuzzy fractional kinetic equations to modelling of the acid hydrolysis reaction. Abstract and Applied Analysis. 2013. article no 610314
53. Anastassiou, GA, and Gal, SG (2001). On a fuzzy trigonometric approximation theorem of Weierstrass-type. Journal of Fuzzy Mathematics. 9, 701-708.
54. Anastassiou, GA (2010). Fuzzy Mathematics: Approximation Theory. Heidelberg, Germany: Springer
55. Chalco-Cano, Y, and Roman-Flores, H (2008). On new solutions of fuzzy differential equations. Chaos, Solitons & Fractals. 38, 112-119. https://doi.org/10.1016/j.chaos.2006.10.043
56. Puri, ML, and Ralescu, DA (1986). Fuzzy random variables. Journal of Mathematical Analysis and Applications. 114, 409-422. https://doi.org/10.1016/0022-247X(86)90093-4
57. Pirzada, UM, and Vakaskar, DC (2017). Existence of Hukuhara differentiability of fuzzy-valued functions. Journal of the Indian Mathematical Society. 84, 27-46. https://doi.org/10.18311/jims/2017/5824
58. Armand, A, Allahviranloo, T, and Gouyandeh, Z (2018). Some fundamental results on fuzzy calculus. Iranian Journal of Fuzzy Systems. 15, 27-46. https://doi.org/10.22111/ijfs.2018.3948
59. Wasques, VF, Esmi, E, Barros, LC, and Sussner, P (2020). The generalized fuzzy derivative is interactive. Information Sciences. 519, 93-109. https://doi.org/10.1016/j.ins.2020.01.042
60. Ma, Y, Li, H, and Zhang, S (2020). Solving two-dimensional fuzzy Fredholm integral equations via sinc collocation method. Advances in Difference Equations. 2020. article no. 290
61. Kaleva, O (1987). Fuzzy differential equations. Fuzzy Sets and Systems. 24, 301-317. https://doi.org/10.1016/0165-0114(87)90029-7
62. Bede, B (2013). Mathematics of Fuzzy Sets and Fuzzy Logic. Heidelberg, Germany: Springer
63. You, C, and Zhang, R (2018). Existence and uniqueness theorems for complex fuzzy differential equation. Journal of Intelligent & Fuzzy Systems. 34, 2213-2222. https://doi.org/10.3233/JIFS-171231

Moa’ath N. Oqielat received his Ph.D. in applied mathematics from Queensland University of Technology (Australia) in 2010. He then began work at Al-Balqa Applied University as assistant professor of applied mathematics and was promoted to associate professor in 2020. His research interests focus on modeling and simulations, fractional calculus, and numerical analysis. Dr. Oqielat has written over 35 research papers.

E-mail: moaathoqily@bau.edu.jo

Ahmad El-Ajou received his Ph.D.degree in mathematics from the University of Jordan (Jordan) in 2009. He began work at Al-Balqa Applied university in 2011 as an assistant professor of applied mathematics. In 2015, he became an associate professor, and in 2020 became a professor, at the same university. His research interests focus on numerical analysis, fractional differential and integral equations, fuzzy differential and integral equations, and simulations and modeling.

E-mail: ajou44@bau.edu.jo

Zeyad Al-Zhour received his Ph.D. degree in applied mathematics from University Putra Malaysia (Malaysia) in 2007. He then started working at Zarqa University (2007–2010) as an assistant professor, and then at Imam Abdulrahman Bin Faisal University (2010-present), and was promoted to Associate Professor in 2018. Dr. Al-Zhours has supervised several MSc. and Ph.D. theses, published more than 60 articles, and completed several projects. He has received many awards (e.g., Distinguished Scientist in Applied Mathematics-VIRA, 2018) and letters of appreciation during his work from many universities and companies. His research interests are in matrix operators and algebraic systems, fractional-order systems, numerical methods, mathematical engineering, and mathematical physics.

Tareq Eriqat received his Ph.D. degree from Yarmouk University (Jordan) in 2021. He began work at the Ministry of Education in 2014 as a mathematics teacher. His research interests are on the numerical solution to fractional differential and integral equations, and numerical analysis.

E-mail: tareq.eriqat92@gmail.com

Mohammed Al-Smadi received his Ph.D. degree from the University of Jordan, Amman in 2011. He then began working as an assistant professor of applied mathematics at Qassim University (2011–2012) and then at Tafila Technical University (2012–2013). Since Sep. 2013, He has worked as an assistant professor of applied mathematics at Al-Balqa Applied University and was promoted to associate professor in 2017. His research interests include applied mathematics, numerical analysis, and fractional calculus.

### Article

#### Original Article

International Journal of Fuzzy Logic and Intelligent Systems 2022; 22(1): 23-47

Published online March 25, 2022 https://doi.org/10.5391/IJFIS.2022.22.1.23

## A New Approach to Solving Fuzzy Quadratic Riccati Differential Equations

1Department of Mathematics, Faculty of Science, Al-Balqa Applied University, Salt, Jordan
2Department of Basic Engineering Sciences, College of Engineering, Imam Abdulrahman Bin Faisal University, Dammam, Saudi Arabia
3Department of Applied Science, Ajloun College, Al-Balqa Applied University, Ajloun, Jordan
4Nonlinear Dynamics Research Centre (NDRC), Ajman University, Ajman, UAE

Received: June 4, 2021; Revised: September 5, 2021; Accepted: November 1, 2021

This is an Open Access article distributed under the terms of the Creative Commons Attribution Non-Commercial License (http://creativecommons.org/licenses/by-nc/3.0/) which permits unrestricted noncommercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

### Abstract

In this paper, we present in detail the power series solutions to fuzzy quadratic Riccati differential equations (QRDEs) along with a suitable fuzzy constant through an interactive derivative, more specifically, the Hukuhara-strongly generalized differentiability (H-SGD) based on our new technique. This technique is called the Laplace residual power series (LRPS) method, and it mainly depends on a new expansion and the combination of the Laplace transform technique with the residual power series method. To validate the accuracy of our proposed algorithm, numerous examples were examined numerically and graphically, and we compared the results of the optimal homotopy asymptotic (OHA), multiagent neural network (MNN), and fourth-order Runge-Kutta (RK-4) methods with the LRPS method at γ = 1.

Keywords: Fuzzy-valued function, Strongly generalized differentiability, Quadratic Riccati differential equation, Laplace residual power series method, Laplace and inverse transforms

### 1. Introduction

Fuzzy set theory is an important topic in the study and modeling of many real-life problems related to numerous physical phenomena and dynamical processes, such as astronomy, quantum mechanics, chromodynamics, quantum optics, electronic mechanisms, electronic-controlled systems, the modeling of anomalous hydraulic diffusion systems, and population dynamics models [114]. Moreover, several problems in applied mathematics that are associated with biology, medicine, physics, and technology can be modeled using fuzzy differential equations [16,1517].

Fuzzy set theory has been explored since the 1920s; however, the term fuzzy derivative was introduced by Chang and Zadeh [15]. Subsequently, Dubois and Prade [16] proposed the fuzzy differential calculus concept. Many researchers have recently presented other results and notations for fuzzy mapping [1821]. Moreover, Bede and his colleagues [22,23] defined strongly generalized differentiability (SGD) as a fuzzy-valued function (FVF). In general, it is not straightforward to obtain the exact solution for these types of problems because of the difficulties involved; therefore, reliable numerical techniques are needed, including a reliable numerical algorithm for handling fuzzy integral equations of the second kind in the Hilbert spaces [4], a residual power series (RPS) for solving uncertain Riccati differential equations

[24], and a fuzzy Picard method for solving fuzzy quadratic Riccati and fuzzy Painlevl equations [25]. Other examples include a novel extended approach for obtaining numerical solutions to fuzzy conformable fractional differential equations [22], a reproducing kernel algorithm for solving second-order two-point fuzzy boundary value problems and second-order fuzzy Volterra integro-differential equations [3,18,26], a Runge-Kutta method applied to the SI epidemiological model [27], a homotopy analysis method (HAM) for solving fuzzy initial value problems (FIVPs) [28], and bidimensional IVP with interactive fuzzy numbers [14]. Additional solutions are a logistic model of population growth through an interactive derivative [29], an HIV viral dynamics model for individuals under antiretroviral treatment [17], a numerical approach in the Hilbert space to solve a fuzzy Atangana-Baleanu fractional hybrid system [30], a novel technique for solving fuzzy fractional delay differential equations [31], and a homotopy perturbation Sumudu transform method to find the analytical fuzzy solution of nonlinear fuzzy integro-differential equations [32].

The existence, uniqueness, and controllability of fuzzy solutions for such nonlinear fuzzy (fractional) differential equations have been studied and discussed in [3339]. Kwun et al. [3336] and Lee et al. [37] respectively studied the existence and uniqueness of fuzzy (periodic) solutions for the nonlinear fuzzy integro-differential equations on $ENn$, semi linear fuzzy integro-differential equation with an initial nonlocal condition, fuzzy differential equations in an n-dimension fuzzy vector space (EN)n, fuzzy differential equations with three types of forcing terms using a Hukuhara derivative, and semilinear fuzzy integro-differential equations using an integral contractor. For the controllability, Yoon et al. [38] and Lee et al. [39] studied the controllability of fuzzy solutions for semilinear fuzzy integro-differential equations with non-local conditions and a forcing term with memory, and abstract fuzzy differential equations in the credibility space from the perspective of the Liu process, respectively.

Quadratic Riccati differential equations (QRDEs) form a precise nonlinear model for defining a certain type of physical scheme with applications in the diffusion process, optimal control, optical networks, and artificial intelligence, and usually have the following general form [4043] of

$u′(t)=A(t)u2(t)+B(t)u(t)+C(t),t≥0,$

subject to

$u(0)=λ.$

Many numerical and analytical techniques have been employed to obtain QRDEs solutions, including the HAM [40], homotopy perpetuation method (HPM) [41], differential transformation method [42], and Adomian decomposition method [43].

Several authors [24,26,28,4045] have reformulated the QRDE as a fuzzy differential equation along with suitable fuzzy constants under the SGD concept, the fuzzy form of which is as follows [24]:

$u^′(t)=Au^2(t)+Bu^(t)+C,0≤t≤T,$

subject to the fuzzy initial condition (FIC):

$u^(0)=λ^,$

where A, B, and C are constants of t, û(t): [0, T] → ℝ is an unknown fuzzy function of the crisp variable t, and λ̂ ∈ ℝ (the set of all fuzzy numbers). Several techniques have been applied to solve (fractional) differential equations [713,4648], and other techniques have been used to find the solution to an FIVP in Eq. (3), including the HPM [44], Euler method [45], fuzzy Picard technique [25], and RPS method [24].

The main objective of the research described in this study is to extend the application of the LRPS method and provide a power series solution to the FIVP in Eq. (3), as well as study the behavior of the fuzzy solution obtained.

The LRPS method is a new method that combines the Laplace transform technique with the RPS method, which is used to solve the integral and differential equations [713]. The Laplace transform method involves transforming the equation from the space within it into a new space, presenting its solution in the new space, and then returning the solution to the original space and thus obtaining the solution to the equation. However, it is sometimes difficult to offer a solution to an equation in a new space. The RPS method assumes a power series solution to the equation and provides an easy technique for determining the series coefficients. The LRPS method uses the two methods together to provide a series solution to the equation in the new space and uses the principle of the RPS method to determine the series coefficients. However, through a new approach, it is easier and faster than the user in the original technique. In fact, our proposed method is an analytical and efficient method that can be easily applied to solve different types of fuzzy problems.

Two fuzzy numbers are considered non-interactive when the joint possibility distribution involving both is given by the minimum t-norm. Thus, any operation between two fuzzy numbers, aside from the one obtained by the minimum t-norm, must be considered to be interactive. Following this reasoning, when we use the Hukuhara difference (H-difference) for two fuzzy numbers, we assume that they are interactive. Thus, all derivatives obtained by a fuzzy arithmetic operation that differs from the standard operation must be considered interactive. In this study, we find analytical and numerical solutions for FIVPs along with a suitable fuzzy constant through an interactive derivative, more specifically, the Hukuhara-strongly generalized differentiability (H-SGD) based on our new LRPS technique. The solution to the FIVP obtained here using H-SGD is more suitable than the classical derivative because the latter uses an interactive derivative obtained by the H-difference, whereas the standard sum is used to compute from the non-interactive operation, that is, the sum operation is based on the minimum t-norm. This is because the arithmetical operations compatible with the H-difference are not known. Thus, it makes no sense to combine interactive and non-interactive H-SGD. Furthermore, with the H-SGD theory, it is possible to extend and generalize our analysis and results in this study to all fuzzy interactive arithmetic operations, and not only to the H-difference operation in future studies. In addition, in the present paper, we introduce a new expansion to Taylor’s formula, which is required for creating the LRPS method and formulating the FIVP. Moreover, several examples are solved and examined numerically and graphically using this new technique, and the results of the optimal homotopy asymptotic (OHA), multiagent neural network (MNN), and fourth-order Runge-Kutta (RK-4) methods are compared with those of the LRPS method at γ = 1. For the calculations and results, Maple 2018 was used to conduct the computations. The results show that the proposed algorithm is efficient, accurate, and robust.

This paper consists of five sections. In Section 2, we study some basic concepts of fuzzy theory, establish a new expansion to Taylor’s formula, and formulate the FIVP. In Section 3, the LRPS algorithm is created and successfully tested to obtain analytical and numerical solutions of the nonlinear FIVP. In Section 4, we test the accuracy of our proposed algorithm (LRPS) by solving three important examples numerically and graphically, and comparing the results of the OHA, MNN, and RK-4 methods with the LRPS method at γ = 1. Finally, the conclusions are presented in the last section..

### 2. Basic Concepts, A New Formula and Drafting

This section reviews some theories and concepts necessary for our study and is divided into three subsections: the first relates to fuzzy calculus, the second relates to a new expansion needed to construct the LRPS solution, and we formulate the FIVP in third.

### 2.1 Basic Concepts and Notations on Fuzzy Theory

In this subsection, we study the most important definitions and notations related to fuzzy analysis and fuzzy calculus [3,4,12, 15,16,1823,4863] used throughout this paper. The provided ℝ denotes the set of all the fuzzy numbers.

Definition 2.1 [19,48,59,60]

Letting v: ℝ → [0, 1], v is then called a fuzzy number if the following conditions hold:

(i) v is normal (i.e., ∃ t0 ∈ ℝ such that v(t0) = 1),

(ii) v is convex (i.e., ∀ t, s ∈ ℝ, and 0 ≤ λ ≤ 1, and it holds that

$v(λt+(1-λ)s)≥min{v(t),v(s)}),$

(iii) v is upper semi-continuous (i.e., {t ∈ ℝ: v(t) ≥ r} is closed for all r ∈ ℝ), and

(iv) v is the bounded support (i.e., $[v]0=supp(v)¯$ is a compact subset of ℝ).

Definition 2.2 [14,17,19,20,53]

Let v ∈ ℝ and γ ∈ [0, 1], the γ-cut (or level) of v is then the crisp set

$[v]γ={t∈ℝ:v(t)≥γ}.$

Therefore, if v ∈ ℝℱ, then γ-levels are closed intervals in ℝ and are denoted by

$[v]γ=[v1(γ),v2(γ)]=[v1γ,v2γ];$

$v1γ=v1(γ)=min{r:r∈[v]γ},v2γ=v2(γ)=max{r:r∈[v]γ}.$

In addition, if v ∈ ℝ, the length of the γ-level set of v is defined as Len[v]γ = v2γv1γ for all γ ∈ [0, 1]. Note that if γ = 0, then Len[v]0 = Diam(v), where Diam(v) is the diameter of a number fuzzy v [14,17]. However, there are several common forms of a fuzzy number v, one of which is the triangular form defined by an ordered triple v = (μ1, μ2, μ3) ∈ ℝ3 with μ1< μ2< μ3 and whose γ-level is

$[v]γ=[μ1+(μ2-μ1)γ,μ3-(μ3-μ2)γ],γ∈[0,1].$

Note that a fuzzy set v of ℝ is said to be a fuzzy number if [v]γ is a non-empty, bounded, and closed interval of ℝ for every γ ∈ [0, 1]. The fuzzy numbers satisfy the following property: vω if and only if [v]γ ⊆ [ω]γ for all γ ∈ [0, 1] [14,27,29,60].

Some conditions must be satisfied by two functions: v1γ, v2γ: [0, 1] → ℝ, such that [v1γ, v2γ] can be a parametric form of a fuzzy number v for all γ ∈ [0, 1]. These conditions are presented in the following theorem.

Theorem 2.1 [57]

Suppose that v1γ, v2γ: [0, 1] → ℝ satisfy the following three conditions:

(i) v1γ is a bounded monotonic non-decreasing left-continuous function on γ ∈ (0, 1],

(ii) v2γ is a bounded monotonic non-increasing left-continuous function on γ ∈ (0, 1],

(iii) v1γ and v2γ are right-continuous at γ = 0, and

(iv) v11v21.

Subsequently, v: ℝ → [0, 1] defined by v(t) = sup{γ: v1γtv2γ} is a fuzzy number with a γ-level set: [v]γ = [v1γ, v2γ]. Conversely, if v is a fuzzy number with [v]γ = [v1γ, v2γ], then the functions v1γ, v2γ: [0, 1] → ℝ satisfy the conditions (i)–(iv).

Definition 2.3 [1016,1821,53]

Let v = [v1γ, v2γ], ω = [ω1γ, ω2γ] ∈ ℝ, μ ∈ ℝ \ {0}, and γ ∈ [0, 1]. If v and ω are non-interactive (i.e., the joint possibility distribution is given by the minimum t-norm), then the following operations for all γ ∈ [0, 1] in the levels are given by the following:

(i) standard sum,

$[v+ω]γ=[v1γ+ω1γ,v2γ+ω2γ];$

(ii) standard differences,

$[v-ω]γ=[v1γ-ω2γ,v2γ-ω1γ];$

(iii) scalar multiplication,

$[μv]γ=μ[v]γ=[min{μv1γ,μv2γ},max{μv1γ,μv2γ}];$

(iv) multiplication,

$[v.ω]γ=[min{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ}, max{v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ}],$

(v) quotient,

$[min {v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ},max {v1γω1γ,v1γω2γ,v2γω1γ,v2γω2γ}], and$

(vi) equality, the two fuzzy numbers v and ω are equal if [v]γ = [ω]γ, that is, v1γ = ω1γ and v2γ = ω2γ.

Definition 2.4 [1021,53]

Letting v = [v1γ, v2γ] and ω = [ω1γ, ω2γ] be fuzzy numbers, and γ ∈ [0, 1], then z is called the H-difference of v and ω and is denoted by vω if there exists an element z ∈ ℝ such that v = ω + z for all γ ∈ [0, 1], and is given in levels by the following:

$[v⊖ω]γ=[v]γ-[ω]γ=[v1γ-ω1γ,v2γ-ω2γ].$

This definition was modified by Bede and Stefanini [23] to the generalized H-difference (gH-difference), as follows:

$v⊖gHω=z⇔v=ω+z or ω=v+(-1)z,$

with γ-levels:

$[v⊖gHω]γ=[min{v1γ-ω1γ,v2γ-ω2γ},max{v1γ-ω1γ,v2γ-ω2γ}].$

Very recently, Wasques et al. [17,27,29,59] generalized the Hukuhara and H-differences and used interactivities for all arithmetic operations and not only for the difference operation. We summarize this work as follows.

In general, a fuzzy subset v of ℝn is described by its membership function ρv: ℝn → [0, 1], where ρv(t) means the degree to which t belongs to v. The r-levels of fuzzy subset v are classical subsets defined as follows:

$[v]γ={t∈ℝn:ρv(t)≥γ},∀γ∈(0,1],$

and

$[v]0={t∈ℝn:ρv(t)>0}.$

Note that the fuzzy subset v of ℝ is a fuzzy number if its γ-levels are closed and nonempty intervals of ℝ and the support of v, supp(v) = {t ∈ ℝ: ρv(t) > 0}, is limited. The family of fuzzy subsets of ℝn with nonempty compact and convex γ-levels is denoted by $ℝℱn$, whereas the family of fuzzy numbers is denoted by ℝ.

Let v and ω ∈ ℝ and ℐ ∈ F(ℝ2). Then, the fuzzy relation ℐ is a joint possibility distribution of v and ω if

$maxt ρℐ(s,t)=ρv(t) and maxs ρℐ(s,t)=ρω(s),$

for all t and s ∈ ℝ.

$maxt ρℐ(s,t)=ρv(t) and maxs ρℐ(s,t)=ρω(s),$

for all t and s ∈ ℝ.

Here, v and ω are the marginal possibility distributions of ℐ. The fuzzy numbers v and ω are said to be non-interactive if and only if their joint possibility distribution ℐ is given by

$ρℐ(s,t)=min{ ρv(t),ρω(s)},$

for all t and s ∈ ℝ.

Otherwise, the fuzzy numbers are said to be interactive.

Definition 2.5 [17]

The Pompeiu-Hausdorff distance $d∞:ℝℱn×ℝℱn→ℝ+∪{0}$ is defined by

$d∞(v,ω)=supγ∈[0,1] dH([v]γ,[ω]γ),$

where dH is the Pompeiu-Hausdorff distance for the compact subsets of ℝn. If v and ω are fuzzy numbers, that is, v and ω ∈ ℝ, then Eq. (20) becomes

$d(v,ω)=supγ∈[0,1] max{∣v1γ-ω1γ∣,∣v2γ-ω2γ∣},$

and gives the Hausdorff distance between [v]γ and [ω]γ.

It is provided that (ℝ, d) is a complete metric space, which processes the following properties [53]:

(i) d(v + z, ω + z) = d(v, ω) for all v, z, ω ∈ ℝ,

(ii) d(μv, μω) = |μ|d(v, ω) for all v, ω ∈ ℝ and μ ∈ ℝ \ {0}, and

(iii) d(v + z, ω + e) ≤ d(v, ω) + d(z, e) for all v, z, ω, e ∈ ℝ.

Note that the metric space d in ℝ has many nice and interesting properties for an addition operation and scalar multiplication, some of which can be found in [19,5356,61].

Definition 2.6 [58,59]

Let [a, b] ⊆ ℝ be a vector space and ℝ be a set of fuzzy numbers. Then, the function û: [a, b] → ℝ is called an FVF on [a, b]. Corresponding to such a function û and γ ∈ [0, 1], we denote two real-valued functions û1γ(t) and û2γ(t) for all t ∈ [a, b] as

$[u^(t)]γ=[u^1γ(t),u^2γ(t)],$

which is called the γ-level representation of an FVF.

Definition 2.7 [20,53]

The complete metric structure on ℝ is given by the Hausdorff distance mapping d: × → ℝ+ ∪ {0} such that

$d(v,ω)=supt∈[0,1] max{∣v1(t)-ω1(t)∣,∣v2(t)-ω2(t)∣}.$
Definition 2.8 [53,55]

Let û1 and û2: [a, b] → ℝ be two FVFs. Then, the uniform distance between û1 and û2 is defined as

$d*(u^1(t),u^2(t))=supt∈[a,b]d(u^1(t),u^2(t)).$
Definition 2.9 [53,54,59]

Let û: [a, b] → ℝ be an FVF and t0 ∈ [a, b]. If ∀ ɛ > 0, ∃ δ > 0 such that d*(û (t), L) < ɛ, whenever |tt0| < δ, then we state that L ∈ ℝ is limit of û at t0, which is denoted by: limtt0û(t) = L. In addition, the FVF û is said to be continuous for t0 ∈ [a, b] if limtt0û(t) = û(t0). Indeed, û is continuous on [a, b] if it is continuous ∀ t ∈ [a, b]. In addition, if û: [a, b] → ℝ is a fuzzy continuous function, then we have [5356] the following:

$d(u^(t),0)=supγ∈[0,1] max{∣u^1γ(t)∣,∣u^2γ(t)∣},∀t∈[a,b].$
Definition 2.10 [20]

Let û: [a, b] → ℝ be an FVF. Then, for a fixed point t0 ∈ [a, b], û is called an SGD at t0 if there exists an element û′(t0) ∈ ℝ such that either of the following hold:

(i) the H-differences û (t0 +ɛ) ⊝ û (t0) and û (t0) ⊝ û (t0ɛ) exist, and

$u^′(t0)=D(1)u^(t0)=limɛ→0+u^(t0+ɛ)⊖u^(t0)ɛ=limɛ→0+u^(t0)⊖u^(t0-ɛ)ɛ,$

or all ɛ > 0 are sufficiently near 0 and the limits are in a metric d.

(ii) The H-differences û (t0) ⊝ û (t0+ɛ) and û (t0ɛ) ⊝ û (t0) exist, and

$u^′(t0)=D(2)u^(t0)=limɛ→0+u^(t0)⊖u^(t0+ɛ)-ɛ=limɛ→0+u^(t0-ɛ)⊖u^(t0)-ɛ,$

for all ɛ > 0 sufficiently near 0 and the limits are in a metric d.

Remark 2.1 [20]

(i) If û is differentiable for all t0 ∈ (a, b), then û is an SGD on (a, b).

(ii) The limits in Eqs. (26) and (27)) are considered in (ℝ, d), and at the edge-point of (a, b), we consider one direction derivative and recalling that as in Eq. (15).

Theorem 2.2 [22,53]

Let û: [a, b] → ℝ be an FVF, and let [û(t)]γ = [û1γ (t), û2γ (t)], ∀t ∈ [0, 1], and D(1)û(t) or D(2)û(t) exists. Then,

(i) if û is (1)-differentiable, then û1γ (t) and û2γ (t) are differentiable functions and

$[u^′(t)]γ=[u^1γ′(t), u^2γ′(t)],$

and (ii) if û is (2)-differentiable, then û1γ (t) and û2γ (t) are differentiable functions and

$[u^′(t)]γ=[u^2γ′(t), u^1γ′(t)].$
Theorem 2.3 [20,22]

Let û: [a, b] → ℝ be an FVF such that û is (1)-differentiable or û is (1)-differentiable at t0 ∈ (a, b). Subsequently, û is continuous at t0.

Definition 2.11 [17,62]

A strongly measurable and limited integrable fuzzy function is known as an integrable function. The fuzzy integral of an FVF û: [a, b] → ℝ with [û(t)]γ = [û1γ (t), û2γ (t)] is defined as follows:

$[Ju^(t)]γ=[∫abu^(t)dt]γ=∫ab[u^(t)]γdt=∫ab[u^1γ(t),u^2γ(t)]dt,$

for all γ ∈ [0, 1]. Provided Eq. (30), a fuzzy number is defined, and the integral exists in ℝ.

Similarly, the fuzzy integral of an FVF û: [a, b] → ℝ with [û(t)]γ = [û1γ (t), û2γ (t)] is defined with reference point t0 as

$[Jt0u^(t)]γ=[∫t0tu^(z)dz]γ=∫t0t[u^(z)]γdz=∫t0t[u^1γ(z),u^2γ(z)]dz,$

for all γ ∈ [0, 1]. Provided Eq. (31) define a fuzzy number and the integral exists in ℝ.

Theorem 2.4 [20]

Let û: [a, b] → ℝ be a FVF. Thus, the following are working:

(i) If û(t) is (1)-differentiable, then

$D(1)(Jt0u^) (t)=u^(t),$

and

$Jt0(D(1)u^(t))=u^(t)⊖u^(t0).$

(ii) If f is û(t) is (2)-differentiable, then.

$D(2)(Jt0u^) (t)-u^(t)=0,$

and

$u^(t0)=u^(t)-Jt0(D(2)u^(t)).$

### 2.2 A New Form to Taylor’s Formula

In this subsection, we present a new form of Taylor’s formula, which we will mainly use in our proposed method for solving FIVPs.

Theorem 2.5

Let u(t) be a piecewise continuous function on [0,∞) and be of exponential order. Suppose that the function U(s) = ℒ[u(t)] (Laplace transform of u(t)) can be expanded as follows:

$U(s)=∑n=0∞cnsn+1,s>0.$

Then cn = u(n)(0).

Proof

We assume that U(s) is a function that can be represented by the expansion in Eq. (34). First, note that if we multiply Eq. (34) by s and taking the limit for both sides as s → ∞, each term in the expansion after the first vanishes; thus, we obtain c0 = lims→∞ sU(s) = u(0), and Eq. (44) becomes

$U(s)=u(0)s+∑n=1∞cnsn+1,s>0.$

For the other aspect as well, multiplying Eq. (35) by s2 leads to the following expansion of c1:

$c1=s2U(s)-su(0)-c2s-c3s2-c4s3-⋯,s>0.$

It is clear that

$c1=lims→∞(s2U(s)-su(0))=lims→∞(s(sU(s)-u(0))=lims→∞s(ℒ[u′(t)] (s))=u′(0).$

If we multiply Eq. (35) by s3 and taking the limit as s → ∞, we obtain

$c2=lims→∞(s3U(s)-s2u(0)-su′(0))=lims→∞s(s2U(s)-su(0)-u′(0))=lims→∞s(ℒ[u′′(t)] (s))=u′′(0).$

Now, we can see the patterns and find the general formula for cn. However, if we continue to multiply Eq. (35) by sn+1 and compute the limit of the expansion of cn as s → ∞, we can obtain cn = u(n)(0), n = 0, 1, 2, · · · . This completes the proof of Theorem 2.5.

### Remark 2.2

The inverse Laplace transform of the expansion in Theorem 2.5 has the following form:

$u(t)=∑n=0∞u(n)(0)n!tn,t≥0,$

which is the classical Taylor’s formula.

The following theorem explains and determines the conditions for convergence of the new form of Taylor’s formula, which is presented in Theorem 2.5, as mentioned above.

### Theorem 2.6

Let U(s) = ℒ[u(t)](s). If |sℒ[u(n+1)(t)]|M, for 0 < sd, where M and d are positive numbers, then the remaining Rn(s) of the new form of Taylor’s formula in Eq. (34) satisfies the following inequality:

$∣Rn(s)∣ ≤Msn+2,0
Proof

First, suppose that ℒ[u(j)(t)](s) is defined as 0 < sd for j = 0, 1, 2, · · ·, n + 1. As given, assume also that

$∣sℒ[u(n+1)(t)]∣ ≤M,0

From the definition of the reminder $Rn(s)=U(s)-∑j=0nu(j)(0)sj+1$, one can obtain the following:

$sn+2Rn(s)=sn+2U(s)-∑j=0nsn+1-j(u(j)) (0)=s(sn+1U(s)-∑j=0nsn-j(u(j)) (0))=sℒ[u(n+1)(t)].$

It follows from Eqs. (41) and (42) that |sn+2Rn(s)|M. Hence,

$-M≤sn+2Rn(s)≤M,0

Thus, by reformulating Eq. (43), we can find the inequality $∣Rn(s)∣ ≤Msn+2$, which completes the proof.

### 2.3 Drafting of FIVP

In this subsection, we state the conditions and assumptions required to formulate the FIVP, as in Eq. (3). First, we point out that the FIVP is expressed in Eq. (3) yields a unique solution [25,45,63]. Second, we assume the γ-level depiction of û′(t), û(t), û2(t) and û(0) as [$u^1γ′(t), u^2γ′(t)$], [û1γ (t), û2γ (t)], [$u^1γ2(t), u^2γ2(t)$] and [û0,1γ (t), û0,2γ (t)], respectively. Thus, the FIVP, as in Eq. (3) can be reformulated as

$[u^′(t)]γ=A[u^2(t)]γ+B[u^(t)]γ+C,t>0,$

with the FIC

$[u^(0)]γ=[u^0]γ.$

### Definition 2.12 [1]

Let û: [a, b] → ℝ, where D(1)û(t) or D(2)û(t) exist. If û and D(1)û(t) satisfy FIVP in Eq. (3), we state that û is (1)-solution for the FIVP in Eq. (3). Otherwise, if û and D(2)û(t) satisfy FIVP, as in Eq. (3), we state that û is (2)-solution for the FIVP, as in Eq. (3).

Finally, to construct the LRPS solution for the target equation in a γ-level depiction, we alter the FIVP, as in Eqs. (44) and (45), into crisp ordinary differential equations (ODEs) systems and assume the following two cases to create the fuzzy solution û(t) for this FIVP regarding to the type of differentiability, where û(t) is either (1)-differentiable or (2)-differentiable.

Case 1

If û(t) is (1)-differentiable, then the IVP in Eq. (40)

$u1γ′(t)=Au^1γ2(t)+Bu^1γ(t)+C,u2γ′(t)=Au^2γ2(t)+Bu^2γ(t)+C,$

with FICs,

$u^1γ(0)=u^0,1γ,u^2γ(0)=u^0.2γ.$

We then have the following procedures:

(i) Construct the LRPS solution for the system in Eqs. (46) and (47).

(ii) Ensure that the solutions [û1γ (t), û2γ (t)] and [$u^1γ′(t), u^2γ′(t)$] are valid γ-level sets, ∀ γ ∈ [0, 1].

(iii) Then, [û1γ (t), û2γ (t)] are the (1)-solutions û in a γ-level depiction.

Case 2

If û(t) is (2)-differentiable, then the FIVP, as in Eqs. (44) and (45), can be transformed into the following ODEs system:

$u^1γ′(t)=Au^2γ2(t)+Bu^2γ(t)+C,u^2γ′(t)=Au^1γ2(t)+Bu^1γ(t)+C,$

with FICs,

$u^1γ(0)=u^0,1γ,u^2γ(0)=u^0.2γ.$

In addition, we have the following procedures:

(i) Construct the LRPS solution for systems (48) and (49).

(ii) Ensure that the solutions [û1γ (t), û2γ (t)] and [$u^1γ′(t), u^2γ′(t)$] are valid γ-level sets, ∀ γ ∈ [0, 1].

(iii) Then, [û2γ (t), û1γ (t)] are the (2)-solutions û in a γ-level depiction.

### 3. LRPS Solution for FIVP

In this section, we present the (1)-solution for the FIVP, as in Eqs. (44) and (45), using the LRPS method. In addition, the same process can be pursued each time û(t) is (2)-differentiable to create the (2)-solution for the FIVP, as in Eqs. (44) and (45). To achieve this, we assumed that û(t) is (1)-differentiable. First, the technique involves applying a Laplace transformation to both sides of the system in Eq. (46); hence,

$ℒ[u^1γ′(t)]=Aℒ[u^1γ2(t)]+Bℒ[u^1γ(t)]+ℒ[C],ℒ[u^2γ′(t)=Aℒ[u^2γ2(t)+Bℒ[u^2γ(t)]+ℒ[C].$

Using the properties of the Laplace transform and the initial conditions in Eq. (47), we can rewrite the system in Eq. (50) into the following algebraic system:

$U^1γ(s)=u^0,1γs+Asℒ[(ℒ-1[U^1γ(s)])2]+BsU^1γ(s)+Cs2,U^2γ(s)=u^0,2γs+Asℒ[(ℒ-1[U^2γ(s)])2]+BsU^2γ(s)+Cs2,$

where Û1γ (s) = ℒ[û1γ (t)], and Û2γ (s) = ℒ[û2γ (t)].

Next, the LRPS technique represents the solution as an infinite series, and the series solution to Eq. (51) is given by the following expansions:

$U^1γ(s)=∑n=0∞ann!sn+1,U^2γ(s)=∑n=0∞bnn!sn+1.$

Because

$lims→∞sU^1γ(s)=u^1γ(0)=u^0,1γ,$

and

$lims→∞sU^2γ(s)=u^2γ(0)=u^0,2γ,$

this leads to a0 = û0,1γ and b0 = û0,2γ.

Therefore, the expansions in Eq. (52) can be written as

$U^1γ(s)=u^0,1γs+∑n=1∞ann!sn+1,U^2γ(s)=u^0,2γs+∑n=1∞bnn!sn+1.$

Consequently, the kth-truncated series of Û1γ (s) and Û2γ (s) can be given by

$U^k,1γ(s)=u^0,1γs+∑n=1kann!sn+1,U^k,2γ(s)=u^0,2γs+∑n=1kbnn!sn+1.$

Prior to applying the LRPS technique to determine the values of the coefficients an and bn in the expansions in Eq. (54), we must define the Laplace residual function for Eq. (51) as

$LRes1γ(s)=U^1γ(s)-u^0,1γs-Asℒ[(ℒ-1[U^1γ(s)])2]-BsU^1γ(s)-Cs2s,LRes2γ(s)=U^2γ(s)-u^0,2γs-Asℒ[(ℒ-1[U^2γ(s)])2]-BsU^2γ(s)-Cs2,$

and the kth-Laplace residual function as follows:

$LResk,1γ(s)=U^k,1γ(s)-u^0,1γs-Asℒ[(ℒ-1[U^k,1γ(s)])2]-BsU^k,1γ(s)-Cs2,LResk,2γ(s)=U^k,2γ(s)-u^0,2γs-Asℒ[(ℒ-1[U^k,2γ(s)])2]-BsU^k,2γ(s)-Cs2.$

It is clear that $limj→∞ LResj,iγ(s)=LRes∞,iγ(s)=LResiγ(s)$ for each s > 0 and i = 1, 2.

Because $lims→∞sLResiγ(s)=0$, it is clear that

$lims→∞sLResk,iγ(s)=0,$

for i = 1, 2 and j = 1, 2, 3, · · ·, k, and thus $lims→∞(sj+1LRes∞,iγ(s))=lims→∞(sj+1LResj,iγ(s))=0, i=1,2,j=1,2,3,⋯,k.$

For the other aspects as well, to determine the value of the first unknown coefficients a1and b1 in Eq. (54), we should substitute Û1,1γ (s) = û0,1γ/s + a1/s2 and Û1,2γ (s) = û0,2γ/s + b1/s2 into the kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), at k = 1 in Eq. (56) to obtain the following:

$LRes1,1γ(s)= U^1,1γ(s)-u^0,1γs-Asℒ[(ℒ-1[U^1,1γ(s)])2] -BsU^1,1γ(s)-Cs2= u^0,1γs+a1s2-u^0,1γs -Asℒ[(ℒ-1[u^0,1γs+a1s2])2] -Bs(u^0,1γs+a1s2)-Cs2= a1s2-As(u^0,1γ2s+2u^0,1γa1s2+2a12s3) -Bs(u^0,1γs+a1s2)-Cs2,LRes1,2γ(s)= U^1,2γ(s)-u^0,2γs-Asℒ[(ℒ-1[U^1,2γ(s)])2] -BsU^1,2γ(s)-Cs2= u^0,2γs+b1s2-u^0,2γs -Asℒ[(ℒ-1[u^0,2γs+b1s2])2] -Bs(u^0,2γs+b1s2)-Cs2= b1s2-As(u^0,2γ2s+2u^0,2γa1s2+2b12s3) -Bs(u^0,2γs+b1s2)-Cs2.$

Now, multiplying both sides of Eq. (56) when k = 1 by s2 gives

$s2LRes1,1γ(s)=a1-A(u^0,1γ2+2u^0,1γa1s+2a12s2)-B(u^0,1γ+a1s)-C,s2LRes1,2γ(s)=b1-A(u^0,2γ2+2u^0,2γa1s+2b12s2)-B(u^0,2γ+b1s)-C.$

Computing the limits of both sides of Eq. (59) as s→∞gives the following expressions:

$lims→∞s2LRes1,1γ(s)=a1-Au^0,1γ2-Bu^0,1γ-C,lims→∞s2LRes1,2γ(s)=b1-Au^0,2γ2-Bu^0,2γ-C.$

Considering j = 1 in Eq. (57) using the expressions in Eq. (60) and solving the resulting algebraic system for a1 and b1, it is easy to obtain $a1=Au^0,1γ2+Bu^0,1γ+C$ and $b1=Au^0,2γ2+Bu^0,2γ+C$. Hence, we obtain the following:

$U^1,1γ(s)=u^0,1γs+(Au^0,1γ2+Bu^0,1γ+C)s2,U^1,2γ(s)=u^0,2γs+(Au^0,2γ2+Bu^0,2γ+C)s2.$

Similarly, to determine the values of the second unknown coefficients a2 and b2, we substitute $U^2,1γ(s)=a0/s+(Aa02+Ba0+C)/s2+2a2/s3$ and $U^2,2γ(s)=b0/s+(Ab02+Bb0+C)/s2+2b2/s3$ into LRes2,1γ (s) and LRes2,2γ (s) in Eq. (56) to obtain the following discretized form:

$LRes2,1γ(s)= a0s+Aa02+Ba0+Cs2+2a2s3-a0s -Asℒ[(ℒ-1[a0s+Aa02+Ba0+Cs2+2a2s3])2] -Bs(a0s+Aa02+Ba0+Cs2+2a2s3) -Cs2,LRes2,2γ(s)= b0s+Ab02+Bb0+Cs2+2b2s3-b0s -Asℒ[(ℒ-1[b0s+Ab02+Bb0+Cs2+2b2s3])2] -Bs(b0s+Ab02+Bb0+Cs2+2b2s3) -Cs2.$

By multiplying both sides of Eq. (56) for k = 2 by s3, we obtain

$s3LRes2,1γ(s)= 2a2-A(2(Aa03+Ba02+Ca0) +4a0a2s+2(Aa02+Ba0+C)2s +12(Aa02a2+Ba0a2+Ca2)s2+24a22s3) -B(Aa02+Ba0+C+2a2s),s3LRes2,2γ(s)= 2b2-A(2(Ab03+Bb02+Cb0) +4b0b2s+2(Ab02+Bb0+C)2s +12(Ab02b2+Bb0b2+Cb2)s2+24b22s3) -B(Ab02+Bb0+C+2b2s).$

Using the fact that $lims→∞sj+1LResj,iγ(s)=0$ in Eq. (57) for j = 2, i = 1, 2, we obtain the following:

$lims→∞s3LRes2,1γ(s) =2a2-A(2(Aa03+Ba02+Ca0)) -B(Aa02+Ba0+C) =2a2-2Aa0(Aa02+Ba0+C) -B(Aa02+Ba0+C)=0,lims→∞s3LRes2,2γ(s) =2b2-A(2(Ab03+Bb02+Cb0)) -B(Ab02+Bb0+C) =2b2-2Ab0(Ab02+Bb0+C) -B(Ab02+Bb0+C)=0.$

By solving the resultant algebraic system in Eq. (64) for a2 and b2, we obtain $a2=Aa0a1+12Ba1$ and $b2=Ab0b1+12Bb1$. Thus,

$U^1,1γ(s)=a0s+(Aa02+Ba0+C)s2+2(Aa0a1+12Ba1)s3,U^1,2γ(s)=b0s+(Ab02+Bb0+C)s2+2(Ab0b1+12Bb1)s3.$

For k = 3, if we substitute Û3,1γ (s) and Û3,2γ (s) into the Laplace residual functions LRes3,1γ (s) and LRes3,2γ (s) of Eq. (65) and the utilize the facts $lims→∞s4LRes3,1γ(s)=0$ and lims→∞ s4LRes3,2γ (s) = 0, then the third coefficients a3 and b3 are given by the following:

$a3=13A(2a0a2+a12)+13Ba2,b3=13A(2b0b2+b12)+13Bb2.$

Therefore, the third LRPS estimation can also be given.

By continuing the same process until random coefficients order k = n, and using the fact that

$lims→∞sn+1LResn,1γ(s)=lims→∞sn+1LResn,2γ(s)=0,$

subsequently an and bn (unknown coefficients) can be achieved. Nevertheless, more precise solutions can be obtained using more iterations. Likewise, if û(t) is (2)-differentiable, the (2)-solution for the FIVP, as in Eqs. (44) and (45), can be obtained.

Finally, the exact or approximate series solution obtained must be transformed into the original space by using the inverse Laplace transform to obtain the series solution to the original FIVP, as in Eqs. (44) and (45).

### 4. Numerical Examples

To validate the accuracy of our proposed algorithm, numerous examples are examined numerically in this section. For calculations and results, Maple 2018 is used to perform the computations.

### Example 4.1

Consider the following FIVP:

$u^′(t)=2u^(t)-u^2(t)+1,t>0,$

with FIC,

$[u^(0)]γ=[γ-1,1-γ],γ∈[0,1].$

For γ = 1, the exact solution of the FIVP, as in Eqs. (67) and (68), can be located as follows:

$u(t)=1+2 tanh (2t+12log(2-12+1)).$

Using Definition 2.10, the FIVP, as in Eqs. (67) and (68), can be represented by the set of ODEs corresponding to their parametric forms as

$u^1γ′(t)=2u^1γ(t)-u^1γ2(t)+1,u^2γ′(t)=2u^2γ(t)-u^2γ2(t)+1,$

subject to

$u^1γ(0)=γ-1,u^2γ(0)=1-γ.$

Applying the LRPS method as below, taking the Laplace transform of both sides of Eq. (70) gives

$sℒ[u^1γ(t)]-u^1γ(0)=2ℒ[u^1γ(t)]-ℒ[u^1γ2(t)]+ℒ[1],sℒ[u^2γ(t)]-u^2γ(0)=2ℒ[u^2γ(t)]-ℒ[u^2γ2(t)]+ℒ[1],$

and Eq. (71) then leads to

$U^1γ(s)=γ-1s+2sU^1γ(s)-1sℒ[(ℒ-1[U^1γ(s)])2]+1s2,U^2γ(s)=1-γs+2sU^2γ(s)-1sℒ[(ℒ-1[U^2γ(s)])2]+1s2,$

where Û1γ (s) = ℒ[û1γ (t)] and Û2γ (s) = ℒ[û2γ (t)].

Now, we assume that an infinite series solution of the algebraic system in Eq. (73) has the following form:

$U^1γ(s)=∑n=0∞ann!sn+1,U^2γ(s)=∑n=0∞bnn!sn+1$

In addition, the kth-truncated series of Eq. (74) is expressed as follows:

$U^k,1γ(s)=u^0,1γs+∑n=1kann!sn+1,U^k,2γ(s)=u^0,2γs+∑n=1kbnn!sn+1.$

Depending on the initial data, $lims→∞sU^1γ(s)=u^1γ(0)=γ-1$ and $lims→∞sU^1γ(s)=u^2γ(0)=1-γ$. Therefore, the series solution to Eq. (73) becomes

$U^k,1γ(s)=γ-1s+∑n=1kann!sn+1,U^k,2γ(s)=1-γs+∫n=1kbnn!sn+1.$

The kth Laplace residual functions, LResk,1γ (s) and LResk,2γ (s) of (73), are as follows:

$LResk,1γ(s)=U^k,1γ(s)-γ-1s-2sU^k,1γ(s)+1sℒ[(ℒ-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-1-γs-2sU^k,2γ(s)+1sℒ[(ℒ-1[U^k,2γ(s)])2]-1s2.$

To determine the first unknown coefficients, a1 and b1, in the expansion (76), we substitute the first-truncated series Û1,1γ (s) = (γ – 1)/s + a1/s2 and Û1,2γ (s) = (1 – γ)/s + b1/s2 into the first Laplace residual functions, LRes1,1γ (s) and LRes1,2γ (s) of Eq. (77), and obtain

$LRes1,1γ(s)= a1s2-2s(γ-1s+a1s2) +1s((γ-1)2s+2(γ-1)a1s2+2a12s3) -1s2,LRes1,2γ(s)= b1s2-2s(1-γs+b1s2) +1s((1-γ)2s+2(1-γ)b1s2+2b12s3) -1s2.$

Now, by remultiplying both sides of Eq. (77) by sk+1, for k = 1, we obtain

$s2LRes1,1γ(s)= a1-2(γ-1+a1s) +((γ-1)2+2(γ-1)a1s+2a12s2) -1,s2LRes1,2γ(s)= b1-2(1-γ+b1s) +((1-γ)2+2(1-γ)a1s+2b12s2) -1.$

Using the fact that $lims→∞sk+1LResk(s)=0$ for k = 1 in Eq. (79) yields a1 = −(γ – 1)2 + 2(γ – 1) + 1 and b1 = −(1 – γ)2 +2(1–γ)+1. Hence, the first-approximate LRPS solution to Eq. (73) can be expressed as

$U^1,1γ(s)=(γ-1)s+(-2+4γ-γ2)s2,U^1,2γ(s)=(1-γ)s+(2-γ2)s2.$

Similarly, to determine the values of the second unknown coefficients a2 and b2, we substitute the second-truncated series Û2,1γ (s) = (γ – 1)/s + (−2 + 4γγ2)/s2 + 2a2/s3, and Û2,2γ (s) = (1–γ)/s+(2–γ2)/s2 +2b2/s3 into the second-Laplace residual functions LRes2,1γ (s) and LRes2,2γ (s) of Eq. (77) to obtain the following discretized form:

$LRes2,1γ(s)= (-2+2γ-γ2)s2+2a2s3 -2s(γ-1s+(-2+4γ-γ2)s2+2a2s3) +1s(24a22s5+(4γ-4s3+-12γ2+48γ-24s4)a2 +γ2-2γ+1s+-2γ3+10γ2-12γ+4s2 +2γ4-16γ3+40γ2-32γ+8s3)-1s2,LRes2,2γ(s)= (2-γ2)s2+2b2s3 -2s(1-γs+(2-γ2)s2+2b2s3) +1s(24b22s5+(-4γ+4s3+-12γ2+24s4)b2 +γ2-2γ+1s+2γ3-2γ2-4γ+4s2 +2γ4-8γ2-32γ+8s3)-1s2.$

Now, re-multiplying s3 by LRes2,1γ (s) and LRes2,2γ (s) in Eq. (81) yields.

$s3LRes2,1γ(s)= 24a22s3 +(2+4γ-8s+-12γ2+48γ-24s2)a2 -2γ3+12γ2-20γ+8 +2γ4-16γ3+40γ2-32γ+8s,s3LRes2,2γ(s)= 24b22s3 +(2-4γs+-12γ2+24s2)b2 +2γ3-4γ+2γ4-8γ2+8s.$

According to the fact: $lims→∞sk+1LResk(s)=0$ for k = 2, then we have the following algebraic equations:

$lims→∞s3LRes2,1γ(s)=2a2-2γ3+12γ2-20γ+8=0,lims→∞s3LRes2,2γ(s)=2b2+2γ3-4γ=0.$

Solving Eq. (83) for a2 and b2, the following discretized forms are obtained:

$a2=γ3-6γ2+10γ-4,b2=-γ3+2γ.$

Hence, the second-approximate LRPS solution to Eq. (73) is expressed as

$U^2,1γ(s)=(γ-1)s+(-2+4γ-γ2)s2+2(γ3-6γ2+10γ-4)s3,U^2,2γ(s)=(1-γ)s+(2-γ2)s2+2(-γ3+2γ)s3.$

The kth-LRPS approximation can be obtained sequentially by using the same procedure. For further documentation, we will report the first two coefficients of ak and bk and compute the new component forms for future use in the succeeding approximation to obtain the following:

$a1=-2+4γ-γ2,a2=γ3-6γ2+10γ-4,a3=13(-20+64γ-64γ2+24γ3-3γ4),a4=13(-32+128γ-180γ2+110γ3-30γ4+3γ5),a5=115(-256+1232γ-2228γ2+1920γ3-840γ4+180γ5-15γ6),a6=145(-1232+6920γ-15288γ2+17108γ3-10500γ4+3570γ5-630γ6+45γ7),a7=1315(-13840+88832γ-231872γ2+319872γ3-255024γ4+120960γ5-33600γ6+5040γ7-315γ8),$

and

$b1=2-γ2,b2=-γ3+2γ,b3=13(-4+8γ2-3γ4),b4=13(-8γ+10γ3-3γ5),b5=115(16-68γ2+60γ4-15γ6),b6=145(136γ-308γ3+210γ5-45γ7),b7=1315(-272+1984γ2-3024γ4+1680γ6-315γ8).$

The LRPS solution to Eq. (73) in an infinite series is as follows:

$U^1γ(s)=(γ-1)s+(-2+4γ-γ2)s2+2(γ3-6γ2+10γ-4)s3+⋯,U^2γ(s)=(1-γ)s+(2-γ2)s2+2(-γ3+2γ)s3+⋯.$

By applying the inverse Laplace transform to Eq. (87), we determine the LRPS solution of Eqs. (70) and (71) as follows:

$u^1γ(t)=(γ-1)+(-2+4γ-γ2)t+(γ3-6γ2+10γ-4)t2+⋯,u^2γ(t)=(1-γ)+(2-γ2)t+(-γ3+2γ)t2+⋯.$

As a special case, when γ = 1, the LRPS solutions to Eqs. (87) and (88) have a general pattern form, which coincides with the exact solution in terms of the infinite series

$u(t)=t+t2+13t3-12t4-115t5-745t6+53315t7+71315t8+⋯.$

Table 1 shows the numerical results of the exact and approximate solutions using the LRPS method for Example 4.1 at γ = 1. Moreover, the relative and absolute errors at some chosen grid points ti in [0,1] with a step-size of 0.1 are also given in Table 1 for n = 51 and γ = 1. The results indicate that the LRPS approximate solutions are almost as good as the exact solutions. However, a more accurate solution can be obtained using more iterations.

Table 2 shows a numerical comparison between the tenth-LRPS solution and the OHA method [2], the MNN method [3], and the RK-4 method [27]. We concluded that the results obtained using the LRPS method were similar to those obtained using other methods.

Figure 1 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8 with diverse values of (t ∈ {0.0, 0.25, 0.5, 0.75}). Figure 2 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

### Example 4.2

$u^′(t)=2u^(t)-u^2(t)+1,t>0,$

with FIC,

$[u^(0)]γ=[0.1+0.1γ,0.3-0.1γ],γ∈[0,1].$

For γ = 1, the initial condition in Eq. (91) becomes u(0) = 0.2. Therefore, the exact solution to Eq. (90) can be expressed as

$u(t)=1+2 tanh (2t+12log(52-452+4)).$

According to Definition 2.10, the FIVPs, as in Eqs. (90) and (91), can be condensed to the set of ODEs corresponding to their parametric forms as follows:

$u^1γ′(t)=2u^1γ(t)-u^1γ2(t)+1,u^2γ′(t)=2u^2γ(t)-u^2γ2(t)+1.$

subject to

$u^1γ(0)=0.1+0.1γ,u^2γ(0)=0.3-0.1γ.$

According to the procedure of the LRPS algorithm presented in Section 3, the Laplace transform of Eq. (90) can be expressed as

$U^1γ(s)=0.1+0.1γs+2sU^1γ(s)-1sℒ[(ℒ-1[U^1γ(s)])2]+1s2,U^2γ(s)=0.3-0.1γs+2sU^2γ(s)-1sℒ[(ℒ-1[U^2γ(s)])2]+1s2,$

where Û1γ (s) = ℒ[û1γ (t)], and Û2γ (s) = ℒ[û2γ (t)].

Depending on the initial conditions in Eq. (94), the proposed LRPS solutions Û1γ (t) and Û2γ (t) of the algebraic system in Eq. (95) can be given by

$U^k,1γ(s)=0.1+0.1γs+∑n=1kann!sn+1,U^k,2γ(s)=0.3-0.1γs+∑n=1kbnn!sn+1,$

and the kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), as

$LResk,1γ(s)=U^k,1γ(s)-0.1+0.1γs-2sU^k,1γ(s)+1sℒ[(ℒ-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-0.3-0.1γs-2sU^k,2γ(s)+1sℒ[(ℒ-1[U^k,2γ(s)])2]-1s2.$

By utilizing lims→∞(sk+1LResk,1γ (s))=0 and lims→∞(sk+1LResk,2γ (s)) = 0, for k = 1, 2, · · ·, the first few terms ak and bk are

$a1=1102(119+18γ-γ2),a2=1103(1071+43γ-27γ2+γ3),a3=13×104(5117-5652γ-658γ2+108γ3-3γ4),a4=13×105(-168147-64585γ+5130γ2+1430γ3-135γ4+3γ5),a5=13×106(-1537123+11682γ+151955γ2+540γ3-2445γ4+162γ5-3γ6),a6=19×107(695079+18187783γ+2825739γ2-719285γ3-40635γ4+11109γ5-567γ6+9γ7),$

and

$b1=1102(151-14γ-γ2),b2=1103(1057+53γ-21γ2-γ3),b3=13×104(-8003+7084γ-82γ2-84γ3-3γ4),b4=13×105(-267421+30985γ+10710γ2-470γ3-105γ4-3γ5),b5=13×106(-935747-560126γ+108755γ2+13020γ3-1005γ4-126γ5-3γ6),b6=19×107(42289513-20342887γ-1706523γ2+685685γ3+38955γ4-5061γ5-441γ6-9γ7).$

We can express the LRPS solution of Eq. (95) in an infinite series as follows:

$U^1γ(s)=(0.1+0.1γ)s+1102(119+18γ-γ2)s2+2103(1071+43γ-27γ2+γ3)s3+⋯,U^2γ(s)=(0.3-0.1γ)s+1102(151-4γ-γ2)s2+2103(1057+53γ-21γ2-γ3)s3+⋯.$

Applying the inverse Laplace transform to Eq. (99), the LRPS solutions of Eqs. (93) and (94) are then given by

$u^1γ(t)=(0.1+0.1γ)+1102(119+18γ-γ2)t+1103(1071+43γ-27γ2+γ3)t2+⋯,u^2γ(t)=(0.3-0.1γ)+1102(151-4γ-γ2)t+1103(1057+53γ-27γ2-γ3)t2+⋯.$

As a special case, when γ = 1, the LRPS solution of Eqs. (90) and (91) have a general pattern form, which coincides with the exact solution in terms of infinite series:

$u(t)= 15+3425t+136125t2-681875t3-70729365t4 -2148846875t5+163744703125t6+1145963224609375t7 +16217728123046875t8+⋯.$

Tables 3 and 4 present the numerical results of the fuzzy solution upper and lower bounds of the system in Eqs. (93) and (94) at the selected grid points with γ = 1 and n = 51. We observed that for γ = 1, the LRPS upper- and lower-bound solutions are the same and coincide with the crisp solution.

Figure 3 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8with diverse values of (t ∈ {0.25, 0.5, 0.75, 1}), where the center (crisp) solution is represented by the midline for γ = 1. Figure 4 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

### Example 4.3

Consider the following FIVP:

$u^′(t)=-u^2(t)+1,t>0,$

with FIC,

$[u^(0)]γ=[γ-1,1-γ],γ∈[0,1].$

For γ = 1 and the crisp initial condition u(0) = 0, the exact solution of (102) is located as

$u(t)=(e2t-1) (e2t+1)-1.$

The FIVP, as in Eqs. (102), and (103), can be reduced to the set of ODEs corresponding to their parametric forms as

$u^1γ′(t)=-u^1γ2(t)+1,u^2γ′(t)=-u^2γ2(t)+1,$

subject to

$u^1γ(0)=γ-1,u^2γ(0)=1-γ.$

Similar to the steps used in previous examples, the Laplace transform of Eq. (102) can be expressed as follows:

$U^1γ(s)=γ-1s-1sℒ[(ℒ-1[U^1γ(s)])2]+1s2,U^2γ(s)=1-γs-1sℒ[(ℒ-1[U^2γ(s)])2]+1s2.$

According to the initial data in Eq. (106), the LRPS solutions, Û1γ (t) and Û2γ (t) of the system in Eq. (107) can be represented as

$U^k,1γ(s)=γ-1s+∑n=1kann!sn+1,U^k,2γ(s)=1-γs+∑n=1kbnn!sn+1.$

The kth-Laplace residual functions, LResk,1γ (s) and LResk,2γ (s), can be defined as follows:

$LResk,1γ(s)=U^k,1γ(s)-γ-1s+1sℒ[(ℒ-1[U^k,1γ(s)])2]-1s2,LResk,2γ(s)=U^k,2γ(s)-1-γs+1sℒ[(ℒ-1[U^k,2γ(s)])2]-1s2.$

By utilizing the following facts: $lims→∞(sk+1LResk,1γ(s))=0$ and lims→∞(sk+1LResk,2γ (s)) = 0, for k = 1, 2, · · ·, the first few coefficients ak and bk are

$a1=(2γ-γ2),a2=(-2γ+3γ2-γ3),a3=13(4γ-14γ2+12γ3-3γ4),a4=13(-2γ+15γ2-25γ3+15γ4+3γ5),a5=115(4γ-62γ2+180γ3-195γ4+90γ5-15γ6),a6=145(-4γ+126γ2-602γ3+1050γ4-840γ5+315γ6-45γ7),$

and

$b1=(2γ-γ2),b2=(2γ-3γ2+γ3),b3=13(4γ-14γ2+12γ3-3γ4),b4=13(2γ-15γ2+25γ3-15γ4+3γ5),b5=115(4γ-62γ2+180γ3-195γ4+90γ5-15γ6),b6=145(4γ-126γ2+602γ3-1050γ4+840γ5-315γ6+45γ7).$

We can express the LRPS solution of Eq. (107) as in an infinite series by

$U^1γ(s)=(γ-1)s+(2γ-γ2)s2+2(-2γ+3γ2-γ3)s3+⋯,U^2γ(s)=(1-γ)s+(2γ-γ2)s2+2(2γ-3γ2+γ3)s3+⋯.$

Applying the inverse Laplace transform to Eq. (111), we then obtain the LRPS solution of Eqs. (105) and (106) as

$u^1γ(t)=(γ-1)+(2γ-γ2)t+(-2γ+3γ2-γ3)t2+⋯,u^2γ(t)=(1-γ)+(2γ-γ2)t+(2γ-3γ2+γ3)t2+⋯.$

As a special case, when γ = 1, the LRPS solutions to Eqs. (102) and (103) have a general pattern form, which coincides with the exact solution in terms of infinite series:

$u(t)=t-13t3+215t5-17315t7+622835t9+⋯.$

Table 5 presents the error analysis of the proposed algorithm for the system in Eqs. (105) and (106) at some nodes t in [0, 1] with a step size of 0.16 for γ = 1 and n = 51. Table 6 presents a numerical comparison between the tenth-LRPS solution and RK-4, OHAM, and MNN given in the literature. We observed that the numerical results of LRPS were as good as those obtained using other techniques.

Figure 5 shows the upper and lower bounds of the triangular fuzzy LRPS solution at n = 8 with diverse values of t (t ∈ {0.25, 0.5, 0.75, 1}), where the center (crisp) solution is represented by the midline for γ = 1. Figure 6 shows the tenth-LRPS approximate solution surface plot for all t ∈ [0, 1] and γ ∈ [0, 1]. The blue color corresponds to the upper bounds of the tenth-LRPS fuzzy solution, whereas the yellow color corresponds to the lower bounds of the solution.

### 5. Conclusion

Many numerical and analytical methods have been applied to solve fuzzy differential equations, some of which have advantages over others. Some of them are accurate and effective; however, they require mathematical operations that can be difficult, long, or sometimes failed. Others are simple and fast; however, they may not provide precise solutions. Our new method, called LRPS, is characterized by its accuracy, speed, and simplicity in finding exact and accurate approximate solutions of linear and nonlinear fuzzy differential equations. In this study, we succeeded in providing accurate approximate and sometimes exact solutions to FIVPs using the newly proposed technique. The proposed expansion of our study enabled us to obtain a serial solution for equations in the Laplace transform space. The LRPS method provides an easy and fast technique for determining the proposed series coefficients as a solution to the equation. Unlike the traditional RPS method, which specifies series coefficients, the derivative must be calculated each time, whereas LRPS only requires the concept of the limit at infinity in determining series coefficients, which distinguishes this method. It is worth noting that the LRPS method can be applied to solve other types of fuzzy ordinary and partial differential equations for an integer or fractional order, such as fuzzy Poisson, KdV, Schrödinger, and apple equations. Furthermore, based on the H-SGD theory, it is possible to extend and generalize our analysis and results to all fuzzy interactive arithmetic operations mentioned in Section 2 in future studies.

### Fig 1.

Figure 1.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.1. Solution plots for (a) t = 0, (b) t = 0.25, (c) t = 0.5, and (d) t = 0.75.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

### Fig 2.

Figure 2.

The 3-dim plot for Example 4.1. Blue and yellow are the lower and upper bounds, respectively, of the tenth-LRPS fuzzy solution.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

### Fig 3.

Figure 3.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.2. Solution plots for (a) t = 0.25, (b) t = 0.5, (c) t = 0.75, and (d) t = 1.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

### Fig 4.

Figure 4.

3-dim plot for Example 4.2. Blue and yellow are the lower and upper bounds, respectively, of the tenth-LRPS fuzzy solution.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

### Fig 5.

Figure 5.

Triangular fuzzy solution plots for the eighth-LRPS solutions of Example 4.3. Solution plots for (a) t = 0, (b) t = 0.25, (c) t = 0.5, and (d) t = 0.75.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

### Fig 6.

Figure 6.

3-dim plot for Example 4.3. Blue and yellow are lower and upper bounds, respectively, of tenth-LRPS fuzzy solution.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 23-47https://doi.org/10.5391/IJFIS.2022.22.1.23

Numerical results of the LRPS solutions for Example 4.1 at γ = 1.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.1102951969160.1102951969171.52656 × 10−161.38406 × 10−15
0.20.2419767996210.2419767996215.55112 × 10−172.29407 × 10−16
0.30.3951048486600.3951048486602.77556 × 10−167.02486 × 10−16
0.40.5678121662930.5678121662933.33067 × 10−165.86579 × 10−16
0.50.7560143934310.7560143934313.33067 × 10−164.40556 × 10−11
0.60.9535662164720.9535662164721.11022 × 10−161.16429 × 10−16
0.71.1529489669791.1529489669805.32907 × 10−144.62212 × 10−14
0.81.3463636553681.3463636554225.39830 × 10−114.00954 × 10−11
0.91.5269113132801.5269113368462.35658 × 10−81.54336 × 10−8)

Approximate solutions of Example 4.1 by various methods at γ = 1.

tiExact10th LRPSOHAMNNRK-4
0.10.1102950.1102950.1103280.1102950.100000
0.20.2419770.2419770.2422730.2419760.219000
0.30.3951050.3951050.3961750.3950890.358004
0.40.5678120.5678120.5702310.5676600.516788
0.50.7560140.7560140.759550.7551340.693439
0.60.9535660.9535660.9550490.9499640.884041
0.71.1529491.1529491.1424441.1414231.082696
0.81.3463641.3463581.3005691.3157231.282012
0.91.5269111.5268141.4004441.4565451.474059

Numerical results of lower bound LRPS-solution of Example 4.2.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.346763995050.3467639950591.11022 × 10−163.20167 × 10−16
0.20.513897275840.5138972758411.11022 × 10−162.16040 × 10−16
0.30.697990872510.6979908725142.22045 × 10−163.18120 × 10−16
0.40.893472160240.8934721602443.33067 × 10−163.72778 × 10−16
0.51.093134993061.0931349930602.22045 × 10−162.03126 × 10−16
0.61.289137171941.2891371719472.22045 × 10−161.72243 × 10−16
0.71.474194376961.4741943769631.05871 × 10−127.18161 × 10−13
0.81.642601362551.6426016335771.02182 × 10−96.22073 × 10−10
0.91.790797918331.7907983475064.29173 × 10−72.39655 × 10−7

Numerical results of upper bound LRPS-solution of Example 4.2.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.10.346763995050.3467639950591.66533 × 10−164.80250 × 10−16
0.20.513897275840.5138972758411.11022 × 10−162.16040 × 10−16
0.30.697990872510.6979908725142.22045 × 10−163.18120 × 10−16
0.40.893472160240.8934721602443.33067 × 10−163.72778 × 10−16
0.51.093134993061.0931349930602.22045 × 10−162.03126 × 10−16
0.61.289137171941.2891371719472.22045 × 10−161.72243 × 10−16
0.71.474194376961.4741943769631.05871 × 10−127.18161 × 10−13
0.81.642601362551.6426016335771.02182 × 10−96.22073 × 10−10
0.91.790797918331.7907983475064.29173 × 10−72.39655 × 10−7

Numerical results of the LRPS solutions of Example 4.3 at γ = 1.

tiExact solutionApproximate solutionAbsolute errorRelative error
0.160.1586485042970.1586485042975.55112 × 10−173.49900 × 10−16
0.3203095069212130.3095069212135.55112 × 10−171.79354 × 10−16
0.480.4462436102490.4462436102495.55112 × 10−171.24397 × 10−16
0.640.5648995528460.5648995528460.000000.00000
0.800.6640367702680.6640367702683.33067 × 10−165.01579 × 10−16
0.960.7442768673620.7442768673584.29645 × 10−125.77265 × 10−12

Numerical comparison of Example 4.3 at γ = 1.

tiExact10th-LRPSOHAMMNNRK-4
0.10.0996680.0996680.0996680.0996680.100000
0.20.1973750.1973750.1973760.1973750.199000
0.30.2913130.2913130.2913150.2913130.295040
0.40.3799490.3799490.3799490.3799490.386335
0.50.4621170.4621210.4620920.4621210.471410
0.60.5370500.5370780.5369100.5370770.549187
0.70.6043680.6045140.6038150.6045130.619026
0.80.6640370.6646410.6622450.6646400.680707
0.90.7162980.7183920.7112870.7183900.734371

### References

1. Oberguggenberger, M, and Pittschmann, S (1999). Differential equations with fuzzy parameters. Mathematical and Computer Modelling of Dynamical Systems. 5, 181-202. https://doi.org/10.1076/mcmd.5.3.181.3683
2. Gumah, G, Naser, MFM, Al-Smadi, M, Al-Omari, SKQ, and Baleanu, D (2020). Numerical solutions of hybrid fuzzy differential equations in a Hilbert space. Applied Numerical Mathematics. 151, 402-412. https://doi.org/10.1016/j.apnum.2020.01.008
3. Alaroud, M, Al-Smadi, M, Rozita Ahmad, R, and Salma Din, UK (2019). An analytical numerical method for solving fuzzy fractional Volterra integro-differential equations. Symmetry. 11. article no. 205
4. Al-Smadi, M (2019). Reliable numerical algorithm for handling fuzzy integral equations of second kind in Hilbert spaces. Filomat. 33, 583-597. https://doi.org/10.2298/FIL1902583A
5. Hasan, S, El-Ajou, A, Hadid, S, Al-Smadi, M, and Momani, S (2020). Atangana-Baleanu fractional framework of reproducing kernel technique in solving fractional population dynamics system. Chaos, Solitons & Fractals. 133. article no. 109624
6. Ahmadian, A, Salahshour, S, and Chan, CS (2017). Fractional differential systems: a fuzzy solution based on operational matrix of shifted Chebyshev polynomials and its applications. IEEE Transactions on Fuzzy Systems. 25, 218-236. https://doi.org/10.1109/TFUZZ.2016.2554156
7. El-Ajou, A, Oqielat, MAN, Ogilat, O, Al-Smadi, M, Momani, S, and Alsaedi, A (2019). Mathematical model for simulating the movement of water droplet on artificial leaf surface. Frontiers in Physics. 7. article no. 132
8. El-Ajou, A, Oqielat, MN, Al-Zhour, Z, and Momani, S (2019). Analytical numerical solutions of the fractional multi-pantograph system: two attractive methods and comparisons. Results in Physics. 14. article no. 102500
9. El-Ajou, A, Al-Zhour, Z, Oqielat, MA, Momani, S, and Hayat, T (2019). Series solutions of nonlinear conformable fractional KdV-Burgers equation with some applications. The European Physical Journal Plus. 134. article no. 402
10. El-Ajou, A, Oqielat, MAN, Al-Zhour, Z, Kumar, S, and Momani, S (2019). Solitary solutions for time-fractional nonlinear dispersive PDEs in the sense of conformable fractional derivative. Chaos: An Interdisciplinary Journal of Nonlinear Science. 29. article no. 093102
11. Oqielat, MAN, El-Ajou, A, Al-Zhour, Z, Alkhasawneh, R, and Alrabaiah, H (2020). Series solutions for nonlinear time-fractional Schrödinger equations: comparisons between conformable and Caputo derivatives. Alexandria Engineering Journal. 59, 2101-2114. https://doi.org/10.1016/j.aej.2020.01.023
12. El-Ajou, A (2020). Taylor’s expansion for fractional matrix functions: theory and applications. Journal of Mathematics and Computer Science. 21, 1-17. https://doi.org/10.22436/jmcs.021.01.01
13. Shqair, M, El-Ajou, A, and Nairat, M (2019). Analytical solution for multi-energy groups of neutron diffusion equations by a residual power series method. Mathematics. 7. article no. 633
14. Santo Pedro, F, de Barros, LC, and Esmi, E (2019). Population growth model via interactive fuzzy differential equation. Information Sciences. 481, 160-173. https://doi.org/10.1016/j.ins.2018.12.076
15. Chang, SSL, and Zadeh, LA (1972). On fuzzy mapping and control. IEEE Transactions on Systems, Man, and Cybernetics. 2, 30-34. https://doi.org/10.1109/TSMC.1972.5408553
16. Dubois, D, and Prade, H (1982). Towards fuzzy differential calculus. Part 1: integration of fuzzy mappings. Fuzzy Sets and Systems. 8, 1-17. https://doi.org/10.1016/0165-0114(82)90025-2
17. Wasques, V, Laiate, B, Santo Pedro, F, Esmi, E, and Barros, LCD 2020. Interactive fuzzy fractional differential equation: application on HIV dynamics. Information Processing and Management of Uncertainty in Knowledge-Based Systems, Array, pp.198-211. https://doi.org/10.1007/978-3-030-50153-2_15
18. Gumah, GN, Naser, MF, Al-Smadi, M, and Al-Omari, SK (2018). Application of reproducing kernel Hilbert space method for solving second-order fuzzy Volterra integro-differential equations. Advances in Difference Equations. 2018. article no. 475
19. Goetschel, R, and Voxman, W (1986). Elementary fuzzy calculus. Fuzzy Sets and Systems. 18, 31-43. https://doi.org/10.1016/0165-0114(86)90026-6
20. Arqub, OA, and Al-Smadi, M (2020). Fuzzy conformable fractional differential equations: novel extended approach and new numerical solutions. Soft Computing. 24, 12501-12522. https://doi.org/10.1007/s00500-020-04687-0
21. Seikkala, S (1987). On the fuzzy initial value problem. Fuzzy Sets and Systems. 24, 319-330. https://doi.org/10.1016/0165-0114(87)90030-3
22. Bede, B, and Gal, SG (2004). Almost periodic fuzzy-number-valued functions. Fuzzy Sets and Systems. 147, 385-403. https://doi.org/10.1016/j.fss.2003.08.004
23. Bede, B, and Stefanini, L (2013). Generalized differentiability of fuzzy-valued functions. Fuzzy Sets and Systems. 230, 119-141. https://doi.org/10.1016/j.fss.2012.10.003
24. Alshammari, M, Al-Smadi, M, Alshammari, S, Arqub, OA, Hashim, I, and Alias, M (2020). An attractive analytic-numeric approach for the solutions of uncertain Riccati differential equations using residual power series. Applied Mathematics and Information Sciences. 14, 177-190. https://doi.org/10.18576/amis/140202
25. Behzadi, SS, Vahdani, B, Allahviranloo, T, and Abbasbandy, S (2016). Application of fuzzy Picard method for solving fuzzy quadratic Riccati and fuzzy Painlevé I equations. Applied Mathematical Modelling. 40, 8125-8137. https://doi.org/10.1016/j.apm.2016.05.003
26. Arqub, OA, Al-Smadi, M, Momani, S, and Hayat, T (2017). Application of reproducing kernel algorithm for solving second-order, two-point fuzzy boundary value problems. Soft Computing. 21, 7191-7206. https://doi.org/10.1007/s00500-016-2262-3
27. Wasques, VF, Esmi, E, Barros, LC, and Bede, B 2019. Comparison between numerical solutions of fuzzy initial-value problems via interactive and standard arithmetics. Fuzzy Techniques: Theory and Applications, Array, pp.704-715. https://doi.org/10.1007/978-3-030-21920-8_62
28. Arqub, OA, El-Ajou, A, Momani, S, and Shawagfeh, N (2013). Analytical solutions of fuzzy initial value problems by HAM. Applied Mathematics & Information Sciences. 7. article no. 1903
29. Wasques, VF, Esmi, E, Barros, LC, and Sussner, P (2018). Numerical solutions for bidimensional initial value problem with interactive fuzzy numbers. Fuzzy Information Processing. Cham, Switzerland: Springer, pp. 84-95 https://doi.org/10.1007/978-3-319-95312-0_8
30. Hasan, S, Al-Smadi, M, El-Ajou, A, Momani, S, Hadid, S, and Al-Zhour, Z (2021). Numerical approach in the Hilbert space to solve a fuzzy Atangana-Baleanu fractional hybrid system. Chaos, Solitons & Fractals. 143. article no. 110506
31. Padmapriya, V, Kaliyappan, M, and Manivannan, A (2020). Numerical solutions of fuzzy fractional delay differential equations. International Journal of Fuzzy Logic and intelligent Systems. 20, 247-254. https://doi.org/10.5391/IJFIS.2020.20.3.247
32. Ahmad, J, Iqbal, A, and Ul Hassan, QM (2021). Study of nonlinear fuzzy integro-differential equations using mathematical methods and applications. International Journal of Fuzzy Logic and Intelligent Systems. 21, 76-85. https://doi.org/10.5391/IJFIS.2021.21.1.76
33. Kwun, YC, Han, CW, Kim, SY, and Park, JS (2004). Existence and uniqueness of fuzzy solutions for the nonlinear fuzzy integro-differential equation on E. International Journal of Fuzzy Logic And Intelligent Systems. 4, 40-44. https://doi.org/10.5391/IJFIS.2004.4.1.040
34. Kwun, YC, Park, JS, Kim, SY, and Park, JH (2006). Existence and uniqueness of solutions for the semilinear fuzzy integrodifferential equations with nonlocal conditions and forcing term with memory. International Journal of Fuzzy Logic and Intelligent Systems. 6, 288-292. https://doi.org/10.5391/IJFIS.2006.6.4.288
35. Kwun, YC, Kim, WH, Nakagiri, SI, and Park, JH (2009). Existence and uniqueness of solutions for the fuzzy differential equations in n-dimension fuzzy vector space. International Journal of Fuzzy Logic and Intelligent Systems. 9, 16-19. https://doi.org/10.5391/IJFIS.2009.9.1.016
36. Kwun, YC, Kim, JS, Hwang, JS, and Park, JH (2010). Existence of periodic solutions for fuzzy differential equations. International Journal of Fuzzy Logic and Intelligent Systems. 10, 184-193. https://doi.org/10.5391/IJFIS.2010.10.3.184
37. Lee, BY, Kwun, YC, Ahn, YC, and Park, JH (2009). The existence and uniqueness of fuzzy solutions for semilinear fuzzy integrodifferential equations using integral contractor. International Journal of Fuzzy Logic and Intelligent Systems. 9, 339-342.
38. Yoon, JH, Kwun, YC, Park, JS, and Park, JH (2007). Controllability for the semilinear fuzzy integrodifferential equations with nonlocal conditions and forcing term with memory. International Journal of Fuzzy Logic and Intelligent Systems. 7, 34-40. https://doi.org/10.5391/IJFIS.2007.7.1.034
39. Lee, BY, Youm, HE, and Kim, JS (2014). Exact controllability for fuzzy differential equations in credibility space. International Journal of Fuzzy Logic and Intelligent Systems. 14, 145-153. https://doi.org/10.5391/IJFIS.2014.14.2.145
40. El-Ajour, A, and Alawneh, A (2008). Modified Homotopy Analysis Method: Application to Linear and Nonlinear Ordinary Differential Equations of Fractional Order. Amman, Jordan: University of Jordan
41. Odibat, Z, and Momani, S (2008). Modified homotopy perturbation method: application to quadratic Riccati differential equation of fractional order. Chaos, Solitons & Fractals. 36, 167-174. https://doi.org/10.1016/j.chaos.2006.06.041
42. Biazar, J, and Eslami, M (2010). Differential transform method for quadratic Riccati differential equation. International Journal of Nonlinear Science. 9, 444-447.
43. Abbasbandy, S (2004). Solving Riccati differential equation using Adomian decomposition method. Appl Math Computers. 157, 54-63.
44. Tapaswini, S, and Chakraverty, S (2013). Approximate solution of fuzzy quadratic Riccati differential equation. Coupled Systems of Mechanics. 23, 255-269. https://doi.org/10.12989/csm.2013.2.3.255
45. Hooshangian, L (2016). Numerical solution for fuzzy Riccati equation under generalized derivation. International Journal of Applications of Fuzzy Sets and Artificial Intelligence. 6, 133-143.
46. Gnitchogna, R, and Atangana, A (2018). New two step Laplace Adam-Bashforth method for integer a noninteger order partial differential equations. Numerical Methods for Partial Differential Equations. 34, 1739-1758. https://doi.org/10.1002/num.22216
47. Atangana, A, and Doungmo Goufo, EF (2014). Extension of matched asymptotic method to fractional boundary layers problems. Mathematical Problems in Engineering. 2014. article no. 107535
48. Bede, B, and Gal, SG (2005). Generalizations of the differentiability of fuzzy-number-valued functions with applications to fuzzy differential equations. Fuzzy Sets and Systems. 151, 581-599. https://doi.org/10.1016/j.fss.2004.08.001
49. Diamond, P, and Kloeden, PE (1994). Metric Spaces of Fuzzy Sets: Theory and Applications. Singapore: World Scientific Publishing
50. Dubois, D, and Prade, H (1982). Towards fuzzy differential calculus. Part 3: Differentiation. Fuzzy Sets and Systems. 8, 225-233. https://doi.org/10.1016/s0165-0114(82)80001-8
51. Zimmermann, HJ (1992). Fuzzy Set Theory and its Applications. Boston, MA: Kluwer Academic Publisher
52. Ghaemi, F, Yunus, R, Ahmadian, A, Salahshour, S, Suleiman, M, and Saleh, SF (2013). Application of fuzzy fractional kinetic equations to modelling of the acid hydrolysis reaction. Abstract and Applied Analysis. 2013. article no 610314
53. Anastassiou, GA, and Gal, SG (2001). On a fuzzy trigonometric approximation theorem of Weierstrass-type. Journal of Fuzzy Mathematics. 9, 701-708.
54. Anastassiou, GA (2010). Fuzzy Mathematics: Approximation Theory. Heidelberg, Germany: Springer
55. Chalco-Cano, Y, and Roman-Flores, H (2008). On new solutions of fuzzy differential equations. Chaos, Solitons & Fractals. 38, 112-119. https://doi.org/10.1016/j.chaos.2006.10.043
56. Puri, ML, and Ralescu, DA (1986). Fuzzy random variables. Journal of Mathematical Analysis and Applications. 114, 409-422. https://doi.org/10.1016/0022-247X(86)90093-4
57. Pirzada, UM, and Vakaskar, DC (2017). Existence of Hukuhara differentiability of fuzzy-valued functions. Journal of the Indian Mathematical Society. 84, 27-46. https://doi.org/10.18311/jims/2017/5824
58. Armand, A, Allahviranloo, T, and Gouyandeh, Z (2018). Some fundamental results on fuzzy calculus. Iranian Journal of Fuzzy Systems. 15, 27-46. https://doi.org/10.22111/ijfs.2018.3948
59. Wasques, VF, Esmi, E, Barros, LC, and Sussner, P (2020). The generalized fuzzy derivative is interactive. Information Sciences. 519, 93-109. https://doi.org/10.1016/j.ins.2020.01.042
60. Ma, Y, Li, H, and Zhang, S (2020). Solving two-dimensional fuzzy Fredholm integral equations via sinc collocation method. Advances in Difference Equations. 2020. article no. 290
61. Kaleva, O (1987). Fuzzy differential equations. Fuzzy Sets and Systems. 24, 301-317. https://doi.org/10.1016/0165-0114(87)90029-7
62. Bede, B (2013). Mathematics of Fuzzy Sets and Fuzzy Logic. Heidelberg, Germany: Springer
63. You, C, and Zhang, R (2018). Existence and uniqueness theorems for complex fuzzy differential equation. Journal of Intelligent & Fuzzy Systems. 34, 2213-2222. https://doi.org/10.3233/JIFS-171231