cf_566_div2_E

此文更新于2019.6.13

题意：
    f(x)=c^(2x-6)f(x-1)f(x-2)f(x-3)
    输入f(1) f(2) f(3) n c
    输出f(n)
数据范围：
    f(1)<1e9
    f(2)<1e9
    f(3)<1e9
    c<1e9
    n<1e18

法1：
令g(x)=f(x)*c^x
则g(x)=g(x-1)g(x-2)g(x-3)
令h(x)=lg(5,g(x))
则h(x)=h(x-1)+h(x-2)+g(x-3)
于是成了bsgs算法+矩阵快速幂

#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
const ll mod=1e9+7;

struct Sarray{
    static const ll mod=1e9+6;
    static const ll LEN=3;
    ll len,data[LEN][LEN];

    Sarray(ll len,ll flag):len(len){
        for(ll i=0;i<len;i++){
            for(ll j=0;j<len;j++)data[i][j]=0;
            data[i][i]=flag;
        }
    }

    Sarray operator *(const Sarray&a){
        Sarray tem(a.len,0);
        for(ll i=0;i<len;i++){
            for(ll j=0;j<len;j++){
                for(ll k=0;k<len;k++){
                    tem.data[i][j]=(tem.data[i][j]+data[i][k]*a.data[k][j])%mod;
                }
            }
        }
        return tem;
    }

    Sarray operator +(const Sarray&a){
        Sarray tem(a.len,0);
        for(ll i=0;i<len;i++){
            for(ll j=0;j<len;j++){
                tem.data[i][j]=(data[i][j]+a.data[i][j])%mod;
            }
        }
        return tem;
    }
};

Sarray qpow(Sarray a,ll b){//会更改a，不能按引用传递
    Sarray tem(a.len,1);
    while(b){
        if(b&1)tem=a*tem;
        a=a*a;
        b>>=1;
    }
    return tem;
}

ll qpow(ll a,ll b){
    ll ret=1;
    while(b){
        if(b&1) ret=ret*a%mod;
        a=a*a%mod;
        b>>=1;
    }
    return ret;
}

ll oria[3][3]={
        0,0,0,
        0,0,0,
        1,0,0
};

ll orib[3][3]={
        0,0,0,
        1,0,0,
        0,0,0
};

ll oric[3][3]={
        1,0,0,
        0,0,0,
        0,0,0
};

ll transs[3][3]={
        1,1,1,
        1,0,0,
        0,1,0
};

ll get(ll n,ll beg[3][3]){
    Sarray orn(3,1),trans(3,1);
    memcpy(orn.data,beg,sizeof(orn.data));
    memcpy(trans.data,transs,sizeof(trans.data));
    if(n==1) return beg[2][0];
    if(n==2) return beg[1][0];
    return (qpow(trans,n-3)*orn).data[0][0];
}

int  main(){
    ll n,f1,f2,f3,c;
    cin>>n>>f1>>f2>>f3>>c;
    ll g1=f1*c%mod;
    ll g2=f2*c%mod*c%mod;
    ll g3=f3*c%mod*c%mod*c%mod;

    ll A=get(n,oria);
    ll B=get(n,orib);
    ll C=get(n,oric);

    ll gn=qpow(g1,A)*qpow(g2,B)%mod*qpow(g3,C)%mod;
    ll t=qpow(c,n);
    cout<<gn*qpow(t,mod-2)%mod<<endl;

}

法2：
令g(x)=f(x)*c^x
则g(x)=g(x-1)g(x-2)g(x-3)
可以肯定出g(x)=g(1)^a*g(2)^b*g(3)*c
于是a,b,c都是x的函数
且a(x)=a(x-1)+a(x-2)+a(x-3) bc同理

#include<bits/stdc++.h>
using namespace std;

typedef long long ll;
const int LEN=3;

void sarray_cpy(int a[][LEN],int b[][LEN],int n){
    for(int i=0;i<n;i++){// a/b可以为同一个数组
        for(int j=0;j<n;j++) b[i][j]=a[i][j];
    }
}

void sarray_mul(int a[][LEN],int b[][LEN],int ret[][LEN],int n,int mod){
    static int c[LEN][LEN];// a/b/ret可以为同一个数组
    for(int i=0;i<n;i++) {
        for(int j=0;j<n;j++){
            c[i][j]=0;
            for(int k=0;k<n;k++){
                c[i][j]=(c[i][j]+1ll*a[i][k]*b[k][j])%mod;
            }
        }
    }
    sarray_cpy(c,ret,n);
}

