Article Search
닫기 International Journal of Fuzzy Logic and Intelligent Systems 2022; 22(2): 144-154

Published online June 25, 2022

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

© The Korean Institute of Intelligent Systems

## Solutions of a System of Fuzzy Fractional Differential in Terms of the Matrix Mittag-Leffler Functions

1Vellore Institute of Technology, Chennai Campus, India
2New Prince Shri Bhavani Arts and Sciences College, Chennai, India
3Division of Mathematics, School of Advanced Sciences, Vellore Institute of Technology, Chennai Campus, India

Correspondence to :

Received: August 18, 2021; Revised: August 18, 2021; Accepted: April 28, 2022

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.

The matrix Mittag-Leffler functions play a crucial role in several applications related to systems with fractional dynamics. These functions represent a generalization for fractional-order systems of the matrix exponential function involved in integer-order systems. Computational techniques for evaluating the matrix Mittag-Leffler functions are therefore of particular importance. In this study, a fuzzy system of differential equations with fractional derivatives was solved in terms of the matrix Mittag-Leffler functions. The matrix Mittag-Leffler function was evaluated based on the Jordan canonical form and the minimal polynomial of the matrix.

Keywords: Mittag-Leffler function, Fuzzy calculus, Fractional calculus, System of fuzzy fractional differential equations, Fractional differential equations

Over the past few centuries, fractional differential (FDEs) have been successfully implemented in several mathematical models in the fields of biological science, chemistry, physics, and engineering . The behavior of a number of physical systems can be properly described using a fractional-order system. Atanackovic and Stankovic  studied an FDE system that arises during lateral motion. The uniqueness, existence, and stability results for a system of FDEs have been presented in [3,4]. In the research community, several analytical and numerical techniques have been employed to determine solutions to a system of FDEs, such as the Laplace transform method [5,6], Adomian decomposition method , and exponential integrators .

A salient feature of the response of an FDE system is that the matrix Mittag-Leffler function serves as the basis of the solution, instead of the matrix exponential function as in the case of ordinary differential Owing to the difficulty in acquiring solutions in the form of fractional systems, numerical methods are often used to solve such fractional systems. Moret and Novati  applied the Krylov subspace method to evaluate the matrix Mittag-Leffler function. Matychyn and Onyshchenko  solved the Bagley-Tovik based on the Jordan canonical form.

Garrappa and N. Popolizio  investigated the computation of the Mittag-Leffler function and evaluated a three-parameter Mittag-Leffler function using the inverse Laplace transform method; the results have been presented in . Duan and his colleague [13,14] employed different methods, such as the inverse Laplace transform method, Jordan canonical matrix method, and minimal polynomial method, for solving a system of FDEs, wherein the solutions were expressed in terms of the matrix Mittag-Leffler function.

Motivated by [13,14], in this study, we intend to apply the Jordan canonical method and the minimal polynomial method to a fuzzy system with fractional derivatives. In modeling real-world phenomena, the crisp quantities of a system of FDEs may cause imprecision and uncertainty. This uncertainty gave rise to several studies on the fuzzy systems of FDEs. We primarily aim to develop the Jordan canonical and minimal polynomial proposed by Duan  for computing the matrix Mittag-Leffler function in the solutions of a system of fuzzy FDEs.

Recently, several researchers have focused on fuzzy differential with fractional-order derivatives. Agarwal et al.  introduced a fuzzy FDE with a combination of the Hukuhara difference and Riemann-Liouville derivative. The existence of a solution of FDEs with uncertainty, involving the Riemann-Liouville operator, was proved by Arshad and Lupulescu  and Allahviranloo et al. . However, the H-differentiable functions have an increasing length of support in the time variable. To overcome this issue, certain authors have considered the generalized Hukuhara differentiability (gH-differentiability), as in . This derivative can enhance the set of fuzzy solutions and provide further insight into the behavior of fuzzy solutions. Balooch Shahriyar et al.  obtained solutions of a system of fuzzy fractional differential (SFFDEs) based on the method of eigenvalues and eigenvectors. The collocation method for discontinuous piecewise polynomial spaces has also been applied to the SFFDE in .

In our research, we compared our obtained solutions with the solutions obtained by the eigenvalue and eigenvector methods, presented in . Our approach proves to be beneficial because it is non-differentiable, non-integral, and easy to implement using a computer owing to its matrix-based structure.

This remainder of this article is structured as follows: Section 2 provides the basic concepts of fuzzy and fractional calculus. Section 3 describes the solution procedure for fuzzy FDEs. Section 4 presents the computation of the matrix Mittag-Leffler function. In Section 5, several numerical examples are presented. Finally, Section 6 outlines the conclusions.

In this section, we review the basic definitions of fuzzy calculus, fractional calculus, and the matrix Mittag-Leffler function. Consider ℜF to be a collection of all fuzzy numbers on the real number ℜ.

### Definition 2.1 

The membership function for a fuzzy set w is defined by w : ℜ → [0, 1], and w(p) is denoted as the degree of membership of an element p in the fuzzy set w for each p ∈ ℜ. If w is a convex fuzzy set, then for each p, q ∈ ℜ and λ ∈ [0, 1], wehavew(λp+(1−λ)q) ≥ min{w(p), w(q)}. If w is a normal fuzzy set, then there exists a p ∈ ℜ such that w(p) = 1. The support for w is denoted by sup(w) and is defined by sup(w) = cl{p ∈ ℜ : w(p) > 0}. Moreover, w is upper semi-continuous.

### Definition 2.2 [25,26]

If w is a fuzzy number, then it is defined by [w]r = [wr, w̄r], which satisfies the following conditions:

• (1) wr is an increasing function, and r is a decreasing function such that wrr.

• (2) wr and r are bounded and left-hand continuous functions in (0,1].

• (3) wr and r are bounded and right-hand continuous functions at r = 0.

### Definition 2.3 

Let w, y ∈ ℜF. If there exists a z ∈ ℜF such that w = y+z, then z is referred to as the H-difference between w and y and is denoted by wH y. The symbol −H represents the H-difference, and wH y w + (−1)y.

### Definition 2.4 

Let I = [a, b]. For F : I → ℜF and p0I, if F is strongly H-differentiable at p0, then there exists an element F′(p0) ∈ ℜF such that

• (1) H-differences F(p0+ω)−HF(p0) and F(p0)−HF(p0ω) exist for eachω > 0 sufficiently close to 0; additionally, lim

limω0F(p0+ω)-HF(p0)ω=limω0F(p0)-HF(p0-ω)ω=F(p0).

• (2) H-differences F(p0) −H F(p0 + ω) and F(p0ω) −H F(p0) exist for eachω > 0 sufficiently close to 0; additionally,

limω0F(p0)-HF(p0+ω)-ω=limω0F(p0-ω)-HF(p0)-ω=F(p0).

### Definition 2.5 

Let F : I → ℜF be a fuzzy function, and let F(p; r) = [F(p; r), (p; r)] for r ∈ [0, 1].

• (1) F is assumed to be (1)-type differentiable if F(p; r) and (p; r) are differentiable functions, and F(p;r)=[F_(p;r),F¯(p;r)]

• (2) F is assumed to be (2)-type differentiable if F(p; r) and (p; r) are differentiable functions, and F(p;r)=[F¯(p;r),F_(p;r)]

### Definition 2.6 [30,31]

Let F : I → ℜF and FCF (I) ∩ LF (I). The fuzzy Riemann-Liouville fractional integral of order α is described as follows:

Jα(F(p))=1Γ(α)apF(v)(p-v)1-αdv,0<α1.

### Definition 2.7 [31,32]

Let F : I → ℜF and FCF (I) ∩ LF (I). Then, F is Caputo’s H-differentiable at p if

(DCαF)(p)=1Γ(1-α)apF(v)(p-v)αdv.0<α1.

In addition, F is Caputo [(1)-type - α] differentiable when F is (1)-type differentiable, and F is Caputo [(2)-type - α] differentiable when F is (2)-type differentiable, where CF (I) and LF (I) represent the space of all continuous and Lebesgue integrable fuzzy valued functions on I, respectively.

### Definition 2.8

The Mittag–Leffler function is crucial for representing the solution of FDEs. The standard definition of the Mittag-Leffler function with one parameter and two parameters [33,34] is given as

