Numerical stability of Grünwald-Letnikov method for time fractional delay differential equations

This paper is concerned with the numerical stability of time fractional delay differential equations (F-DDEs) based on Gr\"{u}nwald-Letnikov (GL) approximation (also called fraction backward Euler scheme) for the Caputo fractional derivative, in particular, the numerical stability region and the Mittag-Leffler stability. Using the boundary locus technique, we first derive the exact expression of the numerically stability region in the parameter plane, and show that the fractional backward Euler scheme based on GL scheme is not $\tau(0)$-stable, which is different from the backward Euler scheme for integer DDE models. Secondly, we also prove the numerical Mittag-Leffler stability for the numerical solutions provided that the parameters fall into the numerical stability region, by employing the singularity analysis of generating function. Our results show that the numerical solutions of F-DDEs are completely different from the classical integer order DDEs, both in terms of $\tau(0)$-stabililty and the long-time decay rate.