本蒟蒻看到一道数学题，就顺手切了。感觉单单对这一题而言，部分评论区的大佬过于复杂了
【分析】
先讲讲拉格朗日插值法：
对于给定的 ((n+1)) 个点，我们可以确定唯一的一个 至多(n)次 函数 (f(x))
这里简单解释一下为什么是唯一的至多 (n) 次：
若给定 ((n+1)) 个点，理论上可以确定一个 (n) 次函数 (f(x)) 。但 (f(x)) 可以前几项系数为 (0) 。因此，本蒟蒻按照数学上的定义，认为它是至多 (n) 次的。
若给定 ((n+1)) 个点，肯定还能给出其他的 (k) 个点使得函数 (f(x)) 变为至多 ((n+k)) 次。由于其他 (k) 个点未知，故原本的这些点能确定无数个更高次数的函数。
拉格朗日插值法的优秀就在于，它能不需要用高斯消元法，就能求出这个至多 (n) 次函数
按照朴素思路，我们是构造一个矩阵，然后高斯消元法 (O(n^3)) ，肯定是不够的 （P.S.因为本蒟蒻实在太菜，不懂有没有更优的算法优化高斯消元法）。
拉格朗日插值法的思路很简单，按照题目样例一给的三个点：
1 4
2 9
3 16
我们假设有以下 (3) 个函数：
(G_1(x)=(x-2)(x-3))
(G_2(x)=(x-1)(x-3))
(G_3(x)=(x-1)(x-2))
因此，显然有
(G_1(2)=G_1(3)=0 , G_2(1)=G_2(3)=0 , G_3(1)=G_3(2)=0)
所以我们需要的原函数可以写为：
(f(x)=A_1G_1(x)+A_2G_2(x)+A_3G_3(x))
那就很显然有
(egin{cases}f(1)=A_1G_1(1)+A_2 imes0+A_3 imes 0=A_1G_1(1)\f(2)=A_1 imes 0+A_2G_2(2)+A_3 imes 0=A_2G_2(2)\f(3)=A_1 imes 0+A_2 imes 0+A_3G_3(3)=A_3G_3(3)end{cases})
再根据题意和计算：
(egin{cases}f(1)=4,G_1(1)=2\f(2)=9,G_2(2)=-1\f(3)=16,G_3(3)=2end{cases})
倒代回去可以解得：
(f(x)=2(x-2)(x-3)-9(x-1)(x-3)+8(x-1)(x-2)=(x+1)^2)
当然，你会发现，这样是在代入 (k) 只需要 (O(n)) ，但这么做并没有更简单。
当然啦 ，但是我们求的不是 (f(x)) ，只是 (f(k)) 啊！
所以我们先提取一下上面的思路：
给定 (n) 个点 ((x_1,y_1),(x_2,y_2),(x_3,y_3)dots(x_n,y_n))
我们先设定一个 (displaystyle G_i(x)=prod_{j=1}^{i-1}(x-x_j)cdotprod_{j=i+1}^n(x-x_j))
即 (displaystyle G_i(x)=prod_{j=1&j eq i}^n(x-x_j))
再另 (displaystyle g(x)=prod_{i=1}^n(x-x_i))
因此 (forall iin[1,n]igcap Z,g(x_i)=0)
所以有 (G_i(k)(k-x_i)=g(k))
且 (forall j eq i,G_i(x_j)=0)
由上述式子设 (displaystyle f(x)=sum_{i=1}^nA_iG_i(x))
倒代入点的坐标得：
(displaystyle y_i=f(x_i)=sum_{i=1}^nA_iG_i(x_i)=A_iG_i(x_i))
因此 (A_i=y_i[G_i(x_i)]^{-1})
所以代入原式得：
(displaystyle f(x)=sum_{i=1}^ny_i[G_i(x_i)]^{-1}G_i(x))
我们分类讨论一下：
-
若 (exists iin [1,n]igcap Z,k=x_i) ， (f(k)=f(x_i)=y_i) ，算都不用算
-
若 (forall iin[1,n]igcap Z,k eq x_i)
(g(k) eq0)
由 (G_i(k)(k-x_i)=g(k)) 得
(G_i(k)=g(k)(k-x_i)^{-1})
代入原式：
(displaystyle f(k)=sum_{i=1}^ny_i[G_i(x_i)]^{-1}G_i(k)=sum_{i=1}^ny_i[G_i(x_i)]^{-1}g(k)(k-x_i)^{-1})
看起来更复杂了
但是，我们可以注意到，(g(k)) 对于 (i) 而言，是定值
还有 ([G_i(x_i)]^{-1}) 与 ((k-x_i)^{-1}) 肯定是可以合并的
所以我们按这个方法，整理原式得：
(displaystyle f(k)=g(k)cdotsum_{i=1}^ny_i[G_i(x_i)(k-x_i)]^{-1})
另 (Inv(x)) 表示数 (x) 在模 (998244353) 下的逆元
就有： (displaystyle f(k)=g(k)cdotsum_{i=1}^ny_icdot Inv[G_i(x_i)(k-x_i)])
至于 (g(k)) 是一个 (O(n)) 扫过去的事
(G_i(x_i)) 也是一个 (O(n)) 扫过去的事
乘一下，用 exgcd 或者 欧拉定理加快速幂 可以一个 (O(log n)) 求出逆元
所以求和符号内的复杂度为 (O(n)+O(log n)=O(n))
套个求和符号就是 (O(n imes n)=O(n^2))
复杂度上肯定是过得了了
总结一下，就是：
(f(k)=egin{cases}y_i,k=x_i\g(k)sum_{i=1}^ny_iInv[G_i(x_i)(k-x_i)]end{cases})
那求逆元本蒟蒻是用了非递归的 exgcd ，众位大犇可以直接看本蒟蒻的代码
【代码】
那本蒟蒻就放 我码风极丑的 代码了
#include<cstdio>
using namespace std;
#define f(a,b,c,d) for(register int a=b,c=d;a<=c;a++)
#define g(a,b,c,d) for(register int a=b,c=d;a>=c;a--)
typedef int i32;
typedef long long i64;
const i64 Mod=998244353;
inline char gc(){
static char s[1<<20|1],*p1=s,*p2=s;
return (p1==p2)&&(p2=(p1=s)+fread(s,1,1<<20,stdin),p1==p2)?EOF:*(p1++);
}
inline i64 read(){
register i64 ans=0;register bool neg=0;register char c=gc();
while(c<48||c>57) neg^=!(c^'-'),c=gc();
while(c>=48&&c<=57) ans=(ans<<3)+(ans<<1)+(c^48),c=gc();
return neg?-ans:ans;
}
inline void input(i64 &lld_X,i64 &lld_Y,i64 &lld_Bas,i64 lld_K){
lld_X=read();
lld_Y=read();
lld_K-=lld_X;
if(lld_K<0) lld_K+=Mod;
if(lld_K<0) lld_K+=Mod;
lld_Bas*=lld_K;
lld_Bas%=Mod;
}
inline i64 Inv(i64 lld_A){
register i64 lld_Tmp[10000]={0},lld_Cur=1;
register i64 lld_X=1,lld_Y=0,lld_B=Mod;
lld_Tmp[0]=lld_A/lld_B;
while(lld_B^=lld_A^=lld_B^=lld_A%=lld_B)
lld_Tmp[lld_Cur++]=lld_A/lld_B;
while(lld_Cur--)
lld_Y^=lld_X^=lld_Y^=lld_X-=lld_Tmp[lld_Cur]*lld_Y;
lld_X%=Mod;
if(lld_X<0) lld_X+=Mod;
return lld_X;
}
inline i64 calc(i64 lld_X[],i64 lld_Y,i32 d_Pos,i32 d_N,i64 lld_K){
i64 lld_Ans=1;
f(i,1,I,d_Pos-1)
lld_Ans=lld_Ans*(lld_X[d_Pos]-lld_X[i])%Mod;
f(i,d_Pos+1,I,d_N)
lld_Ans=lld_Ans*(lld_X[d_Pos]-lld_X[i])%Mod;
lld_Ans=lld_Ans*(lld_K-lld_X[d_Pos])%Mod;
if(lld_Ans<0) lld_Ans+=Mod;
return Inv(lld_Ans)*lld_Y%Mod;
}
inline void print(i64 lld_Ans){
char s[20]={0};
fwrite(s,1,sprintf(s,"%lld",lld_Ans),stdout);
}
i32 main(){
i32 d_N=read();
i64 lld_K=read(),lld_X[2048]={0},lld_Y[2048]={0},lld_Bas=1,lld_Ans=0;
f(i,1,I,d_N){
input(lld_X[i],lld_Y[i],lld_Bas,lld_K);
if(lld_X[i]==lld_K){
print(lld_Y[i]);
return 0;
}
}
f(i,1,I,d_N){
lld_Ans+=calc(lld_X,lld_Y[i],i,d_N,lld_K);
if(lld_Ans>=Mod) lld_Ans-=Mod;
}
lld_Ans=lld_Ans*lld_Bas%Mod;
print(lld_Ans);
return 0;
}
最后安利一下 本蒟蒻的博客