[洛谷]P5401 [CTS2019] 珍珠 题解

[洛谷]P5401 [CTS2019] 珍珠 题解

题意概述

D 种珍珠,每种有无限颗,现在等概率的从 D 种珍珠中抽 n 次珍珠,每次抽 1 个珍珠,记第 i 种珍珠最后一共抽到 cnti 个,求有要多少抽法使得 i=1ncnti2m 的方案数。两种方案不同当且仅当存在一轮抽取,一种方案抽到第 i 种珍珠,另一种方案抽到第 j 种珍珠, ij

n1e9,m1e9,D1e5

题目分析

容易将题意转化,转化为求 ans=i=1n[cnti1mod2]n2×m ,记 mx=min(D,n2m) ,那么容易构造 EGF

ans=k=0mx(Dk)×n!×[xn](i0x2i+1)k×(i0x2i)Dk=n!×[xn]k=0mx(Dk)×(exex2)k×(ex+ex2)Dk=n!2D×[xn]eDx×k=0mx(Dk)×(e2x1)k×(e2x+1)Dk=n!2D×[xn]eDx×k=0mx(Dk)×i=0k(1)ki(ki)e2ix×j=0Dk(Dkj)e2jx

这样做的话也是可以往下做的,但是因为我数学太弱了,并没有做出来(但下文会继续写到如何按照这种方法做)。

可以发现这样做不好做的困难主要在于对于 F(x) 每一项的系数出现了两个不固定的 EGF 相乘,所以说导致无法出现一个单一的和式(即将 (e2x1)k(e2x+1)Dk 二项式展开之后会出现 i=0k×j=0Dk),无法直接取出 EGFxn 项系数,所以其实可以考虑改变一下求解的方式。

我们可以尝试构造一个多项式 F(x) ,使得 F(x)xi 的系数为 ans=i 的方案数。

然后设 f(x)xi 项系数表示 ans 至少为 i 的方案数,那么这样就有:

Fi(x)=j=iD(1)ji(ji)fj(x)fi(x)=(Di)×n!×[xn](exex2)i×(ex)Di

这样求解 F(x) 只需要求出 f(x) 然后卷积即可,而求 f(x) 的难度也较低:

fi(x)=(Di)×n!2i[xn]eix(e2x1)i×e(Di)x=(Di)×n!2i[xn]e(D2i)x×(e2x1)i=(Di)×n!2i[xn]j=0i(1)ij(ij)e(D2i)x×e2jx=D!×n!2i×i!×(Di)![xn]j=0i(1)iji!j!(ij)!e(D2(ij))x=D!2i×(Di)!j=0i(1)ij×1j!(ij)!×(D2(ij))n=D!2i×(Di)!j=0i1j!×(1)ij×(D2(ij))n(ij)!

对于求 F(x) 的卷积就比较套路,令 Gi(x)=FDi(x), gj(x)=fDj(x) ,则:

GDi(x)=j=iD(1)jij!i!(ji)!gDj(x)Gi(x)=j=0i(1)ij(Dj)!(Di)!(ij)!gj(x)=1(Di)!j=0i[(Dj)!×gj(x)]×(1)ij(ij)!

综上,两次卷积就可以求出答案,ans=i=0min(D,n2m)Fi(x) ,时间复杂度 O(Dlog2n) ,下面附上代码:

#include <bits/stdc++.h>
using namespace std;
const int N = 3e5 + 30;
const int M = 3e5;
const int mod = 998244353;
const int pr = 3;
const int ig = 332748118;

inline int add(int x, int y) {
	x += y;
	return x >= mod ? x - mod : x;
}

inline int del(int x, int y) {
	x -= y;
	return x < 0 ? x + mod : x;
}

inline int read() {
	char ch = getchar(); int x = 0;
	while(!isdigit(ch)) ch = getchar();
	while(isdigit(ch)) {x = (x << 3) + (x << 1) + ch - 48; ch = getchar();}
	return x;
}

