拉格朗日插值Python代码 拉格朗日插值程序框图

拉格朗日插值：给定 \(k+1\) 个点对 \((x_i,y_i)\) (\(x_i\)各不相同)能够唯一确定一个最高次为 \(k\) 次的多项式，那么如何进行构造，来求该多项式呢？我们先以经过 \((x_1,1),(x_2,0),(x_3,0)\) 这三个点的4次多项式为例：那么我们可以进行构造设 \(f(X)=\frac{(X-x_2)*(X-x_3)}{(x_1-x_2)*(x_1-x_3)}\)，我们带入这三个点发现显然该多项式是我们所求的多项式。那么如果我们将该三个点扩展到 \((x_1,y_1),(x_2,y_2),(x_3,y_3)\) 呢？
显然：

• \(f_1(X)=\frac{(X-x_2)*(X-x_3)}{(x_1-x_2)*(x_1-x_3)}*y_1\)经过\((x_1,y_1)\)，且\(f_1(x_2) = 0,f_1(x_3) = 0\)
• \(f_2(X)=\frac{(X-x_1)*(X-x_3)}{(x_2-x_1)*(x_2-x_3)}*y_2\)经过\((x_2,y_2)\)，且\(f_2(x_1) = 0,f_2(x_3) = 0\)
• \(f_3(X)=\frac{(X-x_1)*(X-x_2)}{(x_3-x_1)*(x_3-x_2)}*y_3\)经过\((x_3,y_3)\)，且\(f_3(x_1) = 0,f_3(x_2) = 0\)

所以我们最后构造的多项式为\(F(X) = f_1(X) + f_2(X) + f_3(X)\)
所以扩展到 \(n\) 项即\(F(X) = \sum_{i=1}^{i=n}f_i(X)\)
这样就完成了拉格朗日插值；模板：
\(Code:\)

/* -*- encoding: utf-8 -*-
'''
@File    :   拉格朗日模板.cpp
@Time    :   2021/06/28 11:21:28
@Author  :   puddle_jumper 
@Version :   1.0
@Contact :   [e-mail]
'''

# here put the import lib*/
#include<set>
#include<iostream>
#include<cstring>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<map>
#include<algorithm>
#include<vector>
#include<queue>
#define ch() getchar()
#define pc(x) putchar(x)
#include<stack>
#include<unordered_map>
#define rep(i,a,b) for(auto i=a;i<=b;++i)
#define bep(i,a,b) for(auto i=a;i>=b;--i)
#define lowbit(x) x&(-x)
#define ll long long
#define ull unsigned long long
#define pb push_back
#define mp make_pair
#define PI acos(-1)
using namespace std;
template<typename T>void read(T&x){
	static char c;
	static int f;
	for(c=ch(),f=1; c<'0'||c>'9'; c=ch())if(c=='-')f=-f;
	for(x=0; c>='0'&&c<='9'; c=ch())x=x*10+(c&15);
	x*=f;
}
template<typename T>void write(T x){
	static char q[65];
	int cnt=0;
	if(x<0)pc('-'),x=-x;
	q[++cnt]=x%10,x/=10;
	while(x)
		q[++cnt]=x%10,x/=10;
	while(cnt)pc(q[cnt--]+'0');
}

const int N = 2e3+10;
const ll mod = 998244353;
int n;ll k;
ll x[N],y[N];

ll qpow(ll a,ll b){
	ll res = 1ll;
	while(b){
		if(b&1){
			res *= a;res%=mod;
		}
		b>>=1;
		a*=a;a%=mod;
	}
	return res%mod;
}
ll inv(ll x){
	return qpow(x,mod-2);
}
void solve(){
	read(n);read(k);
	rep(i,1,n)read(x[i]),read(y[i]);
	ll ans = 0ll;
	rep(i,1,n){
		ll res = y[i];
		rep(j,1,n){
			if(i == j)continue;
			res *= (k-x[j]);
			res %= mod;
			res *= inv(x[i] - x[j]);
			res %= mod;
		}
		ans += (res%mod + mod)%mod;
		ans %= mod;
	}
	write(ans);pc('\n');
}

signed main(){solve();return 0; }