HDU 4686 Arc of Dream 矩阵快速幂,线性同余 难度:1

http://acm.hdu.edu.cn/showproblem.php?pid=4686

当看到n为小于64位整数的数字时,就应该有个感觉,acm范畴内这应该是道矩阵快速幂

Ai,Bi的递推式题目已经给出,

Ai*Bi=Ax*Bx*(Ai-1*Bi-1)+Ax*By*Ai-1+Bx*Ay*Bi-1+Ay*By

AoD(n)=AoD(n-1)+AiBi

构造向量I{AoD(i-1),Ai*Bi,Ai,Bi,1}

初始向量为I0={0,A0*B0,A0,B0,1}

构造矩阵A{

1,0,0,0,0,

1,Ax*Bx,0,0,0,
0,Ax*By,Ax,0,0,

0,Bx*Ay,0,Bx,0,

0,Ay*By,Ay,By,1,

}

则I0*(A)^n则为包含结果的向量,此时I[0]即为解答

注意:

1.据说有n=0,直接输出0

2.long long的使用

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
typedef long long ll;
const ll mod = 1e9+7;
ll n,ao,bo,ax,bx,ay,by;
void init(ll** ans,ll** base){
        fill(ans[0],ans[0]+5,0);
        fill(base[0],base[0]+25,0);
        ans[0][1]=ao*bo%mod;
        ans[0][2]=ao;
        ans[0][3]=bo;
        ans[0][4]=1;
        base[0][0]=1;
        base[1][0]=1;base[1][1]=ax*bx%mod;
        base[2][1]=ax*by%mod;base[2][2]=ax;
        base[3][1]=bx*ay%mod;base[3][3]=bx;
        base[4][1]=ay*by%mod;base[4][2]=ay;base[4][3]=by;base[4][4]=1;
}
void print(ll **t ,int m,int r){
        for(int i=0;i<m;i++){
                for(int j=0;j<r;j++){
                        printf("%I64d ",t[i][j]);
                }
                puts("");
        }
}

void multi(ll** t,ll** b,ll ** tmp,int m,int r,int s){
        fill(tmp[0],tmp[0]+m*s,0);
        for(int i=0;i<m;i++){
                for(int j=0;j<s;j++){
                        for(int k=0;k<r;k++){
                                tmp[i][j]+=(t[i][k]*b[k][j])%mod;
                                tmp[i][j]%=mod;
                        }
                }
        }
}

void qpow(ll** ans,ll ** base,ll**tmp,ll time,int r,int m){
        while(time>0){
                if(time&1){
                        multi(ans,base,tmp,r,m,m);
                        copy(tmp[0],tmp[0]+r*m,ans[0]);
                        time--;
                }
                multi(base,base,tmp,m,m,m);
                copy(tmp[0],tmp[0]+m*m,base[0]);
                time>>=1;
        }
}

int main(){
        ll ** ans = new ll*[1],**tmp=new ll*[5],** base =new ll*[5];
        ans[0]=new ll[5],tmp[0]=new ll[25],base[0]=new ll[25];
        for(int i=1;i<5;i++){
                tmp[i]=tmp[0]+5*i;
                base[i]=base[0]+5*i;
        }
        while(scanf("%I64d",&n)==1){
                scanf("%I64d%I64d%I64d",&ao,&ax,&ay);
                ao%=mod;ax%=mod;ay%=mod;
                scanf("%I64d%I64d%I64d",&bo,&bx,&by);
                bo%=mod;bx%=mod;by%=mod;
                if(n==0)puts("0");
                else {
                        init(ans,base);
                        qpow(ans,base,tmp,n,1,5);
                        printf("%I64d\n",ans[0][0]);
                }

        }
        delete ans;delete base;delete tmp;
        return 0;
}

  

posted @ 2015-02-26 17:21  雪溯  阅读(201)  评论(0编辑  收藏  举报