Solution -「LOJ #150」挑战多项式 ||「模板」多项式全家桶

(mathcal{Description})

  Link.

  给定 (n) 次多项式 (F(x)),在模 (998244353) 意义下求

[G(x)equivleft{left[1+lnleft(2+F(x)-F(0)-exp int frac{1}{sqrt{F(t)}} ext dt ight) ight]^k ight}'pmod{x^n} ]

其中保证 (F(0)) 是模数的二次剩余,开根取模意义下较小常数项值。

  (nle10^5)

(mathcal{Solution})

  先这样再那样,嗯。

  涉及到多项式 (exp)(ln),求逆,开根,求幂,积分求导,具体推导可以看 这里 哟~

(mathcal{Code})

/*~Rainybunny~*/

#include <cstdio>
#include <cassert>
#include <cstdlib>
#include <algorithm>

#define rep( i, l, r ) for ( int i = l, rep##i = r; i <= rep##i; ++i )
#define per( i, r, l ) for ( int i = r, per##i = l; i >= per##i; --i )

const int MAXL = 1 << 18, MOD = 998244353, INV2 = MOD + 1 >> 1;

inline void subeq( int& a, const int b ) { ( a -= b ) < 0 && ( a += MOD ); }
inline int sub( int a, const int b ) { return ( a -= b ) < 0 ? a + MOD : a; }
inline int mul( const long long a, const int b ) { return int( a * b % MOD ); }
inline int add( int a, const int b ) { return ( a += b ) < MOD ? a : a - MOD; }
inline void addeq( int& a, const int b ) { ( a += b ) >= MOD && ( a -= MOD ); }
inline int mpow( int a, int b ) {
    int ret = 1;
    for ( ; b; a = mul( a, a ), b >>= 1 ) ret = mul( ret, b & 1 ? a : 1 );
    return ret;
}

namespace QResidue {

/*
 * PolyOper::polySqrt() need this.
 * */

int SQRI;

struct Complex {
	int x, y;
	Complex(): x( 0 ), y( 0 ) {}
	Complex( const int tx, const int ty ): x( tx ), y( ty ) {}
	inline Complex operator * ( const Complex& c ) const {
		return Complex( add( mul( x, c.x ), mul( SQRI, mul( y, c.y ) ) ),
			add( mul( x, c.y ), mul( y, c.x ) ) );
	}
};

inline bool isResidue( const int x ) { return mpow( x, MOD - 1 >> 1 ) == 1; }

inline Complex cpow( Complex a, int b ) {
	Complex ret( 1, 0 );
	for ( ; b; a = a * a, b >>= 1 ) if ( b & 1 ) ret = ret * a;
	return ret;
}

inline int residue( const int n ) {
	if ( !n ) return 0;
	if ( !isResidue( n ) ) return -1;
	srand( 20120712 ); int a = rand() % MOD;
	while ( !a || isResidue( sub( mul( a, a ), n ) ) ) a = rand() % MOD;
	SQRI = sub( mul( a, a ), n );
	int ret = cpow( Complex( a, 1 ), MOD + 1 >> 1 ).x;
    return ret < MOD - ret ? ret : MOD - ret;
}

} // namespace QResidue.

namespace PolyOper {

/* *
 * - call init() before any operation;
 * - for `poly(n,u,w)', w can't be u;
 * - for `poly(n,u,w)', it will only access u[0..n-1] and w[0..n-1];
 * - for `poly(n,u,w)', there's no need to clear w before calling them.
 * - remember that `n' is always 2^k.
 * */

const int MG = 3;
int inv[MAXL + 5], omega[20][MAXL];

inline void init() {
    inv[1] = 1;
    rep ( i, 2, MAXL + 3 ) inv[i] = mul( MOD - MOD / i, inv[MOD % i] );
    rep ( i, 0, 17 ) {
        int* oi = omega[i];
        oi[0] = 1, oi[1] = mpow( MG, MOD - 1 >> i >> 1 );
        rep ( j, 2, ( 1 << i ) - 1 ) oi[j] = mul( oi[j - 1], oi[1] );
    }
}

inline void ntt( const int n, int* u, const int type ) {
    static int rev[MAXL]; rev[0] = 0;
    int lgn = 1; for ( ; 1 << lgn < n; ++lgn );
    rep ( i, 1, n - 1 ) rev[i] = rev[i >> 1] >> 1 | ( i & 1 ) << lgn >> 1;
    rep ( i, 1, n - 1 ) if ( i < rev[i] ) {
        u[i] ^= u[rev[i]] ^= u[i] ^= u[rev[i]];
    }
    for ( int i = 0, stp = 1; stp < n; ++i, stp <<= 1 ) {
        int* oi = omega[i];
        for ( int j = 0; j < n; j += stp << 1 ) {
            rep ( k, j, j + stp - 1 ) {
                int ev = u[k], ov = mul( oi[k - j], u[k + stp] );
                u[k] = add( ev, ov ), u[k + stp] = sub( ev, ov );
            }
        }
    }
    if ( !~type ) {
        int ivn = MOD - ( MOD - 1 ) / n;
        rep ( i, 0, n - 1 ) u[i] = mul( u[i], ivn );
        std::reverse( u + 1, u + n );
    }
}

inline void polyInv( const int n, const int* u, int* w ) {
    static int tmp[2][MAXL];
    if ( n == 1 ) return void( w[0] = mpow( u[0], MOD - 2 ) );
    polyInv( n >> 1, u, w );
    rep ( i, 0, ( n >> 1 ) - 1 ) tmp[0][i] = w[i], tmp[1][i] = u[i];
    rep ( i, n >> 1, n - 1 ) tmp[0][i] = 0, tmp[1][i] = u[i];
    rep ( i, n, ( n << 1 ) - 1 ) tmp[0][i] = tmp[1][i] = 0;
    ntt( n << 1, tmp[0], 1 ), ntt( n << 1, tmp[1], 1 );
    rep ( i, 0, ( n << 1 ) - 1 ) {
        tmp[0][i] = mul( tmp[0][i], sub( 2, mul( tmp[0][i], tmp[1][i] ) ) );
    }
    ntt( n << 1, tmp[0], -1 );
    rep ( i, 0, n - 1 ) w[i] = tmp[0][i];
}

inline void polyInt( const int n, const int* u, int* w ) {
    w[0] = 0;
    rep ( i, 1, n ) w[i] = mul( u[i - 1], inv[i] );
}

inline void polyDer( const int n, const int* u, int* w ) {
    w[n - 1] = 0;
    rep ( i, 0, n - 2 ) w[i] = mul( u[i + 1], i + 1 );
}

inline void polyLn( const int n, const int* u, int* w ) {
    static int tmp[2][MAXL]; assert( u[0] == 1 );
    polyDer( n, u, tmp[0] ), polyInv( n, u, tmp[1] );
    rep ( i, n, ( n << 1 ) - 1 ) tmp[0][i] = tmp[1][i] = 0;
    ntt( n << 1, tmp[0], 1 ), ntt( n << 1, tmp[1], 1 );
    rep ( i, 0, ( n << 1 ) - 1 ) tmp[0][i] = mul( tmp[0][i], tmp[1][i] );
    ntt( n << 1, tmp[0], -1 ), polyInt( n - 1, tmp[0], w );
}

inline void polyExp( const int n, const int* u, int* w ) {
    static int tmp[2][MAXL];
    if ( n == 1 ) return assert( !u[0] ), w[0] = 1, void();
    polyExp( n >> 1, u, w );
    rep ( i, 0, ( n >> 1 ) - 1 ) tmp[0][i] = w[i], tmp[1][i] = 0;
    rep ( i, n >> 1, ( n << 1 ) - 1 ) tmp[0][i] = tmp[1][i] = 0;
    polyLn( n, tmp[0], tmp[1] ), tmp[1][0] = sub( tmp[1][0], 1 );
    rep ( i, 0, n - 1 ) tmp[1][i] = sub( u[i], tmp[1][i] );
    ntt( n << 1, tmp[0], 1 ), ntt( n << 1, tmp[1], 1 );
    rep ( i, 0, ( n << 1 ) - 1 ) tmp[0][i] = mul( tmp[0][i], tmp[1][i] );
    ntt( n << 1, tmp[0], -1 );
    rep ( i, 0, n - 1 ) w[i] = tmp[0][i];
}

inline void polySqrt( const int n, const int* u, int* w ) {
    static int tmp[2][MAXL];
    if ( n == 1 ) return assert( ~( w[0] = QResidue::residue( u[0] ) ) );
    polySqrt( n >> 1, u, w ); rep ( i, n >> 1, n - 1 ) w[i] = 0;
    polyInv( n, w, tmp[0] );
    rep ( i, n, ( n << 1 ) - 1 ) tmp[0][i] = tmp[1][i] = 0;
    rep ( i, 0, n - 1 ) tmp[1][i] = u[i];
    ntt( n << 1, tmp[0], 1 ), ntt( n << 1, tmp[1], 1 );
    rep ( i, 0, ( n << 1 ) - 1 ) tmp[0][i] = mul( tmp[0][i], tmp[1][i] );
    ntt( n << 1, tmp[0], -1 );
    rep ( i, 0, n - 1 ) w[i] = mul( add( w[i], tmp[0][i] ), INV2 );
}

} // namespace PolyOper.

/*** templates are above. ***/

int n, k, tmpF[MAXL + 5], F[MAXL + 5], G[MAXL + 5];

int main() {
    PolyOper::init();
    scanf( "%d %d", &n, &k ), ++n;
    rep ( i, 0, n - 1 ) scanf( "%d", &F[i] ), tmpF[i] = F[i];
    int len = 1; for ( ; len < n; len <<= 1 );

    PolyOper::polySqrt( len, F, G );
    PolyOper::polyInv( len, G, F );
    PolyOper::polyInt( len, F, G );
    PolyOper::polyExp( len, G, F );
    rep ( i, 0, len - 1 ) F[i] = sub( tmpF[i], F[i] );
    F[0] = sub( add( F[0], 2 ), tmpF[0] );
    PolyOper::polyLn( len, F, G );
    G[0] = add( G[0], 1 );
    PolyOper::polyLn( len, G, F );
    rep ( i, 0, len - 1 ) F[i] = mul( F[i], k );
    PolyOper::polyExp( len, F, G );
    PolyOper::polyDer( len, G, F );

    rep ( i, 0, n - 2 ) printf( "%d%c", F[i], i + 2 < n ? ' ' : '
' );
    return 0;
}

原文地址:https://www.cnblogs.com/rainybunny/p/15005747.html