拉格朗日插值

$拉格朗日插值法可构造一个经过任意n个点(a_1,b_1),(a_2,b_2),(a_3,b_3),\dots ,(a_n,b_n)的函数f,设f=f_1+f_2+\cdots +f_n,需要让f_1,f_2,\dots ,f_n均经过这n个点,且当x=a_i时,f_i=b_i;当x=a_j(j\ne i)时,f_i=0.$

$那么显然,由f_1,f_2,\dots ,f_n相加得到的的f,恰好经过了a_1,a_2,\dots ,a_n.$
$f_1,f_2,\dots ,f_n被称为拉格朗日基函数,表达式为:f_i=b_i{\prod }_{j\ne i}\frac{x-x_j}{x_i-x_j}$

$那么,由拉格朗日插值法得到的函数f的表达式为:{\sum }_{i=1}^n{b}_i{\prod }_{j\ne i}\frac{x-x_j}{x_i-x_j},时间复杂度为O(n^2).$

$若在程序实现中,每一次计算x_i-x_j值时都求一次其的乘法逆元,时间复杂度将乘上logn,我们选择先计算所有x_i-x_j的乘积,再求逆元.$

另外,此处求逆元必须选择快速幂方法,否则你将陷入深深的迷惑中.为什么呢?在此题中,可能会需要求负数的逆元,但根据辗转相除法求负数和-1的最大公约数时,求得的最大公约数为-1(在cmath库函数__gcd中计算也会得到-1,可见其运用的是辗转相除法).可是在数学定义上,gcd(a,b)=gcd(|a|,|b|),也就是说负数和-1的最大公约数为1.贝祖等式ax+by=gcd(a,b)的gcd(a,b)显然是数学定义上的.那么便不能得到ax+by=1,或者说我们算的其实是ax+by=-1,完全错误了.(我在这个问题上困扰了接近1个小时TAT)

代码:

#include <bits/stdc++.h>
#define ll long long

using namespace std;
ll n,k,x[2005],y[2005],mod=998244353,ans;//998244353是个题目常用取模质数 

ll fpm(ll x,ll power,ll mod) 
{
    x %= mod;
    ll ans = 1;
    for (; power; power >>= 1,(x*=x)%=mod)
    	if(power&1) (ans*=x)%=mod;
    return ans;
}

ll inv(ll x)
{
	return fpm(x,mod-2,mod);
}

int main()
{
	cin>>n>>k;
	for(int i=1;i<=n;i++)
	{
		cin>>x[i]>>y[i];
	}
	for(int i=1;i<=n;i++)
	{
		ll ans_up=y[i],ans_down=1;
		for(int j=1;j<=n;j++)
		{
			if(j==i) continue;
			ans_up=(ans_up*(k-x[j]))%mod;
			ans_down=(ans_down*(x[i]-x[j]))%mod;
		}
		ans=(ans+(ans_up*inv(ans_down)))%mod;
	}
	cout<<(ans+mod)%mod<<endl;
	
	return 0;
}