BZOJ-3557: [Ctsc2014]随机数（FFT+多项式快速除法+快速幂）

2019年4月21日 252次阅读 来源: AmadeusChan

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

这道题整整搞了两天常数啊！！！都快哭死了，最后手贱把FFT改成了指针版，然后居然快了差不多一倍，然后就AC了！（好感动）

PS：诸位刷BZOJ的神犇们对不起。。。我又把OJ卡了好几个160s+了 QAQ：

《BZOJ-3557: [Ctsc2014]随机数（FFT+多项式快速除法+快速幂）》 f7246b600c338744074646a5530fd9f9d62aa06c.jpg.png

这道题还是挺神的QAQ（ORZ了好久VFK、Picks、CLJ等神牛的题解才搞懂。。。我果然还是太弱了555）：

对于第一问，把M0,X都看成多项式，那么答案Mk就是：

Mk(x)=(x^k M0) mod ( x^m + X(x) ) ，那么多项式取模+FFT+快速幂 O(m log m log k)即可。

第二问，由于保证X是好的，有由于X会跑出一个环，这个环包含了所有(0,2m)的数，事实上就是生成的元素关于多项式乘法+取模形成了一个循环群，大小为2m-1，那么我们可以得到x^m=x(mod (X(x)+x^m)),然后我们已知

M(k2^l)=M0 * x^( k2^l )，对M0求一次模(xm+X(x))意义下的逆元，M0是上述群的元素，所以逆元就是M0(2m-2)，然后乘到左边，可以得到x(k2l)，然后xk也是群中的元素，那么xk=(xk)(2m)=(xk)(2(m-l)2^l)=( x^( k*2^l ) ) (2(m-l))，然后快速幂求出x^k，乘上M0就可以得到要求的Mk了。

（上述运算均在模(x^m+X(x))意义下进行）

代码：

#include <cstdio>

#include <cstring>

#include <cmath>

#include <cstdlib>

 

#define REP( i , l , r ) for ( int i = l ; i <= r ; ++ i )

#define rep( i , x ) for ( int i = 0 ; i < x ; ++ i )

#define DOWN( i , r , l ) for ( int i = r ; i >= l ; -- i )

 

const int P = ( 479 << 21 ) ^ 1 , G = 3 ;

 

const int maxD = 21 ;

const int maxN = 1 << maxD ;

 

int MAXN ;

 

typedef long long ll ;

 

inline int Max( int x , int y ) {

    return x > y ? x : y ;

}

 

inline int Min( int x , int y ) {

    return x < y ? x : y ;

}

 

inline void Swap( int &x , int &y ) {

    int t = x ; x = y ; y = t ;

}

 

int tra[ maxN ] , A[ maxN ] , w[ 2 ][ maxD ][ maxN >> 1 ] , IN[ maxN + 1 ] ;

 

#define mul( x , y ) ( ( ll ) x * y % P )

#define plu( x , y ) ( ( x + y ) % P )

#define sub( x , y ) ( ( x - y + P ) % P )

 

inline int power( int x , int cnt ) {

    int y = 1 ;

    for ( ; cnt ; cnt >>= 1 ) {

        if ( cnt & 1 ) y = mul( y , x ) ;

        x = mul( x , x ) ;

    }

    return y ;

}

 

inline void Init_fft(  ) {

    for ( int i = 1 , cnt = 0 ; i < MAXN ; i <<= 1 , ++ cnt ) {

        w[ 1 ][ cnt ][ 0 ] = 1 , w[ 1 ][ cnt ][ 1 ] = power( G , ( P ^ 1 ) / ( i << 1 ) ) ;

        for ( int j = 2 ; j < i ; ++ j ) w[ 1 ][ cnt ][ j ] = mul( w[ 1 ][ cnt ][ j - 1 ] , w[ 1 ][ cnt ][ 1 ] ) ;

        rep( j , i ) w[ 0 ][ cnt ][ j ] = power( w[ 1 ][ cnt ][ j ] , P - 2 ) ;

    }

    for ( int i = 1 ; i <= MAXN ; i <<= 1 ) {

        IN[ i ] = power( i , P - 2 ) ;

    }

}

 