Eα(z)=k=0zkΓ(αk+1),α>0,zC,Eα,β(z)=k=0zkΓ(αk+β),α>0,β>0,zC,

where Γ denotes a gamma function Γ(y)=0ty-1e-tdt.

Consider ACnXn. Then, the matrix Mittag-Leffler function  is given as

Eα,β(A)=k=0AkΓ(αk+β),α>0,β>0.

By setting β = 1 in Eq. (1), we obtain the one-parameter matrix Mittag-Leffler function as Eα(A). Moreover, by setting α = β = 1, the matrix Mittag-Leffler function becomes a matrix exponential, i.e.,

E1,1(A)=k=0Akk!=eA.

We consider the following SFFDE:

DtαX˜(t)=AX˜(t)+f˜(t),t>0,

where , A denotes the n-th order crisp matrix, (t) = (1(t), 2(t), …, n(t))T denotes a continuous function, and the elements of (t) are fuzzy numbers. The fuzzy initial condition (0) = (1(0), 2(0), …, n(0))T is specified, and Dtα is the fuzzy Caputo’s fractional derivative, where 0 < α ≤ 1.

By operating on both sides of Eq. (2), we obtain

X˜(t)=X˜(0)+AJαX˜(t)+Jαf˜(t),

where

X˜(0)=[X_(0;r),X¯(0;r)]=((x_1(0;r),x¯1(0;r)),(x_2(0;r),x¯2(0;r)),,(x_n(0;r),x¯n(0;r)))T.

The parametric form of Eq. (3) is given as

X_(t;r)=X_(0;r)+AJαX_(t;r)+Jαf_(t;r),X¯(t;r)=X¯(0;r)+AJαX¯(t;r)+Jαf¯(t;r).

Let us denote the kth approximate solution by ϕk (t; r) and ϕ̄k(t; r), where the zeroth approximate solution is given by

ϕ_0(t;r)=X_(0;r);ϕ¯0(t;r)=X¯(0;r).

For k ≥ 1, we derive the following recursive formula:

ϕ_k(t;r)=X_(0;r)+AJαϕ_k-1(t;r)+Jαf_(t;r),ϕ¯k(t;r)=X¯(0;r)+AJαϕ¯k-1(t;r)+Jαf¯(t;r).

Using Eq. (4), we can calculate the following:

ϕ_1(t;r)=X_(0;r)+AtαΓ(α+1)X_(0;r)+Jαf_(t;r),ϕ¯1(t;r)=X¯(0;r)+AtαΓ(α+1)X¯(0;r)+Jαf¯(t;r).

Similarly,

ϕ_2(t;r)=X_(0;r)+AtαΓ(α+1)X_(0;r)+A2t2αΓ(2α+1)X_(0;r)+Jαf_(t;r)+AJ2αf_(t;r),ϕ¯2(t;r)=X¯(0;r)+AtαΓ(α+1)X¯(0;r)+A2t2αΓ(2α+1)X¯(0;r)+Jαf¯(t;r)+AJ2αf¯(t;r),ϕ_k(t;r)=i=0kAitiαΓ(iα+1)X_(0;r)+i=0kAiJ(i+1)αf_(t;r),ϕ¯k(t;r)=i=0kAitiαΓ(iα+1)X¯(0;r)+i=0kAiJ(i+1)αf¯(t;r).

Then, the formula for the fractional integral is

J(i+1)αf_(t;r)=tiα+α-1Γ(iα+α)*f_(t;r),

and

J(i+1)αf¯(t;r)=tiα+α-1Γ(iα+α)*f¯(t;r).

By applying the limit k → ∞ for ϕk (t; r) and ϕ̃k(t; r), we obtain the solution in a series form:

X_(t;r)=i=0AitiαΓ(iα+1)X_(0;r)+i=0Aitiα+α-1Γ(iα+α)*f_(t;r),X¯(t;r)=i=0AitiαΓ(iα+1)X¯(0;r)+i=0Aitiα+α-1Γ(iα+α)*f¯(t;r).

Using the matrix Mittag-Leffler function, we can express the solution as follows:

X_(t;r)=Eα(Atα)X_(0;r)+tα-1Eα,α(Atα)*f_(t;r),X¯(t;r)=Eα(Atα)X¯(0;r)+tα-1Eα,α(Atα)*f¯(t;r).

where the convolution is defined as tα-1Eα,α(Atα)*f(t;r)=0tτα-1Eα,α(Aτα)f(t-τ;r)dτ. Then, the solution of Eq. (2) is of the form

X˜(t)=Eα(Atα)X˜(0)+tα-1Eα,α(Atα)*f˜(t).

### Theorem 3.1

If A has real distinct eigenvalues λi and the corresponding eigenvector ζi for i = 1, 2, …, n, the solution of the fuzzy fractional system (2) is a fuzzy number X˜(t)=i=1nc˜iEα(λitα)ζi, where (t) = 0.

Proof

Consider X˜(t)=i=1nc˜iEα(λitα)ζi to be a solution to the fuzzy fractional system (2).

Suppose = [1, 2, …, n], where i are fuzzy numbers, and [C˜]r=i=1n[c˜i]r. Then,

X_(t;r)=min{i=1nciEα(λitα)ζi/ci[C˜]r}=i=1nciEα(λitα)ζi_,X¯(t;r)=i=1nciEα(λitα)ζi_by Definition 7in ,X¯(t;r)=max{i=1nciEα(λitα)ζi/ci[C˜]r}=i=1nciEα(λitα)ζi¯=i=1nciEα(λitα)ζi¯.

Differentiating Eqs. (8) and (9), we obtain

DtαX_(t;r)=i=1nciλiEα(λitα)ζi_,

and

DtαX¯(t;r)=i=1nciλiEα(λitα)ζi¯.

Because λi and ζi are the eigenvalues and corresponding eigenvectors of matrix A, respectively, i=λζi for i = 1, 2, …, n. Then, we have

DtαX_(t;r)=Ai=1nciEα(λitα)ζi_=AX_(t;r),DtαX¯(t;r)=Ai=1nciEα(λitα)ζi¯=AX¯(t;r).

i.e.,

(DtαX_(t;r),DtαX¯(t;r))=A(X_(t;r),X¯(t;r)),DtαX˜(t)=AX˜(t).

This completes the proof.

### Theorem 3.2

If A has repeated eigenvalues λi with multiplicity, and k and ζi, i = 1, 2, …, n are the corresponding eigenvectors, then the solution of the fuzzy fractional system (2) is a fuzzy number.

X˜(t)=i=1kc˜i(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)+i=k+1nc˜iEα(λitα)ζi,

where (t) = 0.

Proof

Consider Eq. (10) as the solution of the fuzzy fractional system (2). Suppose = [1, 2, …, n], where i are fuzzy numbers, and [C˜]r=i=1n[c˜i]r. Then,

X_(t;r)=min{i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)+i=k+1nciEα(λitα)ζi/ci[C˜]r}=i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)+i=k+1nciEα(λitα)ζi_=i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)_+i=k+1nciEα(λitα)ζi_,X¯(t;r)=max{i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)+i=k+1nciEα(λitα)ζi/ci[C˜]r}=i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)+i=k+1nciEα(λitα)ζi¯=i=1kci(Eα(k)(λitα)(tα(A-λiI)ζi)+Eα(λitα)ζi)¯+i=k+1nciEα(λitα)ζi¯.

Differentiating Eqs. (11) and (12), we obtain

DtαX_(t;r)=i=1kci(λiEα(k)(λitα)(tα(A-λiI)ζi)+Eα(k)(λitα)(A-λiI)ζi+λiEα(λitα)ζi)_+i=k+1nciλiEα(λitα)ζi_,DtαX¯(t;r)=i=1kci(λiEα(k)(λitα)(tα(A-λiI)ζi)+Eα(k)(λitα)(A-λiI)ζi+λiEα(λitα)ζi)¯+i=k+1nciλiEα(λitα)ζi¯

Because λi and ζi are eigenvalues and the corresponding eigenvectors of matrix A, i = λζi for i = 1, 2, …, n. Then, we have

DtaX_(t;r)=AX_(t;r),DtaX¯(t;r)=AX¯(t;r),

i.e.,

(DtaX_(t;r),DtaX¯(t;r))=A(X_(t;r),X¯(t;r)),DtaX˜(t)=AX˜(t).

This completes the proof.

