多项式插值学习记录

拉格朗日插值公式:∑yi*(x-xj)/(xi-xj)    i=0..k,j=0..k且j!=i

我们可以这么记忆:这个公式的第i项可以满足:把xi带进去等于yi,把其他的x带进去=0

优点:逼格高,好用

牛顿插值:定义阶差f(x0)=y0,f(x0,x1)=(f(x0)-f(x1))/(x0-x1),f(x0,x1,x2)=(f(x0,x1)-f(x1,x2))/(x0-x2) .........

然后f(x)=f(x0)+f(x0,x1)*(x-x0)+f(x0,x1,x2)*(x-x0)*(x-x1)+.......

我们可以这么记忆:每个阶差就用两个小一点的区间求差....然后除以横坐标的差.....就和算斜率差不多...

优点:支持在线,动态插点只要O(n)

每种插值给出一道例题:

1.Codeforces Educational Round 7 F. The sum of the k-th powers

题意:给定n,k 求∑i^k 对10^9+7取模 ( 1<=i<=n )     n<=10^9  k<=10^6

以前见过另一道一模一样的题 但那题n,k是long long范围,模数只有10^6,因为i^k=(i+mod)^k,所以可以考虑直接暴力模数范围内的值。

这道题模数很大,但是k比较小,所以我们考虑一个多项式做法。

令f(x)=题目所求的∑i^k  那么f(x+1)-f(x)=(x+1)^k,显然是一个k次多项式,

那么根据结论:k+1次多项式差分后是k次多项式,我们可以确定f(x)是一个k+1次多项式。

所以可以直接插值,f(x)=∑yi*∏(n-j)/∏(i-j) i=0..k+1,j!=i

只要暴力算出前k+1个的值,然后求和就可以了。其中两个连乘可以ok计算,yi用快速幂logk计算,逆元的复杂度是log mod所以总复杂度大概是klogmod+klogk

#include<iostream>
#include<cstdio>
#define ll long long
#define mod 1000000007
using namespace std;
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-') f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=x*10+ch-'0'; ch=getchar();}
    return x*f;
}

ll pow(ll x,int p)
{
    ll sum=1;
    for(ll i=x;p;p>>=1,i=(i*i)%mod)if(p&1)
        sum=sum*i%mod;
    return sum;
}

int n;
ll ans=0,K,y=0;

void solve()
{
    ll sum=0;
    for(int i=1;i<=n;i++)
        sum=(sum+pow(i,K))%mod;
    printf("%d",(int)sum);
}

int main()
{
    n=read();K=read();
    if(n<=K+1){solve();return 0;}
    ll sum1=n,sum2=1;
    for(int i=1;i<=K+1;i++)
        sum1=sum1*(n-i)%mod,sum2=sum2*(mod-i)%mod;
    cout<<sum2<<endl;
    for(ll i=1;i<=K+1;i++)
    {
        y=(y+pow(i,K))%mod;
        ll sum3=sum1*pow(n-i,mod-2)%mod;
        sum2=sum2*i%mod;
        sum2=sum2*pow(mod-K+i-2,mod-2)%mod;
     //   cout<<i<<" "<<y<<" "<<sum3<<" "<<sum2<<endl;
        ans=(ans+y*sum3%mod*pow(sum2,mod-2))%mod;
    }
    printf("%d\n",(int)ans);
    return 0;
}

 2.hdu 5275  Dylans loves polynomial

题意:给定n个点的坐标xy和q个询问,每次询问给出L,R,x要求求出x在 用L-R的点构成的L-R次的多项式 的值。n,q<=3000,x<=250000

题解:牛顿插值,预处理出所有阶差(线性筛逆元 r[i]=-(mod/i)*rev[mod%i]),然后直接套公式qaq

复杂度maxX+n^2

#include<iostream>
#include<cstdio>
#define ll long long
#define mod 1000000007
using namespace std;
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(ch<'0'||ch>'9'){if(ch=='-') f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=x*10+ch-'0'; ch=getchar();}
    return x*f;
}
ll f[3005][3005];
int n,m;
ll x[3005],y[3005];
ll rev[250005];

ll getrev(ll x)
{
    return x<0?-rev[-x]:rev[x];
}

int main()
{
    rev[1]=1;
    for(int i=2;i<=250000;i++)
        rev[i]=(mod-mod/i)*rev[mod%i]%mod;
    n=read();
    for(int i=1;i<=n;i++)
    {x[i]=read();y[i]=read();}
    for(int i=1;i<=n;i++)f[i][i]=y[i];
    for(int k=1;k<n;k++)
        for(int i=1;i+k<=n;i++)
        {
            int j=i+k;
            f[i][j]=((f[i][j-1]-f[i+1][j])*(getrev(x[i]-x[j]))%mod+mod)%mod;
        }
    m=read();
    for(int i=1;i<=m;i++)
    {
        int l=read(),r=read();
        ll q=read(),sum=1,ans=0;
        for(int j=l;j<=r;j++)
        {
            ans=(ans+f[l][j]*sum)%mod;
            sum=sum*(q-x[j]+mod)%mod;
        }
        printf("%lld\n",ans);
    }
    return 0;
}

 

posted @ 2017-03-06 21:42  FallDream  阅读(278)  评论(0编辑  收藏  举报