一
稍微转化一下,就是找所有和原树差距不超过k的不同构树的个数
一个挺trick的想法是:
由于矩阵树定理的行列式的值是把邻接矩阵数值看做边权的图的所有生成树的边权乘积之和
那么如果把不存在于原树中的边的边权设为x,做矩阵树定理得到n-1次的多项式第i次项系数就是选择新选择i个边的方案数!
带着x不好做,x=1~n带入,然后插值即可
O(n^4)
二
开始碾标算了:
还是可以树形DP,经典的树形DP套路难办的原因是不知道干掉的子树接在哪里
所以我们干脆先不管接在哪里,先都砍断
一个公式:
https://rqy.moe/Solutions/bzoj5475/
证明导数第二步有点生成函数的思想了
对于拆成k个连通块,方案数可以直接算
设f[i][j][k],以i为根子树保留j条边,和i相连的连通块sz为k的方案数。枚举子树保留边数和sz就可以愉快转移了
但是会算重
如果直接乘上n^(k-2)可能会把断的边连回去,不止保留了j个,并且同一个方案可能因为保留的不同的j个计算多次!
不慌不慌
一切都是因为实际上可能连回去保留了不止j个,我们钦定了j个,其他随便选的方案数是f,
枚举具体保留了几个,可以二项式反演
$f[t]=sum_{k=t}^{n} C(k,t)*ans[k]$
ans表示恰好保留了k个
没了
O(n^3)
三
进一步优化!
$Pi a_i$有组合意义,就是对于k个连通块,我们每个连通块恰好选择一个点的方案数!
f[i][j][0/1],表示,,,,,i这个连通块内有没有选择点,总方案数
愉快转移,断边必须y选择。
代码
:
#include<bits/stdc++.h> #define reg register int #define il inline #define numb (ch^'0') using namespace std; typedef long long ll; il void rd(int &x){ char ch;x=0;bool fl=false; while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true); for(x=numb;isdigit(ch=getchar());x=x*10+numb); (fl==true)&&(x=-x); } namespace Miracle{ const int N=55; const int mod=998244353; int n,k; int qm(int x,int y){ int ret=1;while(y){ if(y&1) ret=(ll)ret*x%mod;x=(ll)x*x%mod; y>>=1; }return ret; } struct node{ int nxt,to; }e[2*N]; int hd[N],cnt; int ad(int x,int y){ return x+y>=mod?x+y-mod:x+y; } void add(int x,int y){ e[++cnt].nxt=hd[x]; e[cnt].to=y; hd[x]=cnt; } int sz[N]; int f[N][N][2]; int C[N][N]; ll ans[N],dp[N]; void dfs(int x,int fa){ sz[x]=1; f[x][0][1]=f[x][0][0]=1; int g[2][N][2]; memset(g[0],0,sizeof g[0]); //memset(g,0,sizeof g); g[0][0][1]=g[0][0][0]=1; int tmp=0; for(reg i=hd[x];i;i=e[i].nxt){ int y=e[i].to; if(y==fa) continue; dfs(y,x); tmp^=1; //for(reg j=0;j<sz[x]+sz[y];++j) g[tmp][j][0]=g[tmp][j][1]=0; memset(g[tmp],0,sizeof g[tmp]); for(reg j=0;j<sz[x];++j){ for(reg k=0;k<sz[y];++k){ g[tmp][j+k][0]=ad(g[tmp][j+k][0],(ll)g[tmp^1][j][0]*f[y][k][1]%mod); g[tmp][j+k][1]=ad(g[tmp][j+k][1],(ll)g[tmp^1][j][1]*f[y][k][1]%mod); g[tmp][j+k+1][0]=ad(g[tmp][j+k+1][0],(ll)g[tmp^1][j][0]*f[y][k][0]%mod); g[tmp][j+k+1][1]=ad(g[tmp][j+k+1][1],ad((ll)g[tmp^1][j][1]*f[y][k][0]%mod,(ll)g[tmp^1][j][0]*f[y][k][1]%mod)); } } sz[x]+=sz[y]; } //cout<<" f "<<x<<endl; for(reg j=0;j<sz[x];++j) { f[x][j][0]=g[tmp][j][0],f[x][j][1]=g[tmp][j][1]; //cout<<j<<" n0 "<<" : "<<f[x][j][0]<<" n1 "<<f[x][j][1]<<endl; } } int main(){ rd(n);rd(k); if(k==0){ puts("1");return 0; } int x; for(reg i=2;i<=n;++i){ rd(x);++x; add(i,x);add(x,i); } dfs(1,0); for(reg j=0;j<n;++j) { if(j!=n-1) dp[j]=(ll)f[1][j][1]*qm(n,n-j-2); else dp[j]=1; } ans[n-1]=dp[n-1]; C[0][0]=1; for(reg i=1;i<=n;++i){ C[i][0]=1; for(reg j=1;j<=n;++j){ C[i][j]=(C[i-1][j-1]+C[i-1][j])%mod; } } ll op=ans[n-1]; for(reg i=n-2;i>=n-k-1;--i){ // cout<<" ii "<<i<<" : "<<dp[i]<<endl; for(reg j=i+1;j<n;++j){ dp[i]=(dp[i]+mod-(ll)C[j][i]*ans[j]%mod)%mod; } ans[i]=dp[i]; op=(op+dp[i])%mod; } printf("%lld",op); return 0; } } signed main(){ Miracle::main(); return 0; } /* Author: *Miracle* Date: 2019/3/4 19:02:40 */
思路:
矩阵树:新加k条新边,为了控制边数,所以用x表示边权得到方案数,插值也很漂亮
DP:直接DP难写,受prufer序列的推论启发,考虑对于连通块直接统计,然后二项式反演去重,
组合意义的优化,类似于:[NOI2009]管道取珠