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.
*/

View Code

  1 #include <iostream>
  2 #include <cmath>
  3 #define maxn 3
  4 
  5 using namespace std;
  6 
  7 struct Mat
  8 {
  9     int m, n;
 10     int d[maxn][maxn];
 11     void init(int m, int n) 
 12     {
 13         this->m = m;
 14         this->n = n;
 15         memset(d, 0, sizeof(d));
 16     }
 17     void initE(int size)                              //生成单位阵
 18     {
 19         m = n = size;
 20         for(int i = 0; i < n; i ++)
 21         {
 22             for(int j = 0; j < n; j ++)
 23             {
 24                 d[i][j] = i==j;
 25             }
 26         }
 27     }
 28     Mat operator * (const Mat & mat) const
 29     {
 30         static Mat res;
 31         res.init(m, mat.n);
 32         for(int i = 0; i < res.m; i ++)
 33         {
 34             for(int k = 0; k < n; k ++)
 35             {
 36                 if(d[i][k]==0)    continue;
 37                 for(int j = 0; j < res.n; j ++)
 38                 {
 39                     res.d[i][j] += d[i][k] * mat.d[k][j];
 40                 }
 41             }
 42         }
 43         return res;
 44     }
 45     Mat mul_mod(const Mat & mat, int mod) const
 46     {
 47         static Mat res;
 48         res.init(m, mat.n);
 49         for(int i = 0; i < res.m; i ++)
 50         {
 51             for(int k = 0; k < n; k ++)
 52             {
 53                 if(d[i][k]==0)    continue;
 54                 for(int j = 0; j < res.n; j ++)
 55                 {
 56                     res.d[i][j]=(res.d[i][j]+d[i][k]*mat.d[k][j]) % mod;
 57                 }
 58             }
 59         }
 60         return res;
 61     }
 62     void pow_mod(int k, int mod)                      //this = this^k % mod;
 63     {
 64         static Mat a;
 65         a = *this;
 66         for(this->initE(n); k; k>>=1, a=a.mul_mod(a, mod))
 67             if(k&1)    *this=this->mul_mod(a, mod);
 68     }
 69     void view_mat()
 70     {
 71         for(int i = 0; i < m; i++)
 72         {
 73             printf("%d", d[i][0]);
 74             for(int j = 1; j < n; j++)
 75                 printf(" %d", d[i][j]);
 76             printf("\n");
 77         }
 78     }
 79 };
 80 
 81 int main(void)
 82 {
 83 #ifndef ONLINE_JUDGE
 84     freopen("in.txt", "r", stdin);
 85 #endif
 86     int cas;
 87     int mod = 100000024;
 88     scanf("%d", &cas);
 89     while(cas--)
 90     {
 91         int n;
 92         scanf("%d", &n);
 93         Mat a;
 94         a.d[0][0] = a.d[1][1] = 5;
 95         a.d[0][1] = 12, a.d[1][0] = 2;
 96         a.m = a.n = 2;
 97         a.pow_mod(n - 1, mod);
 98         int tmp = (a.d[0][0] * 5 + a.d[0][1] * 2) % mod;
 99         printf("%d\n", ((2 * tmp - 1) % mod + mod) % mod);
100     }
101     return 0;
102 }

posted on 2012-08-03 01:34 cchun 阅读(170) 评论(0) 编辑收藏举报