### 4. Calculation of the Matrix Mittag-Leffler Function

The solution of the SFFDEs consists of the matrix Mittag-Leffler function Eα(Atα). In this section, we obtain this matrix Mittag-Leffler function using the Jordan canonical method and the minimal polynomial method.

### 4.1 Matrix of the Jordan Canonical Method

Consider that the matrix of the Jordan canonical for A with an order of n×n is J. J = diag(J1, J2, …, Jk) and A = PJP−1, where

Ji=[λi1λi11λi]ni×ni,i=1,2,,k,

and i=1kni=n. Thus, we have

Eα(Atα)=PEα(Jtα)P-1=Pdiag(Eα(J1tα),,Eα(Jktα))P-1,

where Eα(Jitα) are upper triangular matrices.

Eα(Jitα)=[Eα(λitα)tαEα(λitα)t(ni-1)αEα(ni-1)(λitα)(ni-1)!Eα(λitα)t(ni-2)αEα(ni-2)(λitα)(ni-2)!tαEα(λitα)Eα(λitα)],

where Eαk(λitα)=dkEα(z)dzkz=λitα

Consider that the minimal polynomial of the nth order matrix A is m(λ) = (λλ1)n1 (λλ2)n2 … (λλk)nk, where i=1kni=n.

According to the theory of matrix analysis , we can express Eα(Atα) as a matrix polynomial of degree n − 1, i.e.,

Eα(Atα)=a0(t)I+a1(t)A++am-1(t)Am-1.

Eq. (15) holds if and only if Eα(λtα) and the polynomial ψ(λ, t) = a0(t)+a1(t)λ+· · ·+am−1(t)λm−1 are consistent in the spectrum of matrix A; that is,

djψ(λ,t)dλjλ=λi=tjαEα(j)(λitα),i=1,2,k,j=0,1,ni-1.

We can obtain ai(t) by solving Eq. (16).

### 5.1 Example

We consider the following SFFDE:

DtaX˜(t)=AX˜(t),   A=(2-14-3).

With an initial condition,

X˜(0)=(x1_(0)x2_(0)x1¯(0)x2¯(0))=(r+12+r3-r5-r),

where λ1 = 1 and λ2 = −2 denote the eigenvalues of the matrix A.

By applying the Jordan canonical method, we have A = PJP−1, where

P=(1114)   and   J=(100-2),Eα(Atα)=PEα(Atα)P-1=P(Eα(tα)00Eα(-2tα))P-1=13(4-14-1)Eα(tα)+13(-11-44)Eα(-2tα).

We apply the minimal polynomial method as follows:

Next, we consider the minimal polynomial of A as m(λ) = (λ − 1)(λ + 2). Further, we consider Eα(Atα) = a0(t)I2 + a1(t)A, where I2 represents a second-order identity matrix. Thus,

a0(t)+a1(t)=Eα(tα),a0(t)-2a1(t)=Eα(-2tα).

By solving the abovementioned system, we obtain

a0(t)=23Eα(tα)+13Eα(-2tα),a1(t)=13Eα(tα)-13Eα(-2tα).

In addition, the matrix Mittag-Leffler function has the form

Eα(Atα)=[23Eα(tα)+13Eα(-2tα)]I2+[13Eα(tα)-13Eα(-2tα)]A=13(4-14-1)Eα(tα)+13(-11-44)Eα(-2tα).

The results presented in Eqs. (18) and (19) are the same.

The solution of SFFDE (17) parameterized by the order α is

X˜(t)=Eα(Atα)X˜(0).

Then,

X_(t;r)=(r+23)(11)Eα(tα)+13(14)Eα(-2tα),X¯(t;r)=(73-r)(11)Eα(tα)+23(14)Eα(-2tα).

The above results are consistent with those obtained by the eigenvalue and eigenvector methods in [23, Eg. 1].

### 5.2 Example

Consider the following SFFDE:

DtαX˜(t)=AX˜(t),A=(11121-10-11).

With an initial condition,

X˜(0)=(x1_(0)x2_(0)x3_(0)x1¯(0)x2¯(0)x3¯(0))=(0.75+0.25r1.5+0.5r3.75+0.25r1.25-0.25r2.5-0.5r4.25-0.25r).

Subsequently, λ1 = −1 and λ2 = λ3 = 2 denote the eigenvalues of matrix A. By applying the Jordan canonical method, we obtain A = PJP−1, where

P=(-3014102-11),J=(-100021002),Eα(Atα)=PEα(Atα)P-1=P(Eα(-tα)000Eα(2tα)tαEα(2tα)00Eα(2tα))P-1=19(3-3-3-444-222)Eα(-tα)+19(63345-42-27)Eα(2tα)+19(000633-6-3-3)tαEα(2tα).

Further, we apply the minimal polynomial method as follows:

The minimal polynomial of A is m(λ) = (λ + 1)(λ − 2)2 and Eα(Atα) = a0(t)I3 + a1(t)A + a2(t)A2. Thus,

a0(t)-a1(t)+a2(t)=Eα(-tα),a0(t)+2a1(t)+4a2(t)=Eα(2tα),a1(t)+4a2(t)=tαEα(2tα).

We now solve the system to obtain

a0(t)=49Eα(-tα)+59Eα(2tα)-69tαEα(2tα),a1(t)=-49Eα(-tα)+49Eα(2tα)-39tαEα(2tα),a2(t)=19Eα(-tα)-19Eα(2tα)+39tαEα(2tα).

Then, the matrix Mittag-Leffler function assumes the form

Eα(Atα)=[49Eα(-tα)+59Eα(2tα)-69tαEα(2tα)]I3+[-49Eα(-tα)+49Eα(2tα)-39tαEα(2tα)]A+[19Eα(-tα)-19Eα(2tα)+39tαEα(2tα)]A2=19(3-3-3-444222)Eα(-tα)+19(63345-42-27)Eα(2tα)+19(000633-6-3-3)tαEα(2tα).

The results presented in Eqs. (21) and (22) are the same. The solution of SFFDE (20) parameterized by the order α is

X˜(t)=Eα(Atα)X˜(0).

Then,

X_(t;r)=19(-13.5-1.5r18+2r9+r)Eα(-tα)+19(20.25+3.75r-4.5+2.5r24.75+1.25r)Eα(2tα)+19(020.25+3.75r-20.25-3.75r)tαEα(2tα),X¯(t;r)=19(-16.5+1.5r22-2r11-r)Eα(-tα)+19(27.75-3.75r0.5-2.5r27.25-1.25r)Eα(2tα)+19(027.75-3.75r-27.75+3.75r)tαEα(2tα).

The solutions presented above can be rewritten in the following form:

X_(t;r)=19(4.5+0.5r)   (-342)Eα(-tα)+19(-4.5+2.5r)   (01-1)Eα(2tα)+19(20.25+3.75r)   [(01-1)tαEα(2tα)+(101)Eα(2tα)],X¯(t;r)=19(5.5-0.5r)   (-342)Eα(-tα)+19(0.5-2.5r)   (01-1)Eα(2tα)+19(27.75-3.75r)   [(01-1)tαEα(2tα)+(101)Eα(2tα)].

The above results are consistent with those obtained by the eigenvalue and eigenvector methods in [23, Eg. 2].

Figure 1(a) and 1(b) present the solution curves of Example 1 for different values of α. Figure 2(a)–2(c) present the solution curves of Example 2 for different values of α, wherein the purple, blue, and red lines represent α = 0.6, α = 0.8, and α = 1.0, respectively. In the examples above, the results indicate that the solutions obtained by the proposed method are similar to those obtained by the eigenvalue and eigenvector method in .

This study primarily aimed to develop a solution for a system of fuzzy FDEs in terms of the matrix Mittag-Leffler functions via two methods: the Jordan canonical method and the minimal polynomial method. The effectiveness, ability, and simplicity of the proposed techniques were demonstrated using numerical examples. The numerical and graphical results revealed that the solutions obtained by the proposed methods were consistent with the solutions obtained by the eigenvalue and eigenvector method. The solutions were plotted using MATLAB. Fig. 1.

Solutions 1(t) and 2(t) of Example 1 for different values of (a) α = 0.6, (b) α = 0.8, (c) α = 1.0 at t = 1. Fig. 2.