int qpow(int a, int b) {
	int res = 1;
	while(b) {
		if(b & 1) res = 1ll * res * a % mod;
		a = 1ll * a * a % mod;
		b >>= 1;
	}
	return res;
}

struct node {
	int n, mx = 1, f[N], g[N];
	void NTT(int *a, int len) {
		int x, y, g, pw;
		for(int j = len >> 1; j >= 1; j >>= 1) {
			g = qpow(pr, (mod - 1) / (j << 1));
			for(int i = 0; i < len; i += (j << 1)) {
				pw = 1;
				for(int k = 0; k < j; ++k, pw = 1ll * pw * g % mod) {
					x = a[i + k]; y = a[i + j + k];
					a[i + k] = add(x, y);
					a[i + j + k] = 1ll * pw * del(x, y) % mod;
				}
			}
		}
	}
	
	void INTT(int *a, int len) {
		int x, y, g, pw;
		for(int j = 1; j < len; j <<= 1) {
			g = qpow(ig, (mod - 1) / (j << 1));
			for(int i = 0; i < len; i += (j << 1)) {
				pw = 1;
				for(int k = 0; k < j; ++k, pw = 1ll * pw * g % mod) {
					x = a[i + k]; y = 1ll * pw * a[i + j + k] % mod;
					a[i + k] = add(x, y);
					a[i + j + k] = del(x, y);
				}
			}
		}
		int Inv = qpow(len, mod - 2);
		for(int i = 0; i < len; ++i) a[i] = 1ll * a[i] * Inv % mod;
	}
	
	void mul(int *a, int *b, int len) {
		NTT(a, len); NTT(b, len);
		for(int i = 0; i < len; ++i) a[i] = 1ll * a[i] * b[i] % mod;
		INTT(a, len);
	}
	
	void init() {while(mx <= n + n) mx <<= 1;}
	void clear() {for(int i = n + 1; i < mx; ++i) f[i] = g[i] = 0;}
}F;
int D, n, m, k, fac[N], inv[N], ip[N], ans;

int main() {
	D = read(); n = read(); m = read(); k = min(n - 2 * m, D);
	if(k < 0) return puts("0"), 0;
	fac[0] = ip[0] = 1;
	for(int i = 1; i <= D; ++i) ip[i] = 499122177ll * ip[i - 1] % mod;
	for(int i = 1; i <= D; ++i) fac[i] = 1ll * fac[i - 1] * i % mod;
	inv[D] = qpow(fac[D], mod - 2);
	for(int i = D; i >= 1; --i) inv[i - 1] = 1ll * inv[i] * i % mod;
	F.n = D; F.init();
	for(int i = 0; i <= F.n; ++i) F.g[i] = inv[i];
	for(int i = 0; i <= F.n; ++i) F.f[i] = 1ll * ((i & 1) ? mod - 1 : 1) * qpow(del(D, 2 * i), n) % mod * inv[i] % mod;
	F.mul(F.f, F.g, F.mx); F.clear();
	for(int i = 0; i <= F.n; ++i) F.f[i] = 1ll * fac[D] * ip[i] % mod * F.f[i] % mod * fac[i] % mod * inv[D - i] % mod;
	reverse(F.f, F.f + F.n + 1);
	for(int i = 0; i <= F.n; ++i) F.g[i] = 1ll * ((i & 1) ? mod - 1 : 1) * inv[i] % mod;
	F.mul(F.f, F.g, F.mx); F.clear(); reverse(F.f, F.f + F.n + 1);
	for(int i = 0; i <= k; ++i) F.f[i] = 1ll * F.f[i] * inv[i] % mod, ans = add(ans, F.f[i]);
	printf("%d\n", ans);
	return 0;
}

接下来再说另一种做法,也就是本文最开始就提到的做法:

ans=n!2D×[xn]eDx×k=0mx(Dk)×(e2x1)k×(e2x+1)Dk

