HDU 4686 Arc of Dream

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

题意：fn = ∑aibi  ( 0<=i<= n-1 )

a 0 = A0 

a i = a i-1*AX+AY 

b 0 = B0 

b i = b i-1*BX+BY

给出n，求fn对1,000,000,007的模

思路：写出相邻项的转移式，

f(n+1) = f(n) + a(n) * b(n)

a(n) * b(n) = ( a(n-1) * AX + AY ) * ( b(n-1) * BX + BY ) = AX*BX * a(n-1) * b(n-1)  +  AX*BY * a(n-1)  +  AY*BX * b(n-1)  +  AY*BY

AX，AY，BX，BY都是常数，不用管。

可以弄一个5*5的转移矩阵 

ai*bi   fi   ai   bi   1 * {   AX*BX1 0 0 0

   01 0 0 0

       AX*BY0 AX 0 0

       AY*BX0 0 BX 0

       AY*BY0 AY BY 1   }

#include <cstdio>#include <cmath>#include <cstring>#include <string>#include <cstdlib>#include <iostream>#include <algorithm>#include <stack>#include <map>#include <set>#include <vector>#include <sstream>#include <queue>#include <utility>using namespace std;#define rep(i,j,k) for (int i=j;i<=k;i++)#define Rrep(i,j,k) for (int i=j;i>=k;i--)#define Clean(x,y) memset(x,y,sizeof(x))#define LL long long#define ULL unsigned long long#define inf 0x7fffffff#define mod 1000000007LL n;LL a0,b0,ax,bx,ay,by;struct node{    LL a[5][5];    void P()    {        Clean(a,0);        a[0][0] = ax*bx % mod;        a[1][0] = ax*by % mod;        a[2][0] = ay*bx % mod;        a[4][0] = ay*by % mod;        a[1][1] = ax % mod;        a[2][2] = bx % mod;        a[4][1] = ay % mod;        a[4][2] = by % mod;        a[0][3] = a[3][3] = a[4][4] = 1;    }    void E()    {        Clean(a,0);        rep(i,0,4) a[i][i] = 1;    }};node multi( node &x , node &y ){    node ans;    rep(i,0,4)        rep(j,0,4)        {            ans.a[i][j] = 0;            rep(k,0,4)                ans.a[i][j] = ( ans.a[i][j] + x.a[i][k] * y.a[k][j] % mod ) % mod;        }    return ans;}int main(){    while( scanf("%I64d",&n) == 1 )    {        scanf("%I64d %I64d %I64d",&a0,&ax,&ay);        scanf("%I64d %I64d %I64d",&b0,&bx,&by);        node ans,temp;        ans.E() , temp.P();        while( n )        {            if( n & 1 ) ans = multi( ans , temp );            temp = multi( temp , temp );            n >>= 1;        }        LL S = 0;        LL x[] = { a0%mod*b0%mod , a0%mod , b0%mod , 0 , 1 };        rep(i,0,4)            S = ( S + x[i] * ans.a[i][3] % mod ) % mod;        printf("%I64d\n",S);    }    return 0;}