拉格朗日插值

P4781 【模板】拉格朗日插值

广为人知的是，\((n+1)\) 个点可以唯一确定一个 \(n\) 次多项式，也就是多项式的点值表示法。
如果我们给出这 \(n+1\) 个点，想要求函数上的另外一个值，就需要用到拉格朗日插值法。

具体的，我们并没有将这个函数的表达式具体的求出来，而是直接求出我们想要的值是多少。

事实上，我们有

\[f(x)=\sum_{i=1}^n y_i\times \prod_{j\ne i} \frac{x-x_j}{x_i-x_j} \]

也就是说我们可以直接把整个函数都给算出来。
如果我们只是算某个 \(x\) 坐标上的值的话，那我们直接将这个值带入即可。具体的，我们设我们要求 \(x\) 坐标为 \(k\) 的值，有

\[f(x)=\sum_{i=1}^n y_i\times \prod_{j\ne i} \frac{k-x_j}{x_i-x_j} \]

照着这个东西算即可。可以发现时间复杂度 \(O(n^2)\)。

code

点击查看代码

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=2e3+7,p=998244353;
int ksm(int x,int k){
	int res=1ll;x%=p;
	while(k){
		if(k&1)res=res*x%p;
		x=x*x%p;k>>=1;
	}
	return res;
}
int inv(int x){return ksm(x,p-2);}
int x[N],y[N];
signed main(){
	ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);
	int n,k,ans=0;cin>>n>>k;
	for(int i=1;i<=n;i++)cin>>x[i]>>y[i];
	for(int i=1;i<=n;i++){
		int s1=1ll,s2=1ll;
		for(int j=1;j<=n;j++){
			if(i==j)continue;
			s1=s1*(k-x[j]+p)%p;
			s2=s2*(x[i]-x[j]+p)%p;
		}
		ans=(ans+y[i]*s1%p*inv(s2)%p)%p;
	}
	cout<<ans;return 0;
}

对于连续的横坐标的优化

在解决实际问题的时候，我们给定的点的横坐标很多时候都是连续的，在这种情况下我们可以将其复杂度优化到 \(O(n)\)。

我们发现，当 \(x_i\) 连续的时候，上面的式子中 \(\prod\) 后面的东西就比较有规律了。我们可以预处理出 \(g_{1,j}=\prod_{i=1}^j k-x_i\) 与 \(g_{2,j}=\prod_{i=j}^n k-x_i\)，这样就可以快速计算分子。
分母的话是好求的，就是 \((-1)^{(n-i)}(i-1)!(n-i)!\)，然后要预处理一个阶乘逆元。
然后就做完了，需要一些精细实现。复杂度 \(O(n)\)。

code

需要注意的是那个符号问题。同时为了保证复杂度需要注意一些预处理的东西。
由于没有横坐标连续的模板题，因此就写个了和普通模板题一样的格式。大概是对的。

点击查看代码

#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=2e6+7,p=998244353;
int x[N],y[N],pre[N],ni[N],s1[2][N];
int ksm(int x,int k){
	int res=1;x%=p;
	while(k){
		if(k&1)res=res*x%p;
		x=x*x%p;k>>=1;
	}
	return res;
}
int inv(int x){return ksm(x,p-2);}
void init(int n,int k){
	pre[0]=1;for(int i=1;i<=n;i++)pre[i]=pre[i-1]*i%p;
	ni[n]=inv(pre[n]);for(int i=n-1;i>=0;i--)ni[i]=ni[i+1]*(i+1)%p;
	s1[0][0]=s1[1][n+1]=1;
	for(int i=1;i<=n;i++)s1[0][i]=s1[0][i-1]*(k-x[i]+p)%p;
	for(int i=n;i>=1;i--)s1[1][i]=s1[1][i+1]*(k-x[i]+p)%p;
}
signed main(){
	ios::sync_with_stdio(false),cin.tie(0),cout.tie(0);
	int n,k,s2=1,ans=0;cin>>n>>k;for(int i=1;i<=n;i++)cin>>x[i]>>y[i];
	init(n,k);
	for(int i=1;i<=n;i++){
		s2=ni[i-1]*ni[n-i]%p;if((n-i)&1)s2=p-s2;//注意这里的符号问题。 
		ans=(ans+y[i]*(s1[0][i-1]*s1[1][i+1]%p)%p*s2%p)%p;
	}
	cout<<ans;return 0;
}