其实主要就是对于 k=0mx(Dk)×(e2x1)k×(e2x+1)Dk 无法化简,展开。现在我们设 G(x)=k=0mx(Dk)(x1)k×(x+1)Dk ,那么只要我们知道了 G(x) 的各项系数,这个式子就是可做的, ans 可化为:

ans=n!2D×[xn]eDx×G(e2x)=n!2D×[xn]i=0Dgi×(e2x)i×eDx=n!2D×[xn]i=0Dgi×e(2iD)x=n!2D×i=0Dgi×(2iD)nn!=12D×i=0Dgi×(2iD)n

这一步可以 O(D+DlogD×logn)=O(DlogDn) 通过欧拉筛筛 n 次幂求解。

所以现在的问题就是如何求 G(x) 的系数。这一步可以考虑求导。

G(x)=k=0mxDk!(Dk)!×(k×(x1)k1×(x+1)Dk+(Dk)×(x1)k×(x+1)Dk1)=k=0mxD×(D1)!k!(Dk)!×k(x1)k1(x+1)Dk+D×(D1)!k!(Dk)!(Dk)(x1)k(x+1)Dk1=k=0mxD×((D1)!(k1)!(Dk)!×(x1)k1(x+1)Dk+(D1)!k!(Dk1)!(x1)k(x+1)Dk1)=k=0mxD×((D1k1)×(x1)k1(x+1)Dk+(D1Dk1)(x1)k(x+1)Dk1)

恰好 (D1k1)+(D1Dk1)=(Dk) ,所以尝试将 G(x) 变一下形:

G(x)=k=0mx((D1k1)(x1)k×(x+1)Dk+(D1Dk1)(x1)k×(x+1)Dk)

这样 G(x)G(x) 形式就相近了,然后可以作差之类的找到 G(x)G(x) 之间的系数关系,也即找到了 G(x)xi 项和 xi1 项之间的关系,所以可以凑一下式子, G(x) 差一个常数 DG(x) 次数低了 1 ,于是考虑用 D×G(x)x×G(x)

DG(x)xG(x)=k=0mxD×((D1k1)(x1)k1(x+1)Dk+(D1Dk1)(x1)k(x+1)Dk1)=D×(k=0mx(D1k1)(x1)k1(x+1)Dk+k=0mx(D1k)(x1)k(x+1)Dk1)=D×(k=1mx+1(D1k1)(x1)k1(x+1)Dkk=0mx(D1k1)(x1)k1(x+1)Dk)=D×(D1mx)(x1)mx(x+1)Dmx1

这样就会变得简洁很多,并且由于 [xi](DG(x)xG(x))=(Dn)[xi]G(x) ,所以可以通过求 D×(D1mx)(x1)mx(x+1)Dmx1 的系数来得出 G(x) 的系数,这个可以直接卷积,也可以再次求导降低复杂度 ,下面是求系数方法:

H(x)=(x1)a(x+1)b ,则 :

H(x)=a(x1)a1(x+1)b+b(x1)a(x+1)b1H(x)H(x)=ax1+bx+1x2H(x)H(x)=(a+b)xH(x)+(ab)H(x)(i1)hi1(i+1)hi+1=(a+b)hi1+(ab)hihi+1=1i+1((a+bi+1)hi1+(ab)hi)

这样,令 a=mx,b=Dmx1 然后就可以递推出 D×(D1mx)(x1)mx(x+1)Dmx1 的系数,但有几个要注意的地方:

1.根据 H(x)=(x1)a(x+1)b 这个式子可知 h0=(1)a,h1=(1)a×b+(1)a1×a

2.观察递推式可以发现 hD 通过递推得到为 0 ,但显然 hD 不为 0 ,所以要单独特殊计算,即 hD=i=0a(Di)

递推是 O(D) 的,所以总的复杂度为 O(DlogDn) 。代码如下:

#include <bits/stdc++.h>
using namespace std;
const int N = 3e5 + 30;
const int M = 3e5;
const int mod = 998244353;
const int pr = 3;
const int ig = 332748118;

