The Riemann theta functions are multivariate complex-analytic functions with deep links to the theory of abelian varieties over the complex numbers. For applications in algorithmic number theory, one would like to efficiently evaluate them at a given point (a period matrix) to very large precision with certified error bounds. In this talk, I will present a new evaluation algorithm for theta functions in uniform quasi-linear time in the required precision, in any fixed dimension (joint work with Noam D. Elkies). I will also discuss its implementation in the FLINT number theory library.