解题报告（三）多项式求值与插值（拉格朗日插值）（ACM / OI）
整理的算法模板合集： ACM模板
点我看算法全家桶系列！！！
实际上是一个全新的精炼模板整合计划
繁凡出品的全新系列：解题报告系列 —— 超高质量算法题单，配套我写的超高质量的题解和代码，题目难度不一定按照题号排序，我会在每道题后面加上题目难度指数（1∼51 \sim 51∼5），以模板题难度 111 为基准。
这样大家在学习算法的时候就可以执行这样的流程：
%
阅读【学习笔记】 / 【算法全家桶】学习算法 ⇒\Rightarrow⇒ 阅读相应算法的【解题报告】获得高质量题单 ⇒\Rightarrow⇒ 根据一句话题解的提示尝试自己解决问题 ⇒\Rightarrow⇒ 点开详细题解链接学习巩固（好耶）
%
要是26个英文字母用完了我就接上24个希腊字母，我就不信50道题不够我刷的hhh%
解题报告系列合集：【解题报告系列】超高质量题单 + 题解（ICPC / CCPC / NOIP / NOI / CF / AT / NC / P / BZOJ）
本题单前置知识：【学习笔记】多项式全家桶（包含全套证明）
拉格朗日插值法
有拉格朗日插值公式：
f(x)=∑i=1nyi∏j≠ix−xjxi−xjf(x)=\sum_{i=1}^ny_i\prod_{j\neq i}\dfrac {x-x_j} {x_i-x_j}f(x)=i=1∑nyij=i∏xi−xjx−xj
A. P4781 【模板】拉格朗日插值
Weblink
https://www.luogu.com.cn/problem/P4781
Problem
我们将 nnn个点带入拉格朗日插值公式计算 f(k)f(k)f(k) 即可。时间复杂度O(n2)O(n^2)O(n2)。
注意本题还要求逆元，为了防止求逆元的时间复杂度影响整体的时间复杂度，所以我们分别计算出分子和分母，再将分子乘进分母的逆元，累加进最后的答案，时间复杂度的瓶颈就不会在求逆元上，总体的时间复杂度为O(n2)O(n^2)O(n2)。
Code
本题需要累乘，会爆int，记得开long long
#include <cstdio>
#include <iostream>
#include <cstring>
#include <algorithm>using namespace std;const int N = 500007;
typedef long long ll;
const int mod = 998244353;
ll n, m, k;struct Point
{ll x, y;
}A[N];ll qpow(ll a, ll b, ll c)
{ll res = 1;while(b) {if(b & 1) res = res * a % c;a = a * a % c;b >>= 1;}return res;
}ll inv(ll x) {return qpow(x, mod - 2, mod);}int main()
{scanf("%lld%lld", &n, &k);for(int i = 1; i <= n; ++ i) {scanf("%lld%lld", &A[i].x, &A[i].y);}ll ans = 0;for(int i = 1; i <= n; ++ i) {ll s1 = A[i].y % mod;ll s2 = 1ll;for(int j = 1; j <= n; ++ j) {if(i != j) {s1 = s1 * (k - A[j].x) % mod;s2 = s2 * (A[i].x - A[j].x) % mod;}}ans += s1 * inv(s2) % mod;}printf("%lld\n", (ans % mod + mod) % mod);return 0;
}
我们使用拉格朗日插值公式，对于一个数 kkk ，我们很容易在 O(n2)O(n^2)O(n2) 时间求得 F(k)F(k)F(k) 的数值，如果 xix_ixi 是连续的，我们甚至可以利用预处理在 O(n)O(n)O(n) 时间内得到 F(k)F(k)F(k) 的数值。
但是如果 xix_ixi 不连续，又有多组查询，就需要得到这个多项式的系数以保证求一个函数值的时间为 O(n)O(n)O(n) 。
重心拉格朗日插值法
考虑对拉格朗日插值公式进行优化：
f(x)=∑i=1nyi∏j≠ix−xjxi−xjf(x)=\sum_{i=1}^n y_i\prod_{j\neq i}\dfrac {x-x_j} {x_i-x_j}f(x)=i=1∑nyij=i∏xi−xjx−xj
设 h=∏i=1nx−xih=\prod_{i=1}^n x-x_ih=i=1∏nx−xi
带入得：
=h∑i=1n∏j≠iyi(xi−xj)(x−xi)=h\sum_{i=1}^n \prod_{j\neq i}\dfrac {y_i} {(x_i-x_j)(x-x_i)}=hi=1∑nj=i∏(xi−xj)(x−xi)yi
设
ti=∏j≠iyixi−xjt_i=\prod_{j\neq i} \dfrac {y_i} {x_i-x_j}ti=j=i∏xi−xjyi
带入得：
=h∑i=1ntix−xi=h\sum_{i=1}^n \dfrac {t_i} {x-x_i}=hi=1∑nx−xiti
在每次加入一个新点时，计算出它的 tit_iti ，并且更新别的点的 tit_iti ，时间复杂度 O(n)O(n)O(n)，加入 nnn 个点就是 n2n^2n2 。
这样我们就可以 O(n)O(n)O(n) 求 hhh ，然后再 O(n)O(n)O(n) 求出 ∑\sum∑，总时间复杂度为 O(n)O(n)O(n)。
B. 拉格朗日插值法求系数
P4781 【模板】拉格朗日插值
Solution
首先观察式子，发现对于每个 iii , 上面部分是 (x−x1)×(x−x2)…×(x−xn)(x−xi)\cfrac{(x-x_1) \times (x-x_2)… \times (x-x_n)}{(x-x_i)}(x−xi)(x−x1)×(x−x2)…×(x−xn), 下面那部分和 yiy_iyi 都是常数。
所以可以 O(n2)O(n^2)O(n2) 处理出(x−x1)×(x−x2)…×(x−xn)(x-x_1) \times (x-x_2)… \times (x-x_n)(x−x1)×(x−x2)…×(x−xn) 这个 nnn 次多项式，然后通过模拟长除法，O(n)O(n)O(n)时间内可以得到(x−x1)×(x−x2)…×(x−xn)(x−xi)\cfrac{(x-x_1) \times (x-x_2)… \times (x-x_n)}{(x-x_i)}(x−xi)(x−x1)×(x−x2)…×(x−xn)，然后常系数直接 O(n)O(n)O(n) 暴力算出来，就得到了 xix_ixi 对应的多项式，最后把所有多项式加起来就得到了最终的系数。
Code
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 5007;ll mod = 998244353;ll qpow(ll a, ll b) {ll res = 1;while (b) {if (b & 1)res = res * a % mod;a = a * a % mod, b >>= 1;}return res;
}
ll a[N], b[N], c[N], temp[N];
ll x[N], y[N];
int n;
void mul(ll *f, int len, ll t) { //len为多项式的次数+1，函数让多项式f变成f*(x+t)for (int i = len; i > 0; --i)temp[i] = f[i], f[i] = f[i - 1];temp[0] = f[0], f[0] = 0;for (int i = 0; i <= len; ++i)f[i] = (f[i] + t * temp[i]) % mod;
}
void dev(ll *f, ll *r, ll t) { //f是被除多项式的系数，r保存f除以x+t的结果for (int i = 0; i <= n; ++i)temp[i] = f[i];for (int i = n; i > 0; --i) {r[i - 1] = temp[i];temp[i - 1] = (temp[i - 1] - t * temp[i]) % mod;}return;
}
void lglr() {memset(a, 0, sizeof a);b[1] = 1, b[0] = -x[1];for (int i = 2; i <= n; ++i) {mul(b, i, -x[i]);}//预处理(x-x1)*(x-x2)...*(x-xn)for (int i = 1; i <= n; ++i) {ll fz = 1;for (int j = 1; j <= n; ++j) {if (j == i)continue;fz = fz * (x[i] - x[j]) % mod;}fz = qpow(fz, mod - 2);fz = fz * y[i] % mod; //得到多项式系数dev(b, c, -x[i]);//得到多项式，保存在b数组for (int j = 0; j < n; ++j)a[j] = (a[j] + fz * c[j]) % mod;}
}
int main() {ll k;cin >> n >> k;for (int i = 1; i <= n; ++i)scanf("%lld%lld", &x[i], &y[i]);lglr();ll ans = 0;ll res = 1;for (int i = 0; i < n; ++i) {ans = (ans + res * a[i]) % mod;res = res * k % mod;}ans = (ans + mod) % mod;cout << ans << endl;
}
C. 2019ICPC 南昌邀请赛 B. Polynomial
Problem
n≤1000,m≤2000,1≤L≤R≤9999990n\le 1000, m\le 2000, 1\le L\le R\le 9999990n≤1000,m≤2000,1≤L≤R≤9999990
Solution
给定了 fif_ifi，i∈0∼ni\in 0\sim ni∈0∼n， nnn 次多项式给定 n+1n+1n+1 个值，并且还都是连续的，显然可以使用拉格朗日插值 O(n)O(n)O(n) 计算多项式 f(x)f(x)f(x)。题目中有 mmm 次询问，每次询问的是一个区间的取值，显然可以用前缀和来维护。
设 S(x)S(x)S(x) 是 fif_ifi 的前缀和，则答案为 S(R)−S(L−1)S(R)-S(L-1)S(R)−S(L−1)。
而我们需要预处理出 S(1)∼S(9999990)S(1)\sim S(9999990)S(1)∼S(9999990)，直接插值的话时间复杂度 O(n×9999990)O(n\times9999990)O(n×9999990)。
我们知道 nnn 次多项式的前缀和是 n+1n+1n+1 次的多项式，也就意味着 S(x)S(x)S(x) 需要用 n+2n+2n+2 个点来通过拉格朗日插值求解。但是我们只有 n+1n+1n+1 个点，我们可以利用拉格朗日插值法求出 f(n+1)f(n+1)f(n+1)，这样就有了 n+2n+2n+2 个点。我们只需要对 S(x)S(x)S(x) 进行插值即可。
Time
O(T×m×n)O(T\times m\times n)O(T×m×n)
Code
#include <bits/stdc++.h>using namespace std;
const int N = 5007, mod = 9999991;#define int long longint n, m, t;
int inv[mod + 7], infact[mod + 7];
int sum[N], a[N];void init(int n)
{inv[1] = 1;for(int i = 2; i <= n; ++ i)inv[i] = (mod - mod / i) * inv[mod % i] % mod;infact[0] = 1;for(int i = 1; i <= n; ++ i) infact[i] = infact[i - 1] * inv[i] % mod;
}int lagrange(int x, int *a, int n)
{int res = 0;int p = 1;for(int i = 0; i <= n; ++ i)p = p * (x - i) % mod;for(int i = 0; i <= n; ++ i) {int f = (n - i) & 1 ? -1 : 1;res = (res + mod + a[i] * f * p % mod * inv[x - i] % mod * infact[i] % mod * infact[n - i] % mod) % mod;}return res;
}signed main()
{init(mod);scanf("%lld", &t);while(t -- ) {scanf("%lld%lld", &n, &m);for(int i = 0; i <= n; ++ i) {scanf("%lld", &a[i]);a[i] %= mod;}a[n + 1] = lagrange(n + 1, a, n); sum[0] = a[0];for(int i = 1; i <= n + 1; ++ i) sum[i] = (sum[i - 1] + a[i]) % mod;while(m -- ) {int l, r;scanf("%lld%lld", &l, &r); if(r <= n + 1) {printf("%lld\n", (sum[r] - sum[l - 1] + mod) % mod);}else if(l - 1 <= n + 1) {printf("%lld\n", (lagrange(r, sum, n + 1) - sum[l - 1] + mod) % mod);}else printf("%lld\n", (lagrange(r, sum, n + 1) - lagrange(l - 1, sum, n + 1) + mod) % mod);}}return 0;
}
D. CF622F The Sum of the k-th Powers（自然数k次幂的和）
Problem
求 ∑i=1nik\sum_{i=1}^n i^k∑i=1nik ，n≤109,k≤106n \leq 10^9,k \leq 10^6n≤109,k≤106
Solution
n=k+2n=k+2n=k+2
我们知道对于一个 nnn 次多项式 f(x)f(x)f(x) ，若有一个数列 f(1),f(2),f(3),...f(n)f(1),f(2),f(3),...f(n)f(1),f(2),f(3),...f(n)，多项式差分为 f(2)−f(1),f(3)−f(2),...f(2)-f(1),f(3)-f(2),...f(2)−f(1),f(3)−f(2),...，显然有结论：多项式差分是关于 iii 的 n−1n-1n−1 次多项式。
题中所给的数列自然数k次幂的和，即 f(i)=∑j=1ijkf(i)=\sum_{j=1}^i j^kf(i)=∑j=1ijk ，进行多项式差分后得到：Δf(i)=(i+1)k\Delta f(i)=(i+1)^kΔf(i)=(i+1)k ，是一个关于 iii 的 kkk 次多项式，可以得出结论：原多项式 f(n)f(n)f(n) 是一个关于 nnn 的 k+1k+1k+1 次多项式。
即：自然数k次幂的和是一个 k+1k+1k+1 次多项式。
所以我们只需要预处理出 k+2k+2k+2 个点的值 f(i)f(i)f(i) ，就可以使用拉格朗日插值法求出 f(n)f(n)f(n)。显然我们可以预处理出 0∼k+20\sim k+20∼k+2 的 f(i)f(i)f(i) 的值，这样得到的是 xxx 取值连续的序列，我们就可以是哟个拉格朗日插值 O(n)O(n)O(n) 计算出 f(n)f(n)f(n) 的值，即为答案。
Time
O(k)O(k)O(k)
Code
#include <bits/stdc++.h>
using namespace std;const int N = 1100007, mod = 1e9 + 7;int read(){int s = 0, ne = 1; char c = getchar();while(c < '0' || c > '9') {if(c == '-') ne = -1; c = getchar();}while(c >= '0' && c <= '9') s = (s << 1) + (s << 3) + c - '0', c = getchar();return s * ne;
}int qpow(int a, int b)
{int res = 1;while(b) {if(b & 1) res = 1ll * res * a % mod;a = 1ll * a * a % mod; b >>= 1;} return res;
}int n, k;
int s[N], pre[N], suf[N], infact[N], fact[N], ans;void init(int k)
{infact[0] = fact[0] = 1;for(int i = 1; i <= k + 10; ++ i)fact[i] = 1ll * fact[i - 1] * i % mod;infact[k + 10] = qpow(fact[k + 10], mod - 2);for(int i = k + 9; i >= 0; -- i) infact[i] = 1ll * infact[i + 1] * (i + 1) % mod;infact[0] = 1;}// i = 1 to n, ∑ i^kvoid lagrange(int n, int k)
{ s[0] = 0;for(int i = 1; i <= k + 2; ++ i)s[i] = (s[i - 1] + qpow(i, k)) % mod;if(n <= k + 2) {printf("%d", s[n]);return ;}pre[0] = 1;for(int i = 1; i <= k + 2; ++ i)pre[i] = 1ll * pre[i - 1] * ((n - i + mod) % mod) % mod;suf[k + 3] = 1;for(int i = k + 2; i; -- i)suf[i] = 1ll * suf[i + 1] * ((n - i + mod) % mod) % mod;for(int i = 1; i <= k + 2; ++ i) {s[i] = 1ll * s[i] * pre[i - 1] % mod * suf[i + 1] % mod * infact[i - 1] % mod * infact[k + 2 - i] % mod;if((k + 2 - i) & 1) ans = (1ll * ans - s[i] + mod) % mod;else ans = (1ll * ans + s[i]) % mod;}printf("%d\n", ans);
}int main()
{scanf("%d%d", &n, &k);init(k);lagrange(n, k);return 0;
}
E. 牛客挑战赛36 D. 排名估算（ “概率论全家桶”，好题，拉格朗日插值求自然数 k 次幂之和）
Weblink
https://ac.nowcoder.com/acm/contest/3782/D
Problem
期末考后，小 C 登录了学校的成绩查询系统，却发现自己的排名被屏蔽了。为了知道自己的排名，小 C 使用了系统中的“好友伴学”功能。每次，系统会在除了小 C 之外的所有考生中随机抽取一名，然后返回 Ta 的排名比小 C 高还是低。这次考试有 nnn 个人参加，小 C 总共使用的 mmm 次 “好友伴学” 功能，却没有一次抽中排名比自己高的人。请问小C在这次考试中的期望排名是多少?
2≤n≤1011,0≤m≤50002≤n≤10^{11},0≤m≤50002≤n≤1011,0≤m≤5000
假设小C的的期望排名化为最简分数后是 p/qp/qp/q，输出 p∗q−1mod998244353p*q^{-1}\mod 998244353p∗q−1mod998244353。
Solution
设 pip_ipi 表示排名为第 n−in -in−i 时，mmm 次没抽中排名比自己高的人（排名高指排名在自己前面hhh）的概率。
显然排名只有 1∼n1\sim n1∼n，则对于 0≤i≤n−10\le i\le n-10≤i≤n−1，pi=(in−1)mp_{i}=(\cfrac{i}{n-1})^{m}pi=(n−1i)m（排名 n−in-in−i，则比自己排名高的一共有 n−in-in−i 个人，没有抽中显然是抽到了自己后面的 iii 排名低的人）
令事件 AAA 为 mmm 次没抽中，事件 BiB_iBi 为排名为 n−in-in−i 。
根据乘法公式：P(ABi)=P(A∣Bi)×P(Bi)=P(Bi∣A)×P(A)P(AB_i)=P(A|B_i)\times P(B_i)=P(B_i|A)\times P(A)P(ABi)=P(A∣Bi)×P(Bi)=P(Bi∣A)×P(A)，
我们首先计算 P(A)P(A)P(A)，根据全概率公式：
P(A)=∑i=0n−1P(Bi)×P(A∣Bi)P(A)=\sum\limits_{i=0}^{n-1}P(B_i)\times P(A|B_i)P(A)=i=0∑n−1P(Bi)×P(A∣Bi)
即排名为 iii 的前提下事件AAA ，mmm 次没抽中的概率之和。
显然
P(Bi)=1n，P(A∣Bi)=pi=(in−1)mP(B_i)=\frac{1}{n}，P(A|B_i)=p_i=(\frac{i}{n-1})^{m}P(Bi)=n1，P(A∣Bi)=pi=(n−1i)m
则
P(A)=∑i=0n−1(in−1)m×1nP(A)=\sum_{i=0}^{n-1}(\frac{i}{n-1})^{m} \times \frac{1}{n}P(A)=i=0∑n−1(n−1i)m×n1
代入贝叶斯公式：
P(Bi∣A)=P(Bi)×P(A∣Bi)∑j=0n−1P(Bj)×P(A∣Bj)=pi×1n∑j=0n−1pj×1n\begin{aligned}P(B_i|A)&=\frac{\displaystyle P(B_i)\times P(A|B_i)}{\displaystyle\sum\limits_{j=0}^{n-1}P(B_j)\times P(A|B_j)}&\\&=\cfrac{p_i\times \dfrac{1}{n}}{\displaystyle\sum_{j=0}^{n-1}p_j\times \frac{1}{n}}\end{aligned} P(Bi∣A)=j=0∑n−1P(Bj)×P(A∣Bj)P(Bi)×P(A∣Bi)=j=0∑n−1pj×n1pi×n1
期望 E(x)=P×xE(x)=P\times xE(x)=P×x
xxx 是排名，ppp 为排名为 xxx 时抽 mmm 次抽不到的概率。
显然答案为：
ans=∑x=1nE(x)=∑x=1nP(Bi∣A)×(n−i)=∑i=0n−1pi×1n×(n−i)∑j=0n−1pj×1n=∑i=0n−1pi×(n−i)∑i=0n−1pi=n−∑i=0n−1im+1∑i=0n−1im\begin{aligned}ans & = \sum_{x=1}^{n}E(x)& \\ & = \sum_{x=1}^{n} P(B_i|A) \times (n - i)&\\&=\sum_{i=0}^{n-1}\cfrac{p_i\times \dfrac{1}{n}\times (n-i)}{\displaystyle\sum_{j=0}^{n-1}p_j\times \frac{1}{n}}&\\&=\displaystyle \frac{\displaystyle \sum_{i=0}^{n-1}p_{i}\times(n-i)}{\displaystyle \sum_{i=0}^{n-1}p_{i}}&\\&=n-\cfrac{\displaystyle\sum_{i=0}^{n-1}i^{m+1}}{\displaystyle\sum_{i=0}^{n-1}i^m}\end{aligned} ans=x=1∑nE(x)=x=1∑nP(Bi∣A)×(n−i)=i=0∑n−1j=0∑n−1pj×n1pi×n1×(n−i)=i=0∑n−1pii=0∑n−1pi×(n−i)=n−i=0∑n−1imi=0∑n−1im+1
显然答案就是一个自然数 kkk 次幂之和，我们使用拉格朗日插值 O(n)O(n)O(n) 计算即可。
特判一下 m=0m=0m=0 的情况，即不抽，那么期望就是：
E(x)=P×x=1n×∑i=1ni=1n×n(n+1)2=n+12E(x)=P\times x=\displaystyle\cfrac{1}{n}\times\sum\limits_{i=1}^{n}i=\cfrac{1}{n}\times \cfrac{n(n+1)}{2}=\cfrac{n+1}{2}E(x)=P×x=n1×i=1∑ni=n1×2n(n+1)=2n+1
当然不需要化简为最简分数，因为除法转换成逆元乘起来是等价的。
Code
// Problem: 排名估算
// Contest: NowCoder
// URL: https://ac.nowcoder.com/acm/contest/3782/D
// Memory Limit: 524288 MB
// Time Limit: 2000 ms
//
// Powered by CP Editor (https://cpeditor.org)#include <bits/stdc++.h>
#define int long long
using namespace std;const int N = 5007, M = 5007, mod = 998244353;int read(){int s = 0, ne = 1; char c = getchar();while(c < '0' || c > '9') {if(c == '-') ne = -1; c = getchar();}while(c >= '0' && c <= '9') s = (s << 1) + (s << 3) + c - '0', c = getchar();return s * ne;
}int qpow(int a, int b)
{int res = 1;while(b) {if(b & 1) res = 1ll * res * a % mod;a = 1ll * a * a % mod; b >>= 1;} return res;
}int s[N], pre[N], suf[N], infact[N], fact[N];void init(int k)
{ infact[0] = fact[0] = 1;for(int i = 1; i <= k + 10; ++ i)fact[i] = 1ll * fact[i - 1] * i % mod;infact[k + 10] = qpow(fact[k + 10], mod - 2);for(int i = k + 9; i >= 0; -- i) infact[i] = 1ll * infact[i + 1] * (i + 1) % mod;infact[0] = 1;
}int lagrange(int n, int k)
{ int ans = 0; s[0] = 0;for(int i = 1; i <= k + 2; ++ i)s[i] = (s[i - 1] + qpow(i, k)) % mod;if(n <= k + 1) return s[n]; pre[0] = 1;for(int i = 1; i <= k + 2; ++ i)pre[i] = 1ll * pre[i - 1] * ((n - i + mod) % mod) % mod;suf[k + 3] = 1;for(int i = k + 2; i; -- i)suf[i] = 1ll * suf[i + 1] * ((n - i + mod) % mod) % mod;for(int i = 0; i <= k + 2; ++ i) {s[i] = 1ll * s[i] * pre[i - 1] % mod * suf[i + 1] % mod * infact[i - 1] % mod * infact[k + 2 - i] % mod;if((k + 2 - i) & 1) ans = (1ll * ans - s[i] + mod) % mod;else ans = (1ll * ans + s[i]) % mod;}return (ans + mod) % mod;
}int n, m;signed main()
{init(M - 5);scanf("%lld%lld", &n, &m);if(m == 0) {printf("%lld\n", (n + 1) % mod * qpow(2, mod - 2) % mod);return 0;}int up = lagrange(n - 1, m + 1); int down = lagrange(n - 1, m); down = qpow(down, mod - 2);printf("%lld\n", (n - up * down % mod + mod) % mod);
}
F. P4593 [TJOI2018]教科书般的亵渎
Weblink
https://www.luogu.com.cn/problem/P4593
Problem
Solution
炉石6年老玩家报道！
看好了，我将为大家表演一波教科书般的亵渎（&3/-+7%&567%*19&%…+/*-+…_+%*……&&*）(●ˇ∀ˇ●) ————> 扭了，扭了（扭曲虚空）
老师的术士打的还不够多
显然我们需要使用 k=m+1k=m+1k=m+1 张亵渎，那么每使用一张亵渎，获得的贡献就是一段自然数k次幂的和，我们减掉中间缺掉的随从的贡献就行了。
Time
O(m2)O(m^2)O(m2)
Code
G. P4463 [集训队互测 2012] calc（拉格朗日优化DP）
Weblink
https://www.luogu.com.cn/problem/P4463
Problem
k≤109,n≤500k≤10 ^9 ,n≤500k≤109,n≤500，p≤109p \le 10^9p≤109，并且 ppp 为素数，p>k>n+1p>k>n+1p>k>n+1。
Solution
显然对于一种取值的合法序列，这个序列不管怎么排列，合法序列的值都一样的，我们先考虑暴力计算，设 dp(i,j)dp(i,j)dp(i,j) 表示前 iii 个数取值域范围 [1,j][1,j][1,j] 的所有取值不同的合法序列的值之和。直接转移很不方便，我们可以只考虑递增的序列，即我们仅需讨论第 iii 个数取还是不取 jjj
即：
dp[i][j]=j∗dp[i−1][j−1]+dp[i][j−1]dp[i][j] = j * dp[i - 1][j - 1] + dp[i][j - 1]dp[i][j]=j∗dp[i−1][j−1]+dp[i][j−1]
显然答案就是所有取值不同的合法序列的值之和乘上排列的方案数 n!n!n!。答案就是 dp[n][k]dp[n][k]dp[n][k]，但是 k≤1e9k\le 1e9k≤1e9，考虑优化。
一个DP的递推式可以看作是一个多项式，多项式 fn(i)f_n(i)fn(i) 就是 dp[n][i]dp[n][i]dp[n][i]，那么答案就是 dp[n][k]=fn(k)dp[n][k] = f_{n}(k)dp[n][k]=fn(k)
代入递推式得：
fi(j)−fi(j−1)=j∗fi−1(j−1)f_{i}(j)-f_{i}(j - 1)=j*f_{i - 1}(j - 1)fi(j)−fi(j−1)=j∗fi−1(j−1)
设 fi(j)f_i(j)fi(j) 是 g(n)g(n)g(n) 次多项式
前面是一个差分的形式，显然有结论：
两个 nnn 次多项式的差分是一个 n−1n-1n−1 次多项式
两个 nnn 次多项式的前缀和是一个 n+1n+1n+1 次多项式
自己代入展开算一下就知道了
则 fi(j)−fi(j−1)f_i(j)-f_i(j - 1)fi(j)−fi(j−1) 是 g(n)−1g(n)-1g(n)−1 次多项式，j∗fi−1(j−1)j*f_{i - 1}(j - 1)j∗fi−1(j−1) 是一个 g(n−1)+1g(n-1)+1g(n−1)+1 次多项式（乘上了一个 jjj 嘛），即：
g(n)−1=g(n−1)+1，g(n)=g(n−1)+2g(n)-1=g(n-1)+1，g(n)=g(n-1)+2g(n)−1=g(n−1)+1，g(n)=g(n−1)+2
显然 g(0)=0g(0)=0g(0)=0，则 g(n)=2×ng(n)=2\times ng(n)=2×n。也就意味着我们只需要计算出 fff 的前 2×n+12\times n+12×n+1 项的值，就可以唯一确定一个 2×n2\times n2×n 项的多项式，并且因为我们得到的 2×n+12\times n+12×n+1 项的值中的 xxx 还是连续的，也就意味着我们可以用拉格朗日插值法 ，在 O(2∗n)/O(n2)O(2*n) / O(n^2)O(2∗n)/O(n2) 的复杂度下直接求出 fn(k)f_n(k)fn(k)，既是所求的答案。由于需要 O(n2)O(n^2)O(n2) 预处理dpdpdp 数组，所以复杂度为 O(n2)O(n^2)O(n2)。
当然本题还有生成函数的做法，利用多项式科技可以做到 O(nlogn)O(nlogn)O(nlogn) ：P5850 calc加强版。
Time
O(n2)O(n^2)O(n2)
Code
#include <bits/stdc++.h>using namespace std;
#define int long long
const int N = 5007;int n, m, k, mod;
int dp[N][N];int qpow(int a, int b)
{int res = 1;while(b) {if(b & 1) res = res * a % mod;a = a * a % mod;b >>= 1;}return res;
}int inv(int x)
{return qpow(x, mod - 2);
}signed main()
{scanf("%lld%lld%lld", &k, &n, &mod);int m = 2 * n + 1;for(int i = 0; i <= m; ++ i)dp[0][i] = 1;for(int i = 1; i <= n; ++ i) {for(int j = 1; j <= m; ++ j) {dp[i][j] = (dp[i][j - 1] + dp[i - 1][j - 1] * j % mod) % mod;}}int ans = 0, fact = 1;for(int i = 1; i <= n; ++ i)fact = fact * i % mod;if(k <= m) {ans = dp[n][k];printf("%lld\n", ans * fact % mod);return 0;}for(int i = 1; i <= m; ++ i) {int up = dp[n][i], down = 1;for(int j = 1; j <= m; ++ j) {if(i != j) {up = up * (k - j + mod) % mod; down = down * (i - j + mod) % mod;}}ans = (ans + up * inv(down) % mod) % mod;} printf("%lld\n", ans * fact % mod);
}
H. CF995F Cowmpany Cowmpensation（拉格朗日优化DP）
Weblink
https://www.luogu.com.cn/problem/CF995F
Problem
%给定n个点m条边的有向无环图，其中没有入度的点被视为源点，没有出度的点被视为汇点。 保证源点和汇点数目相同。 考虑所有把源汇点两两配对，并用两两不相交的路径把它们两两连接起来的所有方案。 如果这个方案中，把源点按标号1到n排序后，得到的对应汇点序列的逆序数对的个数是奇数，那么A给B一块钱，否则B给A一块钱。 问最后A的收益，对大质数取模。 n ≤ 600
树形结构，给每个节点分配工资（[1,d]），子节点不能超过父亲节点的工资，问有多少种分配方案
Solution
Code
杜教的代码！%%%
#include <bits/stdc++.h>
using namespace std;
#define rep(i,a,n) for (int i=a;i<n;i++)
#define per(i,a,n) for (int i=n-1;i>=a;i--)
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(),(x).end()
#define fi first
#define se second
#define SZ(x) ((int)(x).size())
typedef vector<int> VI;
typedef long long ll;
typedef pair<int,int> PII;
const ll mod=1000000007;
ll powmod(ll a,ll b) {ll res=1;a%=mod; assert(b>=0); for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;}
ll gcd(ll a,ll b) { return b?gcd(b,a%b):a;}
// headnamespace polysum {const int D=101000;ll a[D],f[D],g[D],p[D],p1[D],p2[D],b[D],h[D][2],C[D];ll calcn(int d,ll *a,ll n) {if (n<=d) return a[n];p1[0]=p2[0]=1;rep(i,0,d+1) {ll t=(n-i+mod)%mod;p1[i+1]=p1[i]*t%mod;}rep(i,0,d+1) {ll t=(n-d+i+mod)%mod;p2[i+1]=p2[i]*t%mod;}ll ans=0;rep(i,0,d+1) {ll t=g[i]*g[d-i]%mod*p1[i]%mod*p2[d-i]%mod*a[i]%mod;if ((d-i)&1) ans=(ans-t+mod)%mod;else ans=(ans+t)%mod;}return ans;}void init(int M) {f[0]=f[1]=g[0]=g[1]=1;rep(i,2,M+5) f[i]=f[i-1]*i%mod;g[M+4]=powmod(f[M+4],mod-2);per(i,1,M+4) g[i]=g[i+1]*(i+1)%mod;}ll polysum(ll n,ll *a,ll m) { // a[0].. a[m] \sum_{i=0}^{n-1} a[i]a[m+1]=calcn(m,a,m+1);rep(i,1,m+2) a[i]=(a[i-1]+a[i])%mod;return calcn(m+1,a,n-1);}ll qpolysum(ll R,ll n,ll *a,ll m) { // a[0].. a[m] \sum_{i=0}^{n-1} a[i]*R^iif (R==1) return polysum(n,a,m);a[m+1]=calcn(m,a,m+1);ll r=powmod(R,mod-2),p3=0,p4=0,c,ans;h[0][0]=0;h[0][1]=1;rep(i,1,m+2) {h[i][0]=(h[i-1][0]+a[i-1])*r%mod;h[i][1]=h[i-1][1]*r%mod;}rep(i,0,m+2) {ll t=g[i]*g[m+1-i]%mod;if (i&1) p3=((p3-h[i][0]*t)%mod+mod)%mod,p4=((p4-h[i][1]*t)%mod+mod)%mod;else p3=(p3+h[i][0]*t)%mod,p4=(p4+h[i][1]*t)%mod;}c=powmod(p4,mod-2)*(mod-p3)%mod;rep(i,0,m+2) h[i][0]=(h[i][0]+h[i][1]*c)%mod;rep(i,0,m+2) C[i]=h[i][0];ans=(calcn(m,C,n)*powmod(R,n)-c)%mod;if (ans<0) ans+=mod;return ans;}
}const int N=3010;
int n,d,p,dp[N][N];
VI s[N];
ll pres[N];
void dfs(int u) {rep(i,0,n+1) dp[u][i]=1;for (auto v:s[u]) {dfs(v);rep(i,0,n+1) pres[i]=dp[v][i];rep(i,1,n+1) pres[i]=(pres[i]+pres[i-1])%mod;rep(i,0,n+1) dp[u][i]=dp[u][i]*pres[i]%mod;}
}int main() {scanf("%d%d",&n,&d); --d;polysum::init(3456);rep(i,2,n+1) {scanf("%d",&p);s[p].pb(i);}dfs(1);rep(i,0,n+1) pres[i]=dp[1][i];rep(i,1,n+1) pres[i]=(pres[i]+pres[i-1])%mod;printf("%lld\n",polysum::calcn(n,pres,d));
}
I. （2018牛客多校（一）F）Sum of Maximum（组合数学+拉格朗日插值）
http://tokitsukaze.live/2018/07/19/2018niuke1.F/
%[https://blog.csdn.net/qq_42819598/article/details/95497117?utm_medium=distribute.pc_relevant.none-task-blog-2%7Edefault%7EBlogCommendFromMachineLearnPai2%7Edefault-3.control&depth_1-utm_source=distribute.pc_relevant.none-task-blog-2%7Edefault%7EBlogCommendFromMachineLearnPai2%7Edefault-3.control](https://blog.csdn.net/qq_42819598/article/details/95497117?utm_medium=distribute.pc_relevant.none-task-blog-2~default~BlogCommendFromMachineLearnPai2~default-3.control&depth_1-utm_source=distribute.pc_relevant.none-task-blog-2~default~BlogCommendFromMachineLearnPai2~default-3.control)
多项式求值与插值
A. P5050 【模板】多项式多点求值
Weblink
https://www.luogu.com.cn/problem/P5050
Problem
Solution
Time
O(nlog2n)O(nlog^2n)O(nlog2n)
Code
B. P5158 【模板】多项式快速插值
Weblink
https://www.luogu.com.cn/problem/P5158
Problem
Solution
好长时间没用笔写过字了，怎么这么丑…
Time
O(nlog2n)O(nlog^2n)O(nlog2n)
Code
%https://blog.csdn.net/a_forever_dream/article/details/112547349
%https://blog.csdn.net/C20190102/article/details/106693455
解题报告（三）多项式求值与插值（拉格朗日插值）（ACM / OI）相关推荐
- 【解题报告系列】超高质量题单 + 题解（ACM / OI）超高质量题解
整理的算法模板合集: ACM模板 点我看算法全家桶系列!!! 实际上是一个全新的精炼模板整合计划 繁凡出品的全新系列:解题报告系列 -- 超高质量算法题单,配套我新写的超高质量的题解和代码,题目难度不 ...
- matlab多项式的求值,多项式求值的MATLAB实现
公茂果老师的课件中,给出了四种多项式求值的算法,下面给出代码示例: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %E-mail: [email protected] ...
- 《C#零基础入门之百识百例》（五十二）封装介绍 -- 二维多项式求值
C#零基础入门 面向对象 -- 封装介绍 -- 二维多项式求值 前言 一,封装概念 二,封装属性 三,实例练习 -- 二维多项式求值 3.1 题目描述 3.2 问题分析 3.3 参考代码 前言 本文属 ...
- 【洛谷新手村解题报告三 字符串/递归前】C++语言，一题多解，思路和WA反思
[洛谷新手村解题报告三] 简单字符串 过程函数与递归 首先下一部分!字符串,这个前面两题新手跳过吧 简单字符串 第三题 统计单词数 [2/2] 给定一个单词,请你输出它在给定的文章中出现的次数和第一次 ...
- 6-2 多项式求值 (15 分)
6-2 多项式求值 (15 分) 本题要求实现一个函数,计算阶数为n,系数为a[0] ... a[n]的多项式f(x)=∑i=0n(a[i]×xi) 在x点的值. 函数接口定义: double f( ...
- C语言实现一维多项式求值
/* * plyv.h * * Created on: Oct 14, 2010 * Author: jenson */ #ifndef PLYV_H_ #define PLYV ...
- 1002: 简单多项式求值 ZZULIOJ
1002: 简单多项式求值 题目描述 对用户输入的任一整数,输出以下多项式的值. y=2x^2+x+8 输入 输入整数x的值. 输出 输出一个整数,即多项式的值. 样例输入 Copy 1 样例输出 C ...
- 【C++】用递归函数实现Hermite多项式求值
要求:用递归函数实现Hermite多项式求值.当x>1时,Hermite多项式定义为: 当输入实数x和整数n后,求出Hermite多项式前n项的值. #include <iostream& ...
- Matlab多项式求值（2）（roots和poly的使用）
Matlab多项式求值(2)(roots和poly的使用) #Matlab学习笔记 1.多项式求根 n次多项式具有n个根(包括实根和共轭复根) roots函数可以用于求解多项式的全部根,其调用格式为: ...
最新文章
- 打一场AI竞赛，让你知道我的厉害
- R语言ggplot2可视化：使用ggplot2绘制按时间顺序排列的时间线图（chronological timeline plot）
- Linux怎么让文件按大小排序,linuxsortlinux系统下，按文件的大小进行排序的命令...
- centos7 开启 关闭 NetworkManager
- opensuse13.2 安裝五笔
- python中readlines_python中read() readline()以及readlines()用法
- PCM设备在电网系统中的应用介绍
- mysql数据库ACID实现原理
- win10鼠标灵敏度怎么调_和平精英灵敏度怎么设置才最合适 调最适合自己用的
- oracle 函数_oracle中row_number、rank、dense_rank() 的区别（over函数）
- linux脚本猜数字游戏,101个shell脚本 猜数字游戏代码
- 紧急！你用的是这20家 VPS 吗？它们疑似马上要携款跑路！
- Windows 11 将“扼杀”第三方浏览器？Firefox 绝地反击！
- 2018icpc ecfinal 西安游记
- 通过IIS安装包安装IIS
- 台式计算机上的硬磁盘,如何在台式计算机上安装机械硬盘驱动器？在台式计算机上安装机械硬盘驱动器的详细步骤...
- contiki笔记2-contiki的第一个程序
- 给定一个字符串s，请你找出其中不含有重复字符的最长子串的长度。
- CSS简介代码实现表格(table)只有横线没有竖线的效果，一个简约卡片化精美的数据表格(类似Bootstrap表格)
- Java Import