HDU2256_Problem of Precision_矩阵巧妙解特定高精度

/*
*State: 2256    0MS    260K    2351 B    C++
*题目大意：
*        求(sqrt(2) + sqrt(3))^2n取下界之后mod1024
*        (1 <= n <= 10^9)。
*解题思路：
*        看到了n那么大，肯定是二分求幂，但是精度是个
*        问题。可以(sqrt(2) + sqrt(3))^2n = (5 + 2sqrt(6))^n
*        然后令5+2sqrt(6)^n = an + bn * sqrt(6)
*                           = (5 + 2sqrt(6)) * (an-1 + bn-1 * sqrt(6))
*                           = (5an-1 + 12bn-1) + (2an-1 + 5bn-1)* sqrt(6)
*        所以an = 5an-1 + 12bn-1, bn = 2an-1 + 5bn-1.
*        之后求出an跟bn之后，不能用(int)(an + bn * sqrt(6)) % mod这个求，因为
*        注意mod运算不能支持浮点的运算。
*        可以由an+bn * sqrt(6) + an-bn * sqrt(6) - (an-bn * sqrt(6))
*            = 2 * an - (an-bn * sqrt(6)),
*        由于an-bn *sqrt(6) =  (sqrt(2) - sqrt(3))^2n = (0.1……)^n 无限接近于0.
*        而结果要求取下界，所以是2 * an - 1.
*/

#include <iostream>
#include <cmath>
#define maxn 3

using namespace std;

struct Mat
{
    int m, n;
    int d[maxn][maxn];
    void init(int m, int n) 
    {
        this->m = m;
        this->n = n;
        memset(d, 0, sizeof(d));
    }
    void initE(int size)                              //生成单位阵
    {
        m = n = size;
        for(int i = 0; i < n; i ++)
        {
            for(int j = 0; j < n; j ++)
            {
                d[i][j] = i==j;
            }
        }
    }
    Mat operator * (const Mat & mat) const
    {
        static Mat res;
        res.init(m, mat.n);
        for(int i = 0; i < res.m; i ++)
        {
            for(int k = 0; k < n; k ++)
            {
                if(d[i][k]==0)    continue;
                for(int j = 0; j < res.n; j ++)
                {
                    res.d[i][j] += d[i][k] * mat.d[k][j];
                }
            }
        }
        return res;
    }
    Mat mul_mod(const Mat & mat, int mod) const
    {
        static Mat res;
        res.init(m, mat.n);
        for(int i = 0; i < res.m; i ++)
        {
            for(int k = 0; k < n; k ++)
            {
                if(d[i][k]==0)    continue;
                for(int j = 0; j < res.n; j ++)
                {
                    res.d[i][j]=(res.d[i][j]+d[i][k]*mat.d[k][j]) % mod;
                }
            }
        }
        return res;
    }
    void pow_mod(int k, int mod)                      //this = this^k % mod;
    {
        static Mat a;
        a = *this;
        for(this->initE(n); k; k>>=1, a=a.mul_mod(a, mod))
            if(k&1)    *this=this->mul_mod(a, mod);
    }
    void view_mat()
    {
        for(int i = 0; i < m; i++)
        {
            printf("%d", d[i][0]);
            for(int j = 1; j < n; j++)
                printf(" %d", d[i][j]);
            printf("\n");
        }
    }
};

int main(void)
{
#ifndef ONLINE_JUDGE
    freopen("in.txt", "r", stdin);
#endif
    int cas;
    int mod = 100000024;
    scanf("%d", &cas);
    while(cas--)
    {
        int n;
        scanf("%d", &n);
        Mat a;
        a.d[0][0] = a.d[1][1] = 5;
        a.d[0][1] = 12, a.d[1][0] = 2;
        a.m = a.n = 2;
        a.pow_mod(n - 1, mod);
        int tmp = (a.d[0][0] * 5 + a.d[0][1] * 2) % mod;
        printf("%d\n", ((2 * tmp - 1) % mod + mod) % mod);
    }
    return 0;
}