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*BX 1 0 0 0

   0 1 0 0 0

       AX*BY 0 AX 0 0

       AY*BX 0 0 BX 0

       AY*BY 0 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 1000000007

LL 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;
}