1 算法介绍
众所周知的，\(n\) 个平面上的 \(x\) 坐标不同的点始终能用 \(n-1\) 次多项式来拟合。
假设拟合函数为 \(f\)，给出点的坐标为 \((x_i,y_i)\)，那么有：
把每个点坐标带入式子会发现这是对的，然后把式子展开不难得到这是个 \(n-1\) 次多项式，于是我们可以 \(O(n^2)\) 求某个点的函数值。
嗯，这就是拉插的全部内容了，但我想我们还需要一些实用的拓展🤨：
2 算法拓展
2.1 当保证需要拟合的点坐标连续
标题的意思是保证 \(x_i=i\)。
把此时的式子写出来：
然后我们发现分母预处理阶乘后可以随意求，分子只需要预处理 \(k-i\) 的前后缀积也不难求。
发现这是 \(O(n)\) 的，妙中妙🥰。
2.2 当需要确切的多项式值||更通用的优化方法
当你写代码时需要确切的多项式值时，瞎找个位置暴力模拟多项式乘/加法展开后面的式子即可做到 \(O(n^3)\)，然而这个复杂度连小学生暴写高斯消元都能做到，劣中劣😂。
但我们定睛一看🤔：
主播你是不是能把 \(g(k)=\prod\limits_{j=1}^{n} (k-x_j)\) 提出来啊，这样就不用重复算了。
发现你现在模拟多项式暴算这个式子是 \(O(n^2)\) 的。
又发现把这个求出来后你以后的询问甚至可以直接带入 \(O(n)\) 算，这对于任意情景都能做到 \(O(n^2)\) 预处理，\(O(n)\) 查询，强中强😁。
2.3 重心拉格朗日插值法||更更通用的优化方法
这个时候有人就要问了，主播你写了一坨模拟多项式计算我给个强制在线怎么就死完了啊🤣。
你先别管啥神人会让你强制在线搞这个，但一坨模拟多项式计算的确丑陋。
于是我们还是看这个式子：
这次我们设 \(w(i)=\prod\limits_{i\ne j}\frac{1}{x_i-x_j}\)，这就是重心权。
我们加新点的时候计算他本身的 \(w\) 并更新其他的 \(w\) 即可😤。
这又实现了 \(O(n)\) 加新点，\(O(n)\) 查单点值，注意到对于普通题目他的预处理复杂度仍是 \(O(n^2)\) 的。
对于重心拉格朗日插值，我们有另外的方法能优化常数、缩短码长。
我们考虑对于一个恒等于 \(1\) 的函数 \(h\) 进行重心拉格朗日插值：
然后我们对于原本的重心拉格朗日插值左右同时除以上面的东西：
注意到 \(h(k)=1\)，\(g(k)\) 能抵消，\(\prod\limits_{i\ne j}\frac{1}{x_i-x_j}=w(i)\)，于是化简得：
优点是优化掉了一个 \(g\) 函数，代码能短几行。
2.5 小结
这个时候有人又要问了，主播你给了这么一坨优化方法我该怎么选啊？
没关系我们有美观表格：
| 是否需要 \(x_i=i\) | 在线插入复杂度 | 离线预处理复杂度 | 查询复杂度 | 求原系数复杂度 | |
|---|---|---|---|---|---|
| 朴素拉插 | ❌ | \(O(1)\) | \(O(1)\) | \(O(\textcolor{red}{n^2})\) | \(O(n^3)\) |
| 分别预处理分子分母 | ✔ | \(O(1)\) | \(O(n)\) | \(O(n)\) | \(O(n^3)\) |
| 模拟多项式计算求原系数 | ❌ | \(O(\textcolor{red}{n^2})\) | \(O(n^2)\) | \(O(n)\) | \(O(n^2)\) |
| 重心拉格朗日插值 | ❌ | \(O(n)\) | \(O(n^2)\) | \(O(n)\) | \(O(n^2)\) |
注意到朴素拉插看似插入快但不难发现他只有一种情景会比重心拉插快——只插入不查询😱，所以这是被完爆的。
然后我们的模拟多项式计算求原系数各项数值都不比重心拉插优，但是我给它求原系数复杂度 \(O(n^2)\) 的原因是他的式子 可以 \(O(n^2)\) 做但你还是要写模拟多项式计算，我说这是集其之长了。
所以说手动实现下 2.1 和 2.3 就好了，2.2 和最初始的不用写。
2.1 代码，这里以 CF622F 为例
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int mod=1e9+7;
int n,k,suml[1001000],sumr[1001000],fac[1001000];
int qpow(int x,int y){assert(x>=0);int ans=1;while(y){if(y&1){ans=ans*x%mod; }x=x*x%mod;y>>=1;}return ans;
}
signed main(){cin>>n>>k;suml[0]=sumr[k+3]=fac[0]=1;for(int i=1;i<=k+2;i++){fac[i]=fac[i-1]*i%mod;suml[i]=suml[i-1]*((n-i+mod)%mod)%mod;}for(int i=k+2;i>=1;i--){sumr[i]=sumr[i+1]*((n-i+mod)%mod)%mod;}int sum=0,ans=0;for(int i=1;i<=k+2;i++){sum=(sum+qpow(i,k))%mod;
// cout<<i<<' '<<sum<<" "<<suml[i-1]*sumr[i+1]%mod<<' '<<fac[i-1]*fac[k+2-i]%mod*((k-i)&1?-1:1)<<'\n';ans=(ans+sum*suml[i-1]%mod*sumr[i+1]%mod*qpow((fac[i-1]*fac[k+2-i]%mod*((k-i)&1?-1:1)+mod)%mod,mod-2)%mod)%mod;}cout<<ans<<'\n';return 0;
}
2.3 代码（这里写的第二种）
struct Murasame{
int n,x[10010],y[10010],w[10010];
void operator ()(){n=0;
}
int operator ()(int X){int up=0,dw=0;for(int i=1;i<=n;i++){if(x[i]==X) return y[i];up=(up+w[i]*inv(X-x[i])%mod*y[i]%mod)%mod;dw=(dw+w[i]*inv(X-x[i])%mod)%mod;}return (up*inv(dw)%mod+mod)%mod;
}
void operator ()(int fi,int se){x[++n]=fi;y[ n]=se;
}
void init(){for(int i=1;i<=n;i++){w[i]=1;for(int j=1;j<=n;j++){if(j==i) continue;w[i]=w[i]*inv(x[i]-x[j])%mod;}}
}
}
3 应用
3.1 求 k 次方的前缀和
标题的意思是求 \(\sum\limits_{i=1}^{n} i^k\)。
传送门。
有人要问了：啊主播主播你给的 2.1 的例题我咋不会啊。
那么就引出了拉格朗日插值的一个经典应用了——求 \(\sum\limits_{i=1}^{n} i^k\)。
具体解释是这个东西会是 \(d+1\) 次多项式，然后就可以拉格朗日插值了。
这个东西的严谨证明冗长而晦涩，于是我决定感性速通一手。
首先我认为高阶等差数列的定义小学二年级找规律会讲，不知道的建议重学。
点击重学
给一个 4 阶等差数列的例子，看完不懂定义是这个👎1 1 1 11 2 3 4 51 2 4 7 11 161 2 4 8 15 26 42
1 2 4 8 16 31 57 99
定理：数列 \(a_n\) 是一个 \(k\) 阶等差数列当且仅当数列的通项 \(a_n\) 为一个 \(n\) 的 \(k\) 次多项式。
证明：
考虑转为证明一次差分后通项公式的次数会 \(-1\)，这显然与原命题同真伪。
即：若有 \(k\) 次多项式 \(f(x)\)，其差分函数 \(g(x)=f(x+1)-f(x)\) 中的 \(x^k\) 项系数为 \(0\)。
令 \(f\) 中 \(x^k\) 项的系数为 \(a\)，现在只考虑 \(x^k\) 项有：\(a\times (x+1)^k-a\times x^k\) 注意到 \((x+1)^k\) 中的 \(x^k\) 项系数为 \(1\)，故最后 \(g\) 中的 \(x^k\) 项全部抵消。
看完证明后我们有一个牛的：对于任意多项式通项函数，有差分后次数 \(-1\)，前缀和后次数 \(+1\)。
那么回到原问题中，发现我们求的是 \(f(x)=x^k\) 的前缀和，于是这就是一个 \(k+1\) 次多项式，直接拉插即可。
3.2 优化 dp
传送门。
首先你把这个 \(a\) 改成无序集合，最后将答案乘 \(n!\) 即可。
然后有我奶奶都会的 dp，设 \(f(i,j)\) 为在 \([1,j]\) 中选 \(i\) 个数的权值和，有转移：
这是 \(O(nk)\) 的，我们发现这个状态 \(i\) 小 \(j\) 大，考虑搞成关于 \(j\) 的多项式然后拉格朗日插值求。
设 \(f(i,j)\) 为 \(g(i)\) 次多项式，整理式子得：
左边是个差分的样子，次数应该是 \(g(i)-1\)，右边因为又乘了个 \(j\) 所以次数应该是 \(g(i-1)+1\)，得 \(g(i)-1=g(i-1)+1\)，即 \(g(i)=g(i-1)+2\)，故我们只需要知道前 \(2i\) 个点就可以拉插搞出答案，复杂度优化为 \(O(n^2)\)。