The Gamma function (over $\mathbb{R}$) allows us to express the factorial function as an Euler integral: $$ n! = \int_0^\infty x^n \, e^{-x} \, dx $$ What happens if we use the trapezoid rule on the improper integral on the right hand side? For general functions it says: $$ \int_a^b f(x) \, dx \approx (b-a) \left[ \frac{f(a)+ f(b)}{2} \right] $$ We can partition $\mathbb{R}_{\geq 0} $ as the union of segments $[n,n+1]$ for over all possible integers $n \geq 0$. Therefore: $$ a! \approx \sum_{m \geq 0} 1 \cdot \frac{m^a \, e^{-m}+ m^{a+1}e^{-(m+1)}}{2} = - \frac{1}{2} + \sum_{m \geq 0} m^a e^{-m}$$
How bad is the trapezoid rule in this kind of estimate? The error term is usually:
$$ \left|\; \int_a^b f(x) \, dx \;- \;(b-a) \left[ \frac{f(a)+ f(b)}{2} \right] \; \right| \leq \frac{(b-a)^3}{12}f''(\xi) $$
So for each term we are losing $\frac{1}{12} \times \text{constant}$. Can we understand this error term better?
Related: An Elementary Proof of Error Estimates for the Trapezoidal Rule D. Cruz-Uribe and C. J. Neugebauer