inline void FFT( int *a , int N , bool flag ) {

    tra[ 0 ] = 0 ;

    for ( int i = 1 , j = N >> 1 ; i < N ; i <<= 1 , j >>= 1 ) for ( int k = 0 ; k < i ; ++ k ) tra[ i + k ] = j ;

    for ( int i = 1 ; i < N ; i <<= 1 ) for ( int j = 0 ; j < i ; ++ j ) tra[ i + j ] |= tra[ j ] ;

    int *pt = A , *px , *pw ;

    rep( i , N ) *( pt ++ ) = a[ tra[ i ] ] ;

    int x , y ;

    for ( int i = 1 , cnt = 0 ; i < N ; i <<= 1 , ++ cnt ) {

        for ( int k = 0 ; k < N ; k += ( i << 1 ) ) {

            pt = A + k , px = A + k + i , pw = w[ flag ][ cnt ] ;

            for ( int j = 0 ; j < i ; ++ j , ++ pt , ++ px , ++ pw ) {

                if ( *px == 0 ) {

                    if ( *pt == 0 ) continue ;

                    *px = *pt ; continue ;

                }

                x = *pt , y = ( ll ) *px * *pw % P ;

                *pt = plu( x , y ) , *px = sub( x , y ) ;

            }

        }

    }

    if ( ! flag ) {

        y = IN[ N ] ;

        rep( i , N ) a[ i ] = mul( A[ i ] , y ) ;

    } else rep( i , N ) a[ i ] = A[ i ] ;

}

 

struct pol {

    int a[ maxN ] , n ;

} ;

 

inline void trans( pol &x , int N , int a[] ) {

    REP( i , 0 , x.n ) a[ i ] = x.a[ i ] ;

    for ( int i = x.n + 1 ; i < N ; ++ i ) a[ i ] = 0 ;

}

 

inline void Copy( pol &x , pol &y ) {

    x.n = y.n ;

    REP( i , 0 , x.n ) x.a[ i ] = y.a[ i ] ;

}

 

int a[ maxN ] , b[ maxN ] , c[ maxN ] ;

 

inline void Mul( pol &x , pol &y ) {

    int m = Max( x.n , y.n ) + 1 , N ; for ( N = 1 ; N < m ; N <<= 1 ) ; N <<= 1 ;

    trans( x , N , a ) , trans( y , N , b ) ;

    FFT( a , N , true ) , FFT( b , N , true ) ;

    rep( i , N ) c[ i ] = mul( a[ i ] , b[ i ] ) ;

    FFT( c , N , false ) ;

    x.n = x.n + y.n ;

    REP( i , 0 , x.n ) x.a[ i ] = c[ i ] & 1 ;

}

 

inline void Inc( pol &x , pol &y ) {

    if ( x.n < y.n ) {

        REP( i , 0 , x.n ) x.a[ i ] ^= y.a[ i ] ;

        REP( i , ( x.n + 1 ) , y.n ) x.a[ i ] = y.a[ i ] ;

    } else REP( i , 0 , y.n ) x.a[ i ] ^= y.a[ i ] ;

}

 

inline void Dec( pol &x , pol &y ) {

    if ( x.n < y.n ) {

        REP( i , 0 , x.n ) x.a[ i ] ^= y.a[ i ] ;

        REP( i , ( x.n + 1 ) , y.n ) x.a[ i ] = y.a[ i ] ;

    } else REP( i , 0 , y.n ) x.a[ i ] ^= y.a[ i ] ;

}

 

int ic[ maxD ] ;

 

inline void Inv( pol &x , pol &y , int t ) {

    int cnt = 0 , N ;

    for ( ; ; ) {

        ic[ cnt ++ ] = t ;

        if ( t == 1 ) break ;

        t = ( t + 1 ) >> 1 ;

    }

    y.n = 0 , y.a[ 0 ] = 1 ;

    DOWN( i , ( cnt - 2 ) , 0 ) {

        t = ic[ i ] ; for ( N = 1 ; N < t ; N <<= 1 ) ; N <<= 1 ;

        rep( j , t ) a[ j ] = x.a[ j ] ;

        for ( int j = t ; j < N ; ++ j ) a[ j ] = 0 ;

        trans( y , N , b ) ;

        FFT( a , N , true ) , FFT( b , N , true ) ;

        rep( j , N ) c[ j ] = mul( a[ j ] , b[ j ] ) ;

        FFT( c , N , false ) ;

        for ( int j = t ; j < N ; ++ j ) c[ j ] = 0 ;

        rep( j , t ) c[ j ] &= 1 ;

        FFT( c , N , true ) ;

        rep( j , N ) c[ j ] = mul( c[ j ] , b[ j ] ) ;

        FFT( c , N , false ) ;

        y.n = t - 1 ;

        rep( j , t ) y.a[ j ] = c[ j ] & 1 ;

    }

}

 

inline void Res( pol &x ) {

    for ( int i = 0 , j = x.n >> 1 ; i <= j ; ++ i ) Swap( x.a[ i ] , x.a[ x.n - i ] ) ;

}

 

pol pa ;

 

