题目大意 : 给出 (n, k)，求 (sum_{i=1}^{n} i^k) 对 (10^9 +7) 取模的值
(1 leq n leq 10^9, 0 leq k leq 10^6)
解题思路 :
考虑 (k) 比较大，且模数不太好 (NTT)，于是考虑插值做法， 学习了一波插值之后，发现就是个板子题 (雾)
拉格朗日插值：
对于 (n) 个点值表达，可以确定一个 (n - 1) 次的多项式，考虑已知点值怎么构造这个多项式
如果对于第 (i) 个点值，能够构造出这么一个多项式 (f_i(x)) ，将每一个点值表达 (j) 带进去的结果如下
[f_i (x) = egin{cases} y_i , & ext{if $ i= j$} \ 0. \, & ext{otherwise} end{cases} ]那么最终的多项式 (F(x) = sum_{i=1}^{n} f_i(x))
通过贺资料发现当 $f_i(x) = y_i imes prod_{j=1, j eq i}^{n} frac{x-x_j}{x_i-x_j} $ 时即可满足条件，不妨证明一下为什么
考虑分子部分，如果带进去的 (x eq x_i) 那么当 (j = i) 的时候分子是 (0)，结果也是 (0)
考虑分母部分，如果 (x = x_i) 那么分子就等于分母，且分子分母都不存在 (0) ，所以结果就是 (y_i)
大力观察题面发现，(k) 次幂的前缀和一定能用一个 (k + 1) 次的多项式来表示，不妨递推出 (x = 1 .. (k+2)) 时的点值表示直接做插值
那么第 (n) 项的答案就等于 $F(n) = sum y_i imes prod_{j=1, j eq i}^{k+2} frac{n-x_j}{x_i-x_j} $
因为 (x_i = i) ，所以 (F(n) = sum y_i imes prod_{j=1, j eq i}^{k+2} frac{n-j}{i-j})， 后面部分也就是 (frac{(n-1) imes(n-2) imes.. imes(n-k-2)}{(n-i) imes(i-1)! imes(n-i)! imes (-1)^{n-i}})
那么阶乘与处理一下枚举一下 (i) 即可，注意特判一下 (i = n) 的情况不能直接除，应该记一个没有 ((n-i)) 这一项的答案
/*program by mangoyang*/
#include<bits/stdc++.h>
#define inf (0x7f7f7f7f)
#define Max(a, b) ((a) > (b) ? (a) : (b))
#define Min(a, b) ((a) < (b) ? (a) : (b))
typedef long long ll;
using namespace std;
template <class T>
inline void read(T &x){
int f = 0, ch = 0; x = 0;
for(; !isdigit(ch); ch = getchar()) if(ch == '-') f = 1;
for(; isdigit(ch); ch = getchar()) x = x * 10 + ch - 48;
if(f) x = -x;
}
#define int ll
int js[2000005], n, k;
const int Mod = 1000000007;
inline int Pow(int a, int b){
int ans = 1;
for(; b; b >>= 1, a = a * a % Mod)
if(b & 1) ans = ans * a % Mod;
return ans;
}
inline int Inv(int x){ return Pow((x + Mod) % Mod, Mod - 2); }
main(){
int res = 1, res2 = 1; read(n), read(k), js[0] = 1;
for(int i = 1; i <= k + 2; i++){
js[i] = (js[i-1] * i) % Mod;
res = (res * (n - i) + Mod) % Mod;
if(i != n) res2 = (res2 * (n - i) % Mod + Mod) % Mod;
}
int now = 0, ans = 0;
for(int i = 1; i <= k + 2; i++){
now = (now + Pow(i, k)) % Mod;
int A = (i != n) ? res * Inv(n - i) % Mod : res2;
int B = js[i-1] * js[k+2-i] % Mod;
if((k + 2 - i) & 1) B = (-B + Mod) % Mod;
ans = (ans + now * (A * Inv(B) % Mod) % Mod) % Mod;
}
cout << (ans % Mod + Mod) % Mod;
return 0;
}