HDU 4549 矩阵连乘

M斐波那契数列F[n]是一种整数数列，它的定义如下：
F[0] = a
F[1] = b
F[n] = F[n-1] * F[n-2] ( n > 1 )

现在给出a, b, n，你能求出F[n]的值吗？输出F[n]对1000000007取模后的值即可

不难推出 f(n)=a^fib(n-2)*b^fib(n-1)%1000000007，所以通过欧拉定理或者费马小定理降幂，再用快速幂取模相乘即可。

#include <iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define M 1000000007
const int MAX=2;
typedef struct
{
    long long m[MAX][MAX];
} Matrix;
Matrix P=
{
    0,1,
    1,1,
};
Matrix I=
{
    1,0,
    0,1,
};
Matrix matrixmul(Matrix a,Matrix b) //矩阵乘法
{
    int i,j,k;
    Matrix c;
    for (i = 0 ; i < MAX; i++)
        for (j = 0; j < MAX; j++)
        {
            c.m[i][j] = 0;
            for (k=0; k<MAX; k++)
                c.m[i][j]+=((a.m[i][k]%(M-1))*(b.m[k][j]%(M-1)))%(M-1);
            c.m[i][j] %= (M-1);
        }
    return c;
}
Matrix quickpow(long long n)
{
    Matrix m = P, b = I;
    while (n >= 1)
    {
        if (n & 1)
            b = matrixmul(b,m);
        n = n >> 1;
        m = matrixmul(m,m);
    }
    return b;
}
long long modular(long long a,long long b)
{
    long long ret=1;
    while(b)
    {
        if(b&1)
            ret=ret*a%M;
        b>>=1;
        a=a*a%M;
    }
    return ret;
}
int main()
{
    long long a,b,n;
    while(cin>>a>>b>>n)
    {
        a%=M,b%=M;
        if(n==0)
            printf("%I64d\n",a%M);
        else if(n==1)
            printf("%I64d\n",b%M);
        else if(n==2)
            printf("%I64d\n",a*b%M);
        else
        {
            Matrix q;
            q=quickpow(n-2);
            long long f1=(q.m[0][0]%(M-1)+q.m[0][1]%(M-1))%(M-1),f2=(q.m[1][0]%(M-1)+q.m[1][1]%(M-1))%(M-1);
            long long ans=modular(a,f1)*modular(b,f2)%M;
            printf("%I64d\n",ans);
        }
    }
    return 0;
}