FoxHfunctions in fractional diffusion

Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 www.elsevier.com/locate/cam Fox H functions in fractional diffusion Francesco Mainardia,∗ , Gianni Pagninib , R.K. Saxenac a Dipartimento di Fisica, Università di Bologna and INFN, Sezione di Bologna, Via Irnerio 46, I-40126 Bologna, Italy b Istituto per le Scienze dell’Atmosfera e del Clima (ISAC) del CNR, Via Gobetti 101, I-40129 Bologna, Italy c Department of Mathematics and Statistics, Jan Narain Vyas University, Jodhpur 342005, India Received 10 November 2003; received in revised form 18 June 2004 Abstract The H functions, introduced by Fox in 1961, are special functions of a very general nature, which allow one to treat several phenomena including anomalous diffusion in a unified and elegant framework. In this paper we express the fundamental solutions of the Cauchy problem for the space–time fractional diffusion equation in terms of proper Fox H functions, based on their Mellin–Barnes integral representations. We pay attention to the particular cases of space-fractional, time-fractional and neutral-fractional diffusion. © 2004 Elsevier B.V. All rights reserved. MSC: 26A33; 33C20; 33C60; 33E12; 33E20; 33E30; 44A15; 60G18; 60J60 Keywords: Fox H-functions; Mellin–Barnes integrals; Fractional derivatives; Fractional diffusion; Probability distributions 1. Introduction The H functions, introduced by Fox [4] in 1961 as symmetrical Fourier kernels, can be regarded as the extreme generalization of the generalized hypergeometric functions p Fq , beyond the Meijer G functions. Like the Meijer G functions, the Fox H functions turn out to be related to the Mellin–Barnes integrals and to the Mellin transforms, but in a more general way. After Fox, the H functions were carefully investigated by Braaksma [2], who provided their convergent and asymptotic expansions in the complex plane, based on their Mellin–Barnes integral representation. ∗ Corresponding author. Tel.: +39 051 209 1098; fax: +39 051 247244. E-mail address: mainardi@bo.infn.it (F. Mainardi). 0377-0427/$ - see front matter © 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.cam.2004.08.006 322 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 More recently, the H functions, being related to the Mellin transforms, have been recognized to play a fundamental role in the probability theory and in fractional calculus as well as in their applications, including non-Gaussian stochastic processes and phenomena of nonstandard (i.e. anomalous) relaxation and diffusion, see e.g. [1,11–13,15,26–29,32–36]. In section 2, we summarize the essential definitions and notations for the Fox H functions. In section 3, we introduce the partial differential equation of fractional order (both in space and in time), that is intended to generalize in a proper way the standard equation for normal diffusion. We also recall the main results of this generalized equation based on the Fourier–Laplace representation of its fundamental solution, the so-called Green function. Then, in section 4, we provide for the general Green function a representation in terms of Mellin–Barnes integrals and, consequently, in terms of Fox H functions. We then concentrate our attention to the particular but relevant cases of space fractional, time fractional and neutral fractional diffusion for which the corresponding Green functions are clearly interpreted as probability densities. Further properties regarding the Green function in the general cases of space–time fractional diffusion can be extracted from the analysis contained in [21] where, however, the passage from the Mellin–Barnes integrals to the corresponding H-functions is not treated. 2. The Fox H functions According to a standard notation, the Fox H function is defined as
1 m,n Hp,q (z) = Hm,n (s) zs ds, 2
i L p,q (2.1) where L is a suitable path in the complex plane C to be disposed later, zs = exp{s(log |z| + i arg z)}, and Hm,n p,q (s) = A(s) = C(s) = m  j =1 A(s) B(s) , C(s) D(s) (2.2) (bj − j s), q  j =m+1 B(s) = (1 − bj + j s), n  j =1 (1 − aj + j s), D(s) = p  j =n+1 (aj − j s) (2.3) (2.4) with 0
n
p , 1
m
q , {aj , bj } ∈ C , {j , j } ∈ R+ . An empty product, when it occurs, is taken to be one so n = 0 ⇐⇒ B(s) = 1, m = q ⇐⇒ C(s) = 1, n = p ⇐⇒ D(s) = 1. Due to the occurrence of the factor zs in the integrand of (2.1), the H function is, in general, multi-valued, but it can be made one-valued on the Riemann surface of log z by choosing a proper branch. We also note that when the ’s and ’s are equal to 1, we obtain the Meijer’s G-functions Gm,n p,q (z). F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 323 The above integral representation of the H functions, by involving products and ratios of Gamma functions, is known to be of Mellin–Barnes integral type.1 A compact notation is usually adopted for (2.1):     (a ,  )  j j j =1,...,,p m,n m,n Hp,q (z) = Hp,q z  . (2.5)  (bj , j )j =1,...,q Thus, the singular points of the kernel H are the poles of the Gamma functions entering the expressions of A(s) and B(s), that we assume do not coincide. Denoting by P(A) and P(B), the sets of these poles, we write P(A) ∩ P(B) = ∅. The conditions for the existence of the H-functions can be made by inspecting the convergence of integral (2.1), which can depend on the selection of the contour L and on certain relations between the parameters {ai , i } (i = 1, . . . , p) and {bj , j } (j = 1, . . . , q). For the analysis of the general case we refer to the specialized treatises on H functions, e.g. [27,28,35] and, in particular to the paper by Braaksma [2], where an exhaustive discussion on the asymptotic expansions and analytical continuation of these functions is found; see also [12]. In the following we limit ourselves to recall the essential properties of the H functions preferring to later analyse in detail those functions related to fractional diffusion. As it will be shown later, this phenomenon depends on one real independent variable and three parameters; in this case we shall have z = x ∈ R and m
2, n
2, p
3, q
3. The contour L in (2.1) can be chosen as follows: (i) L = L−i∞,+i∞ chosen in a manner to go from −i∞ to +i∞ leaving to the right all the poles of P(A), namely the poles sj,k = (bj + k)/j ; j = 1, 2, . . . , m; k = 0, 1, . . . of the functions  entering A(s), and to left all the poles of P(B), namely the poles sj,l = (aj − 1 − l)/j ; j = 1, 2, . . . , n; l = 0, 1, . . . of the functions  entering B(s). (ii) L = L+∞ is a loop beginning and ending at +∞ and encircling once in the negative direction all the poles of P(A), but none of the poles of P(B). (iii) L = L−∞ is a loop beginning and ending at −∞ and encircling once in the positive direction all the poles of P(B), but none of the poles of P(A). Braaksma has shown that, independently of the choice of L the Mellin–Barnes integral makes sense and defines an analytic function of z in the following two cases:  > 0,  = 0, 0 < |z| < ∞, 0 < |z| < , where where = = q  j =1 p  j =1 j − − j j p  j , (2.6) j =1 q  j j . (2.7) j =1 On account of the following useful and important formula for the H-function        (a ,  )   j j 1,p n,m 1  (1 − bj , j )1,q m,n Hp,q z  , = Hq,p  (bj , j )1,q z  (1 − aj , j )1,p (2.8) 1 As a historical note, we point out that the names refer to the two authors, who in the first 1910s developed the theory of these integrals using them for a complete integration of the hypergeometric differential equation. However, these integrals were first used in 1888 by Pincherle, see e.g. [23]. Recent treatises on Mellin–Barnes integrals are those in [25,30]. 324 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 we can transform the H-function with  < 0 and argument z to one with  > 0 and argument 1/z. This property is suitable to compare the results of the theory of H functions based on (2.1) with zs with the other one with z−s , often used in the literature. Other important properties of the Fox H functions, that can be easily derived from their definition, are included in the list below: (i) The H-function is symmetric in the set of pairs (a1 , 1 ), . . . , (an , n ), (an+1 , n+1 ), . . . , (ap , p ), (b1 , 1 ), . . . , (bm , m ) and (bm+1 , m+1 ), . . . , (bq , q ). (ii) If one of the (aj , j ), j = 1, . . . , n, is equal to one of the (bj , j ), j = m + 1, . . . , q; [or one of the pairs (aj , j ), j = n + 1, . . . , p is equal to one of the (bj , j ), j = 1, . . . , m], then the H-function reduces to one of the lower order, that is, p, q and n [or m] decrease by a unity. Provided n  1 and q > m, we have         (a ,  ) (a ,  ) j j 1,p   j j 2,p m,n−1 m,n Hp,q z z = Hp−1,q−1 , (2.9)  (bj , j )1,q−1 (a1 , 1 )  (bj , j )1,q−1        (a ,  )  (a ,  ) (b ,  )  j j 1,p−1  j j 1,p−1 1 1 m−1,n m,n . (2.10) z z = Hp−1,q−1 Hp,q  (bj , j )2,q  (b1 , 1 ) (bj , j )2,q (iii) z (iv)  m,n Hp,q        (a +  ,  )  (a ,  ) j j 1,p  j  j j 1,p m,n z . = Hp,q z   (bj + j , j )1,q  (bj , j )1,q        1 m,n  (aj , j )1,p c  (aj , cj )1,p m,n z H , = Hp,q z   (bj , j )1,q  (bj , cj )1,q c p,q (2.11) c > 0. (2.12) The convergent and asymptotic expansions (for z → 0 or z → ∞) are mostly obtained by applying the residue theorem in the poles (assumed to be simple) of the Gamma functions entering A(s) or B(s) that are found inside the specially chosen path. In some cases (in particular if n = 0 ⇐⇒ B(s) = 1) we find an exponential asymptotic behaviour. In the presence of a multiple pole s0 of order N the treatment becomes more cumbersome because we need to expand in power series at the pole the product of the involved functions, including zs , and to take the first N terms up to (s − s0 )N −1 inclusive. Then the coefficient of (s − s0 )N −1 is the required residue. Let us consider the case N = 2 (double pole) of interest for the fractional diffusion. Then, the expansions for the Gamma functions are of the type (s) = (s0 ) [1 + (s0 )(s − s0 ) + O((s − s0 )2 )], (s) = s → s0 , s0  = 0, −1, −2, . . . , (−1)k [1 + (k + 1)(s + k) + O((s + k)2 )], (k + 1)(s + k)] s → −k, where k = 0, 1, 2, . . . and (z) denotes the logarithmic derivative of the  function, (z) = d ′ (z) log (z) = , dz (z) F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 325 whereas the expansion of zs yields the logarithmic term zs = zs0 [1 + log z(s − s0 ) + O((s − s0 )2 ], s → s0 . 3. The fractional diffusion equation An interesting way to generalize the classical diffusion equation * 2 *x 2 * u(x, t) = *t u(x, t), −∞ < x < + ∞, t 0 (3.1) is to replace in (3.1) the partial derivatives in space and time by suitable linear integro-differential operators, to be intended as derivatives of noninteger order, that allows the corresponding Green function (see below) to be still interpreted as a probability density evolving in time but with an appropriate similarity law. It turns out that this generalized diffusion equation, that we refer to as space–time fractional diffusion equation, is xD   u(x, t) = t D∗ u(x, t), −∞ < x < + ∞, t  0, (3.2) where the , ,  are real parameters restricted as follows: 0 < 
2, | |
min(, 2 − ), 0 < 
2. (3.3)  Here x D  and t D∗ are integro–differential operators, the Riesz–Feller space-fractional derivative of order  and asymmetry and the Caputo time-fractional derivative of order , respectively. The allowed region for the parameters  and in the plane {, } is called the Feller–Takayasu diamond, see e.g. [8,9,21]. The relevant cases of the space–time fractional diffusion equation (3.2) include, in addition to the standard case of normal diffusion { = 2,  = 1}, the space-fractional diffusion {0 <  < 2,  = 1}, the time-fractional diffusion { = 2, 0 <  < 2} and the neutral-fractional diffusion {0 <  =  < 2}. Let us now resume the essential definitions of the fractional derivatives in (3.2) based on their Fourier and Laplace representations. By denoting the Fourier transform of a sufficiently well-behaved (generalized) function f (x), f ( ) = +∞ F {f (x); } = −∞ e +i x f (x) dx, ∈ R, the Riesz–Feller space-fractional derivative of order  and skewness turns out to be defined by F x D  f (x); (3.4) = − ( )f ( ),  ( ) = | | ei(sign )
/2 , 0 < 
2, | |
min {, 2 − }. (3.5) Thus, we recognize that the Riesz–Feller derivative is required to be the pseudo-differential operator whose symbol − ( ) is the logarithm of the characteristic function of a general Lévy strictly stable probability density with index of stability  and asymmetry parameter (improperly called skewness) according to Feller’s parameterization as revisited by Gorenflo et al., see e.g. [8,9]. 326 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 For = 0 we have a symmetric operator with respect to x, that can be interpreted as /2 d2  . x D0 = − − dx 2 (3.6) This can be formally deduced by writing −| | = −( 2 )/2 . For 0 <  < 2 and | |
min {, 2 − }, the Riesz–Feller derivative can be shown to admit the integral representation in the x domain, 
∞ (1 + ) f (x + ) − f (x)  sin [( + )
/2] d x D f (x) = 1+
0 
∞ f (x − ) − f (x) d . (3.7) + sin[( − )
/2] 1+ 0 By denoting the Laplace transform of a sufficiently well-behaved (generalized) function f (t), f(s) = ∞ L {f (t); s}= 0 e −st f (t) dt, R(s) > af , the Caputo time-fractional derivative of order  (m−1 < 
m, m ∈ N) turns out to be defined through  L{t D∗ f (t); s} = s  f(s) − m−1  s −1−k f (k) (0+ ), k=0 m − 1 < 
m. (3.8) This leads to define, see e.g. [7,31],  (m) ( ) d 1 t f   , m − 1 <  < m,  (m − ) 0 (t − )+1−m (3.9) t D∗ f (t) := m   d f (t),  = m. dt m The reader should observe that the Caputo fractional derivative represents a sort of regularization in the time origin for the classical Riemann–Liouville fractional derivative,2 see e.g. [7,31]. When the diffusion equations (3.1), (3.2) are equipped by the initial and boundary conditions u(x, 0+ ) = (x), u(± ∞, t) = 0, (3.10) +∞ their solution reads u(x, t) = −∞ G( , t) (x − ) d , where G(x, t) denotes the fundamental solution (known as the Green function) corresponding to (x) = (x), the Dirac generalized function.3 It is straightforward to derive from (3.2) the Fourier–Laplace transform of the Green function by taking into account the Fourier transform for the Riesz–Feller space-fractional derivative, see (3.4)–(3.5), and the Laplace transform for the Caputo time-fractional derivative, see (3.8). We have    ( , s) − s −1 ,  ( , s) = s  G − ( )G ,  ,  (3.11) 2 We note that the Caputo fractional derivative coincides with that introduced independently by Djrbashian and Nersesian, which has been adopted by Kochubei [16], for treating initial value problems in the presence of fractional derivatives. 3 We note that when 1 < 
2 to Eq. (3.2) we must add a second initial condition of type u (x, 0+ ) = (x), which implies t two Green functions corresponding to {u(x, 0+ ) = (x), ut (x, 0+ ) = 0} and {u(x, 0+ ) = 0, ut (x, 0+ ) = (x)}. Here we limit ourselves to consider only the first Green function. For the time-fractional diffusion equation the second Green function has been investigated in [22]. F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 327 so that   ( , s) = G ,  s −1 s  +  ( ) (3.12) . By using the known scaling rules for the Fourier and Laplace transforms, we infer without inverting the two transforms, G, (x, t) = t − K, (x/t  ), (3.13)  = /, where the one-variable function K, is the reduced Green function and x/t  is the similarity variable.   (0, s) = 1/s ⇐⇒ G  (0, t) = 1, the normalization property We note from G ,  ,
+∞ −∞ G, (x, t) dx =
+∞ −∞ K, (x) dx = 1, (3.14) and, from  ( ) =  (− ) = −  (− ), the symmetry relation K, (−x) = K−, (x), (3.15) which allows us to restrict our attention to x > 0. When  = 2 ( = 0) and  = 1 the inversion of the Fourier–Laplace transform in (3.12) is trivial: we recover the Gaussian density, evolving in time with variance 2 = 2t, typical of the normal diffusion, 1 G02,1 (x, t) = √ exp(−x 2 /(4t)), 2
t x ∈ R, t > 0, (3.16) which exhibits the similarity law (3.13) with  = 21 . 4. Mellin–Barnes and Fox H representations of the Green function Mainardi et al. [21] have inverted the Fourier–Laplace transform (3.12) of the Green function by passing through the Mellin transform. Here we recall and complement their main results by introducing the representation of the reduced Green function in terms of proper Fox H functions, starting from its general Mellin–Barnes integral representation for x > 0,
+i∞ 1 1 (s/)(1 − s/)(1 − s) K, (x) = x s ds, (4.1) x 2
i −i∞ ([( − )/2]s)(1 − [( − )/2]s)(1 − (/)s) with 0 <  < min(, 1) under the condition | |
2 − . From the Mellin–Barnes representation (4.1) we now derive the representation of K, (x) in terms of a proper H function, taking into account the theory of Fox functions briefly summarized in Section 2. At first we distinguish the two cases (a)  <  and (b)  >  for which the corresponding H functions turn out to be singular in x = 0 and ∞, respectively. Taking x > 0, we get     (0, 1 ) (0, 1) (0, − ) 1 1   2 H 1,2 K, (x) = ,  < , (4.2a)  x 3,3 x  (0, 1 ) (0,  ) (0, − )   2 328 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331     (1, 1 ) (1,  ) (1, − ) 1    2 , H 2,1 x  K, (x) =  (1, 1 ) (1, 1) (1, − ) x 3,3  2  > . (4.2b) When  =  the corresponding H function is singular in z = x = 1. However, the singularity is removable because, surprisingly, the corresponding (reduced) Green function can be expressed (in explicit form) in terms of a (nonnegative) elementary function, that we denote by N (x), as it is shown in [21]. We refer to this case as to neutral-fractional diffusion and the corresponding representation through H functions is redundant. Explicitly we write, for x > 0, Neutral diffusion : 0 <  =  < 2;
min{, 2 − }, K, := N (x) = x −1 sin[(
/2)( − )] .
1 + 2x  cos[(
/2)( − )] + x 2 1 (4.3) As far as we know, this case of fractional diffusion seems not so well treated in the literature. We note that N (x) may be considered the fractional generalization (with skewness) of the well-known (symmetric) Cauchy density. For the other particular cases outlined in Section 3 we have to properly use properties (2.9)–(2.11) in general expressions (4.2a)–(4.2b) in order to obtain the corresponding representations in terms of simpler Fox H functions. Normal diffusion :  = 2,  = 1; = 0. The case of normal (or standard) diffusion is known to be characterized by the Gaussian probability density function. Indeed the reduced Green function reads 1 0 K2,1 (x) := D(x) = √ exp(−x 2 /4), 2
so, for x > 0, we have 1 1 D(x) = 2x 2
i
+i∞ −i∞ x ∈ R,    (1, 1) 1  1,0 x 2 2 . x s ds = H1,1  (0, 1) (1 − s/2) 2 (1 − s) (4.4) (4.5) Space-fractional diffusion: 0 <  < 2,  = 1; | |
min{, 2 − }. In this case the reduced Green function K,1 (x) is known to be the -strictly stable Lévy density that we denote by L (x). Then, for x > 0, we have
+i∞ (s/)(1 − s) 1 1 x s ds (4.6) K,1 (x) := L (x) = x 2
i −i∞ ([( − )/2]s)(1 − [( − )/2]s) with 0 <  < min(, 1), Then, by distinguishing the two cases as in Eqs. (4.2), we obtain: (a) 0 <  < 1; | |
,    − − 1 1,1 1  (1, 1) ( 2 , 2 ) . L (x) = H2,2   x  ( 1 , 1 ) ( 2− , 2− ) (b) 1 <  < 2; | |
2 − ,     (1 − 1 , 1 ) (1 − − , − ) 1 1,1    2 2 L (x) = H2,2 x  .  (0, 1) (1 − − , − )  2 2 (4.7a) (4.7b) 329 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 We note that it was Schneider [33], who first in 1986 has recognized that all stable probability densities can be represented in terms of Fox H functions. Previously, the stable (non-Gaussian) densities were known in general through their (convergent and 1 −1 asymptotic) series representations that in a few particular cases, { = 31 , = −1 3 }, { = 2 , = 2 }, { = 23 , = 0}, { = 23 , = 21 }, were interpreted in terms of known special functions. In his remarkable (but almost entirely neglected) article, Schneider has also pointed out the errors present in the literature for some of the above particular cases. More recently, the representation of the stable densities through Mellin–Barnes integrals has been exhaustively treated by a number of authors as in [36,21]. Time-fractional diffusion :  = 2, 0 <  < 2; = 0. In this case the reduced Green function K,1 (x) is known to be a probability density with stretched exponential tails, that we denote (for historical reasons) by 21 M/2 where M/2 denotes a Wright-type function.4 We thus write 0 K2,  (x) := 1 2 (4.8) M/2 (x), where, for x > 0, 1 1 M/2 (x) = x 2
i
+i∞ −i∞ (1 − s) (1 − s/2) s x ds 1,0 = H1,1     x  (1 − 2 , 2 ) (0, 1)  (4.9) with 0 <  < 1. As a check we note that the simpler H function in (4.5) for the Gaussian density is recovered from (4.7a) in the limit  = 2 and from (4.8)–(4.9) in the limit  = 1. 4 The function M (z) is defined for any order  ∈ (0, 1) and ∀z ∈ C by  M (z) = ∞  n=0 ∞ 1  (−z)n−1 (−z)n = (n) sin(
n). n![−n + (1 − )]
(n − 1)! n=1 It turns out that M (z) is an entire function of order  = 1/(1 − ), which provides a generalization of the Gaussian and of the Airy function. In fact, we obtain 1 M1/2 (z) = √ exp(− z2 /4),
M1/3 (z) = 32/3 Ai(z/31/3 ). M (z) is a special case of the Wright function , (z). Originally, Wright [37–39] introduced and investigated the function , (z) := ∞  n=0 zn , n!(n + )   0, z ∈ C, with the restriction   0, in a series of notes starting from 1933 in the framework of the asymptotic theory of partitions. Only later, in 1940, he [40] considered the case −1 <  < 0. We note that in the handbook of the Bateman Project, see [3], Vol. 3, Ch. 18, presumably for a misprint,  is restricted to be nonnegative. In his first analysis of the time fractional diffusion equation, Mainardi [17], aware of the Bateman project but not of the 1940 paper by Wright, introduced the two (Wright-type) auxiliary functions, F (z) := −,0 (−z) and M (z) := −,1− (−z) with 0 <  < 1, inter-related through F (z) =  z M (z). For detailed information on the Wright-type functions (possibly related to time-fractional diffusion equations), the interested reader may consult e.g. [18–20,5,6,10,14]. 330 F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 5. Conclusions As a conclusive remark we point out that the nonnegativity of the above functions, obtained in the particular cases of neutral, space and time fractional diffusion, are relevant in proving that, in the general case of space–time fractional diffusion, the Green functions are still spatial probability densities evolving in time, provided that 0 < 
2 with 0 < 
1 and 1


2, see [21]. The proof is based on the convolution theorem for the Mellin transforms and provides interesting subordination formulas, see [24]. This fact could also be shown by using the properties of the Fox H functions. Acknowledgements The authors are grateful to the anonymous referees for useful comments and suggestions. Research performed under the auspices of the National Group of Mathematical Physics (G.N.F.M.—I.N.D.A.M.) of Italy. FM acknowledges partial support by the Italian Ministry of University (M.I.U.R) and by the National Institute of Nuclear Physics (INFN). References [1] V.V. Anh, N.N. Leonenko, Spectral analysis of fractional kinetic equations with random data, J. Statist. Phys. 104 (2001) 1349–1387. [2] B.L.J. Braaksma, Asymptotic expansions and analytical continuations for a class of Barnes-integrals, Compositio Math. 15 (1962–1963) 239–341. [3] A. Erdélyi, W. Magnus, F. Oberhettinger, F.G. Tricomi, Higher Transcendental Functions, vol. 1, McGraw-Hill, New York, 1953. [4] C. Fox, The G and H functions as symmetrical Fourier kernels, Trans. Amer. Math. Soc. 98 (1961) 395–429. [5] R. Gorenflo, Yu. Luchko, F. Mainardi, Analytical properties and applications of the Wright function, Fract. Cal. Appl. Anal. 2 (1999) 383–414. [6] R. Gorenflo, Yu. Luchko, F. Mainardi, Wright functions as scale-invariant solutions of the diffusion-wave equation, J. Comput. Appl. Math. 118 (2000) 175–191. [7] R. Gorenflo, F. Mainardi, Fractional calculus: integral and differential equations of fractional order, in: A. Carpinteri, F. Mainardi (Eds.), Fractals and Fractional Calculus in Continuum Mechanics, Springer, Wien, 1997, pp. 223–276 (reprinted with revisions in http://www.fracalmo.org). [8] R. Gorenflo, F. Mainardi, Random walk models for space-fractional diffusion processes, Fract. Cal. Appl. Anal. 1 (2) (1998) 167–191. [9] R. Gorenflo, F. Mainardi, D. Moretti, G. Pagnini, P. Paradisi, Discrete random walk models for space–time fractional diffusion, Chem. Phys. 284 (2002) 521–541. [10] R. Gorenflo, F. Mainardi, H.M. Srivastava, Special functions in fractional relaxation-oscillation and fractional diffusionwave phenomena, in: D. Bainov (Ed.), Proceedings VIII International Colloquium on Differential Equations, Plovdiv 1997, VSP International Science Publ., Utrecht, 1998, pp. 195–202. [11] R. Hilfer, Fractional time evolution, in: R. Hilfer (Ed.), Applications of Fractional Calculus in Physics, World Scientific, Singapore, 2000, pp. 87–130. [12] A.A. Kilbas, M. Saigo, On the H functions, J. Appl. Math. Stochastic Anal. 12 (1999) 191–204. [13] A.A. Kilbas, M. Saigo, H-transforms. Theory and Applications, CRC Press, Boca Raton, FL, 2004. [14] A. A Kilbas, M. Saigo, J.J. Trujillo, On the generalized Wright function, Fract. Cal. Appl. Anal. 5 (2002) 437–460. [15] V. Kiryakova, Generalized Fractional Calculus and Applications, Pitman Research Notes in Mathematics, vol. 301, Longman, Harlow, 1994. [16] A.N. Kochubei, Fractional order diffusion, J. Differential Equations 26 (1990) 485–492. F. Mainardi et al. / Journal of Computational and Applied Mathematics 178 (2005) 321 – 331 331 [17] F. Mainardi, On the initial value problem for the fractional diffusion-wave equation, in: S. Rionero, T. Ruggeri (Eds.), Waves and Stability in Continuous Media, World Scientific, Singapore, 1994, pp. 246–251. [18] F. Mainardi, The fundamental solutions for the fractional diffusion-wave equation, Appl. Math. Lett. 9 (6) (1996) 23–28. [19] F. Mainardi, Fractional relaxation-oscillation and fractional diffusion-wave phenomena, Chaos, Soliton. Fract. 7 (1996) 1461–1477. [20] F. Mainardi, Fractional calculus: some basic problems in continuum and statistical mechanics, in: A. Carpinteri, F. Mainardi (Eds.), Fractals and Fractional Calculus in Continuum Mechanics, Springer, Wien and New York, 1997, pp. 248–291 (reprinted in http://www.fracalmo.org). [21] F. Mainardi, Yu. Luchko, G. Pagnini, The fundamental solution of the space–time fractional diffusion equation, Fract. Cal. Appl. Anal. 4 (2) (2001) 153–192 (reprinted in http://www.fracalmo.org). [22] F. Mainardi, G. Pagnini, The Wright functions as solutions of the time-fractional diffusion equations, Appl. Math. Comput. 141 (2003) 51–62. [23] F. Mainardi, G. Pagnini, Salvatore Pincherle: the pioneer of the Mellin–Barnes integrals, J. Comput. Appl. Math. 153 (2003) 331–342. [24] F. Mainardi, G. Pagnini, R. Gorenflo, Mellin transform and subordination laws in fractional diffusion processes, Fract. Cal. Appl. Anal. 6 (2003) 441–459. [25] O.I. Marichev, Handbook of Integral Transforms of Higher Transcendental Functions, Theory and Algorithmic Tables, Ellis Horwood, Chichester, 1983. [26] A.M. Mathai, A few remarks on the exact distributions of certain multivariate statics; II, in: D.G. Kabe, R.P. Gupta (Eds.), Multivariate Statistical Inference, North-Holland, Amsterdam, 1973, pp. 169–181. [27] A.M. Mathai, R.K. Saxena, Generalized Hypergeometric Functions with Applications in Statistics and Physical Sciences, Lecture Notes in Mathematics, vol. 348, Springer, Berlin, 1973. [28] A.M. Mathai, R.K. Saxena, The H-function with Applications in Statistics and Other Disciplines, Wiley Eastern Ltd., New Delhi, 1978. [29] R. Metzler, J. Klafter, The random walk’s guide to anomalous diffusion: a fractional dynamics approach, Phys. Rep. 339 (2000) 1–77. [30] R.B. Paris, D. Kaminski, Asymptotic and Mellin–Barnes Integrals, Cambridge University Press, Cambridge, 2001. [31] I. Podlubny, Fractional Differential Equations, Academic Press, San Diego, 1999. [32] R.K. Saxena, T.F. Nonnenmacher, Application of the H-function in Markovian and non-Markovian chain models, Fract. Cal. Appl. Anal. 7 (2004) 135–148. [33] W.R. Schneider, Stable distributions: Fox function representation and generalization, in: S. Albeverio, G. Casati, D. Merlini (Eds.), Stochastic Processes in Classical and Quantum Systems, Lecture Notes in Physics, vol. 262, Springer, Berlin, 1986, pp. 497–511. [34] W.R. Schneider, W. Wyss, Fractional diffusion and wave equations, J. Math. Phys. 30 (1989) 134–144. [35] H.M. Srivastava, K.C. Gupta, S.P. Goyal, The H-Functions of One and Two Variables with Applications, South Asian Publishers, New Delhi, 1982. [36] V.V. Uchaikin, V.M. Zolotarev, Chance and Stability. Stable Distributions and their Applications, VSP, Utrecht, 1999. [37] E.M. Wright, On the coefficients of power series having exponential singularities, J. London Math. Soc. 8 (1933) 71–79. [38] E.M. Wright, The asymptotic expansion of the generalized Bessel function, Proc. London Math. Soc. (Ser. II) 38 (1935) 257–270. [39] E.M. Wright, The asymptotic expansion of the generalized hypergeometric function, J. London Math. Soc. 10 (1935) 287–293. [40] E.M. Wright, The generalized Bessel function of order greater than one, Quart. J. Math., Oxford Ser. 11 (1940) 36–48.