题目概述

求 $\sum_{i=1}^{n}i^K$ 。

解题报告

这个东西可以用二项式定理和组合数 $O(n^2)$ 递推来着，但是 $n$ 太tm大了，不过 $K$ 很小，所以要想办法搞成只和 $K$ 有关的式子，需要用到伯努利数，大佬传送门。

伯努利数的递推式：$B_0=1,B_n=-{1\over n+1}\sum_{i=0}^{n-1}{n+1\choose i}B_i$ ，变形一下可以得出 $\sum_{i=0}^{t}{t+1\choose i}B_i=0,t>0$ 。

一个数列 $\{f_n\}$ 的指数生成函数为 $\hat{F}(x)=\sum_{i\ge 0} f_i{x^i\over i!}$ ，那么两个指数生成函数的卷积 $\hat{H}(x)=\hat{F}(x)*\hat{G}(x)$ 对应的系数就是：

$$ h_n=\sum_{i+j=n}{n\choose i}f_ig_j $$

这和 $\sum_{i=0}^{t}{t+1\choose i}B_i=0,t>0$ 怎么这么像啊，我们两边补上 $B_{t+1}$ 得到 $\sum_{i=0}^{t+1}{t+1\choose i}B_i=B_{t+1},t>0$ ，把 $t+1$ 换成 $n$ 那就是 $\sum_{i=0}^{n}{n\choose i}B_i=B_{n},n>1$ ，所以左边是 $\hat{B}(x)*e^x$（ $e^x$ 表示全是 $1$ 的数列的指数生成函数）。

右边是 $\hat{B}(x)$ ，但是由于 $n>1$ 所以左边第二项的系数( $1\over 2$ )不等于右边第二项的系数( $-{1\over 2}$ )，因此要补上一个 $1$ ：

$$ \hat{B}(x)*e^x=\hat{B}(x)+x $$

移项得到 $\hat{B}(x)={x\over e^x-1}$ （因此可以多项式求逆得到伯努利数），这时候来看如何求出 $\sum_{i=0}^{n}i^K$ 。

令 $e^{kx}=\sum_{i\ge 0}k^i{x^i\over i!}$ ，不难得出递推式 $e^{kx}=e^{(k-1)x}*e^x=(e^x)^k$ ，那么只要把 $e^{kx},k\in[0,n-1]$ 全部加起来，第 $K+1$ 项就是 $\sum_{i=0}^{n-1}i^K$ ，又因为等比数列求和（显然多项式也可以等比数列求和，感谢万L群大佬的指导）：

$$ \sum_{k=0}^{n-1}e^{kx}={e^{nx}-1\over e^x-1} $$

记 $\hat{S}_K(x)=\sum_{k=0}^{n-1}e^{kx}$ ，又根据 $\hat{B}(x)={x\over e^x-1}$ ，得到：

$$ \hat{S}_K(x)=\hat{B}(x){e^{nx}-1\over x} $$

因为 $e^{nx}-1\over x$ 其实就是 $e^{nx}$ 去掉常数项之后往左移动了 $1$ 位所以只需要当成常数项为 $0$ 的 $e^{nx}$ 然后改变下下标就行了。

$$ {S_K(x)_K\over K!}=\sum_{i+j=K+1,j>0}{B_i\over i!}{n^j\over j!}\\ S_K(x)_K={1\over K+1}\sum_{i+j=K+1,j>0}{K+1\choose i}B_in^j\\ S_K(x)_K={1\over K+1}\sum_{i=0}^{K}{K+1\choose i}B_in^{K+1-i} $$

这样就可以预处理+每次 $O(K)$ 求出 $\sum_{i=0}^{n-1}i^K$ 啦！由于这道题 $K$ 不大所以可以直接递推伯努利数。

示例程序

#include<cstdio>
using namespace std;
typedef long long LL;
const int maxk=2000,MOD=1e9+7;

int te,K,fac[maxk+5],INV[maxk+5],pre[maxk+5],Be[maxk+5];LL n;

#define C(x,y) ((x)<(y)?0:((LL)fac[x]*pre[y]%MOD*pre[(x)-(y)]%MOD))
inline void AMOD(int &x,int tem) {if ((x+=tem)>=MOD) x-=MOD;}
inline void Make(int n){
    pre[0]=INV[1]=1;for (int i=2;i<=n;i++) INV[i]=(LL)(MOD-MOD/i)*INV[MOD%i]%MOD;
    fac[0]=1;for (int i=1;i<=n;i++) fac[i]=(LL)fac[i-1]*i%MOD,pre[i]=(LL)pre[i-1]*INV[i]%MOD;
    for (int i=Be[0]=1;i<=n;Be[i]=(LL)(MOD-Be[i])*INV[i+1]%MOD,i++)
        for (int j=0;j<i;j++) AMOD(Be[i],(LL)C(i+1,j)*Be[j]%MOD);
}
inline int Pow(int w,int b) {int s=1;for (;b;b>>=1,w=(LL)w*w%MOD) if (b&1) s=(LL)s*w%MOD;return s;}
int main(){
    freopen("program.in","r",stdin);freopen("program.out","w",stdout);
    for (Make(maxk+1),scanf("%d",&te);te;te--){
        scanf("%lld%d",&n,&K);n++;int ans=0;
        for (int i=0;i<=K;i++) AMOD(ans,(LL)C(K+1,i)*Be[i]%MOD*Pow(n%MOD,K+1-i)%MOD);
        ans=(LL)ans*INV[K+1]%MOD;printf("%d\n",ans);
    }return 0;
}