并不对劲的拉格朗日插值

Posted 2022-12-05 xzyf

题目大意

给\(n\)（\(n\leq 2000\)）个点\((x_1,y_1),...,(x_n,y_n)\)(\(x,y\leq 998244353\))，求多项式\(f(x)\)使\(\forall i\in [1,n],f(x_i) mod 998244353=y_i\)。

题解

结论：\(f(x)=\sum\limits_i=1^ny_i\frac\prod\limits_j=1,j\neq i^n(x-x_j)\prod\limits_j=1,j\neq i^n(x_i-x_j)\)
该式可以保证当\(x=x_k\)时，\(i=k\) 时 \(\frac\prod\limits_j=1,j\neq i^n(x-x_j)\prod\limits_j=1,j\neq i^n(x_i-x_j)=1\)，\(i\neq k\) 时 \(\frac\prod\limits_j=1,j\neq i^n(x-x_j)\prod\limits_j=1,j\neq i^n(x_i-x_j)=0\)，符合题目要求。
考虑怎么求\(f(x)\)每一项系数。
\(f(x)=\sum\limits_i=1^n\fracy_i\prod\limits_j=1,j\neq i^n(x_i-x_j)\times(\prod\limits_j=1,j\neq i^n(x-x_j))\)
前一部分可以直接求，对于后一部分\(\prod\limits_j=1,j\neq i^n(x-x_j)\)：
先求出\(g(x)=\prod\limits_j=1^n(x-x_j)\)，即求出\((x-x_1)(x-x_2)(x-x_3)...(x-x_n)\)，\(k\)次项系数相当于在\((-x_1),(-x_2),...,(-x_n)\)中选\(n-k\)个数的所有方案的积的和，这部分可以设\(dp(i,j)\)表示考虑前\(i\)个数，选了\(j\)个数；
接下来考虑如何从\(g(x)\)中去掉\((x-x_i)\)得到\(\prod\limits_j=1j\neq i^n(x-x_j)\)，相当于已知\(g(x),x_i\)，求\(h(x)\)使\(h(x)\times(x-x_i)=g(x)\)，设\(g(x)=\sum\limits_j=0^n a_j\times x^j\)，\(h(x)=\sum\limits_j=0^n-1b_j\times x^j\)，就有\(\sum\limits_j=0^n-1b_j\times x^j+1-\sum\limits_j=0^n-1b_j\times x_j\times x^j=\sum\limits_j=0^n a_j\times x^j\)，该式无论\(x\)取何值恒成立，就可以得出\(x\)对应次项的系数相等，即\(b_0=-\fraca_0x_i\)，\(\forall j\in[1,n],b_j=-\fraca_j-b_j-1x_i\)。
时间复杂度\(\Theta(n^2)\)。

代码

#include<algorithm>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<ctime>
#include<iomanip>
#include<iostream>
#include<map>
#include<queue>
#include<set>
#include<stack>
#include<vector>
#define rep(i,x,y) for(register int i=(x);i<=(y);++i)
#define dwn(i,x,y) for(register int i=(x);i>=(y);--i)
#define view(u,k) for(int k=fir[u];~k;k=nxt[k])
#define LL long long
#define maxn 2007
using namespace std;
int read()

    int x=0,f=1;char ch=getchar();
    while(!isdigit(ch)&&ch!='-')ch=getchar();
    if(ch=='-')f=-1,ch=getchar();
    while(isdigit(ch))x=(x<<1)+(x<<3)+ch-'0',ch=getchar();
    return x*f;

void write(int x)

    if(x==0)putchar('0'),putchar('\n');return;
    int f=0;char ch[20];
    if(x<0)putchar('-'),x=-x;
    while(x)ch[++f]=x%10+'0',x/=10;
    while(f)putchar(ch[f--]);
    putchar('\n');
    return;

const LL mod=998244353;
int qx[maxn],qy[maxn],a[maxn],n,K,f[maxn];
int mo(int x)return x>=mod?x-mod:x;
int mul(int x,int y)int res=1;while(y)if(y&1)res=(LL)res*x%mod;x=(LL)x*x%mod,y>>=1;return res;
int getf(int x)

    int ans=0,now=1;
    rep(i,0,n-1)ans=mo(ans+(LL)now*a[i]%mod),now=(LL)now*x%mod;
    return ans;

int main()

    n=read(),K=read();
    rep(i,1,n)qx[i]=read(),qy[i]=read();
    f[0]=1;
    rep(i,1,n)dwn(j,i,1)
        f[j]=mo(f[j]+(LL)f[j-1]*(mod-qx[i])%mod);
    reverse(f,f+n+1);
    rep(i,1,n)
    
        int num=1,nyx=mul(mod-qx[i],mod-2);
        rep(j,1,n)if(j!=i)num=(LL)num*mo(qx[i]-qx[j]+mod)%mod;
        num=(LL)mul(num,mod-2)*qy[i]%mod;
        int lst=0;
        rep(j,0,n-1)
        
            lst=(LL)mo(f[j]-lst+mod)*nyx%mod,a[j]=mo(a[j]+(LL)lst*num%mod);
        
    
    write(getf(K));
    return 0;

一些感想

终于填上这个坑了。。。