多項式入門——拉格朗日插值

插值用來求解這樣一類問題：給定 \(n\) 點 \((x_i,y_i)\) 求過這些點的多項式。

1 簡介

設 \(f(x)\) 為這個多項式，我們有：

\[f(x)\equiv f(a)\bmod (x-a)\tag{1} \]

這是因為：

\[f(x)-f(a)=(a_0-a_0)+a_1(x-a)+a_2(x^2-a^2)+... \]

而後者顯然有一個根為 \(a\) ，所以 \((1)\) 式得證。

通過把這 \(n\) 個點代入我們可以得到：

\[\begin{cases} f(x)\equiv y_1\bmod x-x_1\\ ...\\ f(x)\equiv y_n\bmod x-x_n \end{cases} \]

顯然模數互質，所以我們考慮用中國剩餘定理來解決這個事情。

令 \(M=\prod_{i=1}^n(x-x_i)\) ，\(m_i=\frac{M}{x-x_i}=\prod_{i\not= j}(x-x_j)\) 。並且在模 \(x-x_i\) 的意義下，我們有：

\[m_i^{-1}=\frac{1}{\prod\limits_{i\not=j}(x_i-x_j)} \]

這是因為我們有：

\[\prod_{i\not =j}(x-x_j)\equiv \prod_{i\not =j}(x-x_j-x+x_i)=\prod_{i\not =j}(x_i-x_j) \]

所以上述結論成立。

所以我們有：

\[f(x)\equiv \sum\limits_{i=1}^ny_im_im_i^{-1}\equiv \sum\limits_{i=1}^ny_i\prod\limits_{j\not=i}\frac{x-x_j}{x_i-x_j} \]

這個東西可以在 \(n^2\)

的時間內求出。

2 例題

連結

直接模擬上面的過程即可。

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N 2010
#define M number
using namespace std;

const int INF=0x3f3f3f3f;
const ll mod=998244353;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

inline ll ksm(ll a,ll b,ll mod){
    ll res=1;
    while(b){
        if(b&1) (res*=a)%=mod;
        a=a*a%mod;
        b>>=1;
    }
    return res;
}

inline ll inv(ll a){
    return ksm(a,mod-2,mod);
}

ll n,k,x[N],y[N],ans;

int main(){
    read(n);read(k);
    for(int i=1;i<=n;i++){
        read(x[i]);read(y[i]);
    }
    for(int i=1;i<=n;i++){
        ll fenzi=1,fenmu=1;
        for(int j=1;j<=n;j++){
            if(j==i) continue;
            fenmu*=(x[i]-x[j]);fenmu%=mod;
            fenzi*=(k-x[j]);fenzi%=mod;
        }
        ans+=y[i]*fenzi%mod*inv(fenmu)%mod;ans%=mod;
    }
    printf("%lld\n",(ans%mod+mod)%mod);
    return 0;
}
 

注意：

1. 需要把分母乘出來再求逆元，這樣複雜度的瓶頸就不會是求逆元。

2. 結尾註意化為正數。