拉格朗日插值学习笔记
应用
众所周知，在平面直角坐标系中，对于任意的 \(n\) 个点，都一定有一个不超过 \(n-1\) 次的函数与之相对应。拉格朗日插值适用于求解这 \(n\) 个点对应的函数。
思路
考虑给定的 \(n\) 个点的坐标表示为 \((x_i,y_i)\)，不难构造出如下函数：
\[f(x)=\sum_{i=1}^{n}y_ig_i(x) \]那么此时只需要构造出符合要求的 \(g_i\) 即可。不难发现，为了使 \(f\) 符合条件，\(g_i\) 应该满足：
\[\forall n\in x,g_i(n)=\left\{\begin{matrix}1&n=x_i\\0&n\ne x_i\end{matrix}\right. \]其中，\(x\) 表示所有 \(x_i\) 构成的集合。不难发现，当 \(n\ne x_i\) 时函数值为 \(0\)，说明该函数 \(g_i(x)\) 一定有因式 \((x-x_j)(i\ne j)\)，因此不难看出 \(g_i(x)\) 应当有因式 \(\prod_{j\ne i}(x-x_j)\)。又因为当 \(n=x_i\) 时，函数值为 \(0\)，不难构造出：
\[g_i(n)=\frac{\prod_{j\ne i}(n-x_j)}{\prod_{j\ne i}(x_i-x_j)} \]如此就构造出了符合条件的函数，此时可以 \(O(n^2)\) 求解。在单次求解的时候可以采用。
不妨进行进一步转换，令 \(t_i=\frac{y_i}{\prod_{j\ne i}(x_i-x_j)}\)，不难看出，\(f(x)=t_i\prod_{j\ne i}(x-x_j)\)。不妨计算出 \(\prod_{i=1}^{n}(x-x_i)\) 后将对应项除掉，接着利用多项式乘法算出每一次项前的系数后，就可以按照多项式加法的方式 \(O(n)\) 求解。这样的方式适用于多组查询或需要求出函数对应项系数的情况。
代码
\(O(n^2)\) 写法
#include <bits/stdc++.h>
using namespace std;
#define ll long long
const int N=2e3+5,M=998244353;
int n,k;
int x[N],y[N];
ll QuickPow(ll a,int b){
ll res=1;
while(b>0){
if(b&1)res=res*a%M;
a=a*a%M;
b>>=1;
}
return res;
}
ll Lagrange(int x[],int y[],int X){
ll res=0;
for(int i=1;i<=n;i++){
ll l_i=1;
for(int j=1;j<=n;j++){
if(i!=j)l_i=l_i*(X-x[j])%M*QuickPow(x[i]-x[j],M-2)%M;
}
res=(res+l_i*y[i]%M)%M;
}
return res;
}
int main(){
scanf("%d%d",&n,&k);
for(int i=1;i<=n;i++)scanf("%d%d",&x[i],&y[i]);
printf("%lld",(Lagrange(x,y,k)+M)%M);
return 0;
}
\(O(n)\) 做法
#include <bits/stdc++.h>
using namespace std;
#define ll long long
const int N=2e3+5,M=998244353;
int n,k;
ll t[N],x[N],y[N],g[N],fs[N],f[N];
ll QuickPow(ll a,int b){
ll res=1;
while(b>0){
if(b&1)res=res*a%M;
a=a*a%M;
b>>=1;
}
return res;
}
void Init(){
for(int i=0;i<n;i++){
t[i]=1;
for(int j=0;j<n;j++){
if(i==j)continue;
t[i]=t[i]*(x[i]-x[j])%M;
}
t[i]=QuickPow(t[i],M-2)*y[i]%M;
}
fs[0]=1;
for(int i=0;i<n;i++){
for(int j=n-1;j>0;j--)fs[j]=(fs[j-1]+fs[j]*(M-x[i])%M)%M;
fs[0]=fs[0]*(M-x[i])%M;
}
for(int i=0;i<n;i++){
ll inv=QuickPow(M-x[i],M-2);
g[0]=fs[0]*inv%M;
for(int j=1;j<n;j++)g[j]=(fs[j]-g[j-1])*inv%M;
for(int j=0;j<n;j++)f[j]=(f[j]+t[i]*g[j]%M)%M;
}
return ;
}
ll Lagrange(ll X){
ll res=0,Pow=1;
for(int i=0;i<n;i++){
res=(res+Pow*f[i]%M)%M;
Pow=Pow*X%M;
}
return res;
}
int main(){
scanf("%d%d",&n,&k);
for(int i=0;i<n;i++)scanf("%lld%lld",&x[i],&y[i]);
Init();
printf("%lld",(Lagrange(k)+M)%M);
return 0;
}