码迷,mamicode.com
首页 > 其他好文 > 详细

HDU4565 So Easy!【矩阵连乘】【推导】

时间:2015-04-26 00:01:46      阅读:171      评论:0      收藏:0      [点我收藏+]

标签:

题目链接:

http://acm.hdu.edu.cn/showproblem.php?pid=4565


题目大意:

定义Sn = [( a + √b )^n] % m,[x]表示x向上取整,比如[3.14] = 4。给你a b n m的值,

求Sn是多少。


思路:

这道题很经典,因为(a-1)^2< b < a^2,所以0 < |a-√(b)| < 1,所以

Sn = [( a + √b )^n] % m = ( [( a + √b )^n]  + [( a - √b))^n]  ) % m。

右边其实是一个整数,如果将右边二项式展开的话,除了相互抵消的部分,剩下的部分全为

数。这个式子设An = (a + √b)^n,Bn = (a - √b)^n,Cn = [An],Sn = Cn % m。

先来求Cn的递推公式。

设Cn = pC(n-1) + qC(n-2),特征方程为x^2 = p*x + q,从[( a + √b )^n]  + [( a - √b))^n] 

可以看出特征根分别为a + √b和a - √b。则p = 2*a,q = b - a^2。

则Cn = 2*a*C(n-1) + (b-a^2)C(n-2)。然后开始构造矩阵

[   Cn   ]   =  [ 2*a    b-a^2 ] * [C(n-1)]    

[C(n-1)]       [    1            0   ]   [C(n-2)]


[   Cn   ]   =  [ 2*a    b-a^2 ] ^(n-2)   *   [C2]    

[C(n-1)]       [    1            0   ]                   [C1]

由上边Sn公式可知,C2 = 2*(a^2 + b),C1 = 2*a。

然后用矩阵快速幂求解得到Cn,然后输出就可以了。


AC代码:

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#include<cmath>
#define LL __int64
using namespace std;
const int MAXN = 2;
struct Matrix
{
    LL m[MAXN][MAXN];
};
LL A,B,N,M,ret,mo;
Matrix I =
{
    1,0,
    0,1
};
Matrix Multi(Matrix a, Matrix b,LL mo)
{
    Matrix c;
    for(int i = 0; i < MAXN; ++i)
    {
        for(int j = 0; j < MAXN; ++j)
        {
            c.m[i][j] = 0;
            for(int k = 0; k < MAXN; ++k)
            {
                c.m[i][j] += a.m[i][k] * b.m[k][j];
            }
            c.m[i][j] %= mo;
        }
    }
    return c;
}

Matrix Power(Matrix a, LL k)
{
    Matrix ans = I, p = a;
    while(k)
    {
        if(k&1)
            ans = Multi(ans,p,mo);
        p = Multi(p,p,mo);
        k >>= 1;
    }
    return ans;
}

int main()
{

    Matrix a,ans;
    while(cin >> A >> B >> N >> mo)
    {
        a.m[0][0] = 2*A%mo;
        a.m[0][1] = ((B%mo-A*A%mo)%mo + mo) % mo;
        a.m[1][0] = 1;
        a.m[1][1] = 0;
        if(N == 1)
            cout << 2*A % mo << endl;
        else
        {
            ans = Power(a,N-2);
            ret = (ans.m[0][0]%mo*2*(A*A%mo+B%mo)%mo + 2*A%mo*ans.m[0][1]%mo)%mo;
            ret %= mo;
            cout << ret << endl;
        }
    }

    return 0;
}



HDU4565 So Easy!【矩阵连乘】【推导】

标签:

原文地址:http://blog.csdn.net/lianai911/article/details/45276793

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!