Solutions 1(t), 2(t), and 3(t) of Example 2 for different values of (a) α = 0.6, (b) α = 0.8, (c) α = 1.0 at t = 1.

1. Torvik, PJ, and Bagley, RL (1984). On the appearance of the fractional derivative in the behavior of real materials. Journal of Applied Mechanics. 51, 294-298. https://doi.org/10.1115/1.3167615
2. Atanackovic, TM, and Stankovic, B (2004). On a system of differential equations with fractional derivatives arising in rod theory. Journal of Physics A: Mathematical and General. 37, 1241-1250. https://doi.org/10.1088/0305-4470/37/4/012
3. Daftardar-Gejji, V, and Babakhani, A (2004). Analysis of a system of fractional differential equations. Journal of Mathematical Analysis and Applications. 293, 511-522. https://doi.org/10.1016/j.jmaa.2004.01.013
4. Deng, W, Li, C, and Lu, J (2007). Stability analysis of linear fractional differential system with multiple time delays. Nonlinear Dynamics. 48, 409-416. https://doi.org/10.1007/s11071-006-9094-0
5. Charef, A, and Boucherma, D (2011). Analytical solution of the linear fractional system of commensurate order. Computers & Mathematics with Applications. 62, 4415-4428. https://doi.org/10.1016/j.camwa.2011.10.017
6. Duan, JS, Fu, SZ, and Wang, Z (2013). Solution of linear system of fractional differential equations. Pacific Journal of Applied Mathematics. 5, 93-106.
7. Daftardar-Gejji, V, and Jafari, H (2005). Adomian decomposition: a tool for solving a system of fractional differential equations. Journal of Mathematical Analysis and Applications. 301, 508-518. https://doi.org/10.1016/j.jmaa.2004.07.039
8. Garrappa, R (2013). Exponential integrators for time–fractional partial differential equations. The European Physical Journal Special Topics. 222, 1915-1927. https://doi.org/10.1140/epjst/e2013-01973-1
9. Moret, I, and Novati, P (2011). On the convergence of Krylov subspace methods for matrix Mittag–Leffler functions. SIAM Journal on Numerical Analysis. 49, 2144-2164. https://doi.org/10.1137/080738374
10. Matychyn, I, and Onyshchenko, V (2018). Matrix Mittag Leffler function in fractional systems and its computation. Bulletin of the Polish Academy of Sciences: Technical Sciences. 66, 495-500. https://doi.org/10.24425/124266
11. Garrappa, R, and Popolizio, N (2018). Computing the matrix Mittag-Leffler function with applications to fractional calculus. Journal of Scientific Computing. 77, 129-153. https://doi.org/10.1007/s10915-018-0699-5
12. Garrappa, R (2015). Numerical evaluation of two and three parameter Mittag-Leffler functions. SIAM Journal on Numerical Analysis. 53, 1350-1369. https://doi.org/10.1137/140971191
13. Duan, J (2018). System of linear fractional differential equations and the Mittag-Leffler functions with matrix variable. Journal of Physics: Conference Series. 1053. article no. 012032
14. Duan, J, and Chen, L (2018). Solution of fractional differential equation systems and computation of matrix Mittag–Leffler functions. Symmetry. 10. article no. 503
15. Agarwal, RP, Lakshmikantham, V, and Nieto, JJ (2010). On the concept of solution for fractional differential equations with uncertainty. Nonlinear Analysis: Theory, Methods & Applications. 72, 2859-2862. https://doi.org/10.1016/j.na.2009.11.029
16. Arshad, S, and Lupulescu, V (2011). On the fractional differential equations with uncertainty. Nonlinear Analysis: Theory, Methods & Applications. 74, 3685-3693. https://doi.org/10.1016/j.na.2011.02.048
17. Allahviranloo, T, Salahshour, S, and Abbasbandy, S (2012). Explicit solutions of fractional differential equations with uncertainty. Soft Computing. 16, 297-302. https://doi.org/10.1007/s00500-011-0743-y
18. Allahviranloo, T, Armand, A, and Gouyandeh, Z (2014). Fuzzy fractional differential equations under generalized fuzzy Caputo derivative. Journal of Intelligent & Fuzzy Systems. 26, 1481-1490. https://doi.org/10.3233/IFS-130831
19. Van Hoa, N (2015). Fuzzy fractional functional differential equations under Caputo gH-differentiability. Communications in Nonlinear Science and Numerical Simulation. 22, 1134-1157. https://doi.org/10.1016/j.cnsns.2014.08.006
20. Van Hoa, N (2015). Fuzzy fractional functional integral and differential equations. Fuzzy Sets and Systems. 280, 58-90. https://doi.org/10.1016/j.fss.2015.01.009
21. Long, VH, Son, NTK, and Tam, HTT (2015). Global existence of solutions to fuzzy partial hyperbolic functional differential equations with generalized Hukuhara derivatives. Journal of Intelligent & Fuzzy Systems. 29, 939-954. https//doi.org/10.3233/IFS-151623
22. Long, HV, Son, NTK, and Tam, HTT (2017). The solvability of fuzzy fractional partial differential equations under Caputo gH-differentiability. Fuzzy Sets and Systems. 309, 35-63. https://doi.org/10.1016/j.fss.2016.06.018
23. Balooch Shahriyar, MR, Ismail, F, Aghabeigi, S, Ahmadian, A, and Salahshour, S (2013). An eigenvalue-eigenvector method for solving a system of fractional differential equations with uncertainty. Mathematical Problems in Engineering. 2013. article no. 579761
24. Alijani, Z, Baleanu, D, Shiri, B, and Wu, GC (2020). Spline collocation methods for systems of fuzzy fractional differential equations. Chaos, Solitons & Fractals. 131. article no. 109510
25. Zimmermann, HJ (1991). Fuzzy Set Theory and its Applications. Boston, MA: Kluwer Academic Publishers
26. Friedman, M, Ma, M, and Kandel, A (1999). Numerical solutions of fuzzy differential and integral equations. Fuzzy Sets and Systems. 106, 35-48. https://doi.org/10.1016/S0165-0114(98)00355-8
27. Ma, M, Friedman, M, and Kandel, A (1999). Numerical solutions of fuzzy differential equations. Fuzzy Sets and Systems. 105, 133-138. https://doi.org/10.1016/S0165-0114(97)00233-9
28. 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
29. 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
30. Abasbandy, S, Allahviranloo, T, Shahryari, MB, and Salahshour, S (2012). Fuzzy local fractional differential equation. International Journal of Industrial Mathematics. 4, 231-245.
31. Allahviranloo, T. (2021) . Fuzzy Fractional Differential Operators and Equations. https://doi.org/10.1007/978-3-030-51272-9
32. Armand, A, Allahviranloo, T, and Gouyandeh, Z (2016). General solution of Basset equation with Caputo generalized Hukuhara derivative. Journal of Applied Analysis and Computation. 6, 119-130. https://doi.org/10.11948/2016010
33. Podlubny, I (1999). Fractional Differential Equations. San Diego, CA: Academic Press
34. Mainardi, F, and Gorenflo, R (2000). On Mittag-Leffler-type functions in fractional evolution processes. Journal of Computational and Applied mathematics. 118, 283-299. https://doi.org/10.1016/S0377-0427(00)00294-6
35. Horn, RA, and Johnson, CR (2012). Matrix Analysis. New York, NY: Cambridge University Press V. Padmapriya is a Ph.D. research scholar studying at the Division of Mathematics, Vellore Institute of Technology-Chennai Campus, India. She is currently working as an assistant professor in the Department of Statistics, New prince Shri Bhavani Arts and Science college, Chennai, India. Her research interests include fuzzy differential and fuzzy fractional differential M. Kaliyappan is working as a professor of Mathematics at the Division of Mathematics, School of Advanced Sciences, VIT Chennai, India. He has over 23 years of experience in teaching and research. His research interests include approximation theory, numerical computing, differential equations, fuzzy differential equations, fractional differential equations, fuzzy fractional differential equations, and optimization theory.

E-mail: kaliyappan.m@vit.ac.in

### Article

#### Original Article International Journal of Fuzzy Logic and Intelligent Systems 2022; 22(2): 144-154

Published online June 25, 2022 https://doi.org/10.5391/IJFIS.2022.22.2.144

## Solutions of a System of Fuzzy Fractional Differential in Terms of the Matrix Mittag-Leffler Functions

