高效求解自然数k次幂和的方法

前言

寒假集训比赛时遇到了一道《WYF的盒子》，题目大意就是求$\sum_{i=m}^n i^k$对$p$取模的结果。
其中有的数据点满足$n-m\le 5000$但$k\le{10}^{12}$，剩下的数据点满足$n,m\le10^{12}$且$k\le 2000$。对于所有数据模数$p\le10^{12}$。
其中$n$和$m$相差很小的点显然可以使用快速幂，但是剩下的点怎么做？换言之，我们需要找到一个时间复杂度之和$k$有关的高效的求自然数幂和的算法。

---

普通公式

运用二项式定理我们可以推出一般的公式。

$$S_k=\frac{1}{k+1}[(n+1)^{k+1}-(n+1)-\sum_{i=1}^{k-1}C_{k+1}^{k-i}S_i]$$
推算过程链接： http://wenku.baidu.com/link?url=F07Cj1SxmyNr8_EriB3UrJU8tn4f3aByvpYKPvrLCoqV2HWjeXG1y6XWbgFQL8mBcP5Z0is1Sqpe0zxbDRmMSsBjt-0mQRRvvHTMqq57gpK
当然我们也可以不使用二项式定理，只研究递推关系，也可以得到相同的公式。
推算过程链接： http://blog.csdn.net/acdreamers/article/details/38929067
这种公式固然优美，时间复杂度为 $\mathrm O(k^2)$，但是不可避免地要做除法，如果除数和模数互质那好办，直接逆元即可。但是这里显然不满足。

---

拉格朗日插值法&牛顿插值法

暑假时HeD大神给我们讲了这两种解法，初中党表示一脸懵逼。
拉格朗日插值法其实不难，拿高中选修数论初步看看就懂了。但是还是避免不了除法。
牛顿插值法只需要差分，不需要做除，时间复杂度$\mathrm O(k^2)$，是很优秀的算法，然而我不会。