2015沈阳网络赛 1002（HDU 5451 矩阵快速幂 + 矩阵循环群）

最新推荐文章于 2021-04-13 22:00:04 发布

Flintx

最新推荐文章于 2021-04-13 22:00:04 发布

阅读量1.8k

点赞数

分类专栏： ACM and novicer 数学-矩阵文章标签：数学矩阵快速幂区域赛

本文链接：https://blog.youkuaiyun.com/qq_15714857/article/details/48583019

版权

ACM and novicer 同时被 2 个专栏收录

105 篇文章

订阅专栏

数学-矩阵

5 篇文章

订阅专栏

本文针对HDU5451题目的求解，介绍了一种利用矩阵快速幂来计算特定形式幂次表达式的取模运算的方法，并通过代码实现了具体的算法流程。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

HDU 5451

题意：

输入 x ( 0 <= x <= 2^32 ) 与 M，求：

$⌊ (5 + 2 6 \sqrt) 2 x + 1 ⌋ % M$ $\lfloor (5+2\sqrt 6)^{2^x+1}\rfloor\%M$

思路：

关于 $\lfloor (5+2\sqrt 6)^{n}\rfloor\%M$ 的题，之前的博客HDU 2256已经提到过，利用矩阵快速幂即可。
这次的难点是在于求矩阵乘法的循环节，首先介绍一个引理：

剩余类 $Z_p$ 上的 n 阶可逆矩阵个数为 :

$N u m n (p) = (p n - 1) (p n - p) \dots (p n - p n - 1)$ $Num_n(p) = ( p^n - 1)(p^n-p)\ldots(p^n-p^{n-1})$

有关这个定理的证明可参考http://www.docin.com/p-773748932.html

对于这一题来说，在剩余类 $Z_M$ 上共有 $(p^2-1)(p^2-p)$ 个2阶可逆矩阵，关于模M便构成了一个 $(p^2-1)(p^2-p)$ 阶的矩阵循环群，可以得到循环节为 $(p^2-1)$ 。
关于找循环节的方法可以参考(http://blog.youkuaiyun.com/acdreamers/article/details/25616461)
那么只要把 $2^x$ 利用快速幂模除 $(p^2-1)$ 再套矩阵快速幂即可求出答案。

代码：

/*
* @author FreeWifi_novicer
* language : C++/C
*/
#include<cstdio>
#include<iostream>
#include<cstring>
#include<cstdlib>
#include<cmath>
#include<algorithm>
#include<string>
#include<map>
#include<set>
#include<vector>
#include<queue>

using namespace std;

#define clr( x , y ) memset(x,y,sizeof(x))
#define cls( x ) memset(x,0,sizeof(x))
#define pr( x ) cout << #x << " = " << x << endl 
#define pri( x ) cout << #x << " = " << x << " " 
#define test( t ) int t ; cin >> t ; int kase = 1 ; while( t-- )
#define out( kase ) printf( "Case #%d: " , kase++ )
#define mp make_pair
#define pb push_back
typedef long long lint;
typedef long long ll;
typedef long long LL;

int mod ;
const int maxn = 4;
struct Matrix{
    int n,m;
    lint a[maxn][maxn];
    Matrix(int n , int m){
        this->n = n;
        this->m = m;
        cls(a);
    }
    Matrix operator * (const Matrix &tmp){
        Matrix res(n,tmp.m);
        for(int i = 0 ; i < n ; i++)
            for(int j = 0 ; j < tmp.m ; j++)
                for(int k = 0 ; k < m ; k++)
                    res.a[i][j] = (res.a[i][j] + (a[i][k] * tmp.a[k][j]) % mod) % mod;
        return res;
    }
};

void Matrix_print(Matrix x){
    for(int i = 0 ; i < x.n ; i++){
        for(int j = 0 ; j < x.m ; j++) cout << x.a[i][j] << ' ';
        cout << endl;
    }
}
Matrix fast_pow(Matrix x ,int n){
    Matrix res(x.n,x.m);
    for(int i = 0 ; i <  x.n ; i++) res.a[i][i] = 1;
    while(n){
        if(n&1)
            res = res * x;
        x = x*x;
        n >>= 1;
    }
    return res;
}
lint quick_pow( lint x , lint n , lint p ){
    lint res = 1 ;
    x %= p ;
    while( n ){
        if( n & 1 )
            res = res * x % p ;
        n >>= 1 ;
        x = x * x % p ;
    }
    return res ;
}
void solve(){
    lint n , x;
    cin >> x >> mod ;
    lint MOD = ( mod - 1 ) * ( mod + 1 ) ;
    n = quick_pow( 2 , x , MOD ) + 1 ;
    if( n == 0 ){
        cout << 1 << endl ;
        return ;
    }
    if(n == 1){
        cout << 9 << endl;
        return;
    }
    Matrix base(2,1);
    Matrix fun(2,2);
    base.a[0][0] = 5;
    base.a[1][0] = 2;
    fun.a[0][0] = 5;
    fun.a[0][1] = 12;
    fun.a[1][0] = 2;
    fun.a[1][1] = 5;
    fun = fast_pow(fun,n-1);
    base = fun * base;
    cout << (2*base.a[0][0] - 1) % mod << endl;
}

int main(){
    //freopen( "input.txt" , "r" , stdin ) ;
    test( t ){
        out( kase ) ;
        solve();
    }
    return 0;
}