1Vellore Institute of Technology, Chennai Campus, India
2New Prince Shri Bhavani Arts and Sciences College, Chennai, India
3Division of Mathematics, School of Advanced Sciences, Vellore Institute of Technology, Chennai Campus, India

Received: August 18, 2021; Revised: August 18, 2021; Accepted: April 28, 2022

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

The matrix Mittag-Leffler functions play a crucial role in several applications related to systems with fractional dynamics. These functions represent a generalization for fractional-order systems of the matrix exponential function involved in integer-order systems. Computational techniques for evaluating the matrix Mittag-Leffler functions are therefore of particular importance. In this study, a fuzzy system of differential equations with fractional derivatives was solved in terms of the matrix Mittag-Leffler functions. The matrix Mittag-Leffler function was evaluated based on the Jordan canonical form and the minimal polynomial of the matrix.

Keywords: Mittag-Leffler function, Fuzzy calculus, Fractional calculus, System of fuzzy fractional differential equations, Fractional differential equations

### 1. Introduction

Over the past few centuries, fractional differential (FDEs) have been successfully implemented in several mathematical models in the fields of biological science, chemistry, physics, and engineering . The behavior of a number of physical systems can be properly described using a fractional-order system. Atanackovic and Stankovic  studied an FDE system that arises during lateral motion. The uniqueness, existence, and stability results for a system of FDEs have been presented in [3,4]. In the research community, several analytical and numerical techniques have been employed to determine solutions to a system of FDEs, such as the Laplace transform method [5,6], Adomian decomposition method , and exponential integrators .

A salient feature of the response of an FDE system is that the matrix Mittag-Leffler function serves as the basis of the solution, instead of the matrix exponential function as in the case of ordinary differential Owing to the difficulty in acquiring solutions in the form of fractional systems, numerical methods are often used to solve such fractional systems. Moret and Novati  applied the Krylov subspace method to evaluate the matrix Mittag-Leffler function. Matychyn and Onyshchenko  solved the Bagley-Tovik based on the Jordan canonical form.

Garrappa and N. Popolizio  investigated the computation of the Mittag-Leffler function and evaluated a three-parameter Mittag-Leffler function using the inverse Laplace transform method; the results have been presented in . Duan and his colleague [13,14] employed different methods, such as the inverse Laplace transform method, Jordan canonical matrix method, and minimal polynomial method, for solving a system of FDEs, wherein the solutions were expressed in terms of the matrix Mittag-Leffler function.

Motivated by [13,14], in this study, we intend to apply the Jordan canonical method and the minimal polynomial method to a fuzzy system with fractional derivatives. In modeling real-world phenomena, the crisp quantities of a system of FDEs may cause imprecision and uncertainty. This uncertainty gave rise to several studies on the fuzzy systems of FDEs. We primarily aim to develop the Jordan canonical and minimal polynomial proposed by Duan  for computing the matrix Mittag-Leffler function in the solutions of a system of fuzzy FDEs.

Recently, several researchers have focused on fuzzy differential with fractional-order derivatives. Agarwal et al.  introduced a fuzzy FDE with a combination of the Hukuhara difference and Riemann-Liouville derivative. The existence of a solution of FDEs with uncertainty, involving the Riemann-Liouville operator, was proved by Arshad and Lupulescu  and Allahviranloo et al. . However, the H-differentiable functions have an increasing length of support in the time variable. To overcome this issue, certain authors have considered the generalized Hukuhara differentiability (gH-differentiability), as in . This derivative can enhance the set of fuzzy solutions and provide further insight into the behavior of fuzzy solutions. Balooch Shahriyar et al.  obtained solutions of a system of fuzzy fractional differential (SFFDEs) based on the method of eigenvalues and eigenvectors. The collocation method for discontinuous piecewise polynomial spaces has also been applied to the SFFDE in .

In our research, we compared our obtained solutions with the solutions obtained by the eigenvalue and eigenvector methods, presented in . Our approach proves to be beneficial because it is non-differentiable, non-integral, and easy to implement using a computer owing to its matrix-based structure.

This remainder of this article is structured as follows: Section 2 provides the basic concepts of fuzzy and fractional calculus. Section 3 describes the solution procedure for fuzzy FDEs. Section 4 presents the computation of the matrix Mittag-Leffler function. In Section 5, several numerical examples are presented. Finally, Section 6 outlines the conclusions.

### 2. Basic Concepts

In this section, we review the basic definitions of fuzzy calculus, fractional calculus, and the matrix Mittag-Leffler function. Consider ℜF to be a collection of all fuzzy numbers on the real number ℜ.

### Definition 2.1 

The membership function for a fuzzy set w is defined by w : ℜ → [0, 1], and w(p) is denoted as the degree of membership of an element p in the fuzzy set w for each p ∈ ℜ. If w is a convex fuzzy set, then for each p, q ∈ ℜ and λ ∈ [0, 1], wehavew(λp+(1−λ)q) ≥ min{w(p), w(q)}. If w is a normal fuzzy set, then there exists a p ∈ ℜ such that w(p) = 1. The support for w is denoted by sup(w) and is defined by sup(w) = cl{p ∈ ℜ : w(p) > 0}. Moreover, w is upper semi-continuous.

### Definition 2.2 [25,26]

If w is a fuzzy number, then it is defined by [w]r = [wr, w̄r], which satisfies the following conditions:

• (1) wr is an increasing function, and r is a decreasing function such that wrr.

• (2) wr and r are bounded and left-hand continuous functions in (0,1].

• (3) wr and r are bounded and right-hand continuous functions at r = 0.

### Definition 2.3 

Let w, y ∈ ℜF. If there exists a z ∈ ℜF such that w = y+z, then z is referred to as the H-difference between w and y and is denoted by wH y. The symbol −H represents the H-difference, and wH y w + (−1)y.

### Definition 2.4 

Let I = [a, b]. For F : I → ℜF and p0I, if F is strongly H-differentiable at p0, then there exists an element F′(p0) ∈ ℜF such that

• (1) H-differences F(p0+ω)−HF(p0) and F(p0)−HF(p0ω) exist for eachω > 0 sufficiently close to 0; additionally, lim

$limω→0F(p0+ω)-HF(p0)ω=limω→0F(p0)-HF(p0-ω)ω=F′(p0).$

• (2) H-differences F(p0) −H F(p0 + ω) and F(p0ω) −H F(p0) exist for eachω > 0 sufficiently close to 0; additionally,

$limω→0F(p0)-HF(p0+ω)-ω=limω→0F(p0-ω)-HF(p0)-ω=F′(p0).$

### Definition 2.5 

Let F : I → ℜF be a fuzzy function, and let F(p; r) = [F(p; r), (p; r)] for r ∈ [0, 1].

• (1) F is assumed to be (1)-type differentiable if F(p; r) and (p; r) are differentiable functions, and $F′(p;r)=[F′_(p;r),F′¯(p;r)]$

• (2) F is assumed to be (2)-type differentiable if F(p; r) and (p; r) are differentiable functions, and $F′(p;r)=[F′¯(p;r),F′_(p;r)]$

### Definition 2.6 [30,31]

Let F : I → ℜF and FCF (I) ∩ LF (I). The fuzzy Riemann-Liouville fractional integral of order α is described as follows:

$Jα(F(p))=1Γ(α)∫apF(v)(p-v)1-αdv, 0<α≤1.$

### Definition 2.7 [31,32]

Let F : I → ℜF and FCF (I) ∩ LF (I). Then, F is Caputo’s H-differentiable at p if

$(DCαF) (p)=1Γ(1-α)∫apF′(v)(p-v)αdv.0<α≤1.$

In addition, F is Caputo [(1)-type - α] differentiable when F is (1)-type differentiable, and F is Caputo [(2)-type - α] differentiable when F is (2)-type differentiable, where CF (I) and LF (I) represent the space of all continuous and Lebesgue integrable fuzzy valued functions on I, respectively.

### Definition 2.8

The Mittag–Leffler function is crucial for representing the solution of FDEs. The standard definition of the Mittag-Leffler function with one parameter and two parameters [33,34] is given as

$Eα(z)=∑k=0∞zkΓ(αk+1), α>0, z∈C,Eα, β(z)=∑k=0∞zkΓ(αk+β), α>0, β>0, z∈C,$

