The best previously known algorithm for evaluating the Riemann zeta function,
ζ
(
σ
+
i
t
)
\zeta (\sigma + it)
, with
σ
\sigma
bounded and
t
t
large to moderate accuracy (within
±
t
−
c
\pm {t^{ - c}}
for some
c
>
0
c > 0
, say) was based on the Riemann-Siegel formula and required on the order of
t
1
/
2
{t^{1/2}}
operations for each value that was computed. New algorithms are presented in this paper which enable one to compute any single value of
ζ
(
σ
+
i
t
)
\zeta (\sigma + it)
with
σ
\sigma
fixed and
T
⩽
t
⩽
T
+
T
1
/
2
T \leqslant t \leqslant T + {T^{1/2}}
to within
±
t
−
c
\pm {t^{ - c}}
in
O
(
t
ε
)
O({t^\varepsilon })
operations on numbers of
O
(
log
t
)
O(\log t)
bits for any
ε
>
0
\varepsilon > 0
, for example, provided a precomputation involving
O
(
T
1
/
2
+
ε
)
O({T^{1/2 + \varepsilon }})
operations and
O
(
T
1
/
2
+
ε
)
O({T^{1/2 + \varepsilon }})
bits of storage is carried out beforehand. These algorithms lead to methods for numerically verifying the Riemann hypothesis for the first
n
n
zeros in what is expected to be
O
(
n
1
+
ε
)
O({n^{1 + \varepsilon }})
operations (as opposed to about
n
3
/
2
{n^{3/2}}
operations for the previous method), as well as improved algorithms for the computation of various arithmetic functions, such as
π
(
x
)
\pi (x)
. The new zeta function algorithms use the fast Fourier transform and a new method for the evaluation of certain rational functions. They can also be applied to the evaluation of
L
L
-functions, Epstein zeta functions, and other Dirichlet series.