LGP4781 [LG TPLT] 拉格朗日插值 学习笔记

LGP4718 [LG TPLT] 拉格朗日插值 学习笔记

Luogu Link

前言

服了。复习拉格朗日插值于是找自己学习笔记；怎么又是未找到。

题意简述

对于 \(n+1\) 个 \(x_i\) 互不相同的点 \((x_i,y_i)\)，经过这 \(n\) 个点可唯一确定一个 \(n\) 次多项式。现在给你这样的 \(n\) 个点和一个整数 \(k\)，直接输出 \(f(k)\) 的值。答案对 \(998244353\) 取模。

\(n\le 10^3\)。

做法解析

拉插可以在 \(O(n^2)\) 时间内解决这个问题。

具体来说，这个唯一确定的多项式长成这个样子：\(f(k)=\sum_{i=0}^n y_i\prod_{i\neq j}\frac{k-x_j}{x_i-x_j}\)。正确性显然：你把 \(x_p\) 带进去，在求和式 \(p=i\) 那一项有，\(\prod\) 后的分母和分子完全相等，所以得到一个 \(y_i\)；而在 \(p\neq i\) 的项，\(\prod\) 后总有 \(j=p\) 的分式导致当前项归 \(0\)。

这道题不需要求每一项的系数，所以直接把 \(k\) 往里面代就行了。

如果要求呢？这里给出一个 \(O(n^2)\) 做法。

首先我们记 \(D(k)=\prod_{i=0}^n k-x_i\)，\(c_i=y_i\prod_{i\neq j}\frac{1}{x_i-x_j}\)，\(E_i(k)=\frac{D(k)}{k-x_i}\)。这样就有：\(f(k)=\sum_{i=0}^n E_i(k)\times c_i\)。

首先 \(D(k)\) 的系数显然可以 \(O(n^2)\) 算。\(E_i\) 的系数和 \(c_i\) 每个 \(O(n)\) 算，总共也是 \(O(n^2)\)。最后 \(f(k)\) 的第 \(i\) 项系数，不妨写作 \(f[k^i]\)，就是 \(\sum_{i=0}^n E[k^i]\times c_i\)。

以防你不知道 \(E_i\) 系数怎么算：长除法。

另外，在横坐标为连续整数时，我们还能做到 \(O(n)\) 插值。为什么呢？原柿 \(f(k)=\sum_{i=0}^n y_i\prod_{i\neq j}\frac{k-x_j}{x_i-x_j}\) 在这里可以写成 \(f(k)=\sum_{i=0}^n y_i\prod_{i\neq j}\frac{k-j}{i-j}\)，进一步地，可以写成

\[f(k)=\sum_{i=0}^n y_i\times \frac{\prod_{j=0}^n(k-j)}{(k-i)\times i!\times (n-i)!\times (-1)^{n-i}} \]

预处理各种东西之后，显然这是可以 \(O(n)\) 算的。应用参考CFP622学习笔记。

代码实现

#include <bits/stdc++.h>
using namespace std;
using namespace obasic;
using namespace omodint;
using mint=m998;
const int MaxN=2e3+5;
int N,K,X,Y;
mint uob[MaxN];
namespace omathe{
    struct spoi{int x,y;};
    mint laintp(spoi a[],int n,mint k){
        mint res=0,uprd=1;
        for(int i=0;i<=n;i++){
            uob[i]=k-a[i].x;
            uprd*=uob[i];
        }
        for(int i=0;i<=n;i++){
            mint tmp=1;
            for(int j=0;j<=n;j++){
                if(i==j)continue;
                tmp*=(a[i].x-a[j].x);
            }
            res+=a[i].y*uprd/uob[i]/tmp;
        }
        return res;
    }
};
using namespace omathe;
spoi A[MaxN];
int main(){
    readis(N,K),N--;
    for(int i=0;i<=N;i++)readis(X,Y),A[i]={X,Y};
    writi(miti(laintp(A,N,K)));
    return 0;
}