I am sure one has to use indefinite sums and the Euler-Maclaurin formula.
$$\sum_{x=0}^{n} f(x)=\sum_{x}f(n+1)-\sum_{x}f(0)=\int_{0}^{n+1}f(t) \ dt -\frac{1}{2}f(n+1)+\left(\sum_{k=1}^{\infty}\frac{B_{2k}}{(2k)!}f^{(2k-1)}(n+1)\right)+\frac{1}{2}f(0)-\left(\sum_{k=1}^{\infty}\frac{B_{2k}}{(2k)!}f^{(2k-1)}(0)\right)$$
Where $\sum\limits_{x}$ is the anti-difference operator and $B_{2k}$ is the sequence of Bernoulli numbers.
For example with $\tan(x)$I tried the following.
F[n_] := NIntegrate[Tan[x], {x, 0, n + 1}] - 1/2*Tan[n + 1] +
Sum[(BernoulliB[2 k] D[Tan[n + 1], 2 k - 1])/(Factorial[2 k]), {k,
1, 1000}] + 1/2*Tan[0] -
Sum[(BernoulliB[2 k] D[Tan[0], 2 k - 1])/(Factorial[2 k]), {k, 1,
1000}]
But this does not work since the integral does not converge fast enough at the asymptotes.
One person told me approximating such a sum works at an interval from $[-\pi/2,\pi/2]$. In such cases how do we extend the rate of convergence.
Can programming be used efficiently in such a process?