inline int add(int x, int y) {
	x += y;
	return x >= mod ? x - mod : x;
}

inline int del(int x, int y) {
	x -= y;
	return x < 0 ? x + mod : x;
}

inline int read() {
	char ch = getchar(); int x = 0;
	while(!isdigit(ch)) ch = getchar();
	while(isdigit(ch)) {x = (x << 3) + (x << 1) + ch - 48; ch = getchar();}
	return x;
}

int qpow(int a, int b) {
	int res = 1;
	while(b) {
		if(b & 1) res = 1ll * res * a % mod;
		a = 1ll * a * a % mod;
		b >>= 1;
	}
	return res;
}

int D, n, m, k, fac[N], inv[N], ans, f[N], g[N], a, b, mem, Inv[N], prime[N], pw[N], opt = 1;
bool vis[N];

inline int C(int x, int y) {return 1ll * fac[x] * inv[y] % mod * inv[x - y] % mod;}

void init() {
	pw[1] = 1; pw[0] = 0;
	for(int i = 2; i <= D; ++i) {
		if(!vis[i]) {
			prime[++prime[0]] = i;
			pw[i] = qpow(i, n);
		}
		for(int j = 1; j <= prime[0] && i * prime[j] <= D; ++j) {
			vis[i * prime[j]] = true;
			pw[i * prime[j]] = 1ll * pw[i] * pw[prime[j]] % mod;
			if(i % prime[j] == 0) break;
		}
	}
}

int main() {
	D = read(); n = read(); m = read(); k = min(n - 2 * m, D); 
	if(n & 1) opt = mod - 1;
	if(k == D) return printf("%d\n", qpow(D, n)), 0;
	if(k < 0) return puts("0"), 0;
	fac[0] = 1; Inv[1] = 1; init();
	for(int i = 1; i <= D; ++i) fac[i] = 1ll * fac[i - 1] * i % mod;
	inv[D] = qpow(fac[D], mod - 2);
	for(int i = D; i >= 1; --i) inv[i - 1] = 1ll * inv[i] * i % mod;
	for(int i = 2; i <= D; ++i) Inv[i] = 1ll * (mod - mod / i) * Inv[mod % i] % mod;
	
	a = D - k - 1; b = k; mem = C(D - 1, k);
	
	f[0] = (b & 1) ? (mod - 1) : 1;
	f[1] = add(1ll * ((b & 1) ? mod - 1 : 1) * a % mod, 1ll * (((b - 1) & 1) ? mod - 1 : 1) * b % mod);
	
	for(int i = 1; i < D; ++i)
		f[i + 1] = 1ll * (mod - 1) * del(add(1ll * a * del(f[i - 1], f[i]) % mod, 1ll * b * add(f[i - 1], f[i]) % mod), 1ll * (i - 1) * f[i - 1] % mod) % mod * Inv[i + 1] % mod;
	for(int i = 0; i < D; ++i)
		g[i] = 1ll * f[i] * D % mod * mem % mod * Inv[D - i] % mod;
	for(int i = 0; i <= k; ++i)
		g[D] = add(g[D], C(D, i));

	for(int i = 0; i <= D; ++i)
		if(2 * i >= D) ans = add(ans, 1ll * g[i] * pw[2 * i - D] % mod) % mod;
		else ans = add(ans, 1ll * g[i] * pw[D - 2 * i] % mod * opt % mod);
	ans = 1ll * ans * qpow(2, mod - 1 - D) % mod;
	printf("%d\n", ans);
	return 0;
}
posted @   _zyc  阅读(52)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 分享一个免费、快速、无限量使用的满血 DeepSeek R1 模型,支持深度思考和联网搜索!
· 基于 Docker 搭建 FRP 内网穿透开源项目(很简单哒)
· ollama系列01:轻松3步本地部署deepseek,普通电脑可用
· 25岁的心里话
· 按钮权限的设计及实现
点击右上角即可分享
微信分享提示