[基本操作]拉格朗日插值
大概就是有 $n$ 对点 $(x_i,y_i)$ 让你构造一个 $n-1$ 次多项式函数过这些点,求这个多项式函数在 $k$ 处的点值
这是一个构造题,我们构造一个函数
$\sum\limits_{i=0}^{n-1} y_i \times \prod\limits_{j=0}^{n-1} \frac{x-x_j}{x_i-x_j}$
然后就没有然后了...
展开求这个多项式函数的每一项是 $O(n^3)$ 的,如果要求点值的话就把 $k$ 代进去就行了
例题比如
bzoj2655 calc
一道很缺德的题
一个序列a1,...,an是合法的,当且仅当:
长度为给定的n。
a1,...,an都是[1,A]中的整数。
a1,...,an互不相等。
一个序列的值定义为它里面所有数的乘积,即a1a2...an。
求所有不同合法序列的值的和。
两个序列不同当且仅当他们任意一位不一样。
输出答案对一个数mod取余的结果。
$n \leq 500,A \leq 1e9,mod是质数$
sol:
首先你可以想到一个 $O(n^2)$ 的 dp:
$f(i,j)$ 表示值域为 $i$,有 $j$ 个数的不同合法序列值的和
$f(0,0) = 1; f(i,j) = f(i-1,j-1)*i*j+f(i-1,j)$
然后发现这是一个 $2n$ 次多项式...(打表找规律
插值即可
#include <bits/stdc++.h> #define int long long #define LL long long #define rep(i, s, t) for (register int i = (s), i##end = (t); i <= i##end; ++i) #define dwn(i, s, t) for (register int i = (s), i##end = (t); i >= i##end; --i) using namespace std; inline int read() { int x = 0, f = 1; char ch; for (ch = getchar(); !isdigit(ch); ch = getchar()) if (ch == '-') f = -f; for (; isdigit(ch); ch = getchar()) x = 10 * x + ch - '0'; return x * f; } int A, n, mod; int f[5100][5100]; inline int ksm(int x, int t) { int res = 1; for(; t; x = 1LL * x * x % mod, t >>= 1) if(t & 1) res = 1LL * x * res % mod; return res; } int nodval[5100]; signed main() { A = read(), n = read(), mod = read(); f[0][0] = 1; rep(i, 1, n+n) { f[i][0] = f[i - 1][0]; rep(j, 1, n) f[i][j] = ((((1LL * f[i-1][j-1] * i) % mod) * (LL)j) % mod + f[i-1][j]) % mod; } if(A <= n+n) cout << f[A][n] << endl; else { int ans = 0; rep(i, 0, 2*n) { int t1 = 1, t2 = 1; rep(j, 0, 2*n) { if(i == j) continue; t1 = 1LL * (A - j) * t1 % mod; t2 = 1LL * (i - j) * t2 % mod; } ans = (ans + (1LL * (f[i][n] * t1) % mod * ksm(t2, mod-2)) % mod) % mod; } cout << (ans%mod + mod)%mod << endl; } }
再比如,自然数幂和
就是 $\sum\limits_{i=1}^n i^k$
然后 $k$ 是 $100000$ 左右,$n$很大, $10^{18}$
sol:
根据观察发现自然数 $k$ 次幂和是 $k+1$ 次多项式(拿一次和二次意会一下)
于是可以插值
这时就要用到一个小技巧, 我们插值如果插 $1,2,...,k$ 这样连续的值,就可以 $O(k)$ 求一个点值
具体的,记 $pref_i = \prod\limits_{j=1}^{i}(n-j)$,$suff_i = \prod\limits_{j=i}^{k}(n-j)$
则 $f_n = \sum\limits_{i=1}^k (-1)^{(k-i)} \times y_i \times \frac{pref_{(i-1)} \times suff_{(i+1)}}{(i-1)! \times (k-i)!}$
然后这题就是令上面式子里的 $k$ 为 $k+2$ 即可
给一个提交地址:51nod1258 序列求和 v4
#include <bits/stdc++.h> #define LL long long #define rep(i, s, t) for (register int i = (s), i##end = (t); i <= i##end; ++i) #define dwn(i, s, t) for (register int i = (s), i##end = (t); i >= i##end; --i) using namespace std; inline LL read() { LL x = 0, f = 1; char ch; for (ch = getchar(); !isdigit(ch); ch = getchar()) if (ch == '-') f = -f; for (; isdigit(ch); ch = getchar()) x = 10 * x + ch - '0'; return x * f; } const int maxn = 500010, mod = 1000000007; LL n; int k, fac[maxn], y[maxn], pre[maxn], suf[maxn]; inline int ksm(int x, int t) { int res = 1; for (; t; x = 1LL * x * x % mod, t >>= 1) if (t & 1) res = 1LL * res * x % mod; return res; } int main() { //freopen("TLE.txt","r",stdin); //freopen("TLEOUT.txt","w",stdout); int T = read(); while(T--) { n = read(), k = read(); int ans = 0; rep(i, 1, k + 2) y[i] = (y[i - 1] + ksm(i, k)) % mod; if (n <= k + 2) cout << y[n] << endl; else { fac[0] = 1; pre[0] = 1; suf[k+2+1] = 1; rep(i, 1, k+2) pre[i] = 1LL * pre[i - 1] * ((n - i) % mod) % mod; dwn(i, k+2, 1) suf[i] = 1LL * suf[i + 1] * ((n - i) % mod) % mod; //rep(i, 1, k+2) pre[i] = 1LL * pre[i - 1] * (n - i) % mod; //dwn(i, k+2, 1) suf[i] = 1LL * suf[i + 1] * (n - i) % mod; rep(i, 1, k+2) fac[i] = 1LL * fac[i - 1] * i % mod; rep(i, 1, k+2) { int opt = ((k+2-i) & 1) ? -1 : 1; int cur = y[i], lx = (1LL * pre[i-1] * suf[i+1]) % mod, fm = ksm((1LL * fac[i-1] * fac[k+2-i]) % mod, mod - 2); cur = 1LL * cur * lx % mod; cur = 1LL * cur * fm % mod; (ans += opt * cur) %= mod; ans = ((ans % mod) + mod) % mod; } /* rep(i, 1, k+2) { int t1 = 1, t2 = 1; rep(j, 1, k+2) { if(i == j) continue; t1 = 1LL * t1 * (n - j) % mod; t2 = 1LL * t2 * (i - j) % mod; } (ans += ((1LL * t1 * y[i]) % mod * ksm(t2, mod - 2)) % mod) %= mod; ans = ((ans % mod) + mod) % mod; }*/ //cout << ans << endl; printf("%d\n",ans); } } }
自然数幂和这个东西还可以扩展
比如 bzoj 3453
给定 k,a,n,d,p
f(i)=1^k+2^k+3^k+......+i^k
g(x)=f(1)+f(2)+f(3)+....+f(x)
求(g(a)+g(a+d)+g(a+2d)+......+g(a+nd))mod p
a,n,d,p 很大,p 是质数,k 不超过 123
sol:
发现 $f_i$ 是 $k+1$ 次多项式,$f_i$ 的前缀和 $g_i$ 就是 $k+3$ 次多项式
下面那个式子,相当于可以看成 $g_i$ 的前缀和(意会),就是 $k+5$ 次多项式
先用拉格朗日插值算出下面那个式子每一项的点值,然后再插一遍
#include <bits/stdc++.h> #define LL long long using namespace std; #define rep(i, s, t) for (register int i = (s), i##end = (t); i <= i##end; ++i) #define dwn(i, s, t) for (register int i = (s), i##end = (t); i >= i##end; --i) inline LL read() { LL x = 0, f = 1; char ch = getchar(); for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -f; for (; isdigit(ch); ch = getchar()) x = 10 * x + ch - '0'; return x * f; } const LL mod = 1234567891; LL T, k, a, n, d; inline LL skr(LL x, LL t) { LL res = 1; for(; t; (x *= x) %= mod, t >>= 1) if(t & 1) res = res * x % mod; return res; } LL cal(LL *f, LL pos, LL len) { if(pos <= len) return f[pos]; LL ans = 0; rep(i, 1, len) { LL t1 = 1, t2 = 1; rep(j, 1, len) { if(j == i) continue; t1 = t1 * (pos - j) % mod; t2 = t2 * (i - j) % mod; } (ans += (((f[i] * t1) % mod * skr(t2, mod - 2)) % mod)) %= mod; } ans = ((ans % mod) + mod) % mod; return ans; } LL f[200], g[200]; int main() { T = read(); while(T--) { k = read(), a = read(), n = read(), d = read(); rep(i, 1, k+3) f[i] = skr(i, k); rep(i, 2, k+3) (f[i] += f[i-1]) %= mod; rep(i, 2, k+3) (f[i] += f[i-1]) %= mod; g[0] = cal(f, a, k+3); rep(i, 1, k+5) g[i] = (cal(f, (a + i * d) % mod, k+3) + g[i-1]) % mod; cout << cal(g, n, k+5) << endl; } }