拉格朗日插值

定义

一般地，若已知 (y=f(x))， 在互不相同 (n+1) 个点：(x_0,x_1,dots x_n) 处的函数值：(y_0,y_1,dots y_n) ，即该函数过点：((x_0,y_0),(x_1,y_1),dots (x_n,y_n)) ，则可以考虑构造一个过这 (n+1) 个点的，最高次幂不超过 (n) 的多项式：

[y=P_n(x) ]

使其满足：

[P_n(x_k)=y_k,k=0,1,dots ,n ]

定理

满足插值条件的，次数不超过 (n) 的多项式是存在而且是唯一的。

分析

可以理解为，将该多项式 (P_n(x)) 拆成 (n+1) 个函数的和的形式：

[P_n(x)=sum_{i=0}^{n}{f_i(x)} ]

其中，第 (i) 个函数 (f_i(x)) 满足除 (f_i(x_i)=y_i) 外，其它位置取值都为 (0)。即：

[f_i(x)=a_i prod_{i eq j}{(x-x_j)} ]

代入 (f(x_i)=y_i)，可以得出：

[a_i=frac{y_i}{prod_{i eq j}{(x_i-x_j)}} ]

即插值公式为：

[egin{align} P_n(x) &=sum_{i=0}^{n}{frac{y_i·prod_{i eq j}{(x-x_j)}}{prod_{i eq j}{(x_i-x_j)}}}\ end{align} ]

另外，还有拉格朗日插值的优化：重心拉格朗日。
https://www.cnblogs.com/ECJTUACM-873284962/p/6833391.html

Polynomial-2019南昌邀请赛

题意

定义(f(x)=a_0+a_1 imes x^1+a_2 imes x^2+dots +a_n imes x^n)，现在已知 (f_{(0)},f_{(1)},f_{(2)},dots ,f_{(n)})，现在要求：

[sum_{i=L}^{R}f_{(i)} mod 9999991 ]

题目链接：https://nanti.jisuanke.com/t/40254

分析

根据拉格朗日插值的公式，可以得出 (f(x)) 的通项：

[f(x)=sum_{i=0}^{n}{frac{y_i}{prod_{i eq j}{(i-j)}}·prod_{i eq j}{(x-j)}} ]

其中常数项部分，可以通过通过阶乘和逆元预处理出来，注意 (-1)。而后面与 (x) 有关的项，可以代入 (x) 通过前缀积和后缀积求出。

但是如果要求一段区间内的函数值之和，可以考虑求出前缀和。因为 (n) 次多项式的前缀和是 (n+1) 次多项式，且前缀和需要 (n+1) 个点来求出。因此，可以利用 (f(x)) 求出 (f(n+1)) ，从而求出前缀和的多项式在各点的取值，进而求出前缀和的多项式，求差即可。

代码

#include<bits/stdc++.h>
//先求Pn的多项式，然后求其前缀和的多项式
using namespace std;
typedef long long ll;
const int mod=9999991;
const int N=1e3+10;
const int M=1e7+5;
ll fact[N],sum[M],inv[N],pre[N],suf[N],f[N];
ll power(ll x,ll y)
{
    ll res=1;
    x%=mod;
    while(y)
    {
        if(y&1) res=res*x%mod;
        x=x*x%mod;
        y>>=1;
    }
    return res;
}
void init()
{
    int maxn=1000;
    fact[0]=1;
    for(int i=1;i<=maxn;i++)
        fact[i]=1LL*fact[i-1]*i%mod;
    inv[maxn]=power(fact[maxn],mod-2);
    for(int i=maxn-1;i>=0;i--)
        inv[i]=1LL*inv[i+1]*(i+1)%mod;
}
ll solve(ll h[],int x,int n)//求出对应的函数值
{
    if(x<=n) return h[x];
    ll res=0;
    pre[0]=x,suf[n]=(x-n);
    for(int j=1;j<=n;j++)
        pre[j]=1LL*pre[j-1]*(x-j)%mod;
    for(int j=n-1;j>=0;j--)
        suf[j]=1LL*suf[j+1]*(x-j)%mod;
    for(int i=0;i<=n;i++)
    {
        ll tmp=h[i];
        tmp=tmp*inv[n-i]%mod*inv[i]%mod;
        int flag=((n-i)%2?-1:1);
        tmp=tmp*flag%mod;
        if(i>0) tmp=tmp*pre[i-1]%mod;
        if(i<n) tmp=tmp*suf[i+1]%mod;
        res=(res+tmp+mod)%mod;
    }
    return res;
}
int main()
{
    int T,L,R,n,m;
    scanf("%d",&T);
    init();
    while(T--)
    {
        scanf("%d%d",&n,&m);
        for(int i=0;i<=n;i++)
            scanf("%lld",&f[i]);
        f[n+1]=solve(f,n+1,n);
        sum[0]=f[0];
        for(int i=1;i<=n+1;i++)
            sum[i]=(sum[i-1]+f[i])%mod;
        while(m--)
        {
            scanf("%d%d",&L,&R);
            printf("%lld
",(solve(sum,R,n+1)-solve(sum,L-1,n+1)+mod)%mod);
        }
    }
    return 0;
}