void sarray_qpow(int aa[][LEN],ll b,int ret[][LEN],int n,int mod){
    static int a[LEN][LEN];// aa ret可以为同一个数组
    sarray_cpy(aa,a,n);
    for(int i=0;i<n;i++){
        for(int j=0;j<n;j++) ret[i][j]=0;
        ret[i][i]=1;
    }
    while(b){
        if(b&1) sarray_mul(ret,a,ret,n,mod);
        sarray_mul(a,a,a,n,mod);
        b>>=1;
    }
}

void sarray_add(int a[][LEN],int b[][LEN],int c[][LEN],int n,int mod){
    for(int i=0;i<n;i++){// a,b,c可以为同一个数组
        for(int j=0;j<n;j++){
            c[i][j]=(a[i][j]+b[i][j])%mod;
        }
    }
}

// a^0 a^1 a^2 a^3 ... a^b
void sarray_sum(int a[][LEN],ll b,int ret[][LEN],int n,int mod){
    static int tmp[LEN][LEN];
    if(b==0) sarray_qpow(a,b,ret,n,mod);
    else{
        ll mid=(b-1)>>1;
        sarray_sum(a,mid,ret,n,mod);
        sarray_qpow(a,mid+1,tmp,n,mod);
        for(int i=0;i<n;i++) tmp[i][i]=(tmp[i][i]+1)%mod;
        sarray_mul(ret,tmp,ret,n,mod);
        if((b&1)==0) {
            sarray_mul(ret,a,ret,n,mod);
            for(int i=0;i<n;i++) ret[i][i]=(ret[i][i]+1)%mod;
        }
    }
}

int qpow(int a,int b,int mod){
    int ret=1;
    while(b){
        if(b&1)ret=1ll*ret*a%mod;
        a=1ll*a*a%mod;
        b>>=1;
    }
    return ret;
}

// a^x === b x=lg(a,b)
int bsgs_lg(int a,int b,int mod){
    map<int,int>mp;
    int sqr=sqrt(mod-1)+1;
    for(int i=0;i<sqr;i++) mp[qpow(a,i,mod)]=i; // baby step
    for(int i=0;i<mod-1;i+=sqr){ // giant step
        int tp=1ll*b*qpow(a,mod-1-i,mod)%mod; // a^(-i)
        if(mp.find(tp)!=mp.end()) return i+mp[tp];
    }
    return -1;// error
}

int main(){
    const int mod=1e9+7;
    ll n;
    int f1,f2,f3,c;
    cin>>n>>f1>>f2>>f3>>c;
    int g1=1ll*f1*c%mod;
    int g2=1ll*f2*c%mod*c%mod;
    int g3=1ll*f3*c%mod*c%mod*c%mod;

    g1=bsgs_lg(5,g1,mod);
    g2=bsgs_lg(5,g2,mod);
    g3=bsgs_lg(5,g3,mod);
    // cout<

    int orn[LEN][LEN]={
            g3,0,0,
            g2,0,0,
            g1,0,0
    };

    int trans[LEN][LEN]={
            1,1,1,
            1,0,0,
            0,1,0
    };

    int ans[LEN][LEN],t;
    if(n==1) t=g1;
    else if(n==2) t=g2;
    else if(n==3) t=g3;
    else {
        int ans1[LEN][LEN];
        int ans2[LEN][LEN];
        sarray_sum(trans,n-3,ans1,3,mod-1);
        sarray_mul(ans1,orn,ans1,3,mod-1);
        sarray_sum(trans,n-4,ans2,3,mod-1);
        sarray_mul(ans2,orn,ans2,3,mod-1);
        t=(ans1[0][0]-ans2[0][0]+mod-1)%(mod-1);

        // for(int i=0;i<3;i++)for(int j=0;j<3;j++)ans[i][j]=(ans1[i][j]-ans2[i][j]+mod)%mod;
        // sarray_mul(ans,orn,ans,3,mod-1);
        // t=ans[0][0];
    }
    int gn=qpow(5,t,mod);
    int inv=qpow(c,mod-1-n%(mod-1),mod);
    cout<<1ll*gn*inv%mod<<endl;

}

此文标签
bsgs 矩阵快速幂