where Γ denotes a gamma function $Γ(y)=∫0∞ty-1e-tdt$.

Consider ACnXn. Then, the matrix Mittag-Leffler function  is given as

$Eα, β(A)=∑k=0∞AkΓ(αk+β), α>0, β>0.$

By setting β = 1 in Eq. (1), we obtain the one-parameter matrix Mittag-Leffler function as Eα(A). Moreover, by setting α = β = 1, the matrix Mittag-Leffler function becomes a matrix exponential, i.e.,

$E1, 1(A)=∑k=0∞Akk!=eA.$

### 3. Solution for a System of Fuzzy FDEs

We consider the following SFFDE:

$DtαX˜(t)=AX˜(t)+f˜(t), t>0,$

where , A denotes the n-th order crisp matrix, (t) = (1(t), 2(t), …, n(t))T denotes a continuous function, and the elements of (t) are fuzzy numbers. The fuzzy initial condition (0) = (1(0), 2(0), …, n(0))T is specified, and $Dtα$ is the fuzzy Caputo’s fractional derivative, where 0 < α ≤ 1.

By operating on both sides of Eq. (2), we obtain

$X˜(t)=X˜(0)+AJαX˜(t)+Jαf˜(t),$

where

$X˜(0)=[X_(0;r),X¯(0;r)]=((x_1(0;r),x¯1(0;r)),(x_2(0;r),x¯2(0;r)), …,(x_n(0;r),x¯n(0;r)))T.$

The parametric form of Eq. (3) is given as

$X_(t;r)=X_(0;r)+AJαX_(t;r)+Jαf_(t;r),X¯(t;r)=X¯(0;r)+AJαX¯(t;r)+Jαf¯(t;r).$

Let us denote the kth approximate solution by ϕk (t; r) and ϕ̄k(t; r), where the zeroth approximate solution is given by

$ϕ_0(t;r)=X_(0;r);ϕ¯0(t;r)=X¯(0;r).$

For k ≥ 1, we derive the following recursive formula:

$ϕ_k(t;r)=X_(0;r)+AJαϕ_k-1(t;r)+Jαf_(t;r),ϕ¯k(t;r)=X¯(0;r)+AJαϕ¯k-1(t;r)+Jαf¯(t;r).$

Using Eq. (4), we can calculate the following:

$ϕ_1(t;r)=X_(0;r)+AtαΓ(α+1)X_(0;r)+Jαf_(t;r),ϕ¯1(t;r)=X¯(0;r)+AtαΓ(α+1)X¯(0;r)+Jαf¯(t;r).$

Similarly,

$ϕ_2(t;r)=X_(0;r)+AtαΓ(α+1)X_(0;r)+A2t2αΓ(2α+1)X_(0;r)+Jαf_(t;r)+AJ2αf_(t;r),ϕ¯2(t;r)=X¯(0;r)+AtαΓ(α+1)X¯(0;r)+A2t2αΓ(2α+1)X¯(0;r)+Jαf¯(t;r)+AJ2αf¯(t;r),⋯ϕ_k(t;r)=∑i=0kAitiαΓ(iα+1)X_(0;r)+∑i=0kAiJ(i+1)αf_(t;r),ϕ¯k(t;r)=∑i=0kAitiαΓ(iα+1)X¯(0;r)+∑i=0kAiJ(i+1)αf¯(t;r).$

Then, the formula for the fractional integral is

$J(i+1)αf_(t;r)=tiα+α-1Γ(iα+α)*f_(t;r),$

and

$J(i+1)αf¯(t;r)=tiα+α-1Γ(iα+α)*f¯(t;r).$

By applying the limit k → ∞ for ϕk (t; r) and ϕ̃k(t; r), we obtain the solution in a series form:

$X_(t;r)=∑i=0∞AitiαΓ(iα+1)X⌣_(0;r)+∑i=0∞Aitiα+α-1Γ(iα+α)*f_(t;r),X¯(t;r)=∑i=0∞AitiαΓ(iα+1)X⌣¯(0;r)+∑i=0∞Aitiα+α-1Γ(iα+α)*f¯(t;r).$

Using the matrix Mittag-Leffler function, we can express the solution as follows:

$X_(t;r)=Eα(Atα)X_(0;r)+tα-1Eα, α(Atα)*f_(t;r),X¯(t;r)=Eα(Atα)X¯(0;r)+tα-1Eα, α(Atα)*f¯(t;r).$

where the convolution is defined as $tα-1Eα, α(Atα)*f(t;r)=∫0tτα-1Eα, α(Aτα)f(t-τ;r)dτ$. Then, the solution of Eq. (2) is of the form

$X˜(t)=Eα(Atα)X˜(0)+tα-1Eα, α(Atα)*f˜(t).$

### Theorem 3.1

If A has real distinct eigenvalues λi and the corresponding eigenvector ζi for i = 1, 2, …, n, the solution of the fuzzy fractional system (2) is a fuzzy number $X˜(t)=∑i=1nc˜iEα(λitα)ζi$, where (t) = 0.

Proof

Consider $X˜(t)=∑i=1nc˜iEα(λitα)ζi$ to be a solution to the fuzzy fractional system (2).

Suppose = [1, 2, …, n], where i are fuzzy numbers, and $[C˜]r=∏i=1n[c˜i]r$. Then,

$X_(t;r) =min{∑i=1nciEα(λitα)ζi/ci∈[C˜]r} =∑i=1nciEα(λitα)ζi_,X¯(t;r) =∑i=1nciEα(λitα)ζi_ by Definition 7 in ,$$X¯(t;r)=max{∑i=1nciEα(λitα)ζi/ci∈[C˜]r}=∑i=1nciEα(λitα)ζi¯=∑i=1nciEα(λitα)ζi¯.$

Differentiating Eqs. (8) and (9), we obtain

$DtαX_(t;r)=∑i=1nciλiEα(λitα)ζi_,$

and

$DtαX¯(t;r)=∑i=1nciλiEα(λitα)ζi¯.$

Because λi and ζi are the eigenvalues and corresponding eigenvectors of matrix A, respectively, i=λζi for i = 1, 2, …, n. Then, we have

$DtαX_(t;r)=A∑i=1nciEα(λitα)ζi_=AX_(t;r),DtαX¯(t;r)=A∑i=1nciEα(λitα)ζi¯=AX¯(t;r).$

i.e.,

$(DtαX_(t;r),DtαX¯(t;r))=A(X_(t;r),X¯(t;r)),DtαX˜(t)=AX˜(t).$

This completes the proof.

### Theorem 3.2

If A has repeated eigenvalues λi with multiplicity, and k and ζi, i = 1, 2, …, n are the corresponding eigenvectors, then the solution of the fuzzy fractional system (2) is a fuzzy number.

$X˜(t)=∑i=1kc˜i(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)+∑i=k+1nc˜iEα(λitα)ζi,$

where (t) = 0.

Proof

Consider Eq. (10) as the solution of the fuzzy fractional system (2). Suppose = [1, 2, …, n], where i are fuzzy numbers, and $[C˜]r=∏i=1n[c˜i]r$. Then,

$X_(t;r)=min{∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)+∑i=k+1nciEα(λitα)ζi/ci∈[C˜]r}=∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)+∑i=k+1nciEα(λitα)ζi_=∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)_+∑i=k+1nciEα(λitα)ζi_,$$X¯(t;r)=max{∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)+∑i=k+1nciEα(λitα)ζi/ci∈[C˜]r}=∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)+∑i=k+1nciEα(λitα)ζi¯=∑i=1kci(Eα(k)(λitα) (tα(A-λiI)ζi)+Eα(λitα)ζi)¯+∑i=k+1nciEα(λitα)ζi¯.$

Differentiating Eqs. (11) and (12), we obtain

$DtαX_(t;r)=∑i=1kci(λiEα(k)(λitα) (tα(A-λiI)ζi)+Eα(k)(λitα) (A-λiI)ζi+λiEα(λitα)ζi)_+∑i=k+1nciλiEα(λitα)ζi_,DtαX¯(t;r)=∑i=1kci(λiEα(k)(λitα) (tα(A-λiI)ζi)+Eα(k)(λitα) (A-λiI)ζi+λiEα(λitα)ζi)¯+∑i=k+1nciλiEα(λitα)ζi¯$

Because λi and ζi are eigenvalues and the corresponding eigenvectors of matrix A, i = λζi for i = 1, 2, …, n. Then, we have

