  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.
