The approximation of integral functionals with respect to a stationary Markov process by a Riemann-sum estimator is studied. Stationarity and the functional calculus of the infinitesimal generator of the process are used to get a better understanding of the estimation error and to prove a general error bound. The presented approach admits general integrands and gives a unifying explanation for different rates obtained in the literature. Several examples demonstrate how the general bound can be related to well-known function spaces.