BZOJ-2432: [Noi2011]兔农(矩阵快速幂)

作者: AmadeusChan | 来源:发表于2019-03-15 09:55 被阅读0次

    题目:http://www.lydsy.com/JudgeOnline/problem.php?id=2432

    我们考虑一下所求数列mod k意义下的情况(样例):

    1 1 2 3 5 0 5 5 3 0 3 3 6 2 ...

    发现每个-1的地方都是0,每两个这种相邻地方的中间是一个类FIB数列,然后我们可以拆成一堆小数列的组合:

    1 1 2 3 5 0

    5 5 3 0

    3 3 6 2 ...我们发现每个小数列只要开头确定了,那么这个数列就确定了,下一个数列的开头也确定了,然后小数列的开头只有k中情况,所以如果小数列数量大于k,一定会出现循环,然后对于一个小数列,就是一个斐波那契数列乘上一个常数a,设Fi a = 1 ( mod k ,那么Fi就是a的乘法逆元,如果Fi不存在,那么就这个小数列会一直持续下去,直接矩阵算掉,否则用Fi来查i,然后我就卡了,最后跑去VFK的BLOG找了个神奇的定理,模G意义下的斐波那契数列循环节不会超过6G,那么就相当好做了先预处理出一定长度的FIB数列,算出Fi最早出现的位置,然后每次查一下就好了,然后一个数列一个数列的跳,如果循环了就把循环用矩阵跳掉,这样似乎就可以O(k log k)了QAQ。

    代码(注意处理小数列的循环的时候,我们维护模p意义下的数要用矩阵快速幂去算,因为这个卡了好久额QAQ):

    #include <cstdio>
    
    #include <algorithm>
    
    #include <cstring>
    
      
    
    using namespace std ;
    
      
    
    #define REP( i , l , r ) for ( int i = l ; i <= r ; ++ i )
    
    #define rep( i , x ) for ( int i = 0 ; i ++ < x ; )
    
    #define DOWN( i , r , l ) for ( int i = r ; i >= l ; -- i )
    
    #define Rep( i , x ) for ( int i = 0 ; i < x ; ++ i )
    
      
    
    const int maxk = 1010000 ;
    
      
    
    typedef long long ll ;
    
      
    
    struct mat {
    
            int n , m , a[ 4 ][ 4 ] ;
    
            mat(  ) {
    
                    n = m = 0 ;
    
                    memset( a , 0 , sizeof( a ) ) ;
    
            }
    
            inline void Init( int _n ) {
    
                    n = m = _n ;
    
                    memset( a , 0 , sizeof( a ) ) ;
    
                    rep( i , n ) a[ i ][ i ] = 1 ;
    
            }
    
    } e ;
    
      
    
    inline mat mul( const mat &x , const mat &y , int mod ) {
    
            mat ret ;
    
            ret.n = x.n , ret.m = y.m ;
    
            rep( i , x.n ) rep( j , y.m ) rep( k , x.m ) {
    
                    ret.a[ i ][ j ] = ( ret.a[ i ][ j ] + ( ll ) x.a[ i ][ k ] * y.a[ k ][ j ] % mod ) % mod ;
    
            }
    
            return ret ;
    
    }
    
      
    
    inline mat power( mat x , ll cnt , int mod ) {
    
            mat ret ; ret.Init( x.n ) ;
    
            for ( ; cnt ; cnt >>= 1 ) {
    
                    if ( cnt & 1 ) ret = mul( ret , x , mod ) ;
    
                    x = mul( x , x , mod ) ;
    
            }
    
            return ret ;
    
    }
    
      
    
    inline void Init_mat(  ) {
    
            e.n = e.m = 2 , e.a[ 1 ][ 1 ] = 0 , e.a[ 1 ][ 2 ] = e.a[ 2 ][ 1 ] = e.a[ 2 ][ 2 ] = 1 ;
    
    }
    
      
    
    inline void getf( int f1 , int f2 , ll cnt , int mod , int &x , int &y ) {
    
            if ( cnt == 1 ) {
    
                    y = f1 ; return ;
    
            }
    
            if ( cnt == 2 ) {
    
                    x = f1 , y = f2 ; return ;      
    
            }
    
            mat rec = power( e , cnt - ll( 2 ) , mod ) , ret ;
    
            ret.n = 2 , ret.m = 1 , ret.a[ 1 ][ 1 ] = f1 % mod , ret.a[ 2 ][ 1 ] = f2 % mod ;
    
            ret = mul( rec , ret , mod ) ;
    
            x = ret.a[ 1 ][ 1 ] , y = ret.a[ 2 ][ 1 ] ;
    
    }
    
      
    
    inline int gcd( int x , int y ) {
    
            if ( x < y ) swap( x , y ) ;
    
            for ( int z ; y ; z = y , y = x % y , x = z ) ;
    
            return x ;
    
    }
    
      
    
    void exgcd( int a , int b , int &x , int &y ) {
    
            if ( ! b ) {
    
                    x = 1 , y = 0 ; return ;
    
            }
    
            int tx , ty ; exgcd( b , a % b , tx , ty ) ;
    
            x = ty , y = tx - ( a / b ) * ty ;
    
    }
    
      
    
    inline int Inv( int a , int mod ) {
    
            if ( gcd( a , mod ) != 1 ) return 0 ;
    
            int x , y ;
    
            if ( a > mod ) exgcd( a , mod , x , y ) ; else exgcd( mod , a , y , x ) ;
    
            if ( x > 0 ) x %= mod ;
    
            if ( x < 0 ) x += mod * ( ( - x ) / mod + ( ( - x ) % mod > 0 ) ) ;
    
            return x ;
    
    }
    
      
    
    int pos[ maxk ] , k , p , f1 , f2 ;
    
    ll n ;
    
      
    
    int last[ maxk ] , cnt = 0 ;
    
    ll sum[ maxk ] ;
    
      
    
    int main(  ) {
    
            Init_mat(  ) ;
    
            memset( pos , 0 , sizeof( pos ) ) ;
    
            scanf( "%lld%d%d" , &n , &k , &p ) ;
    
            f1 = f2 = 1 ;
    
            for ( int t , i = 3 , j = 6 * k ; i <= j ; ++ i ) {
    
                    t = ( f1 + f2 ) % k ;
    
                    f1 = f2 ; f2 = t ;
    
                    if ( ! pos[ t ] ) pos[ t ] = i ;
    
            }
    
            memset( last , 0 , sizeof( last ) ) ;
    
            sum[ cnt = 0 ] = 0 ;
    
            if ( n <= 2 ) {
    
                    printf( "1\n" ) ; return 0 ;
    
            }
    
            for ( int a = 1 , b1 = 1 , b2 = 1 , len , inv , x , y ; ; ) {
    
                    inv = Inv( a , k ) ;
    
                    if ( ! inv || ! pos[ inv ] || n <= pos[ inv ] ) {
    
                            getf( b1 , b2 , n , p , x , y ) ;
    
                            if ( y % k == 1 ) y = ( y - 1 + p ) % p ;
    
                            printf( "%d\n" , y ) ;
    
                            return 0 ;
    
                    }
    
                    len = pos[ inv ] ;
    
                    ++ cnt ;
    
                    sum[ cnt ] = sum[ cnt - 1 ] + ll( len ) , last[ a ] = cnt ;
    
                    getf( a , a , len - 1 , k , x , y ) ; a = y ;
    
                    getf( b1 , b2 , len , p , x , y ) ; y = ( y - 1 + p ) % p ;
    
                    b1 = ( x + y ) % p ; b2 = ( y + b1 ) % p ;
    
                    n -= len ;
    
                    if ( last[ a ] && n > sum[ cnt ] - sum[ last[ a ] - 1 ] ) {
    
                            mat Ex , Ea , Eb , Et ;
    
                            Ea.n = Ea.m = Eb.n = Eb.m = 3 ;
    
                            Ea.a[ 1 ][ 1 ] = 0 , Ea.a[ 1 ][ 2 ] = 1 , Ea.a[ 1 ][ 3 ] = 0 ;
    
                            Ea.a[ 2 ][ 1 ] = 1 , Ea.a[ 2 ][ 2 ] = 1 , Ea.a[ 2 ][ 3 ] = 0 ;
    
                            Ea.a[ 3 ][ 1 ] = 0 , Ea.a[ 3 ][ 2 ] = 0 , Ea.a[ 3 ][ 3 ] = 1 ;
    
                            Eb.a[ 1 ][ 1 ] = 1 , Eb.a[ 1 ][ 2 ] = 1 , Eb.a[ 1 ][ 3 ] = p - 1 ;
    
                            Eb.a[ 2 ][ 1 ] = 1 , Eb.a[ 2 ][ 2 ] = 2 , Eb.a[ 2 ][ 3 ] = p - 2 ;
    
                            Eb.a[ 3 ][ 1 ] = 0 , Eb.a[ 3 ][ 2 ] = 0 , Eb.a[ 3 ][ 3 ] = 1 ;
    
                            Ex.Init( 3 ) ;
    
                            REP( i , last[ a ] , cnt ) {
    
                                    Ex = mul( power( Ea , sum[ i ] - sum[ i - 1 ] - 2 , p ) , Ex , p ) ;
    
                                    Ex = mul( Eb , Ex , p ) ;
    
                            }
    
                            ll sz = sum[ cnt ] - sum[ last[ a ] - 1 ] , t ;
    
                            if ( n % sz ) {
    
                                    t = n / sz ;
    
                                    n %= sz ;
    
                            } else {
    
                                    t = n / sz - 1 ;
    
                                    n = sz ;
    
                            }
    
                            Et.n = 3 , Et.m = 1 , Et.a[ 1 ][ 1 ] = b1 , Et.a[ 2 ][ 1 ] = b2 , Et.a[ 3 ][ 1 ] = 1 ;
    
                            Et = mul( power( Ex , t , p ) , Et , p ) ;
    
                            b1 = Et.a[ 1 ][ 1 ] , b2 = Et.a[ 2 ][ 1 ] ;
    
                    }
    
            }
    
            return 0 ;
    
    }
    

    相关文章

      网友评论

        本文标题:BZOJ-2432: [Noi2011]兔农(矩阵快速幂)

        本文链接:https://www.haomeiwen.com/subject/ylokmqtx.html