HDU 4291 A Short problem 矩阵，多重函数求循环节

题意：

　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 <= 10¹⁸), You should solve for 

g(g(g(n))) mod 10 ⁹ + 7
　　where
g(n) = 3g(n - 1) + g(n - 2)
g(1) = 1
g(0) = 0

题解：

假如最外层要求模10^9+7，则令L0=10^9+7

将最外层内部看做一个整体x，则x在L0下可以求出一个循环节L1

那么x的值对L1取模不影响结果

将次外层内部看做一个整体y，则y在L1下可以求出一个循环节L2

那么y的值对L2取模不影响结果

......

求循环节的代码：

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;

#define lint __int64
const lint mod = 1000000007;
lint L1 = 222222224, L2 = 183120, L3;

void test()
{
    printf("start\n");
    lint gn, gn1 = 1, gn2 = 0;
/*
    for(lint i = 2; ; i++)
    {
        gn = (3 * gn1 + gn2) % mod;
        if(gn == 1 && gn1 == 0)
        {
            printf("循环节L1：%I64d\n", i - 1);
            L1 = i - 1; //222222224
            break;
        }
        gn2 = gn1;
        gn1 = gn;
    }

    gn = 3, gn1 = 1, gn2 = 0;
    for(lint i = 2; ; i++)
    {
        gn = (3 * gn1 + gn2) % L1;
        if(gn == 1 && gn1 == 0)
        {
            printf("循环节L2：%I64d\n", i - 1);
            L2 = i - 1; //183120
            break;
        }
        gn2 = gn1;
        gn1 = gn;
    }
*/
    gn = 3, gn1 = 1, gn2 = 0;
    for(lint i = 2; ; i++)
    {
        gn = (3 * gn1 + gn2) % L2;
        if(gn == 1 && gn1 == 0)
        {
            printf("循环节L3：%I64d\n", i - 1);
            L3 = i - 1; //183120
            break;
        }
        gn2 = gn1;
        gn1 = gn;
    }

    printf("over\n");
}

int main()
{
    test();
}

题解：

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
using namespace std;

#define MAXN 3
#define lint __int64
const lint L0 = 1000000007;
const lint L1 = 222222224;
const lint L2  = 183120;
const lint L3 = 240;
lint gn, mod;

class CMatrix
{
public:
    lint a[MAXN][MAXN], n;
    void init(int n, int flag);
    CMatrix operator+(CMatrix);
    CMatrix operator*(CMatrix);
    CMatrix power(lint);
};

void CMatrix::init(int n, int flag)
{
    this->n = n;
    memset(this->a, 0, sizeof(this->a));
    if(flag == 0) return;
    for(int i = 0; i < n; i++)
        this->a[i][i] = 1;
}

CMatrix CMatrix::operator+(CMatrix tt)
{
    CMatrix ret = *this;
    for(int i = 0; i < n; i++)
        for(int j = 0; j < n; j++)
            ret.a[i][j] = (ret.a[i][j] + tt.a[i][j]) % mod;
    return ret;
}

CMatrix CMatrix::operator*(CMatrix tt)
{
    CMatrix ret;
    ret.init(this->n, 0);
    for(int i = 0; i < n; i++)
    {
        for(int j = 0; j < n; j++)
        {
            for(int k = 0; k < n; k++)
                ret.a[i][j] += this->a[i][k] * tt.a[k][j];
            ret.a[i][j] %= mod;
        }
    }
    return ret;
}


CMatrix CMatrix::power(lint e)
{
    CMatrix ret, tmp;
    ret.init(this->n, 1);
    tmp = *this;
    while(e > 0)
    {
        if(e & 1)
            ret = ret * tmp;
        e >>= 1;
        tmp = tmp * tmp;
    }
    return ret;
}

int main()
{
    while(scanf("%I64d", &gn) != EOF)
    {

        CMatrix x, y;
        x.n = 2;
        x.a[0][0] = 3; x.a[0][1] = 1;
        x.a[1][0] = 1; x.a[1][1] = 0;

        mod = L2;
        y = x.power(gn % L3);
        mod = L1;
        y = x.power(y.a[1][0] % L2);
        mod = L0;
        y = x.power(y.a[1][0] % L1);
        printf("%I64d\n", y.a[1][0]);
    }
    return 0;
}