给定 (n,m,k),计算 (f(n,m)) 的值,模 (10^9+7).
(0 leq N leq 10^{18},0 leq M leq 30,1 leq K leq 2.5 imes 10^6).
易发现答案为 (1^k,2^k,3^k,cdots ,n^k) 的 (m) 阶前缀和的第 (n) 项值。
考虑 (i^k) 对答案的贡献次数,是对 (f_0=1,f_k=0(k>0)) 的数组作 (m) 阶前缀和后的第 ((n-i)) 项,也就是 (f_{n-i})。其值为 (inom{n-i+m-1}{m-1}).
从组合意义上考虑,作一阶前缀和为 (1,1,1,cdots),设后面第 (i) 阶前缀和的 (f_j=g_{i,j}),有递推式 (g_{i,j}=g_{i-1,j}+g_{i,j-1}),恰为只能向右、下走的格点计数。
故贡献的次数为 ((1,0) o (m,n-i)) 的路径数,即为 (inom{n-i+m-1}{m-1})。
所以答案为:
(i^k) 为关于 (i) 的 (k) 次单项式。
(inom{n-i+m-1}{m-1}=frac{(n-i+m-1)!}{(m-1)!(n-i)!}),其中 ((n-i+m-1)!) 为关于 (i) 的 ((n-i+m-1)) 次多项式,((n-i)!) 为关于 (i) 的 ((m-1)) 次多项式,((m-1)!) 是常数。
故 (i^kinom{n-i+m-1}{m-1}) 是关于 (i) 的 ((m+k-1)) 次多项式。
对于 (iin [0,n]),对这个多项式求和即为答案。
其可以写成 (ans=a_0s_0+a_1s_1+a_2s_2+cdots+a_{m+k-1}s_{m+k-1}) 的形式,其中 (s_i=sumlimits_{j=0}^ij^k).
由于 (s_i) 是关于 (i) 的 ((i+1)) 次多项式,故答案为关于 (i) 的 ((m+k)) 次多项式。
设 (l=m+k+1),选出前 (l) 个连续的整数,算出其作 (m) 阶前缀和的答案,拉格朗日插值即可。
时间复杂度 (mathcal{O}(mk)).
#include<iostream>
#include<cstdio>
#include<algorithm>
typedef long long ll;
const ll mod = 1000000007;
template <typename T> T Max(T x, T y) { return x > y ? x : y; }
template <typename T> T Min(T x, T y) { return x < y ? x : y; }
template <typename T> T Add(T x, T y) { return (x + y >= mod) ? (x + y - mod) : (x + y); }
template <typename T> T Mod(T x) { return (x >= mod) ? (x - mod) : (x < 0 ? (x + mod) : x); }
template <typename T> T Mul(T x, T y) { return x * y % mod; }
template <typename T>
T &read(T &r) {
r = 0; bool w = 0; char ch = getchar();
while(ch < '0' || ch > '9') w = ch == '-' ? 1 : 0, ch = getchar();
while(ch >= '0' && ch <= '9') r = (r << 3) + (r <<1) + (ch ^ 48), ch = getchar();
return r = w ? -r : r;
}
ll qpow(ll x, ll y) { ll sumq = 1; while(y) { if(y & 1) sumq = sumq * x % mod; x = x * x % mod; y >>= 1; } return sumq; }
const int N = 2600100;
ll n;
int m, k;
ll a[N], fac[N], inv[N], pre[N], suf[N], ans;
signed main() {
read(n); read(m); read(k); int l = m+k+1;
if(!n) { puts("0"); return 0; }
for(int i = 1; i <= l; ++i) a[i] = qpow(i, k);
for(int j = 1; j <= m; ++j)
for(int i = 1; i <= l; ++i)
a[i] = Add(a[i], a[i-1]);
inv[0] = fac[0] = 1; for(int i = 1; i <= l; ++i) fac[i] = fac[i-1] * i % mod;
inv[l] = qpow(fac[l], mod-2); n %= mod;
for(int i = l-1; i; --i) inv[i] = inv[i+1] * (i+1) % mod;
pre[0] = 1; for(int i = 1; i <= l; ++i) pre[i] = pre[i-1] * (n - i) % mod;
suf[l+1] = 1; for(int i = l; i; --i) suf[i] = suf[i+1] * (n - i) % mod;
for(int i = 1; i <= l; ++i)
ans = Add(ans, Mod(a[i] * pre[i-1] % mod * suf[i+1] % mod * inv[i-1] % mod * inv[l-i] % mod * (((l-i)&1) ? -1 : 1)));
printf("%lld
", ans);
return 0;
}