$DtaX_(t;r)=AX_(t;r),DtaX¯(t;r)=AX¯(t;r),$

i.e.,

$(DtaX_(t;r),DtaX¯(t;r))=A(X_(t;r),X¯(t;r)),DtaX˜(t)=AX˜(t).$

This completes the proof.

### 4. Calculation of the Matrix Mittag-Leffler Function

The solution of the SFFDEs consists of the matrix Mittag-Leffler function Eα(Atα). In this section, we obtain this matrix Mittag-Leffler function using the Jordan canonical method and the minimal polynomial method.

### 4.1 Matrix of the Jordan Canonical Method

Consider that the matrix of the Jordan canonical for A with an order of n×n is J. J = diag(J1, J2, …, Jk) and A = PJP−1, where

$Ji=[λi1 λi1 ⋱⋱ ⋱1 λi]ni×ni, i=1, 2, …, k,$

and $∑i=1kni=n$. Thus, we have

$Eα(Atα)=PEα(Jtα)P-1=Pdiag(Eα(J1tα), …,Eα(Jktα))P-1,$

where Eα(Jitα) are upper triangular matrices.

$Eα(Jitα)=[Eα(λitα)tαEα′(λitα)…t(ni-1)αEα(ni-1)(λitα)(ni-1)! Eα(λitα)…t(ni-2)αEα(ni-2)(λitα)(ni-2)! ⋱⋮ ⋱tαEα′(λitα) Eα(λitα)],$

where $Eαk(λitα)=dkEα(z)dzk∣z=λitα$

### 4.2 Minimal Polynomial Method

Consider that the minimal polynomial of the nth order matrix A is m(λ) = (λλ1)n1 (λλ2)n2 … (λλk)nk, where $∑i=1kni=n$.

According to the theory of matrix analysis , we can express Eα(Atα) as a matrix polynomial of degree n − 1, i.e.,

$Eα(Atα)=a0(t)I+a1(t)A+⋯+am-1(t)Am-1.$

Eq. (15) holds if and only if Eα(λtα) and the polynomial ψ(λ, t) = a0(t)+a1(t)λ+· · ·+am−1(t)λm−1 are consistent in the spectrum of matrix A; that is,

$djψ(λ, t)dλj∣λ=λi=tjαEα(j)(λitα),i=1, 2, …k, j=0, 1, …ni-1.$

We can obtain ai(t) by solving Eq. (16).

### 5.1 Example

We consider the following SFFDE:

$DtaX˜(t)=AX˜(t), A=(2-14-3).$

With an initial condition,

$X˜(0)=(x1_(0)x2_(0)x1¯(0)x2¯(0))=(r+12+r3-r5-r),$

where λ1 = 1 and λ2 = −2 denote the eigenvalues of the matrix A.

By applying the Jordan canonical method, we have A = PJP−1, where

$P=(1114) and J=(100-2),Eα(Atα)=PEα(Atα)P-1=P(Eα(tα)00Eα(-2tα))P-1=13 (4-14-1) Eα(tα)+13 (-11-44) Eα(-2tα).$

We apply the minimal polynomial method as follows:

Next, we consider the minimal polynomial of A as m(λ) = (λ − 1)(λ + 2). Further, we consider Eα(Atα) = a0(t)I2 + a1(t)A, where I2 represents a second-order identity matrix. Thus,

$a0(t)+a1(t)=Eα(tα),a0(t)-2a1(t)=Eα(-2tα).$

By solving the abovementioned system, we obtain

$a0(t)=23Eα(tα)+13Eα(-2tα),a1(t)=13Eα(tα)-13Eα(-2tα).$

In addition, the matrix Mittag-Leffler function has the form

$Eα(Atα)=[23Eα(tα)+13Eα(-2tα)]I2+[13Eα(tα)-13Eα(-2tα)]A=13 (4-14-1) Eα(tα)+13 (-11-44) Eα(-2tα).$

The results presented in Eqs. (18) and (19) are the same.

The solution of SFFDE (17) parameterized by the order α is

$X˜(t)=Eα(Atα)X˜(0).$

Then,

$X_(t;r)=(r+23) (11) Eα(tα)+13 (14) Eα(-2tα),X¯(t;r)=(73-r) (11) Eα(tα)+23 (14) Eα(-2tα).$

The above results are consistent with those obtained by the eigenvalue and eigenvector methods in [23, Eg. 1].

### 5.2 Example

Consider the following SFFDE:

$DtαX˜(t)=AX˜(t), A=(11121-10-11).$

With an initial condition,

$X˜(0)=(x1_(0)x2_(0)x3_(0)x1¯(0)x2¯(0)x3¯(0))=(0.75+0.25r1.5+0.5r3.75+0.25r1.25-0.25r2.5-0.5r4.25-0.25r).$

Subsequently, λ1 = −1 and λ2 = λ3 = 2 denote the eigenvalues of matrix A. By applying the Jordan canonical method, we obtain A = PJP−1, where

$P=(-3014102-11), J=(-100021002),Eα(Atα)=PEα(Atα)P-1=P(Eα(-tα)000Eα(2tα)tαEα′(2tα)00Eα(2tα))P-1=19 (3-3-3-444-222) Eα(-tα)+19 (63345-42-27) Eα(2tα)+19 (000633-6-3-3)tαEα′(2tα).$

Further, we apply the minimal polynomial method as follows:

The minimal polynomial of A is m(λ) = (λ + 1)(λ − 2)2 and Eα(Atα) = a0(t)I3 + a1(t)A + a2(t)A2. Thus,

$a0(t)-a1(t)+a2(t)=Eα(-tα),a0(t)+2a1(t)+4a2(t)=Eα(2tα),a1(t)+4a2(t)=tαEα′(2tα).$

We now solve the system to obtain

$a0(t)=49Eα(-tα)+59Eα(2tα)-69tαEα′(2tα),a1(t)=-49Eα(-tα)+49Eα(2tα)-39tαEα′(2tα),a2(t)=19Eα(-tα)-19Eα(2tα)+39tαEα′(2tα).$

Then, the matrix Mittag-Leffler function assumes the form

$Eα(Atα)= [49Eα(-tα)+59Eα(2tα)-69tαEα′(2tα)]I3 +[-49Eα(-tα)+49Eα(2tα)-39tαEα′(2tα)]A +[19Eα(-tα)-19Eα(2tα)+39tαEα′(2tα)]A2= 19 (3-3-3-444222) Eα(-tα)+19 (63345-42-27) Eα(2tα) +19 (000633-6-3-3)tαEα′(2tα).$

The results presented in Eqs. (21) and (22) are the same. The solution of SFFDE (20) parameterized by the order α is

$X˜(t)=Eα(Atα)X˜(0).$

Then,

$X_(t;r)= 19 (-13.5-1.5r18+2r9+r) Eα(-tα) +19 (20.25+3.75r-4.5+2.5r24.75+1.25r) Eα(2tα) +19 (020.25+3.75r-20.25-3.75r)tαEα′(2tα),X¯(t;r)= 19 (-16.5+1.5r22-2r11-r) Eα(-tα) +19 (27.75-3.75r0.5-2.5r27.25-1.25r) Eα(2tα) +19 (027.75-3.75r-27.75+3.75r)tαEα′(2tα).$

The solutions presented above can be rewritten in the following form:

$X_(t;r)= 19(4.5+0.5r) (-342) Eα(-tα) +19(-4.5+2.5r) (01-1) Eα(2tα) +19(20.25+3.75r) [(01-1)tαEα′(2tα)+(101) Eα(2tα)],X¯(t;r)= 19(5.5-0.5r) (-342) Eα(-tα) +19(0.5-2.5r) (01-1) Eα(2tα) +19(27.75-3.75r) [(01-1)tαEα′(2tα)+(101) Eα(2tα)].$

The above results are consistent with those obtained by the eigenvalue and eigenvector methods in [23, Eg. 2].

Figure 1(a) and 1(b) present the solution curves of Example 1 for different values of α. Figure 2(a)–2(c) present the solution curves of Example 2 for different values of α, wherein the purple, blue, and red lines represent α = 0.6, α = 0.8, and α = 1.0, respectively. In the examples above, the results indicate that the solutions obtained by the proposed method are similar to those obtained by the eigenvalue and eigenvector method in .

