HDU 3306 Another kind of Fibonacci(快速幂矩阵)

#include <cstring>
#include <cstdio>
#include <string>
#include <iostream>
#include <algorithm>
#include <vector>
#include <queue>
using namespace std;
#define MOD 10007
#define LL __int64
LL p[4][4],mat[4][4];
int qmod(int n)
{
    LL c[4][4];
    int i,j,k;
    while(n)
    {
        if(n&1)
        {
            memset(c,0,sizeof(c));
            for(i = 0;i < 4;i ++)
            {
                for(j = 0;j < 4;j ++)
                {
                    for(k = 0;k < 4;k ++)
                    {
                        c[i][j] += mat[i][k]*p[k][j];
                        c[i][j] %= MOD;
                    }
                }
            }
            memcpy(mat,c,sizeof(mat));
        }
        memset(c,0,sizeof(c));
        for(i = 0;i < 4;i ++)
        {
            for(j = 0;j < 4;j ++)
            {
                for(k = 0;k < 4;k ++)
                {
                    c[i][j] += p[i][k]*p[k][j];
                    c[i][j] %= MOD;
                }
            }
        }
        memcpy(p,c,sizeof(p));
        n >>= 1;
    }
    return (mat[0][1] + mat[0][0] + mat[0][2] + mat[0][3])%MOD;
}
int main()
{
    LL x,y,n;
    int i,j;
    while(scanf("%I64d%I64d%I64d",&n,&x,&y)!=EOF)
    {
        memset(p,0,sizeof(p));
        p[0][0] = p[0][1] = 1;
        p[1][1] = x*x;
        p[1][2] = y*y;
        p[1][3] = 2*x*y;
        p[2][1] = 1;
        p[3][1] = x;
        p[3][3] = y;
        for(i = 0;i < 4;i ++)
        {
            for(j = 0;j < 4;j ++)
            {
                if(i == j)
                mat[i][j] = 1;
                else
                mat[i][j] = 0;
            }
        }
        printf("%d\n",qmod(n));
    }
    return 0;
}