R, 104 90 bytes
Since R doesn't have arbitrary precision arithmetic, and other answers already have used the limit definition and Taylor series of e, I use a decimal spigot algorithm to spice things up. I don't usually write loop-based code in R, so there is plenty of room for improvement. (divmod
would be useful).
-7 bytes thanks to pajonk's better looping, and another -7 because I remembered R's colon operator can go backwards.
C=rep(1,254)
d=2
for(i in 1:500){cat(d)
d=0
for(j in 254:2)C[j]=(t=10*C[j]+d)-(d=t%/%j)*j}
Attempt This Online!
Based on The Calculation of e to Many Significant Digits by AHJ Sale (1968). The basic idea is to use an \$m\$th order Taylor polynomial, extract an integer part and fractional part, multiply the fractional part by 10, and repeat on the fractional part. Here's an example extracting the first decimal digit using 4 terms: (the 0 coefficients are coincidence)
\begin{alignat}{4}
\newcommand\pp{\phantom{0}}
e
&\approx 2 + {} && \biggl[ &\frac 1 2 \biggl(\pp1 + {}
&&\frac 1 3 \biggl(\pp1 + {} &\frac 1 4 (\pp1) \biggr) \biggr) \biggr]\\
&= 2 + \frac{1}{10} && \biggl[ &\frac 1 2 \biggl( 10 + {}
&&\frac 1 3 \biggl( 10 + {} &\frac 1 4 (10) \biggr) \biggr) \biggr]\\
&= 2 + \frac{1}{10} && \biggl[ &\frac 1 2 \biggl( 10 + {}
&&\frac 1 3 \biggl( 12 + {} &\frac 1 4 (\pp2) \biggr) \biggr) \biggr]\\
&= 2 + \frac{1}{10} && \biggl[ &\frac 1 2 \biggl( 14 + {}
&&\frac 1 3 \biggl(\pp0 + {} &\frac 1 4 (\pp2) \biggr) \biggr) \biggr]\\
&= 2 + \frac{1}{10} && \biggl[7 + {} &\frac 1 2 \biggl(\pp0 + {}
&&\frac 1 3 \biggl(\pp0 + {} &\frac 1 4 (\pp2) \biggr) \biggr) \biggr]
\end{alignat}
Thanks to Steven B. Segletes for helping with MathJax formatting.
I believe the integers involved never exceed \$10m\$. The fractional part is always less than 1, so the computed digit can't affect the previous digit. \$e\$ is special because the coefficients of the Maclaurin series of \$e^x\$ are all 1. Trig functions for \$\pi\$ would have larger coefficients or produce negative digits, which would require amending previous digits. For example, here's Newton's arctan formula for \$\pi\$:
\begin{align}
\frac \pi 2
&\approx 1 + \frac 1 3 \biggl( 1 + \frac 2 5 \biggl( 1 + \frac 3 7 (1) \biggr) \biggr) \\
&= 1 + \frac 1 3 \biggl( 1! + \frac 1 5 \biggl( 2! + \frac 1 7 (3!) \biggr) \biggr)
\end{align}
There is no clean way to transform \$\frac 3 7 (10)\$ into an integer part while the fractional part is a multiple of \$\frac 3 7\$ and less than 1. The solution is modifying the spigot to maintain held predigits that are released once we are sure the digits are correct. See A Spigot Algorithm for the Digits of π by Rabinowitz and Wagon (1995).
\$m\$ can be calculated beforehand from the error of the finite sum using Stirling's approximation. For \$n+1\$ correct digits, find least \$m\$ satisfying \$m! > 10^{n+1}\$, equivalently
$$\frac 1 2 \ln(2\pi m) + m(\ln m - 1) > (n+1) \ln 10$$
Ungolfed algorithm:
n = 500
m = 254
C = rep(1,m)
cat(2)
for(i in 2:n){
d = 0
for(j in m:2){
t = 10*C[j] + d
d = t %/% j
C[j] = t - d*j
}
cat(d)
}