标签：const 加强版 迭代 int LOJ6538 inline frac 烷基 mod
题目
点这里看题目。
省流版：求结点个数为 \(n\) 的结点儿子数不超过 3 的无标号有根树个数。
对于 \(100\%\) 的数据，满足 \(1\le n\le 10^5\)。
分析
首先，遇到这种问题，不难猜想使用生成函数；并且，其余一大堆树计数的问题，都可以用生成函数解决。
设 \(G(x)\) 为问题的生成函数，也就是 \([x^n]G(x)\) 表示“结点个数为 \(n\) 且结点儿子数不超过 3 的无标号有根树个数”。我们认为 \([x^0]G(x)=1\)，也就是包含空树。
尝试列出一个关于 \(G\) 的方程：先确定一个根，而后考虑根的儿子。由于 \(G\) 包含空树，我们将所有的情况都看作是根有三个儿子。
我们不妨认为儿子之间有顺序，这样很容易计算出方案数为 \(G^3(x)\)；现在，如果两个方案是相同的，就存在一种儿子之间的映射，使得映射后子树对应相同。不难发现这就是在求 \(S_3\) 划分下等价类个数。所以可以联想到 Burnside 引理。最终方程为：
\[G(x)=1+\frac x 6(G^3(x)+3G(x^2)G(x)+2G(x^3)) \]稍微移一下项，就可以得到：
\[\frac x 6G^3(x)+\left(\frac x 2G(x^2)-1\right)G(x)+\frac x 3G(x^3)+1=0 \]这个方程有两种处理方法：一者，分治乘法，比较好想但实现难度未知；二者，我们用牛顿迭代解出这个方程。
这里有一个值得注意的细节，牛顿迭代实际上解的是如下的方程：
\[f(G,x)\equiv 0\pmod {x^n} \]其中，\(f(G,x)\) 同时是关于 \(G,x\) 的多项式方程。其中用到的泰勒级数实际上也是在对 \(G\) 这个变量求偏导。所以最终的迭代过程为：
若：
\[f(G_{t},x)\equiv 0\pmod {x^{\lceil\frac n 2\rceil}}\\ \]则：
\[G_{t+1}\equiv G_t-\frac{f(G_{t},x)}{\frac{\partial f}{\partial G}(G_t,x)}\pmod {x^n} \]
在这个情境下，有 \(f(G,x)=\frac x 6G^3(x)+\left(\frac x 2G(x^2)-1\right)G(x)+\frac x 3G(x^3)+1\)。但是存在 \(G(x^2),G(x^3)\) 这样的项，怎么处理呢？
注意到，在 \(G_t\) 的时候，实际上我们已经得到了 \(G_{t+1}(x^3)\equiv G_{t}(x^3),G_{t+1}(x^2)\equiv G_t(x^2)\pmod {x^n}\)（迭代过程中尤其注意次数）。所以我们可以认为 \(G_t(x^2),G_t(x^3)\) 为常量，求偏导得到的是 \(\frac{\partial f}{\partial G}(G,x)=\frac{x}{2}G^2(x)+\frac{x}{2}G(x^2)-1\)。
之后套用牛顿迭代即可，复杂度为 \(O(n\log n)\)。
小结：
- 注意生成函数在树计数、图计数问题中的应用。
- 注意有无标号的处理方法；一般来说，有标号会倾向于使用 EGF，而无标号会倾向于使用先标号，再去重或者排序等处理手段。
- 学习一下
完整的牛顿迭代法。其中的小细节不能忘了！ - 牛顿迭代的倍增过程尤其注意次数；或者说，注意 \(\pmod {x^n}\)，这是一个很好的辅助条件。
代码
#include <cstdio>
#include <iostream>
#include <algorithm>
#define rep( i, a, b ) for( int i = (a) ; i <= (b) ; i ++ )
#define per( i, a, b ) for( int i = (a) ; i >= (b) ; i -- )
const int mod = 998244353, g = 3, phi = mod - 1;
const int inv2 = ( mod + 1 ) / 2, inv3 = ( mod + 1 ) / 3, inv6 = ( mod + 1 ) / 6;
const int MAXN = ( 1 << 18 ) + 5;
template<typename _T>
void read( _T &x ) {
x = 0; char s = getchar(); bool f = false;
while( s < '0' || '9' < s ) { f = s == '-', s = getchar(); }
while( '0' <= s && s <= '9' ) { x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar(); }
if( f ) x = -x;
}
template<typename _T>
void write( _T x ) {
if( x < 0 ) putchar( '-' ), x = -x;
if( 9 < x ) write( x / 10 );
putchar( x % 10 + '0' );
}
int ans[MAXN];
int N;
inline int Qkpow( int, int );
inline int Inv( const int a ) { return Qkpow( a, mod - 2 ); }
inline int Mul( int x, const int v ) { return 1ll * x * v % mod; }
inline int Sub( int x, const int v ) { return ( x -= v ) < 0 ? x + mod : x; }
inline int Add( int x, const int v ) { return ( x += v ) >= mod ? x - mod : x; }
inline int Qkpow( int base, int indx ) {
int ret = 1;
while( indx ) {
if( indx & 1 ) ret = Mul( ret, base );
base = Mul( base, base ), indx >>= 1;
}
return ret;
}
namespace Basic {
const int L = 18;
int w[MAXN];
void NTTInit( const int n = 1 << L ) {
w[0] = 1, w[1] = Qkpow( g, phi >> L );
rep( i, 2, n - 1 ) w[i] = Mul( w[i - 1], w[1] );
}
void NTT( int *coe, const int n, const int typ ) {
for( int i = 1, j = 0 ; i < n ; i ++ ) {
for( int d = n ; j ^= ( d >>= 1 ), ~ j & d ; );
if( j < i ) std :: swap( coe[j], coe[i] );
}
int *wp, p, k;
for( int s = 1 ; s < n ; s <<= 1 )
for( int i = 0 ; i < n ; i += s << 1 ) {
wp = w, p = ( 1 << L ) / ( s << 1 );
for( int j = 0 ; j < s ; j ++, wp += p )
k = Mul( *wp, coe[i + j + s] ),
coe[i + j + s] = Sub( coe[i + j], k ),
coe[i + j] = Add( coe[i + j], k );
}
if( ~ typ ) return ;
std :: reverse( coe + 1, coe + n );
int inv = Inv( n ); rep( i, 0, n - 1 ) coe[i] = Mul( coe[i], inv );
}
}
namespace PolyInv {
int P[MAXN], U0[MAXN], V[MAXN], U[MAXN];
void Newton( const int n ) {
if( n == 1 ) {
U0[0] = Inv( P[0] );
return ;
}
int h = ( n + 1 ) >> 1, L; Newton( h );
for( L = 1 ; L < ( n << 1 ) ; L <<= 1 );
for( int i = 0 ; i < L ; i ++ ) U[i] = V[i] = 0;
for( int i = 0 ; i < h ; i ++ ) U[i] = U0[i];
for( int i = 0 ; i < n ; i ++ ) V[i] = P[i];
Basic :: NTT( U, L, 1 ), Basic :: NTT( V, L, 1 );
for( int i = 0 ; i < L ; i ++ ) U[i] = Mul( U[i], Sub( 2, Mul( U[i], V[i] ) ) );
Basic :: NTT( U, L, -1 );
for( int i = 0 ; i < n ; i ++ ) U0[i] = U[i];
}
void PolyInv( int *ret, const int *inp, const int n ) {
rep( i, 0, n - 1 ) P[i] = U0[i] = V[i] = U[i] = 0;
rep( i, 0, n - 1 ) P[i] = inp[i];
Newton( n );
rep( i, 0, n - 1 ) ret[i] = U0[i];
}
}
namespace Solve {
int U0[MAXN], U[MAXN], V[MAXN], W[MAXN];
int P[MAXN], Q[MAXN];
void Newton( const int n ) {
if( n == 1 ) { U0[0] = 1; return ; }
int h = ( n + 1 ) >> 1, L; Newton( h );
for( L = 1 ; L < ( n << 1 ) ; L <<= 1 );
for( int i = 0 ; i < L ; i ++ ) P[i] = Q[i] = U[i] = V[i] = W[i] = 0;
for( int i = 0 ; i < h ; i ++ ) {
U[i] = U0[i];
if( i * 2 < n ) V[i * 2] = U0[i];
if( i * 3 < n ) W[i * 3] = U0[i];
}
Basic :: NTT( U, L, 1 );
Basic :: NTT( V, L, 1 );
Basic :: NTT( W, L, 1 );
int *wp = Basic :: w, p = ( 1 << Basic :: L ) / L;
for( int i = 0 ; i < L ; i ++, wp += p ) {
P[i] = Add( 1, Add( Mul( Mul( inv3, *wp ), W[i] ),
Add( Mul( U[i], Sub( Mul( Mul( inv2, *wp ), V[i] ), 1 ) ),
Mul( Mul( inv6, *wp ), Mul( Mul( U[i], U[i] ), U[i] ) ) ) ) );
Q[i] = Sub( Add( Mul( Mul( inv2, *wp ), Mul( U[i], U[i] ) ),
Mul( Mul( inv2, *wp ), V[i] ) ), 1 );
}
Basic :: NTT( P, L, -1 );
Basic :: NTT( Q, L, -1 );
for( int i = n ; i < L ; i ++ ) P[i] = Q[i] = 0;
PolyInv :: PolyInv( Q, Q, n );
Basic :: NTT( P, L, 1 );
Basic :: NTT( Q, L, 1 );
for( int i = 0 ; i < L ; i ++ ) P[i] = Mul( P[i], Q[i] );
Basic :: NTT( P, L, -1 );
for( int i = 0 ; i < n ; i ++ ) U0[i] = Sub( U0[i], P[i] );
for( int i = n ; i < L ; i ++ ) U0[i] = 0;
}
void Solve( int *ret, const int n ) {
Newton( n );
rep( i, 0, n - 1 ) ret[i] = U0[i];
}
}
int main() {
Basic :: NTTInit();
read( N );
Solve :: Solve( ans, N + 1 );
write( ans[N] ), putchar( '\n' );
return 0;
}
标签：const,加强版,迭代,int,LOJ6538,inline,frac,烷基,mod 来源： https://www.cnblogs.com/crashed/p/15925039.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。