  In this work, we extend the fractional linear multistep methods in [C.
Lubich, SIAM J. Math. Anal., 17 (1986), pp.704--719] to the tempered fractional
integral and derivative operators in the sense that the tempered fractional
derivative operator is interpreted in terms of the Hadamard finite-part
integral. We develop two fast methods, Fast Method I and Fast Method II, with
linear complexity to calculate the discrete convolution for the approximation
of the (tempered) fractional operator. Fast Method I is based on a local
approximation for the contour integral that represents the convolution weight.
Fast Method II is based on a globally uniform approximation of the trapezoidal
rule for the integral on the real line. Both methods are efficient, but
numerical experimentation reveals that Fast Method II outperforms Fast Method I
in terms of accuracy, efficiency, and coding simplicity. The memory requirement
and computational cost of Fast Method II are $O(Q)$ and $O(Qn_T)$,
respectively, where $n_T$ is the number of the final time steps and $Q$ is the
number of quadrature points used in the trapezoidal rule. The effectiveness of
the fast methods is verified through a series of numerical examples for
long-time integration, including a numerical study of a fractional
reaction-diffusion model.
