标签：fp frac ary sum 计数 2017 prod 集训 mod
在一个\(s\)个点的图中，存在\(s-n\)条边，使图中形成了\(n\)个连通块，第\(i\)个连通块中有\(a_i\)个点。
现在我们需要再连接\(n-1\)条边，使该图变成一棵树。对一种连边方案，设原图中第\(i\)个连通块连出了\(d_i\)条边，那么这棵树\(T\)的价值为：
\[\mathrm{val}(T) = \left(\prod_{i=1}^{n} {d_i}^m\right)\left(\sum_{i=1}^{n} {d_i}^m\right) \]你的任务是求出所有可能的生成树的价值之和，对\(998244353\)取模。
（可能只有我没读出来题目说连通块内的连边方式不计。）
树和每个点的度数可以联想到\(prufer\)序列。那么设\(c_i\)为第\(i\)个点在\(prufer\)序列中出现的次数，则\(c_i=d_i-1\)。考虑对于一个确定的序列\(c_i\)，它对答案的贡献就是
\[\frac{(n-2)!}{\prod_{i=1}^nc_i!}\times (\prod_{i=1}^na_i^{c_i+1})\times (\prod_{i=1}^n(c_i+1)^m)\times (\sum_{i=1}^n(c_i+1)^m) \]第一项是这个序列对应的有标号无根树个树；第二个是因为每个连通块中的点可以任意分配这个连通块的出边；后面两个是题面定义的价值。
那么有一个暴力思路就是递推，在那之前我们把\((n-2)!\prod_{i=1}^na_i\)看作常数项，只考虑剩余的式子。设
\[f_{n,m}=\sum_{\sum_{c_i}=m}\frac{(\prod_{i=1}^na_i^{c_i})\times (\prod_{i=1}^n(c_i+1)^m)}{\prod_{i=1}^nc_i!}\\ g_{n,m}=\sum_{\sum_{c_i}=m}\frac{(\prod_{i=1}^na_i^{c_i})\times (\prod_{i=1}^n(c_i+1)^m)}{\prod_{i=1}^nc_i!}\times (\sum_{i=1}^n(c_i+1)^m) \]即递推前\(n\)个点的总\(c_i\)为\(m\)的所有情况之和。转移就有：
\[f_{n,m}=\sum_{i=0}^mf_{n-1,m-i}\times \frac{a_i^i\times (i+1)^m}{i!}\\ g_{n,m}=\sum_{i=0}^m(g_{n-1,m-i}+f_{n-1,m-i}\times (i+1)^m)\times \frac{a_i^i\times (i+1)^m}{i!} \]答案就是\(g_{n,n-2}\times \frac{(n-2)!}{\prod_{i=1}^na_i}\)。
现在就可以有\(20\)分的好成绩了。如果用\(NTT\)实现上面的转移就可以有\(40\)分的好成绩。
然后发现我这个式子并不好优化（懒得优化两个式子）。瞟一眼题解之后发现开始那个式子可以化得好看些：
我们还是不管前面的常数项。可以发现每个\(c_i\)对式子的贡献就是\(\frac{a_i^{c_i}(c_i+1)^m}{c_i!}\)或者一个序列中仅有一个\(c_i\)贡献为\(\frac{a_i^{c_i}(c_i+1)^{2m}}{c_i!}\)，那么这个式子就可以由若干个次数表示\(c\)的\(EGF\)乘起来（实际上如果尝试用\(EGF\)推一下那个\(n^3\)递推可以更容易发现这种性质）。乍一看一共有\(n\)个系数不同的\(n\)次多项式，似乎不可做，但是第\(i\)个多项式每一项都有\(a_i\)的若干次方，且与\(x\)次数相同，所以这些多项式都可以写成\(F(a_ix)\)的形式。因此设
\[A(x)=\frac{(i+1)^{2m}x^m}{i!},B(x)=\frac{(i+1)^mx^m}{i!} \]答案就是\(\sum_{i=1}^n\frac{A(a_ix)}{B(a_ix)}\prod_{j=1}^nB(a_jx)\)。这样有什么好处呢？这里补充一下这种trick。
如果式子可以写成\(\sum_{i=1}^nF(a_ix)\)的形式，并且对任意\(m\)都求出了\(\sum_{i=1}^na_i^m\)，那么只要求出\(F(x)\)，式子就可以变成
\[\sum_{m=0}([x^m]F(x))\sum_{i=1}^na_i^m \]因此我们求出\(\frac{A(x)}{B(x)}\)和\(\prod_{i=1}^nB(a_ix)\)即可。但后面这个是\(\prod\)，和前面的\(\sum\)不同，这里就要取个\(\ln\)：
\[\prod_{i=1}^nB(a_ix)=e^{\sum_{i=1}^n\ln{B(a_ix)}} \]那么我们算一下\(\ln{B(x)}\)就可以像上面那样算了。
现在唯一的问题就是怎么对每个\(m\)求出\(\sum_{i=1}^na_i^m\)。类似于自然数幂和的推导，我们写出这个东西的\(OGF\)就有
\[G(x)=\sum_{m=0}(\sum_{i=1}^na_i^m)x^m=\sum_{i=1}^n\sum_{m=0}a_i^mx^m=\sum_{i=1}^n\frac{1}{1-a_ix} \]这就有点像P4705玩游戏这题的技巧，因为
\[x(\ln(1-a_ix))'=-\frac{a_ix}{1-a_ix}=-\frac{1}{1-a_ix}+1 \]那设\(H(x)=\sum_{i=1}^n(\ln(1-a_ix))'\)，那\(G(x)=-xH(x)+n\)。求\(H\)就：
\[H(x)=(\sum_{i=1}^n\ln(1-a_ix))'=(\ln\prod_{i=1}^n(1-a_ix))' \]分治\(NTT\)即可。至此这题就解决了，复杂度瓶颈为最后分治\(NTT\)的\(\mathcal{O}(n\log^2n)\)。
#include<bits/stdc++.h>
#define rg register
#define il inline
#define cn const
#define gc getchar()
#define fp(i, a, b) for(int i = (a), ed = (b); i <= ed; ++i)
#define fb(i, a, b) for(int i = (a), ed = (b); i >= ed; --i)
#define go(u) for(int i = head[u]; ~i; i = e[i].nxt)
using namespace std;
typedef cn int cint;
typedef long long LL;
il void rd(int &x){
x = 0;
rg int f(1); rg char c(gc);
while(c < '0' || '9' < c){if(c == '-')f = -1; c = gc;}
while('0' <= c && c <= '9')x = (x<<1)+(x<<3)+(c^48), c = gc;
x *= f;
}
cint maxn = 30010, mod = 998244353, G = 3, invG = (mod+1)/3;
int n, m, a[maxn], fac[maxn], ifac[maxn], inv[maxn], mul = 1;
int lim, hst, rev[maxn<<2], A[maxn<<2], B[maxn<<2], ln[maxn<<2], invB[maxn<<2];
int H[maxn<<2], E[maxn<<2];
il int fpow(int a, int b, int ans = 1){
for(; b; b >>= 1, a = 1ll*a*a%mod)if(b&1)ans = 1ll*ans*a%mod;
return ans;
}
il void ntt(int *a, cint &typ){
fp(i, 0, lim-1)if(i > rev[i])swap(a[i], a[rev[i]]);
for(rg int md = 1; md < lim; md <<= 1){
rg int len = md<<1, Gn = fpow(typ ? invG : G, (mod-1)/len);
for(rg int l = 0; l < lim; l += len){
for(rg int nw = 0, Pow = 1; nw < md; ++nw, Pow = 1ll*Pow*Gn%mod){
rg int x = a[l+nw], y = 1ll*a[l+nw+md]*Pow%mod;
a[l+nw] = (x+y)%mod, a[l+nw+md] = (x-y+mod)%mod;
}
}
}
if(typ){
rg int inv = fpow(lim, mod-2);
fp(i, 0, lim-1)a[i] = 1ll*a[i]*inv%mod;
}
}
il void init(int n){
lim = 1, hst = 0;
while(lim < n)lim <<= 1, ++hst;
fp(i, 0, lim-1)rev[i] = (rev[i>>1]>>1)|((i&1)<<hst-1);
}
int inv_ary[maxn<<2];
void get_inv(int *a, int *f, int n){
if(n == 1)return f[0] = fpow(a[0], mod-2), void();
get_inv(a, f, n+1>>1), init(2*n-1);
fp(i, 0, n-1)inv_ary[i] = a[i];
ntt(f, 0), ntt(inv_ary, 0);
fp(i, 0, lim-1)f[i] = 1ll*f[i]*(2-1ll*f[i]*inv_ary[i]%mod+mod)%mod;
ntt(f, 1);
fp(i, n, lim-1)f[i] = 0;
fp(i, 0, lim-1)inv_ary[i] = 0;
}
int ln_ary[maxn<<2];
il void get_ln(int *a, int *f, int n){
get_inv(a, ln_ary, n), init(2*n-2);
fp(i, 1, n-1)f[i-1] = 1ll*a[i]*i%mod;
ntt(f, 0), ntt(ln_ary, 0);
fp(i, 0, lim-1)f[i] = 1ll*f[i]*ln_ary[i]%mod;
ntt(f, 1);
fp(i, n-1, lim)f[i] = 0;
fp(i, 0, lim)ln_ary[i] = 0;
fb(i, n-1, 1)f[i] = 1ll*f[i-1]*inv[i]%mod;
f[0] = 0;
}
int exp_ary[maxn<<2];
void get_exp(int *a, int *f, int n){
if(n == 1)return f[0] = 1, void();
get_exp(a, f, n+1>>1), get_ln(f, exp_ary, n), init(2*n-1);
fp(i, 0, n-1)exp_ary[i] = (a[i]-exp_ary[i]+mod)%mod;
if((++exp_ary[0]) == mod)exp_ary[0] = 0;
ntt(f, 0), ntt(exp_ary, 0);
fp(i, 0, lim-1)f[i] = 1ll*f[i]*exp_ary[i]%mod;
ntt(f, 1);
fp(i, n, lim-1)f[i] = 0;
fp(i, 0, lim-1)exp_ary[i] = 0;
}
int div_ary[19][maxn<<2];
void divntt(int d, int l, int r){
if(l == r)return div_ary[d][0] = 1, div_ary[d][1] = mod-a[l], void();
int md = l+r>>1, len = r-l+2;
divntt(d, l, md), divntt(d+1, md+1, r), init(len), ntt(div_ary[d], 0), ntt(div_ary[d+1], 0);
fp(i, 0, lim-1)div_ary[d][i] = 1ll*div_ary[d][i]*div_ary[d+1][i]%mod;
ntt(div_ary[d], 1);
fp(i, len, lim-1)div_ary[d][i] = 0;
fp(i, 0, lim-1)div_ary[d+1][i] = 0;
}
int main(){
// freopen("in", "r", stdin);
rd(n), rd(m);
fp(i, 1, n)rd(a[i]), mul = 1ll*mul*a[i]%mod;
fac[0] = 1; fp(i, 1, n)fac[i] = 1ll*fac[i-1]*i%mod;
ifac[n] = fpow(fac[n], mod-2); fb(i, n, 1)ifac[i-1] = 1ll*ifac[i]*i%mod;
inv[1] = 1; fp(i, 2, n)inv[i] = mod-1ll*(mod/i)*inv[mod%i]%mod;
fp(i, 0, n)A[i] = 1ll*fpow(i+1, m<<1)*ifac[i]%mod;
fp(i, 0, n)B[i] = 1ll*fpow(i+1, m)*ifac[i]%mod;
get_ln(B, ln, n+1), get_inv(B, invB, n+1), init(2*n+1);
// fp(i, 0, n)printf("%d ", B[i]);puts("");
// fp(i, 0, n)printf("%d ", ln[i]);puts("");
// fp(i, 0, n)printf("%d ", invB[i]);puts("");
// fp(i, 0, n)printf("%d ", A[i]);puts("");
// fp(i, 0, n)printf("%d ", invB[i]);puts("");
ntt(invB, 0), ntt(A, 0);
fp(i, 0, lim-1)A[i] = 1ll*A[i]*invB[i]%mod;
ntt(A, 1);
fp(i, n+1, lim)A[i] = 0;
// fp(i, 0, n)printf("%d ", A[i]);puts("");
divntt(0, 1, n), get_ln(div_ary[0], H, n+1);
// fp(i, 0, n)printf("%d ", div_ary[0][i]);puts("");
fp(i, 1, n)H[i] = mod-1ll*H[i]*i%mod;
H[0] = n;
fp(i, 0, n)ln[i] = 1ll*ln[i]*H[i]%mod;
get_exp(ln, E, n+1);
fp(i, 0, n)A[i] = 1ll*A[i]*H[i]%mod;
init(2*n+1), ntt(E, 0), ntt(A, 0);
fp(i, 0, lim-1)A[i] = 1ll*A[i]*E[i]%mod;
ntt(A, 1), printf("%lld\n", 1ll*fac[n-2]*mul%mod*A[n-2]%mod);
return 0;
}
标签：fp,frac,ary,sum,计数,2017,prod,集训,mod 来源： https://www.cnblogs.com/akura/p/14488464.html
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。