HDU 4684 - Arc of Dream（矩阵快速幂）

Problem Description

An Arc of Dream is a curve defined by following function:

where
a ₀ = A0
a _i = a _i-1*AX+AY
b ₀ = B0
b _i = b _i-1*BX+BY
What is the value of AoD(N) modulo 1,000,000,007?

 

Input

There are multiple test cases. Process to the End of File.
Each test case contains 7 nonnegative integers as follows:
N
A0 AX AY
B0 BX BY
N is no more than 10 ¹⁸, and all the other integers are no more than 2×10 ⁹.

 

Output

For each test case, output AoD(N) modulo 1,000,000,007.

 

Sample Input

  
  

   
   1
1 2 3
4 5 6
2
1 2 3
4 5 6
3
1 2 3
4 5 6
  
  

 

Sample Output

  
  

   
   4
134
1902
  
  

                                                              

思路：

题目有一些关系式，变量与前一个变量有关系，又n很大，所以可以判定是矩阵快速幂。

An = axAn-1 + ay;

Bn = bxBn-1 + by;

An*Bn = (axAn-1 + ay)(bxBn-1 + by) (展开可以发现与An-1,Bn-1的关系)；

Sn = Sn-1 + An*Bn.

技巧：将有关系的变量列成一列，根据关系式得到矩阵。

CODE：

#include <iostream>
#include <cstdio>
#include <vector>

using namespace std;
typedef long long ll;
typedef vector<ll> vec;
typedef vector<vec> mat;
const ll M = 1000000007;

mat mul(mat &A, mat &B)
{
    mat C(A.size(), vec(B[0].size()));
    for(int i = 0; i < A.size(); ++i) {
        for(int k = 0; k < B.size(); ++k) {
            for(int j = 0; j < B[0].size(); ++j) {
                C[i][j] = (C[i][j] + A[i][k] * B[k][j]) % M;
            }
        }
    }
    return C;
}
mat pow(mat A, ll n)
{
    mat B(A.size(), vec(A.size()));
    for(int i = 0; i < A.size(); ++i) {
        B[i][i] = 1;
    }
    while(n > 0) {
        if(n & 1) B = mul(B, A);
        A = mul(A, A);
        n >>= 1;
    }
    return B;
}
int main()
{
//freopen("in", "r", stdin);
    ll N, a, b;
    ll ax, ay, bx, by;
    while(~scanf("%I64d", &N)) {
        scanf("%I64d %I64d %I64d", &a, &ax, &ay);
        scanf("%I64d %I64d %I64d", &b, &bx, &by);
        mat A(5, vec(5));
        for(int i = 0; i < 5; ++i) {
            for(int j = 0; j < 5; ++j) {
                A[i][j] = 0;
            }
        }
        A[0][0] = ax; A[0][4] = ay;
        A[1][1] = bx; A[1][4] = by;
        A[2][0] = ax*by%M; A[2][1] = ay*bx%M; A[2][2] = ax*bx%M; A[2][4] = ay*by%M; //a和b相乘的时候要模，不然会爆！！
        A[3][2] = 1; A[3][3] = 1;
        A[4][4] = 1;
        A = pow(A, N);
        ll tmp[] = {a, b, a*b%M, 0, 1}; //这里也要记得模。
        ll ans = 0;
        for(int i = 0; i < 5; ++i) {
            ans = (ans + A[3][i] * tmp[i]) % M;
        }
        printf("%lld\n", ans);
    }
    return 0;
}