inline void Mod( pol &x , pol &y ) {

    if ( x.n < y.n ) return ;

    int n = x.n , m = y.n ;

    Res( x ) , Res( y ) ;

    Inv( y , pa , n - m + 1 ) ;

    x.n = n - m ; Mul( pa , x ) ; x.n = n ; pa.n = n - m ;

    Res( x ) , Res( y ) , Res( pa ) ;

    Mul( pa , y ) ; Dec( x , pa ) ;

    x.n = m - 1 ;

}

 

pol pv , px ;

 

inline void Power( pol &x , pol &mod , ll cnt ) {

    Copy( px , x ) ;

    pv.n = 0 , pv.a[ 0 ] = 1 ;

    for ( ; cnt ; cnt >>= 1 ) {

        if ( cnt & 1 ) {

            Mul( pv , px ) ; Mod( pv , mod ) ;

        }

        if ( cnt >> 1 ) {

            Mul( px , px ) ; Mod( px , mod ) ;

        }

    }

}

 

int d[ maxN ] ;

 

inline void POWER( pol &x , pol &mod , int bit[] , int len ) {

    Copy( px , x ) ; pv.a[ 0 ] = 1 , pv.n = 0 ;

    int N ; for ( N = 1 ; N < mod.n + 1 ; N <<= 1 ) ; N <<= 1 ;

    REP( i , 0 , len ) {

        if ( bit[ i ] || i < len ) {

            trans( px , N , d ) ; FFT( d , N , true ) ;

        }

        if ( bit[ i ] ) {

            Mod( px , mod ) ;

            trans( pv , N , a ) ; FFT( a , N , true ) ;

            rep( j , N ) c[ j ] = mul( a[ j ] , d[ j ] ) ;

            FFT( c , N , false ) ;

            pv.n += px.n ;

            REP( j , 0 , pv.n ) pv.a[ j ] = c[ j ] & 1 ;

            Mod( pv , mod ) ;

        }

        if ( i < len ) {

            rep( j , N ) c[ j ] = mul( d[ j ] , d[ j ] ) ;

            FFT( c , N , false ) ;

            px.n <<= 1 ;

            REP( j , 0 , px.n ) px.a[ j ] = c[ j ] & 1 ;

            Mod( px , mod ) ;

        }

    }

}

 

int ch ;

 

inline void getint( int &t ) {

    for ( ch = getchar(  ) ; ch < '0' || ch > '9' ; ch = getchar(  ) ) ;

    t = ch - '0' ;

}

 

pol M0 , Mx , My , Mz ;

int Mn , bit[ maxN ] ;

 

int main(  ) {

    scanf( "%d" , &Mn ) ;

    for ( MAXN = 1 ; MAXN < Mn + 1 ; MAXN <<= 1 ) ; MAXN <<= 1 ;

    Init_fft(  ) ;

    Mx.n = Mn , M0.n = Mn - 1 ;

    rep( i , Mn ) getint( Mx.a[ i ] ) ; Mx.a[ Mn ] = 1 ;

    rep( i , Mn ) getint( M0.a[ i ] ) ;

    int type ; scanf( "%d" , &type ) ;

    if ( ! type ) {

        My.n = My.a[ 1 ] = 1 , My.a[ 0 ] = 0 ;

        ll k ; scanf( "%lld" , &k ) ;

        Power( My , Mx , k ) ;

        Mul( M0 , pv ) ; Mod( M0 , Mx ) ;

        rep( i , Mn ) if ( M0.a[ i ] ) putchar( '1' ) ; else putchar( '0' ) ;

        putchar( '\n' ) ;

    } else {

        int l ; scanf( "%d" , &l ) ;

        rep( i , Mn ) getint( Mz.a[ i ] ) ;

        Mz.n = Mn - 1 ;

        rep( i , Mn ) bit[ i ] = 1 ; bit[ 0 ] = 0 ;

        POWER( M0 , Mx , bit , Mn - 1 ) ;

        Mul( Mz , pv ) ; Mod( Mz , Mx ) ;

        rep( i , ( Mn - l ) ) bit[ i ] = 0 ; bit[ Mn - l ] = 1 ;

        POWER( Mz , Mx , bit , Mn - l ) ;

        Mul( M0 , pv ) ; Mod( M0 , Mx ) ;

        rep( i , Mn ) if ( M0.a[ i ] ) putchar( '1' ) ; else putchar( '0' ) ;

        putchar( '\n' ) ;

    }

    return 0 ;

}
    原文作者：AmadeusChan
    原文地址: https://www.jianshu.com/p/59e524d6cbc8
    本文转自网络文章，转载此文章仅为分享知识，如有侵权，请联系博主进行删除。
点赞