1. 程式人生 > 其它 >拉格朗日插值模板(P4781)

拉格朗日插值模板(P4781)

題目連結:拉格朗日插值
拉格朗日插值:給定 \(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)\)


這樣就完成了拉格朗日插值;模板:
\(Code:\)

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

# 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; }