ICode9

精准搜索请尝试: 精确搜索
首页 > 其他分享> 文章详细

【模板】拉格朗日插值

2022-04-30 09:31:54  阅读:152  来源: 互联网

标签:拉格朗 函数 limits 插值 ne int prod 模板 mod


link

开始学习数学了……

众所周知拉格朗日插值可以解决如下问题:给定N个点,可以在\(O(N^2)\)的时间内求出过这些点的次数和项数都为N-1的多项式,并且求出该多项式在另一个自变量下的取值。

拉格朗日插值的思想是构造。用一句经典的话来说,假如我们有N个项数和次数都为N-1的函数 \(g_1,g_2\dots g_N\) ,并且函数 \(g_i\) 满足

\[g_i(x_j)=\begin{cases}0&j\ne i\\1&j=i\end{cases}(j\in[1,N]) \]

于是可以发现,假如有函数

\[f=\sum\limits_{i=1}^Ng_i\times y_i \]

那么显然这个函数的次数和项数都符合条件,而且某个\(x_i\)只可能在\(g_i\)时有\(y_i\)的贡献,其它时候都被0掉了,所以这个函数是满足上述条件的,也就是说这样构造出的函数是合法的。

但是那个 \(g\) 函数又应当如何构造呢?根据奇怪的因式定理可以尝试构造函数 \(r\):

\[r_i=\prod\limits_{j=1}^N(x-x_j)(j\ne i) \]

那么这个函数 \(r\) 就可以满足上面那个0的条件。那1的条件呢?直接除以一个一样的式子不就好啦。于是:

\[g_i=\frac{r_i}{r_i}=\frac{\prod\limits_{j=1}^N(x-x_j)}{\prod\limits_{j=1}^N(x_i-x_j)}=\prod\limits_{j=1}^N\frac{x-x_j}{x_i-x_j}(j\ne i) \]

合起来就可以得到最终的结论:

\[f=\sum\limits_{i=1}^Ny_i\times\prod\limits_{j=1}^N\frac{x-x_j}{x_i-x_j}(j\ne i) \]

然后发现可以直接带值进行计算。代码挺好写的,有几个地方需要注意一下:

  • 除法部分确实需要用逆元,但为了少求几次可以考虑全部乘起来之后再一次性搞定。
  • 关于取模,最后答案可能是负数,这就意味着必须要用馍加馍的方法把它强行搞成正数。连Wa三次得出的教训。
#include<cstdio>
#define zczc
#define int long long
const int mod=998244353;
const int N=2010;
inline void read(int &wh){
	wh=0;int f=1;char w=getchar();
	while(w<'0'||w>'9'){if(w=='-')f=-1;w=getchar();}
	while(w>='0'&&w<='9'){wh=wh*10+w-'0';w=getchar();}
	wh*=f;return;
}

int m,w,ans,x[N],y[N];

inline int pow(int s1,int s2){
	if(s2==1)return s1;if(s2==0)return 1;
	int an=pow(s1,s2>>1);
	if(s2&1)return an*an%mod*s1%mod;
	else return an*an%mod;
}

signed main(){
	
	#ifdef zczc
	freopen("in.txt","r",stdin);
	#endif
	
	read(m);read(w);
	for(int i=1;i<=m;i++){read(x[i]);read(y[i]);}
	for(int i=1;i<=m;i++){
		int a=y[i],b=1;
		for(int j=1;j<=m;j++){
			if(i==j)continue;
			a*=w-x[j];a%=mod;
			b*=x[i]-x[j];b%=mod;
		}
		ans+=a*pow(b,mod-2)%mod;ans%=mod;
	}
	printf("%lld",(ans+mod)%mod);
	
	return 0;
}

标签:拉格朗,函数,limits,插值,ne,int,prod,模板,mod
来源: https://www.cnblogs.com/dai-se-can-tian/p/16209462.html

本站声明: 1. iCode9 技术分享网(下文简称本站)提供的所有内容,仅供技术学习、探讨和分享;
2. 关于本站的所有留言、评论、转载及引用,纯属内容发起人的个人观点,与本站观点和立场无关;
3. 关于本站的所有言论和文字,纯属内容发起人的个人观点,与本站观点和立场无关;
4. 本站文章均是网友提供,不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属;如您发现该文章侵犯了您的权益,可联系我们第一时间进行删除;
5. 本站为非盈利性的个人网站,所有内容不会用来进行牟利,也不会利用任何形式的广告来间接获益,纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。

专注分享技术,共同学习,共同进步。侵权联系[81616952@qq.com]

Copyright (C)ICode9.com, All Rights Reserved.

ICode9版权所有