zoukankan      html  css  js  c++  java
  • lagrange 插值法模板

    对于n次多项式

    [Gamma(k) = sum_{i=0}^{n}y_iprod_{j=0,{j} eq {i} }^{n}frac{k-x_j}{x_i-x_j} ]

    const int mod=1e9+7;
    template<class T>
    T qpow(T x,int k){
        T e=1;
        while(k){
            if(k&1)e=1ll*e*x%mod;
            x=1ll*x*x%mod;
            k>>=1;
        }
        return e;
    }
    template<class T>
    T lagrange(int n,int x[],int y[],T k){
        T ans=0;
        for(int i=0;i<=n;++i){
            T fz=1,fm=1;//分子分母
            for(int j=0;j<=n;++j){
                if(i==j)continue;
                fm=1ll*fm*(x[i]-x[j])%mod;
                fz=1ll*fz*(k-x[j])%mod;
            }
            ans=(1ll*ans+1ll*y[i]*fz%mod*qpow(fm,mod-2)%mod)%mod;
        }
        return (ans+mod)%mod;
    }
    

    在x取值连续时的做法

    (x_i)​ 的取值都是连续的,可以把上面的算法优化到O(n)复杂度
    首先把(x_i)换成i,新的式子为

    [f(k) = sum_{i=0}^{n} y_i prod_{i ot = j} frac{k - j}{i - j}f(k) ]

    对于 (prod_{i ot = j} ^{n}frac{k - j}{i - j}​) 可以优化
    对于分子来说,我们维护出关于k的前缀积和后缀积,也就是

    [pre_i = prod_{j = 0}^{i-1} k - j ]

    [suf_i = prod_{j = i+1}^n k - j ]

    对于分母来说,观察发现这其实就是阶乘的形式,我们用fac[i]来表示(i!)
    那么式子就变成$$f(k) = sum_{i=0}^n y_i frac{pre_i * suf_i }{fac[i] * fac[N - i]}$$​分母可能会出现符号问题,也就是说,当N - i为奇数时,分母应该取负号

    const int mod=1e9+7;
    template<class T>
    T qpow(T x,int k){
        T e=1;
        while(k){
            if(k&1)e=1ll*e*x%mod;
            x=1ll*x*x%mod;
            k>>=1;
        }
        return e;
    }
    int lagrange(int n,int x[],int y[],int k){
        int pre[n+1],suf[n+1],fac[n+1],ans=0;
        fac[0]=pre[0]=suf[n]=1;
        for(int i=1; i<=n; ++i)pre[i]=1ll*(k-x[i-1])*pre[i-1]%mod;
        for(int i=n-1;i>=0;++i)suf[i]=1ll*(k-x[i+1])*suf[i+1]%mod;
        for(int i=1;i<=n;++i)fac[i]=1ll*i*fac[i-1]%mod;
        for(int i=0;i<=n;++i){
            int fz=1ll*pre[i]*suf[i]%mod;//分子分母
            int fm=1ll*fac[i]*fac[n-i]%mod;
            if((n-i)&1)fm=(mod-fm)%mod;
            ans=(1ll*ans+1ll*y[i]*fz%mod*qpow(fm,mod-2)%mod)%mod;
        }
        return (ans+mod)%mod;
    }
    

    2019南昌邀请赛B. Polynomial

    题目链接
    Define (f_{(x)} = a_0+a_1 imes x^1 + a_2 imes x^2 +…+ a_n imes x^n)

    the number maybe tremendous huge, so take the module of 9999991.

    For a polynomial, XH doesn't know any (a_i) ​ , but he knows the value of (f_{(0)} , f_{(1)}, f_{(2)}…f_{(n)}) . He wants to calculate (displaystylesum_{i=L}^R f_{(i)}) mod 9999991.

    Input

    The first line contains a number (T ( 1 le T le 5 )) ——the number of test cases.

    For each test case. The first line contains two integers (n (1 le n le 1000)) and (m(1 le m le 2000)).

    The second line contains n+1 integers,representing the value of (f_{(0)} , f_{(1)}, f_{(2)}…f_{(n)})

    The next mm lines, each line contains two integers L, R.

    (( 1 le L le R le 9999990))

    Output

    For each test case, output m lines, each line contains one integer, the value of (isplaystylesum_{i=L}^R f_{(i)}) mod 9999991.

    样例输入

    1 
    3 2
    1 10 49 142
    6 7 
    95000 100000 
    

    样例输出

    2519 
    1895570 
    

    题意就是知道(0, (f_0)) , (1, (f_1)) ... (n, (f_n)) n+1个点,求(S(k)=sum_{i=0}^k f(i))
    由n+1个点可求得 (f(n+1))
    再由(S(k)-S(k-1)=f(k))可知,(S(k)=sum_{i=0}^k f(i))是n+1次多项式
    且我们可以求出这n+2个点,分别是
    (0, (f_0)) , (1, (f_0+f_1)), ..., (n, (f_0+f_0+... +f_n)), (n, (f_0+f_0+... +f_n+f_{n+1}))
    这样,求出(S(k)) 之后
    (S(r)-S(l-1))即为答案。
    代码:

    #include<bits/stdc++.h>
    using namespace std;
    char buf[100000],*L=buf,*R=buf;
    #define gc() L==R&&(R=(L=buf)+fread(buf,1,100000,stdin),L==R)?EOF:*L++
    template<typename T>
    inline void read(T&x) {
        int flag=x=0;
        char ch=gc();
        while (ch<'0'||ch>'9')flag|=ch=='-',ch=gc();
        while (ch>='0'&&ch<='9')x=(x<<3)+(x<<1)+ch-48,ch=gc();
        if(flag)x=-x;
    }
    #define Init(arr,val) memset(arr,val,sizeof(arr))
    const int inf=0x3f3f3f3f,mod=9999991,MAXN=1020;
    typedef long long ll;//鬼畜数据,只能全部用ll,不然一个点都过不了,浪费了了我一下午。。。。
    template<class T>
    T qpow(T x,int k) {
        T e=1;
        while(k) {
            if(k&1)e=1ll*e*x%mod;
            x=1ll*x*x%mod;
            k>>=1;
        }
        return e%mod;
    }
    ll fac[MAXN],ffa[MAXN];
    //fac是阶乘
    //ffa是分母,对于m次询问,都是相同的分母,保存起来,
    //因为S(x) 已经确定,分母里的快速幂只需求一次即可,不然坑爹数据一个点也过不了
    ll lagrange(ll n,ll x[],ll y[],ll k) {
        ll pre[MAXN],suf[MAXN],ans=0;
        fac[0]=pre[0]=suf[n]=1;
        for(int i=1; i<=n; ++i)pre[i]=1ll*(k-x[i-1])*pre[i-1]%mod;
        for(int i=n-1; i>=0; --i)suf[i]=1ll*(k-x[i+1])*suf[i+1]%mod;
        for(int i=0; i<=n; ++i) {
            ll fz=1ll*pre[i]*suf[i]%mod;
            ans=(1ll*ans+1ll*y[i]*fz%mod*ffa[i]%mod)%mod;
        }
        return (ans+mod)%mod;
    }
    
    int main() {
        fac[0]=1;
        for(int i=1; i<MAXN; ++i)fac[i]=1ll*i*fac[i-1]%mod;
        ll n,m,x[MAXN],y[MAXN];
        int t;
        read(t);
        while(t--) {
            read(n),read(m);
            for(int i=0; i<=n; ++i) {
                x[i]=i;
                read(y[i]);
            }
            for(int i=0; i<=n; ++i) {//f(x)
            	ll fm=1ll*fac[i]*fac[n-i]%mod;
             	if((n-i)&1)fm=(mod-fm)%mod;
                ffa[i]=1ll*qpow(fm,mod-2)%mod;
       		}
            y[n+1]=lagrange(n,x,y,n+1);//f(n+1)
            x[n+1]=n+1;
            for(int i=1; i<=n+1; ++i)y[i]+=y[i-1];//前缀和就是S(x)个点的y值
            for(int i=0; i<=n+1; ++i) {//m次查询,分母记录下来,求一次即可
            	ll fm=1ll*fac[i]*fac[n+1-i]%mod;
             	if((n+1-i)&1)fm=(mod-fm)%mod;
                ffa[i]=1ll*qpow(fm,mod-2)%mod;
       		}
            while(m--) {
                ll l,r;
                read(l),read(r);
                printf("%lld
    ",(lagrange(n+1,x,y,r)-lagrange(n+1,x,y,l-1)+mod)%mod);
            }
        }
        return 0;
    }
    
  • 相关阅读:
    根据人脸关键点实现平面三角剖分和最近邻搜索 ( KNN, K=1 ), opencv3.4.2, C++
    KDTree  C++实现
    python 保留小数
    Clion提示:Single-argument constructors must be marked explicitly to avoid unintentional implicit conversions 解法办法
    二叉搜索树的C++ 实现
    排列组合之组合问题 网易深度学习工程师面试题 C++ 使用10方法
    OS X 安装命令行看图工具 chafa 以及其依赖libtool
    leetcode704 C++ 72ms 二分查找
    Deep Interest Network for Click-Through Rate Prediction
    归并排序
  • 原文地址:https://www.cnblogs.com/foursmonth/p/14144986.html
Copyright © 2011-2022 走看看