一道利用 Lagrange 插值法解决的多项式问题。
题源:Codeforces Edu Round 7 F。
题解
不妨记 $p(x)=\sum\limits_{i=1}^x i^k$ ,不难发现,这是一个关于 $x$ 的 $k+1$ 次多项式。
(自然数幂次和公式如下,但发现上述性质只需仔细观察题目。)
$$
S_p(n) = 1^p + 2^p + \cdots+n^p= \frac{1}{p+1} \sum_{k=0}^{p} \binom{p+1}{k} B_k n^{p+1-k}
$$
我们希望能够求出 $p(n)$ 。
考虑 Lagrange 插值法,我们只需要确定 $k+2$ 个点,就能完全确定这个 $k+1$ 次多项式。
我们可以选取 $1,2,\cdots,k+2$ 这连续的 $k+2$ 个点,求出这些点的坐标所需时间复杂度为 $O(k\log k)$。
对于连续点的 Lagrange 插值,可以做如下优化:
$$
p(x)=\sum\limits_{i=1}^{k+2}y_i\prod\limits_{j=1,j\ne i}^{k+2}\frac{x-x_j}{x_i-x_j}=\sum\limits_{i=1}^{k+2}y_i\prod\limits_{j=1,j\ne i}^{k+2}\frac{x-j}{i-j}
$$
连乘项中,分子维护住前缀积 $pre[i]=\prod\limits_{j=1}^{i}(x-j)$,后缀积 $suf[i]=\prod\limits_{j=i}^{k+2}(x-j)$。
分母维护住阶乘前缀积 $pref[i]=\prod\limits_{j=1}^ij$,这样,我们的连乘项可以化成:
$$
\prod\limits_{j=1,j\ne i}^{k+2}\frac{x-j}{i-j}=\frac{pre[i-1]\times suf[i+1]}{pref[i-1]\times (-1)^{k+2-i}\cdot pref[k+2-i]}
$$
这样,就能在 $O(k)$ 时间完成 $p(n)$ 的计算。
总复杂度 $O(k\log k)$,瓶颈在于获取前 $k+2$ 项的初值。代码实现:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83
| #include<bits/stdc++.h> #define int long long using namespace std;
constexpr int mod = 1e9 + 7;
int qmi(int a, int k, int p){ int res = 1; while(k){ if(k & 1) res = res * a % p; k >>= 1; a = a * a % p; } return res; }
int inv(int x){ return qmi(x, mod - 2, mod); }
struct modint{ int v; modint(int x = 0) : v((x % mod + mod) % mod) {}
modint operator+ (const modint &o){ return modint(v + o.v); } modint operator- (const modint &o){ return modint(v - o.v); } modint operator* (const modint &o){ return modint(v * o.v); } modint operator/ (const modint &o){ return modint(v * inv(o.v)); } }; using Z = modint;
void solve(){ int n, k; cin >> n >> k; vector<Z>p(k + 3); p[1] = 1; for(int i = 2; i <= k + 2; i++){ p[i] = p[i - 1] + Z(qmi(i, k, mod)); }
if(n <= k + 2){ cout << p[n].v << '\n'; return; }
vector<Z>pre(k + 3), suf(k + 4); pre[0] = 1; for(int i = 1; i <= k + 2; i++){ pre[i] = pre[i - 1] * Z(n - i); } suf[k + 3] = 1; for(int i = k + 2; i >= 1; i--){ suf[i] = suf[i + 1] * Z(n - i); }
vector<Z>pref(k + 3); pref[0] = 1; for(int i = 1; i <= k + 2; i++){ pref[i] = pref[i - 1] * Z(i); }
Z res(0); for(int i = 1; i <= k + 2; i++){ Z mo = pre[i - 1] * suf[i + 1]; Z de = pref[i - 1] * pref[k + 2 - i] * Z((k + 2 - i) % 2 ? -1 : 1); res = res + p[i] * mo / de; } cout << res.v << '\n'; }
signed main(){ ios::sync_with_stdio(0); cin.tie(0); solve(); return 0; }
|