P1397 [NOI2013]矩阵游戏（递推）

P1397 [NOI2013]矩阵游戏

一波化式子，$f[1][m]=a^{m-1}+bsum_{i=0}^{m-2}a^i$，用快速幂+逆元求等比数列可以做到$logm$

设$v=a^{m-1},k=sum_{i=0}^{m-2}a^i$

那么$f[1][m]=v+bk$

再对纵列化一波式子，$f[i][m]=f[i-1][m]*vc+bk+vd$

如果你直接上个矩乘可以拿到65的好分数

#include<iostream>
#include<cstdio>
#include<cstring>
#define ri register int
using namespace std;
typedef long long ll;
const ll P=1e9+7;
const ll W=1e9;
char q[1000005];
struct bnum{
    ll a[10005],len;
    bnum(){memset(a,0,sizeof(a));len=0;}
    void init(){
        scanf("%s",q); int z=1; len=1;
        for(ri i=strlen(q)-1;i>=0;--i){
            a[len]+=(q[i]-48)*z; z*=10;
            if(z==W) z=1,++len;
        }
        while(!a[len]&&len) --len;
    }
    ll mod(){
        ll re=0;
        for(ri i=1;i<=len;++i) re=(re*W+a[i])%P;
        return re;
    }
    void rem1(){
        --a[1];
        for(ri i=1;a[i]<0;++i) a[i]+=W,--a[i+1];
        while(!a[len]&&len) --len;
    }
    bnum div2(){
        bnum c; c.len=len; ll x=0;
        for(ri i=len;i;--i)
             x=x*W+a[i],c.a[i]=x/2,x%=2;
        while(!c.a[c.len]&&c.len) --c.len;
        return c;
    }
}n,m;
struct mat{
    ll a[2][2];
    mat(){memset(a,0,sizeof(a));}
    mat operator * (const mat &b) const{
        mat c;
        for(ri i=0;i<2;++i)
            for(ri j=0;j<2;++j)
                for(ri k=0;k<2;++k)
                    c.a[i][j]=(c.a[i][j]+a[i][k]*b.a[k][j]%P)%P;
        return c;
    }
}s,p;
ll Pow(ll x,ll y){
    ll re=1;
    for(;y;y>>=1,x=x*x%P) if(y&1) re=re*x%P;
    return re;
}
ll Pow_b(ll x,bnum y){
    ll re=1;
    for(;y.len;y=y.div2(),x=x*x%P) if(y.a[1]&1) re=re*x%P;
    return re; 
}
mat Pow_m(mat x,bnum y){
    mat re; for(ri i=0;i<2;++i) re.a[i][i]=1;
    for(;y.len;y=y.div2(),x=x*x) if(y.a[1]&1) re=re*x;
    return re;
}
int main(){
    ll a,b,c,d,v,k;
    n.init(); m.init(); n.rem1(); m.rem1();
    scanf("%lld%lld%lld%lld",&a,&b,&c,&d);
    v=Pow_b(a,m);
    if(a>1) k=(v-1+P)%P*Pow((a-1+P)%P,P-2)%P;
    else k=m.mod();
    s.a[0][0]=(b*k+v)%P; s.a[0][1]=(b*k%P+v*d%P)%P;
    p.a[0][0]=v*c%P; p.a[1][0]=p.a[1][1]=1;
    p=Pow_m(p,n); s=s*p;
    printf("%lld",s.a[0][0]);
    return 0;
}

观察发现，这个矩乘可以再化：

$f[n][m]=(v+bk)(vc)^{n-1}+(bk+vd)sum_{i=0}^{n-2}(vc)^i$

观察这个式子，复杂度主要在快速幂上，复杂度$O(logn+logm)$

考虑缩小$n,m$

假设存在：$a^n=a^{n-x}\, mod ;  p$，$p$为质数

$ herefore  a^x=1\, mod ;  p$

根据费马小定理，$x=p-1$

$ herefore  a^n=a^{n\, mod\, p-1}\, mod ;  p$

输入$n,m$时记下它们$mod\, p,p-1$的值，代入式子即可

注意等比数列公比$=1$时需要特判

#include<iostream>
#include<cstdio>
#include<cstring>
#define ri register int
using namespace std;
typedef long long ll;
const ll P=1e9+7;
ll a,b,c,d,v,k,n,m,_n,_m,a_,k_,v_;
void init(){
    char c=getchar();
    while(c<'0'||c>'9') c=getchar();
    while('0'<=c&&c<='9'){
        n=(n*10+c-48)%P;
        _n=(_n*10+c-48)%(P-1);
        c=getchar();
    }c=getchar();
    while(c<'0'||c>'9') c=getchar();
    while('0'<=c&&c<='9'){
        m=(m*10+c-48)%P;
        _m=(_m*10+c-48)%(P-1);
        c=getchar();
    }
}
ll Pow(ll x,ll y){
    ll re=1;
    for(;y;y>>=1,x=x*x%P) if(y&1) re=re*x%P;
    return re;
}
int main(){
    init(); scanf("%lld%lld%lld%lld",&a,&b,&c,&d);
    n=(n-1+P)%P; _n=(_n-1+P-1)%(P-1);
    m=(m-1+P)%P; _m=(_m-1+P-1)%(P-1);
    v=Pow(a,_m);
    k=a>1?(v-1+P)*Pow(a-1+P,P-2)%P:m;
    a_=v*c%P; v_=Pow(a_,_n);
    k_=a_>1?(v_-1+P)*Pow(a_-1+P,P-2)%P:n;
    printf("%lld",((v+b*k)%P*v_%P+(b*k+v*d)%P*k_%P)%P);
    return 0;
}