Applied Mathematics and Computation 338 (2018) 733–738
Contents lists available at ScienceDirect
Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc
Short Communication
Some notes on properties of the matrix Mittag-Leffler function Amir Sadeghi a,∗, João R. Cardoso b,c a
Department of Mathematics, Robat Karim Branch, Islamic Azad University, Tehran, Iran Coimbra Polytechnic – ISEC, Coimbra, Portugal c Institute of Systems and Robotics, University of Coimbra, Pólo II, Coimbra 3030-290, Portugal b
a r t i c l e
i n f o
a b s t r a c t
Keywords: Matrix Mittag-Leffler function Matrix exponential Matrix function
We have come across with publications where some properties of the matrix exponential were incorrectly extended to the matrix Mittag-Leffler function, and then used as a key tool to solve certain linear matrix fractional differential equations. The main purpose of these notes is to give some clarifications on properties of the matrix Mittag-Leffler function, by explaining in detail why some identities do not hold and by providing a list of (valid) properties of this function. A sufficient condition to identify very special cases of pairs of matrices satisfying the semigroup property is given as well as examples illustrating the results. © 2018 Elsevier Inc. All rights reserved.
1. Introduction Given A ∈ Cn×n , the matrix Mittag-Leffler (ML) function with two parameters is defined through the convergent series
Eα ,β (A ) =
∞ k=0
Ak
(α k + β )
=
1
(β )
I+
A
(α + β )
+
A2
(2α + β )
+ ··· ,
(1)
where α > 0 and β are assumed to be real numbers (check for instance the recent paper [6] and the references therein). We recall that denotes the well-known Gamma function, which, for any real number x = 0, −1, −2, . . . , has been commonly defined by the convergent improper integral
(x ) =
∞
t x−1 e−t dt.
0
Setting β = 1 in (1), yields the one parameter matrix ML function Eα (A). To simplify, we will often use the notation Eα (A) instead of Eα , 1 (A). For α = β = 1, the matrix ML function is the matrix exponential, that is,
e A = E1,1 ( A ) =
∞ Ak , k! k=0
which satisfies some elegant properties such as the ones in the following lemma, where the symbols and stand, respectively, for the Kronecker product and the Kronecker sum [10, Chapter 4]. ∗
Corresponding author. E-mail addresses:
[email protected] (A. Sadeghi),
[email protected] (J.R. Cardoso).
https://doi.org/10.1016/j.amc.2018.06.037 0 096-30 03/© 2018 Elsevier Inc. All rights reserved.
734
A. Sadeghi, J.R. Cardoso / Applied Mathematics and Computation 338 (2018) 733–738
Lemma 1. Let A, B ∈ Cn×n . Then (i) e(A+B )t = eAt eBt , for all t ∈ R, if and only if AB = BA; (ii) eAB = eA eB . Proof. (i) See, for instance, [8, Section 10.1]. (ii) Follows from (i), by noticing that AI and IB always commute, regardless of A and B being or not commutative. If t = 1 in Lemma 1(i), we can conclude that if the matrices A and B commute then reciprocal may not be true. For instance, if
0 A= 0
0 2π i
0 and B = 0
e A+B
=
eA eB .
Note that, however, the
1 , 2π i
we have eA+B = eA eB = I, and AB = BA. Most of the material in the next sections relies on the general theory of matrix functions. For a detailed account of this theory, we refer the reader to the books [8] and [10, Chapter 6]. For details on the scalar ML function, see, for instance, [7]. The organization of the paper is as follows. In the next section, we provide a discussion about properties that are valid or non valid for the matrix ML function. In Section 3, Lemma 5, we give conditions that pairs of nonzero matrices must satisfy in order to verify the semigroup property. We stress out that these conditions are very restrictive so that those pairs of matrices are exceptional instances. The paper finishes with two examples illustrating the results. 2. Properties of matrix ML function When searching for literature on matrix ML functions, one have noticed that the papers [1–3] assume that the property (i) in Lemma 1, known in the literature as the semigroup property, is also valid for the one parameter matrix ML function Eα , for a general real positive number α . This property in general does not hold, even to ML scalar functions. Indeed, it is shown in [13] that Eα (z1 + z2 ) = Eα (z1 )Eα (z2 ) does not hold in general to the one parameter ML function, unless α = 1 or z1 = 0 or z2 = 0. The two parameter ML function is investigated in [4], where it is shown that the identity Eα ,β (z1 + z2 ) = Eα ,β (z1 )Eα ,β (z2 ) holds only in the following cases: • α = β = 1 (exponential case); • β = 1 ∧ ( z1 = 0 ∨ z2 = 0 ); • β = 2 ∧ ( z1 = 0 ∨ z2 = 0 ). Since the definition of a matrix function f(A) is based on the behavior of the corresponding scalar function f(z) and its derivatives on the eigenvalues of A, it is immediate that the semigroup property does not hold in general, that is, for two given complex commutating matrices A and B, we have in general
Eα ( A + B ) = Eα ( A )Eα ( B ). About the extension of (ii) in Lemma 1 to the ML-function, we will give an example in Section 4 showing that in general
Eα ( A B ) = Eα ( A ) Eα ( B ).
(2)
Hence, the methods proposed in [1–3] to solve some linear matrix fractional differential equations may not be trustworthy. Moreover, other properties of the matrix exponential, such as, eA is nonsingular and (eA )−1 = e−A , are not valid in general to the matrix ML function. Since Eα ,β (z) may be zero for some α ,β and z (see, for instance, [7, Corollary 3.10] and [9]), Eα ,β (A) may be a singular matrix. Note, however, that if α ≈ 1 then
Eα ( A + B ) ≈ Eα ( A )Eα ( B ), Eα ( A B ) ≈ Eα ( A ) Eα ( B ). To understand this claim, suppose that A is a fixed matrix and α varies in the interval [0, ∞]. Let us consider the function f : [0, ∞] −→ Cn×n defined by
f (α ) =
∞ k=0
Ak . (α k + 1 )
Then it is easy to see that f is a continuous function on [0, ∞] and that
lim f (α ) = eA .
α →1
Lemma 2. Let A, B ∈ Cn×n , and α and β are real numbers with α > 0. Let us denote the identity and zero matrix by I and 0, respectively. The following properties hold.
A. Sadeghi, J.R. Cardoso / Applied Mathematics and Computation 338 (2018) 733–738
735
(i) If there exists a non-singular matrix S such that A = SBS−1 , then Eα ,β (A ) = S Eα ,β (B ) S−1 . (ii) Eα ,1 (0 ) = Eα ,2 (0 ) = I. (iii) For each A, there exists a polynomial p(z) with complex coefficients such that Eα ,β (A ) = p(A ). (iv) (v) (vi) (vii)
If A = diag (a11 , . . . , ann ), then Eα ,β (A ) = diag Eα ,β (a11 ), . . . , Eα ,β (ann ) . Am Eα ,β (A ) = Eα ,β (A )Am , for any integer m. B Eα ,β (AB ) = Eα ,β (BA ) B. Eα ,β (A In ) = Eα ,β (A ) I, and Eα ,β (I A ) = I Eα ,β (A );
H
T
(viii) Eα ,β AH = Eα ,β (A ) (XH denotes the conjugate transpose of X) and Eα ,β AT = Eα ,β (A ) . (ix) If A is Hermitian (respectively, symmetric), then Eα ,β (A) is Hermitian (resp., symmetric). (x) Eα ,β (A ) = 21π i C Eα ,β (z )(zI − A )−1 dz (Cauchy integral formula), where C is a closed contour enclosing the spectrum of A. Proof. Since Eα ,β is an entire function (check [6]), the properties above follow easily from the general theory of matrix functions (see [8, Chapter 1] and [10, Chapter 6]). In addition, many specific properties of the scalar ML function can be easily generalized to the matrix case. Some of them are included in the next lemma. However, some care must be taken when properties of scalar functions involving derivatives or integrals are generalized to the matrix scenario. We recall, for instance, that the derivative that has a major role in the matrix case, and whose extension is non-trivial, is the Fréchet derivative. This derivative is the key for investigating the conditioning of matrix functions. See our discussion on this derivative below. Let m ≥ 2 be a positive integer. Given a matrix A ∈ Cn×n with eigenvalues not belonging to the closed negative real axis, there exists a unique matrix X such that X m = A, whose eigenvalues lie on the sector of the complex plane defined by −π /m < arg(z ) < π /m, where arg(z ) denotes the argument of the complex z. This unique matrix X is called the principal mth root of A and will be denoted by A1/m . For background on matrix mth roots see [8] and the references therein. Note that the case m = 2 corresponds to the well-known matrix square root. Lemma 3. Let A ∈ Cn×n , and α and β be real numbers with α > 0. Assume that m ≥ 1 and r are natural numbers. The following identities hold:
mEmα ,β (Am ) =
m −1
Eα ,β e2π ki/m A ;
(3)
k=0
mAr Emα ,β +rα (Am ) =
m −1
e2π ki(m−r )/m Eα ,β e2π ki/m A ; (m > r )
(4)
k=0
Am Eα ,β +mα (A ) = Eα ,β (A ) −
m −1 k=0
Ak ( β ≥ 0 ). (kα + β )
(5)
Moreover, if A has no eigenvalues belonging to the closed negative real axis, then
Eα ,β (A ) =
m−1 1 E mα ,β e2π ki/m A1/m ; m
(6)
k=0
E2α ,β (A ) =
1 Eα ,β A1/2 + Eα ,β −A1/2 . 2
(7)
Proof. Formulae (3) and (4) can be derived from the scalar counterparts included in [14, Section 1.2.5]. Formula (6) comes from (3) by replacing α by α /m and A by A1/m . The scalar version of (5) is due to [15]. In Lemma 3, we have assumed that the eigenvalues of A lie off the closed negative real axis to guarantee that the principal mth root function f (A ) = A1/m is defined in accordance with the standard definition of matrix function (see [8, Definition 1.1]). We recall that if A has an eigenvalue λ on the closed negative real axis, its scalar complex counterpart f (z ) = z1/m is not holomorphic on λ and thus the general properties of matrix functions are not valid for f(A). We believe that some properties in Lemma 3 may be exploited to design efficient algorithms for computing the matrix ML function. For instance, the identity (7) can be used consecutively to reduce the value of the parameter α and to put the spectrum of A close to the identity matrix I. We recall that A1/m → I as m → ∞. Given a map f : Cn×n → Cn×n , the Fréchet derivative of f at X ∈ Cn×n in the direction of ∈ Cn×n is a linear operator Lf (X) that maps the “direction matrix” to Lf (X, ) such that
f (X + ) − f (X ) − L f (X, ) = 0. → 0 lim
736
A. Sadeghi, J.R. Cardoso / Applied Mathematics and Computation 338 (2018) 733–738
The Fréchet derivative of the matrix function f may not exist at X, but if it does it is unique and coincides with the directional (or Gâteaux) derivative of f at X in the direction . Hence the existence of the Fréchet derivative guarantees that, for any ∈ Cn×n ,
L f (X, ) = lim t→0
f (X + t ) − f (X ) . t
(8)
Attending to the properties of the scalar ML function and [8, Thm. 3.8], for any A, ∈ Cn×n the Fréchet derivative of Eα ,β , here denoted by LE (A, ), exists and can be evaluated from (8). Since the matrix ML function can be represented by means of a Taylor series, some calculation using (8) leads to (check [12, Eq. (2.2)])
LE (A, ) =
∞
m −1
1
(mα + β ) m=1
Ak Am−1−k .
k=0
Another possibility of finding the Fréchet derivative is through the computation of the ML function of a block matrix of size 2n × 2n (see [8, Eq. (3.16)]):
Eα ,β
A 0
A
E (A ) = α ,β 0
LE (A, ) , Eα ,β (A )
and reading off the block (1,2). We end this section by proposing an explicit formula for the ML function of a diagonalizable matrix of order 2.
Lemma 4. Let A =
a c
b d
be a diagonalizable matrix of order 2 with eigenvalues λ1 = (a+d2)− and λ2 = (a+d2)+ (not nec-
essarily distinct), where := function of A is
(a − d )2 + 4bc. Set ei := Eα ,β (λi ), i = 1, 2. Provided that and c are not zero, the matrix ML
1 (d − a )(e1 − e2 ) + (e1 + e2 ) Eα ,β (A ) = −2c (e1 − e2 ) 2 where =
−2b(e1 − e2 ) , −(d − a )(e1 − e2 ) + (e1 + e2 )
(9)
(a − d )2 + 4bc.
Proof. Since A is diagonalizable, we can write A = SDS−1 , where
D = diag(λ1 , λ2 ) =
( a + d ) − 2
0
and S =
( a + d )+
0
2
Then the result follows from (i) of Lemma 2.
a−d− 2c
1
a−d+ 2c
1
.
If, in Lemma 4, c or are zero, simpler closed formulae to the ML function of a diagonalizable matrix of order 2 can be easily derived. 3. Pairs of matrices satisfying the semigroup property A trivial pair of matrices satisfying the semigroup property (10) is to take one of the matrices, say A, as the zero matrix and consider B as being arbitrary. If A and B are nonzero, a complete characterization of all pairs of matrices, for which the semigroup property is valid, looks a very a hard task. However, we are able to provide, in next lemma, some special pairs of nonzero matrices A and B satisfying that property. Lemma 5. Let A, B ∈ Cn×n and α > 0. If A and B are nilpotent with index 2 (i.e., A2 = B2 = 0) and AB = BA = 0, then
Eα ( A + B ) = Eα ( A )Eα ( B ).
(10)
Proof. Using the series expansion (1), we have
Eα ( A + B ) = I +
A+B ( A + B )2 + + ··· . (α + 1 ) (2α + 1 )
Since A and B satisfy the conditions mentioned in the statement of the lemma, (11) simplifies to
Eα ( A + B ) = I +
A+B
(α + 1 )
.
Now, a simple calculation, shows that
Eα ( A )Eα ( B ) = I +
A+B
(α + 1 )
and then the result follows.
,
(11)
A. Sadeghi, J.R. Cardoso / Applied Mathematics and Computation 338 (2018) 733–738
737
Notice that all the eigenvalues of a nilpotent matrix are zero and that the Jordan canonical form of a nilpotent matrix with index 2 involves only Jordan blocks of order 1 and 2 (see [11]):
0 J1 (0 ) = 0 or J2 (0 ) = 0
1 . 0
If we restrict the result of Lemma 5 to the case 2 × 2, it is possible to find a characterization of all pairs of commuting nilpotent matrix A and B with index 2. Indeed, a little calculation shows that every 2 × 2 nilpotent matrix with index 2 has the form
A=
±a2 , −ab
ab ∓b2
(12)
for some complex numbers a and b not simultaneously zero. If B is an 2 × 2 nilpotent of index 2 and commutes with A, then B = μA, for some μ ∈ C. Hence, all these pairs A and B involve nonzero matrices and verify (10). We have also investigated if the previous lemma can be extended to pairs of commuting nilpotent matrices of order n ≥ 3 with index 3, but the answer is in general negative. Surprisingly, an analogue version of Lemma 5 to the identity
Eα ( A B ) = Eα ( A ) Eα ( B )
(13)
is of little interest. Indeed, if a condition like A B = 0 is imposed, then A = 0 or B = 0. If A and B are arbitrary matrices of order n, we know that the pairs (0, B) and (A, 0) satisfy (13). However, we are not aware if there is any pair (A, B), with both A and B being nonzero matrices, for which (13) holds. 4. Examples The computations performed in the first example aim at illustrating why (13) may not be valid. They have been carried out in MATLAB, with unit roundoff u ≈ 1.1 × 10−16 . The ML function of the involved 2 × 2 matrices has been computed by formula (9), whereas the scalar ML function and the ML function of the matrix 4 × 4 have been evaluated, respectively, by the methods of [5] and [6]. The results are rounded and just 4 decimal places are shown. Example 1. Consider the following diagonalizable commutating matrices
−3 , 1
2 A= −1
−1 B= −1
−3 . −2
By (9), the matrix ML function evaluated at A and B for α = 0.8 is, respectively,
E0.8 ( A ) =
68.7603 −29.5424
One the other hand,
−88.6273 , 39.2179
⎡
63.9381 ⎢−24.7984 E0.8 ( A ) E0.8 ( B ) = ⎣ −27.4706 10.6545 Since
⎡
−3 0 0 −1
1 ⎢−1 AB=⎣ −1 0 and
−3 0 0 −1
⎡
73.7723 ⎢−31.7450 E0.8 ( A B ) = ⎣ −31.7450 13.6758
E0.8 ( B ) =
−74.3952 39.1397 31.9634 −16.8161
0.9299 −0.3606
−82.4118 31.9634 36.4675 −14.1439
−1.0819 . 0.5692
⎤
95.8903 −50.4484⎥ . −42.4318 ⎦ 22.3236
⎤
0 −3⎥ −3⎦ −1 −95.2350 42.0273 41.0273 −18.0693
−95.2350 41.0273 42.0273 −18.0693
⎤
123.0818 −54.2078⎥ , −54.2078⎦ 23.9580
we observe that E0.8 (AB) = E0.8 (A)E0.8 (B). This examples illustrates the relation (2). Example 2. The following matrices A and B are instead examples of the very special cases of matrices satisfying the conditions of Lemma 5 for which, consequently, as an exceptional case the semigroup property (10) hold:
⎡
−35 ⎢ 70 A=⎣ 14 77
−4 6 4 8
−7 13 4 15
⎤
−11 24 ⎥ , 2 ⎦ 25
⎡
98 ⎢ −112 B=⎣ −140 −182
4 −6 −4 −8
16 −19 −22 −30
⎤
38 −42⎥ . −56⎦ −70
738
A. Sadeghi, J.R. Cardoso / Applied Mathematics and Computation 338 (2018) 733–738
Acknowledgments The work of Amir Sadeghi is supported by Robat Karim branch, Islamic Azad University, Tehran, Iran, and the work of João R. Cardoso by ISR-University of Coimbra (project UID/EEA/0 0 048/2013) funded by “Fundação para a Ciência e a Tecnologia” (FCT). References [1] Z. Al-Zhour, The general (vector) solutions of such linear (coupled) matrix fractional differential equations by using Kronecker structures, Appl. Math. Comput. 232 (2014) 498–510. [2] Z. Al-Zhour, The general solutions of singular and non-singular matrix fractional time-varying descriptor systems with constant coefficient matrices in Caputo sense, Alexandria Eng. J. 55 (2016) 1675–1681. [3] Z. Al-Zuhiri, Z. Al-Zhour, K. Jaber, The exact solutions of such coupled linear matrix fractional differential equations of diagonal unknown matrices by using hadamard product, J. Appl. Math. Phys. 4 (2016) 432–442. [4] S.K. Elagan, On the invalidity of semigroup property for the Mittag-Leffler function with two parameters, J. Egypt. Math. Soc. 24 (2016) 200–203. [5] R. Garrappa, Numerical evaluation of two and three parameter Mittag-Leffler function, SIAM J. Numer. Anal. 53 (2015) 135–169. [6] R. Garrappa, M. Popolizio, Computing the matrix Mittag-Leffler function with applications to fractional calculus, J. Sci. Comput. (2018), doi:10.1007/ s10915- 018- 0699- 5. [7] R. Gorenflo, A. Kilbas, F. Mainardi, S. Rogosin, Mittag-Leffler Functions, Related Topics and Applications, Springer-Verlag, 2014. [8] N.J. Higham, Functions of Matrices: Theory and Computation, Society for Industrial and Applied Mathematics, Philadelphia, 2008. [9] R. Hilfer, H.J. Seybold, Computation of the generalized Mittag-Leffler function and its inverse in the complex plane, Integral Transforms Special Functions 17 (9) (2006) 637–652. [10] R.A. Horn, C.R. Johnson, Topics in Matrix Analysis, Cambridge University Press, Cambridge, 1994. Paperback Edition. [11] R.A. Horn, C.R. Johnson, Matrix Analysis, second ed., Cambridge University Press, 2013. [12] C.S. Kenney, A.J. Laub, Condition estimates for matrix functions, SIAM J. Matrix Anal. Appl. 10 (1989) 191–209. [13] J. Peng, K. Li, A note on property of the Mittag-Leffler function, J. Math. Anal. Appl. 370 (2010) 635–638. [14] I. Podlubny, Fractional Differential Equations, first ed., Academic Press, New York, 1999. [15] R.K. Saxena, Certain properties of generalized Mittag-Leffler function, in: Proceedings of the 3rd Annual Conference of the Society for Special Functions and Their Applications, Chennai, India, 2002, pp. 77–81.