### 6. Conclusion

This study primarily aimed to develop a solution for a system of fuzzy FDEs in terms of the matrix Mittag-Leffler functions via two methods: the Jordan canonical method and the minimal polynomial method. The effectiveness, ability, and simplicity of the proposed techniques were demonstrated using numerical examples. The numerical and graphical results revealed that the solutions obtained by the proposed methods were consistent with the solutions obtained by the eigenvalue and eigenvector method. The solutions were plotted using MATLAB.

### Fig 1. Figure 1.

Solutions 1(t) and 2(t) of Example 1 for different values of (a) α = 0.6, (b) α = 0.8, (c) α = 1.0 at t = 1.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 144-154https://doi.org/10.5391/IJFIS.2022.22.2.144

### Fig 2. Figure 2.

Solutions 1(t), 2(t), and 3(t) of Example 2 for different values of (a) α = 0.6, (b) α = 0.8, (c) α = 1.0 at t = 1.

The International Journal of Fuzzy Logic and Intelligent Systems 2022; 22: 144-154https://doi.org/10.5391/IJFIS.2022.22.2.144

### References

1. Torvik, PJ, and Bagley, RL (1984). On the appearance of the fractional derivative in the behavior of real materials. Journal of Applied Mechanics. 51, 294-298. https://doi.org/10.1115/1.3167615
2. Atanackovic, TM, and Stankovic, B (2004). On a system of differential equations with fractional derivatives arising in rod theory. Journal of Physics A: Mathematical and General. 37, 1241-1250. https://doi.org/10.1088/0305-4470/37/4/012
3. Daftardar-Gejji, V, and Babakhani, A (2004). Analysis of a system of fractional differential equations. Journal of Mathematical Analysis and Applications. 293, 511-522. https://doi.org/10.1016/j.jmaa.2004.01.013
4. Deng, W, Li, C, and Lu, J (2007). Stability analysis of linear fractional differential system with multiple time delays. Nonlinear Dynamics. 48, 409-416. https://doi.org/10.1007/s11071-006-9094-0
5. Charef, A, and Boucherma, D (2011). Analytical solution of the linear fractional system of commensurate order. Computers & Mathematics with Applications. 62, 4415-4428. https://doi.org/10.1016/j.camwa.2011.10.017
6. Duan, JS, Fu, SZ, and Wang, Z (2013). Solution of linear system of fractional differential equations. Pacific Journal of Applied Mathematics. 5, 93-106.
7. Daftardar-Gejji, V, and Jafari, H (2005). Adomian decomposition: a tool for solving a system of fractional differential equations. Journal of Mathematical Analysis and Applications. 301, 508-518. https://doi.org/10.1016/j.jmaa.2004.07.039
8. Garrappa, R (2013). Exponential integrators for time–fractional partial differential equations. The European Physical Journal Special Topics. 222, 1915-1927. https://doi.org/10.1140/epjst/e2013-01973-1
9. Moret, I, and Novati, P (2011). On the convergence of Krylov subspace methods for matrix Mittag–Leffler functions. SIAM Journal on Numerical Analysis. 49, 2144-2164. https://doi.org/10.1137/080738374
10. Matychyn, I, and Onyshchenko, V (2018). Matrix Mittag Leffler function in fractional systems and its computation. Bulletin of the Polish Academy of Sciences: Technical Sciences. 66, 495-500. https://doi.org/10.24425/124266
11. Garrappa, R, and Popolizio, N (2018). Computing the matrix Mittag-Leffler function with applications to fractional calculus. Journal of Scientific Computing. 77, 129-153. https://doi.org/10.1007/s10915-018-0699-5
12. Garrappa, R (2015). Numerical evaluation of two and three parameter Mittag-Leffler functions. SIAM Journal on Numerical Analysis. 53, 1350-1369. https://doi.org/10.1137/140971191
13. Duan, J (2018). System of linear fractional differential equations and the Mittag-Leffler functions with matrix variable. Journal of Physics: Conference Series. 1053. article no. 012032
14. Duan, J, and Chen, L (2018). Solution of fractional differential equation systems and computation of matrix Mittag–Leffler functions. Symmetry. 10. article no. 503
15. Agarwal, RP, Lakshmikantham, V, and Nieto, JJ (2010). On the concept of solution for fractional differential equations with uncertainty. Nonlinear Analysis: Theory, Methods & Applications. 72, 2859-2862. https://doi.org/10.1016/j.na.2009.11.029
16. Arshad, S, and Lupulescu, V (2011). On the fractional differential equations with uncertainty. Nonlinear Analysis: Theory, Methods & Applications. 74, 3685-3693. https://doi.org/10.1016/j.na.2011.02.048
17. Allahviranloo, T, Salahshour, S, and Abbasbandy, S (2012). Explicit solutions of fractional differential equations with uncertainty. Soft Computing. 16, 297-302. https://doi.org/10.1007/s00500-011-0743-y
18. Allahviranloo, T, Armand, A, and Gouyandeh, Z (2014). Fuzzy fractional differential equations under generalized fuzzy Caputo derivative. Journal of Intelligent & Fuzzy Systems. 26, 1481-1490. https://doi.org/10.3233/IFS-130831
19. Van Hoa, N (2015). Fuzzy fractional functional differential equations under Caputo gH-differentiability. Communications in Nonlinear Science and Numerical Simulation. 22, 1134-1157. https://doi.org/10.1016/j.cnsns.2014.08.006
20. Van Hoa, N (2015). Fuzzy fractional functional integral and differential equations. Fuzzy Sets and Systems. 280, 58-90. https://doi.org/10.1016/j.fss.2015.01.009
21. Long, VH, Son, NTK, and Tam, HTT (2015). Global existence of solutions to fuzzy partial hyperbolic functional differential equations with generalized Hukuhara derivatives. Journal of Intelligent & Fuzzy Systems. 29, 939-954. https//doi.org/10.3233/IFS-151623
22. Long, HV, Son, NTK, and Tam, HTT (2017). The solvability of fuzzy fractional partial differential equations under Caputo gH-differentiability. Fuzzy Sets and Systems. 309, 35-63. https://doi.org/10.1016/j.fss.2016.06.018
23. Balooch Shahriyar, MR, Ismail, F, Aghabeigi, S, Ahmadian, A, and Salahshour, S (2013). An eigenvalue-eigenvector method for solving a system of fractional differential equations with uncertainty. Mathematical Problems in Engineering. 2013. article no. 579761
24. Alijani, Z, Baleanu, D, Shiri, B, and Wu, GC (2020). Spline collocation methods for systems of fuzzy fractional differential equations. Chaos, Solitons & Fractals. 131. article no. 109510
25. Zimmermann, HJ (1991). Fuzzy Set Theory and its Applications. Boston, MA: Kluwer Academic Publishers
26. Friedman, M, Ma, M, and Kandel, A (1999). Numerical solutions of fuzzy differential and integral equations. Fuzzy Sets and Systems. 106, 35-48. https://doi.org/10.1016/S0165-0114(98)00355-8
27. Ma, M, Friedman, M, and Kandel, A (1999). Numerical solutions of fuzzy differential equations. Fuzzy Sets and Systems. 105, 133-138. https://doi.org/10.1016/S0165-0114(97)00233-9
28. 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
29. 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
30. Abasbandy, S, Allahviranloo, T, Shahryari, MB, and Salahshour, S (2012). Fuzzy local fractional differential equation. International Journal of Industrial Mathematics. 4, 231-245.
31. Allahviranloo, T. (2021) . Fuzzy Fractional Differential Operators and Equations. https://doi.org/10.1007/978-3-030-51272-9
32. Armand, A, Allahviranloo, T, and Gouyandeh, Z (2016). General solution of Basset equation with Caputo generalized Hukuhara derivative. Journal of Applied Analysis and Computation. 6, 119-130. https://doi.org/10.11948/2016010
33. Podlubny, I (1999). Fractional Differential Equations. San Diego, CA: Academic Press
34. Mainardi, F, and Gorenflo, R (2000). On Mittag-Leffler-type functions in fractional evolution processes. Journal of Computational and Applied mathematics. 118, 283-299. https://doi.org/10.1016/S0377-0427(00)00294-6
35. Horn, RA, and Johnson, CR (2012). Matrix Analysis. New York, NY: Cambridge University Press