CF1349F2

F1-Solution

方便起见给 \(k\) 减 \(1\),考虑答案为:

\[\begin{aligned} &Ans_k=\sum f_{i,k}\binom{n}{i}(n-i)! \end{aligned}\]

对于 \(f_{i,k}\) 考虑通过容斥计算,设 \(g_{i,k}\) 表示长度为 \(i\) 的序列存在 \(k\) 处 < 的方案数。

那么就有:

\[g_{i,k}=\sum f_{i,j}\binom{j}{k} \]

对于 \(g_{i,k}\),考虑其代表了 \(i-k\) 个连通块,此时相当于将 \(i\) 个带标号球分配给 \(i-k\) 个盒子的方案数,又相当于给 \(1\sim i\) 染上 \(i-k\) 种颜色,且每类颜色非空的方案数,即 \(i!\cdot (e^x-1)^{i-k}[x^i]\)

于是:

\[f_{i,k}=\sum_j \binom{j}{k}(-1)^{j-k}g_{i,j} \]

故:

\[Ans_k=\sum_i \binom{n}{i}(n-i)!\sum_j \binom{j}{k}(-1)^{j-k}i!\cdot (e^x-1)^{i-j}[x^i] \]

\[Ans_k=n!\sum_i \sum_j \binom{j}{k}(-1)^{j-k}\cdot (e^x-1)^{i-j}[x^i] \]

\[Ans_k=n!\sum_j \binom{j}{k}(-1)^{j-k}\sum_{i\ge j}(e^x-1)^{i-j}[x^i] \]

设后面那个东西为 \(f_j\),此时前半部分就是差值卷积。

  • \(\frac{n!}{k!}\sum_j j!\frac{(-1)^{j-k}}{(j-k)!}f_j\)

只需要考虑计算 \(f_i\)

为此考虑:

\[f_j=\sum_{i\ge j}(e^x-1)^{i-j}[x^i] \]

并考虑:

\[[x^j]\sum_{i=0}^{n-j}(\frac{e^x-1}{x})^i \]

设为 \(F(x)\),即:

\[[x^j]\frac{1-F(x)^{n-j+1}}{1-F(x)} \]

分成两个 Part:

  • Part1,计算 \(\frac{1}{1-F(x)}\),由于常数项为 \(0\) 无法直接求逆,引入负数幂的定义后表示为 \(x^kF^{}(x)\to x^{-k}F^{-1}(x)\)
  • Part2:

\[[x^j]\frac{F(x)^{n-j+1}}{1-F(x)}\to [x^{n+1}]\frac{(xF(x))^{n-j+1}}{1-F(x)} \]

由于需要对于每个 \(j\) 进行计算,这里考虑加入一个元以提取信息,我们转为求:

\[[x^{n+1}]\sum_j \frac{(xyF(x))^{n-j+1}}{1-F(x)} \]

\[[x^{n+1}][y^{n-j+1}]\sum_j \frac{(xyF(x))^{j}}{1-F(x)}\to \frac{1}{1-F(x)}\times \frac{1}{1-xyF(x)} \]

我们希望将其表述为多项式复合的形式,这样我们才可以通过拉格朗日反演来描述答案,事实上,我们可以设:

\(W(x)=xF(x),H(x)\) 为一个多项式满足 \(H(W(x))=F(x)\),并设 \(G(x)=\frac{1}{1-H(x)}\frac{1}{1-xy}\)

这样我们就有 \(G(W(x))=\frac{1}{1-H(W(x))}\times \frac{1}{1-W(x)y}\)

于是我们所求就是:

