HDU 4291 (矩阵快速幂+循环节)

A Short problem
Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
Total Submission(s): 1026    Accepted Submission(s): 407

Problem Description
According to a research, VIM users tend to have shorter fingers, compared with Emacs users.
Hence they prefer problems short, too. Here is a short one:
Given n (1 <= n <= 1018), You should solve for  g(g(g(n))) mod 109 + 7
where
g(n) = 3g(n - 1) + g(n - 2)
g(1) = 1
g(0) = 0

Input
There are several test cases. For each test case there is an integer n in a single line.
Please process until EOF (End Of File).

Output

　　For each test case, please print a single line with a integer, the corresponding answer to this case.

Sample Input
0
1
2
 
Sample Output
0
1

42837

循环节：对于连续的模操作，总会出现循环的情况。在这道题里，可以首先暴力求出每一层中的循环节的大小。代码如下：

/*
int main()
{
    Int g0=0,g1=1,tem=-1;
    for(Int i=1;;i++)
    {
        tem = (3*g1%mod+g0%mod)%mod;
        g0=g1;
        g1=tem;
        if(g0 == 0 && g1 == 1)
        {
            cout << i << endl;
            break;
        }
    }
    return 0;
}
*/
//1000000007
//222222224
//183120

即，对于外层要模1000000007，求得循环节是222222224，也就是在第二层里每222222224次操作就回出现循环。同理，第三层根据222222224可以求得循环节是183120.

接下来就是进行三层的矩阵快速幂，因为g(n)=3*g(n-1)+g(n-2)，故根据矩阵乘法运算可设2*2矩阵a为3，1；1，0。很明显a*a=10,3;3,1。a[0][0] = g(3),a[0][1]=a[1][0]=g(2),a[1][1]=g(1).

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define mod1 1000000007LL
#define mod2 222222224LL
#define mod3 183120LL
using namespace std;

typedef __int64 Int;

struct Matrix
{
    Int mt[2][2];
};

Matrix cal(Matrix m,Matrix w,Int mod)
{
    Matrix ret;
    for(int i=0; i<=1; i++)
        for(int j=0; j<=1; j++)
            ret.mt[i][j] = 0;
    for(int i=0; i<=1; i++)
    {
        for(int j=0; j<=1; j++)
        {
            for(int k=0; k<=1; k++)
            {
                ret.mt[i][j] += m.mt[i][k]%mod*w.mt[k][j];
                ret.mt[i][j] %= mod;
            }
        }
    }
    return ret;
}

Int fastPower(Matrix a,Int tim,Int mod)
{
    Matrix ret;
    ret.mt[0][0] = 1,ret.mt[0][1]=ret.mt[1][0] = 0,ret.mt[1][1] = 1;
    while(tim)
    {
        if(tim & 1)
            ret = cal(ret,a,mod);
        tim >>= 1;
        a = cal(a,a,mod);
    }
    return ret.mt[0][1];
}

int main()
{
    Int n;
    while(cin >> n)
    {
        if(n < 2)
        {
            printf("%I64d\n",n);
            continue;
        }
        Matrix cst;
        cst.mt[0][0] = 3,cst.mt[0][1]=cst.mt[1][0] = 1,cst.mt[1][1] = 0;
        Int ans = fastPower(cst,fastPower(cst,fastPower(cst,n,mod3),mod2),mod1);
        cout << ans << endl;
    }
}