HDU 6395 Sequence （矩阵快速幂+分块）

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=6395

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

#define debug puts("YES");
#define rep(x,y,z) for(int (x)=(y);(x)<(z);(x)++)
#define ll long long

#define lrt int l,int r,int rt
#define lson l,mid,rt<<1
#define rson mid+1,r,rt<<1|1
#define root l,r,rt
#define mst(a,b) memset((a),(b),sizeof(a))
const int  maxn =3;
const int mod=1e9+7;
const int ub=3;
ll powmod(ll x,ll y){ll t; for(t=1;y;y>>=1,x=x*x%mod) if(y&1) t=t*x%mod; return t;}
ll gcd(ll x,ll y){return y?gcd(y,x%y):x;}
/*
矩阵快速幂+分块
*/

ll a,b,c,d,p,n;

struct mat
{
    ll a[maxn][maxn];
    mat()
    {
        memset(a,0,sizeof(a));
        for(int i=0;i<maxn;i++) a[i][i]=1;
    }
    mat operator*(const mat& y ) const
    {
        mat ret;
        for(int i=0;i<maxn;i++) for(int j=0;j<maxn;j++)
        {
            ret.a[i][j]=0;
            for(int k=0;k<maxn;k++) (ret.a[i][j]+=a[i][k]*y.a[k][j]%mod)%=mod;
        }
        return ret;
    }
};

mat quick_pow(mat x,ll n)
{
    mat ret;
    for(;n;n>>=1,x=x*x) if(n&1) ret=x*ret;
    return ret;
}

int main()
{
    int t;scanf("%d",&t);
    while(t--)
    {
        scanf("%lld%lld%lld%lld%lld%lld",&a,&b,&c,&d,&p,&n);
        if(n==1) printf("%lld\n",a);
        else if(n==2) printf("%lld\n",b);
        else
        {
            mat x,y;memset(x.a,0,sizeof(x.a));
            x.a[0][0]=d,x.a[0][1]=c,x.a[1][0]=1,x.a[2][2]=1;
            int i,j;
            for( i=3;i<=min(n,p);i=j+1)
            {
                j=min(min(n,p),p/(p/i)),x.a[0][2]=p/i;
                y=quick_pow(x,j-i+1)*y;
            }
            x.a[0][2]=0;y=quick_pow(x,max(0LL,n-i+1))*y;
            ll ans=b*y.a[0][0]%mod+a*y.a[0][1]%mod;
            ans%=mod;
            ans=(ans+y.a[0][2])%mod;
            printf("%lld\n",ans);
        }
    }
    return 0;
}