\[[x^{n+1}]G(W(x))=\frac{1}{n+1}G(x)'P(x)^{-(n+1)}[x^{-1}] \]

其中 \(P(x)\) 为 \(W(x)\) 的复合逆,因此,答案即:

\[\frac{1}{n+1}G(x)'\Big(\frac{x}{P(x)}\Big)^{n+1}[x^n] \]

我们先求出 \(G(x)\) 的导数:

\[\frac{y(1-H(x))+(1-xy)H'(x)}{(1-H(x))^2(1-xy)^2} \]

此时答案即:

\[\begin{aligned} &\frac{1}{n+1}G(x)'\Big(\frac{x}{P(x)}\Big)^{n+1}[x^n] \\&\frac{1}{n+1}\frac{y(1-H(x))+(1-xy)H'(x)}{(1-H(x))^2(1-xy)^2}\Big(\frac{x}{P(x)}\Big)^{n+1}[x^n] \end{aligned}\]

另一边,求复合逆一直是一个非常困难的问题,但事实上考虑到 \(H(W(x))=F(x)\),考虑 \(P(x)=\frac{x}{H(x)}\to P(W(x))=\frac{W(x)}{H(W(x))}=x\)

于是我们等价于求:

\[\begin{aligned} \\&[x^n]\frac{1}{n+1}\frac{y(1-H(x))+(1-xy)H'(x)}{(1-H(x))^2(1-xy)^2}\Big(H(x)\Big)^{n+1} \\&[x^n]\frac{1}{n+1}\times \bigg(y\frac{1}{(1-H(x))(1-xy)^2}+\frac{H'(x)}{(1-H(x))^2(1-xy)}\bigg)H(x)^{n+1} \\&[x^n]\frac{1}{n+1}\times \Bigg(\frac{1}{1-H(x)}\sum (i+1)x^iy^{i+1}+\frac{H'(x)}{(1-H(x))^2}\sum x^iy^i\Bigg)H(x)^{n+1} \\&[x^ny^{j}]\frac{1}{n+1}\times \Bigg(j\frac{H(x)^{n+1}x^{j-1}}{1-H(x)}+\frac{x^jH'(x)H(x)^{n+1}}{(1-H(x))^2}\Bigg) \end{aligned}\]

于是答案即:

\[\frac{1}{n+1}\times \Bigg([x^{n-j+1}]j\frac{H(x)^{n+1}}{1-H(x)}+[x^{n-j}]\frac{H'(x)H(x)^{n+1}}{(1-H(x))^2}\Bigg) \]

注意 \(f_j=[x^n][y^{n-j+1}]\) 即:

\[\frac{1}{n+1}\times \Bigg([x^{j}](n-j+1)\frac{H(x)^{n+1}}{1-H(x)}+[x^{j+1}]\frac{H'(x)H(x)^{n+1}}{(1-H(x))^2}\Bigg) \]

最后,我们的瓶颈在于求解 \(H(x)\)

注意到 \(H(x)=\frac{x}{P(x)}\),考虑 \(P(x)\) 为 \(W(x)\) 即 \(e^x-1\) 的复合逆。

于是 \(P(x)=\ln(x+1)\),于是 \(H(x)=\frac{x}{\ln(x+1)}\)

特别的,计算得到的 \(H(x)[x^0]\) 仍为 \(1\),这意味着 \(1-H(x)\) 是需要进行偏移以保证存在逆元的。

最终复杂度 \(\mathcal O(n\log n)\)

注意偏移量。

\(Code:\)

#include<bits/stdc++.h>
using namespace std ;
#define drep( i, s, t ) for( register int i = (t); i >= (s); -- i )
#define rep( i, s, t ) for( register int i = (s); i <= (t); ++ i )
#define re register
#define int long long
//视题目调节 
const int P = 998244353 ; 
const int Gi = 332748118 ;
const int G = 3 ; 
// 基本参数 
const int N = 8e5 + 5 ; 
int gi() {
	char cc = getchar() ; int cn = 0, flus = 1 ;
	while( cc < '0' || cc > '9' ) {  if( cc == '-' ) flus = - flus ; cc = getchar() ; }
	while( cc >= '0' && cc <= '9' )  cn = ( cn * 10 + cc - '0' ) % P, cc = getchar() ;
	return cn * flus ;
}
int fpow( int x, int k ) {
	int ans = 1, base = x ;
	while( k ) {
		if( k & 1 ) ans = ans * base % P ; 
		base = base * base % P, k >>= 1 ;
	} return ans % P ; 
}
namespace Poly {
	int limit, L, Inv, R[N], inv[N], gi[N << 1], igi[N << 1] ; 
	void Init( int x ) {
		limit = 1, L = 0 ; while( limit < x ) limit <<= 1, ++ L ;
		rep( i, 0, limit ) R[i] = ( R[i >> 1] >> 1 ) | ( ( i & 1 ) << ( L - 1 ) ) ; 
		Inv = fpow( limit, P - 2 ) ;
	}
	void _init( int x ) { // max x range
		for( re int i = 0; i <= x; ++ i ) inv[i] = fpow( i, P - 2 ) ;
		Init( x + x ) ; int num = 0 ;
		for( re int k = 1; k < limit; k <<= 1 ) {
			int d = fpow( G, ( P - 1 ) / ( k << 1 ) ) ;
			for( re int j = 0, g = 1; j < k; ++ j, g = g * d % P ) gi[++ num] = g ;
		} num = 0 ;
		for( re int k = 1; k < limit; k <<= 1 ) {
			int d = fpow( Gi, ( P - 1 ) / ( k << 1 ) ) ;
			for( re int j = 0, g = 1; j < k; ++ j, g = g * d % P ) igi[++ num] = g ;
		} 
	}
	void NTT( int *a, int type ) { //NTT多项式乘法 
		for( re int i = 0; i < limit; ++ i ) if( i < R[i] ) swap( a[i], a[R[i]] ) ;
		int num = 0 ;
		for( re int k = 1; k < limit; k <<= 1 ) {
			for( re int i = 0; i < limit; i += ( k << 1 ) )
			for( re int j = i, d = num + 1; j < i + k; ++ j, ++ d ) {
				int Nx = a[j], Ny = a[j + k] * ( ( type == 1 ) ? gi[d] : igi[d] ) % P ;
				a[j] = ( Nx + Ny ) % P, a[j + k] = ( Nx - Ny + P ) % P ; 
			} num += k ;
		} 
		if( type != 1 ) rep( i, 0, limit ) a[i] = a[i] * Inv % P ; 
	}
	// 求逆 
	int H[N], Gx[N] ;
	void PolyInv( int *a, int x ) {
		Gx[0] = fpow( a[0], P - 2 ) ; int lim = 1 ; 
		while( lim < x ) {
			lim <<= 1, Init( lim << 1 ) ; 
			for( re int i = 0; i < lim; ++ i ) H[i] = a[i] ;
			NTT( H, 1 ), NTT( Gx, 1 ) ;
			for( re int i = 0; i < limit; ++ i ) Gx[i] = ( 2ll - Gx[i] * H[i] % P + P ) % P * Gx[i] % P ;
			NTT( Gx, 0 ) ; for( re int i = lim; i < limit; ++ i ) Gx[i] = 0 ;
		}
		for( re int i = 0; i <= x; ++ i ) a[i] = Gx[i] ;
		rep( i, 0, limit ) Gx[i] = H[i] = 0 ;
	} 
	void Deriv( int *a, int x ) { //求导 
		rep( i, 1, x ) a[i - 1] = a[i] * i % P ; a[x] = 0 ; 
	}
	void Integ( int *a, int x ) { //积分 
		drep( i, 1, x ) a[i] = a[i - 1] * inv[i] % P ; a[0] = 0 ;
	}
	int Gt[N], F[N] ;
	void Polyln( int *a, int x ) {
		rep( i, 0, x ) F[i] = Gt[i] = a[i] ;
		PolyInv( Gt, x ), Deriv( F, x ) ; 
		Init( x << 1 ), NTT( Gt, 1 ), NTT( F, 1 ) ;
		rep( i, 0, limit ) F[i] = F[i] * Gt[i] % P ;
		NTT( F, 0 ), Integ( F, x ) ;
		rep( i, 0, x - 1 ) a[i] = F[i] ;
		rep( i, 0, limit ) F[i] = Gt[i] = 0 ;
	}
	int _Gt[N], _F[N], _H[N] ; 
	void Polyexp( int *a, int x ) {
		_Gt[0] = 1 ; int lim = 1 ; 
		while( lim < x ) {
			lim <<= 1 ; rep( i, 0, lim - 1 ) _F[i] = a[i], _H[i] = _Gt[i] ;
			Polyln( _H, lim ), Init( lim << 1 ) ;
			NTT( _Gt, 1 ), NTT( _H, 1 ), NTT( _F, 1 ) ;
			rep( i, 0, limit ) _Gt[i] = ( 1ll - _H[i] + _F[i] + P ) % P * _Gt[i] % P ;
			NTT( _Gt, 0 ) ; rep( i, lim, limit ) _Gt[i] = 0 ;
		}
		rep( i, 0, x - 1 ) a[i] = _Gt[i] ;
		rep( i, 0, limit ) _F[i] = _H[i] = _Gt[i] = 0 ;
	}
	void PolyKSM( int *a, int x, int k ) {
		Polyln( a, x ) ;
		rep( i, 0, x - 1 ) a[i] = a[i] * k % P ;
		Polyexp( a, x ) ;
	}
}
int n, ret, dh[N], h[N], H[N], dH[N], IH[N], IH2[N] ;
int f[N], ans[N], fac[N], inv[N] ; 
signed main()
{
	n = gi(), ret = n + 10, fac[0] = inv[0] = 1 ; 
	Poly::_init(ret * 2) ; 
	rep( i, 1, ret ) fac[i] = fac[i - 1] * i % P ; 
	rep( i, 0, ret ) inv[i] = fpow(fac[i], P - 2) ; 
	h[0] = h[1] = 1, Poly::Polyln(h, ret) ; 
	rep( i, 0, ret ) h[i] = h[i + 1] ; 
	Poly::PolyInv(h, ret) ; 
	rep( i, 0, ret ) IH2[i] = IH[i] = dH[i] = H[i] = h[i] ; 
	Poly::PolyKSM(H, ret, n + 1), Poly::Deriv(dH, ret) ; 
	rep( i, 0, ret ) IH[i] = IH2[i] = (P - h[i]) % P ;
	IH[0] = IH2[0] = 0 ; 
	Poly::Init(ret + ret), Poly::NTT(IH2, 1) ; 
	rep( i, 0, Poly::limit ) IH2[i] = IH2[i] * IH2[i] % P ; 
	Poly::NTT(IH2, 0) ; 
	rep( i, ret + 1, Poly::limit ) IH2[i] = 0 ; 
	rep( i, 0, ret ) IH2[i] = IH2[i + 2], IH[i] = IH[i + 1] ; 
	Poly::PolyInv(IH2, ret), Poly::PolyInv(IH, ret) ; 
	rep( i, 0, ret ) ans[i] = (P - inv[i + 2]) % P ; 
	Poly::PolyInv(ans, ret) ;
	rep( i, 1, n ) f[i] = ans[i + 1] ; 
	//----------------------------------------------
	Poly::Init(ret * 4) ; 
	Poly::NTT(IH2, 1), Poly::NTT(H, 1), 
	Poly::NTT(dH, 1), Poly::NTT(IH, 1) ; 
	rep( i, 0, Poly::limit ) 
		dh[i] = H[i] * IH[i] % P,
		IH2[i] = IH2[i] * dH[i] % P * H[i] % P,
		IH[i] = IH[i] * H[i] % P ; 
	Poly::NTT(IH, 0), Poly::NTT(IH2, 0) ; 
	//---------------------------------------------
	int iv = fpow(n + 1, P - 2) ;
	rep( i, 1, n ) f[i] = (f[i] - iv * ((n - i + 1) * IH[i + 1] % P + IH2[i + 1]) % P + P) % P ; 
	f[0] = n ;
	rep( i, 0, ret ) h[i] = 0 ; 
	rep( i, 0, n ) h[i] = ((n - i) & 1) ? P - inv[n - i] : inv[n - i] ; 
	rep( i, 0, n ) f[i] = f[i] * fac[i] % P ; 
	Poly::Init(n * 2 + 2), Poly::NTT(f, 1), Poly::NTT(h, 1) ;
	rep( i, 0, Poly::limit ) h[i] = f[i] * h[i] % P ; 
	Poly::NTT(h, 0) ; 
	rep( i, 0, n ) f[i] = h[i + n] ; 
	rep( i, 0, n ) f[i] = f[i] * fac[n] % P * inv[i] % P ; 
	rep( i, 1, n ) printf("%lld ", f[i - 1] ) ; 
	return 0 ;
}
上一篇:相关性研究思路及代码实现(MIC-最大信息系数、Relif-F特征选择算法、pearson、spearman、kendall、卡方检验、fisher精确检验、F检验、简单粗暴的分层聚合)


下一篇:web service 部署