1. 程式人生 > >[2018.6.22集訓]admirable-拆系數FFT-多項式相關

[2018.6.22集訓]admirable-拆系數FFT-多項式相關

spa 兩個 pro 處理 -i 同時 需要 n) ble

題目大意

給出一棵樹,現在需要將$k$條不相同的路徑覆蓋到這棵樹上。
定義一種合法的路徑覆蓋方案為,能使得樹上的每條邊的被覆蓋次數$t \in {0,1,k}$的方案。
求合法方案的數量,對$10^9 +9$取模。

$n,k \leq 10^5$。

題解

很容易想到一個計算答案的方法:
被覆蓋$k$次的部分一定是一條鏈,枚舉這條鏈,並對每條枚舉的鏈,計算可行的方案數並累加。

對於每條覆蓋$k$次的鏈,可以發現,如果分別計算出兩個端點處的方案,將其貢獻乘起來,即可得到這條鏈的方案數。
可以發現,對於一個端點處的合法方案,每一條路徑在這一側的端點均在這條鏈的端點的不同子樹中或恰好在這個端點上,因為若有兩條路徑在同一子樹中,會經過同一條邊多次,方案會不合法。

考慮如何計算這個方案數,可以使用多項式乘法解決,構造形如$size[v]*x+1$的多項式,並將每個子樹的多項式乘起來,再取出前$k$項的系數,采用全排列和組合數計算貢獻。

因此,每個點所$u$需要的多項式,為$\prod\limits_{(u,v) \in E}size[v]*x+1$,除掉與當前枚舉的鏈相交的一側的兒子的多項式得到。
乘或除掉一個二項式是$O(n)$的,同時如果預處理,乘除操作的次數將為$O(n)$,因此預處理一個節點處的多項式再枚舉,復雜度$O(n^2)$。

考慮優化。
對於計算初始每個節點處所有兒子乘起來得到的多項式,考慮分治FFT,$O(n\log^2n)$。
考慮除掉二項式的復雜度,可以發現對於兩個$size$相同的多項式,得到的結果相同,而一個節點的子樹的不同的$size$的數量為$O(\sqrt{n})$,於是改為對子樹中每種出現了的$size$預處理,復雜度降至$O(n\sqrt{n})$。

最後,對於一個節點,同一棵子樹內的節點與它組成的鏈,得到的這個點的貢獻相同,因此可以一起算,復雜度降至$O(n)$,具體可見代碼。

於是最後復雜度為$O(n\log^2n+n\sqrt{n})$。
由於模數很惡心,需要使用拆系數FFT,這裏使用了三模數+中國剩余定理合並的版本。

代碼:

#include<map>
#include<cstdio>
#include<algorithm>
using namespace std;

typedef long long ll;
const int N=4e5+9;
const int md=1e9+9;
const
int K=33; inline int read() { int x=0;char ch=getchar(); while(ch<'0' || '9'<ch)ch=getchar(); while('0'<=ch && ch<='9')x=x*10+(ch^48),ch=getchar(); return x; } inline void chk(ll &x){if(x>=md)x-=md;} inline ll qpow(ll a,ll b) { ll ret=1; while(b) { if(b&1)ret=ret*a%md; a=a*a%md;b>>=1; } return ret; } int n,k; int to[N<<1],nxt[N<<1],beg[N],tot; int fa[N],deg[N],siz[N],len[N]; ll f[N],t[N],sumf[N]; map<int,ll> g[N]; ll fac[N],inv[N],ans; inline ll c(ll a,ll b) { if(a<b)return 0; return fac[a]*inv[b]%md*inv[a-b]%md; } inline void init() { fac[0]=1; for(ll i=1;i<N;i++) fac[i]=fac[i-1]*i%md; inv[N-1]=qpow(fac[N-1],md-2); for(ll i=N-1;i>=1;i--) inv[i-1]=inv[i]*i%md; for(ll i=0;i<=k;i++) t[i]=fac[i]*c(k,i)%md; } inline void add(int u,int v) { to[++tot]=v; nxt[tot]=beg[u]; beg[u]=tot; deg[v]++; } inline ll mul(ll *f,ll x,int &l) { f[++l]=0; for(int i=l-1;i>=0;i--) chk(f[i+1]+=f[i]*x%md); } inline ll imul(ll *f,ll x,int &l) { ll invv=qpow(x,md-2); for(int i=l-1;i>=0;i--) { f[i+1]=f[i+1]*invv%md; chk(f[i]=f[i]+md-f[i+1]); } for(int i=0;i<l;i++) f[i]=f[i+1]; l--; } inline ll calc(ll *f,int l) { ll ret=0; for(int i=0;i<=l;i++) chk(ret+=f[i]*t[i]%md); return ret; } namespace ntt { const ll md1=998244353; const ll md2=1004535809; const ll md3=469762049; const ll M=md1*md2; int rev[N]; ll b[K][N],len[K],seg[N],top; const ll muls(ll a,ll b,ll p) { return ((a*b-(ll)((long double)a*b/p)*p)%p+p)%p; } inline void initrev(int n) { for(int i=0;i<n;i++) rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1)); } inline ll qpow(ll a,ll b,ll p) { ll ret=1;a%=p; while(b) { if(b&1)ret=ret*a%p; a=a*a%p;b>>=1; } return ret; } inline void ntt(ll *a,int n,int f,ll md) { for(int i=0;i<n;i++)if(i<rev[i])swap(a[i],a[rev[i]]); for(int h=2;h<=n;h<<=1) { ll w=qpow(3,(md-1)/h,md); if(f)w=qpow(w,md-2,md); for(int j=0;j<n;j+=h) { ll wn=1ll,x,y; for(int k=j;k<j+(h>>1);k++) { x=a[k],y=a[k+(h>>1)]*wn%md;wn=wn*w%md; a[k]=(x+y)%md;a[k+(h>>1)]=(x-y+md)%md; } } } if(f) for(ll i=0,invn=qpow(n,md-2,md);i<n;i++) a[i]=a[i]*invn%md; } inline void cp(ll *a,ll *b,ll n,ll m) { for(int i=0;i<n;i++) b[i]=a[i]; for(int i=n;i<m;i++) b[i]=0; } inline void bfmul(ll *a,int n,ll *b,int m,ll *c) { static ll d[N]; for(int i=0;i<n+m-1;i++) d[i]=0; for(int i=0;i<n;i++) for(int j=0;j<m;j++) chk(d[i+j]+=a[i]*b[j]%md); for(int i=0;i<n+m-1;i++) c[i]=d[i]; } inline void mul(ll *a,int n,ll *b,int m,ll *c) { static ll d[N],e[N],f[N],g[N],l; if(n+m-1<=1000) { bfmul(a,n,b,m,c); return; } for(l=1;l<=n+m;l<<=1);initrev(l); cp(a,d,n,l);cp(b,e,m,l); ntt(d,l,0,md1);ntt(e,l,0,md1); for(int i=0;i<l;i++) f[i]=d[i]*e[i]%md1; ntt(f,l,1,md1); cp(a,d,n,l);cp(b,e,m,l); ntt(d,l,0,md2);ntt(e,l,0,md2); for(int i=0;i<l;i++) g[i]=d[i]*e[i]%md2; ntt(g,l,1,md2); ll inv2=qpow(md2,md1-2,md1); ll inv1=qpow(md1,md2-2,md2); for(int i=0;i<n+m-1;i++) { f[i]=muls(f[i]*md2%M,inv2,M); g[i]=muls(g[i]*md1%M,inv1,M); f[i]=(f[i]+g[i])%M; } cp(a,d,n,l);cp(b,e,m,l); ntt(d,l,0,md3);ntt(e,l,0,md3); for(int i=0;i<l;i++) g[i]=d[i]*e[i]%md3; ntt(g,l,1,md3); for(int i=0;i<n+m-1;i++) { ll k=((g[i]-f[i])%md3+md3)%md3*qpow(M,md3-2,md3)%md3; f[i]=(k%md*(M%md)%md+f[i]%md)%md; } for(int i=0;i<n+m-1;i++) c[i]=f[i]; } inline void work(int l,int r) { if(l==r) { len[++top]=2; b[top][0]=1; b[top][1]=seg[l]; return; } int mid=l+r>>1; work(l,mid); work(mid+1,r); mul(b[top-1],len[top-1],b[top],len[top],b[top-1]); len[top-1]=len[top-1]+len[top]-1; top--; } } inline void dfs_pre(int u) { siz[u]=1;sumf[u]=0;int cnt=0; for(int i=beg[u];i;i=nxt[i]) if(to[i]!=fa[u]) { fa[to[i]]=u; dfs_pre(to[i]); siz[u]+=siz[to[i]]; chk(sumf[u]+=sumf[to[i]]); } len[u]=0;f[0]=1; for(int i=beg[u];i;i=nxt[i]) if(to[i]!=fa[u]) ntt::seg[++cnt]=siz[to[i]]; if(u) ntt::seg[++cnt]=n-siz[u]; ntt::top=0; ntt::work(1,cnt); len[u]=ntt::len[1]-1; for(int i=0;i<=len[u];i++) f[i]=ntt::b[1][i]; for(int i=beg[u];i;i=nxt[i]) if(to[i]!=fa[u] && !g[u].count(siz[to[i]])) { imul(f,siz[to[i]],len[u]); g[u][siz[to[i]]]=calc(f,len[u]); mul(f,siz[to[i]],len[u]); } if(!g[u].count(n-siz[u])) { imul(f,n-siz[u],len[u]); g[u][n-siz[u]]=calc(f,len[u]); mul(f,n-siz[u],len[u]); } chk(sumf[u]+=g[u][n-siz[u]]); } namespace chain { int mina() { ans=c(n,2); for(int i=1;i<=n;i++) { chk(ans+=c(k,1)*(i-1)%md*(n-i)%md); chk(ans+=c(k,1)*c(n-i,2)%md*(1+c(k,1)*(i-1)%md)%md); } printf("%lld\n",ans); return 0; } } inline void dfs_calc(int u,ll sum) { if(fa[u]) chk(ans+=g[u][n-siz[u]]*sum%md); for(int i=beg[u];i;i=nxt[i]) if(to[i]!=fa[u]) { chk(ans+=(g[u][siz[to[i]]]*sumf[to[i]])%md); sum+=sumf[to[i]]; } for(int i=beg[u];i;i=nxt[i]) if(to[i]!=fa[u]) dfs_calc(to[i],(sum-sumf[to[i]]+g[u][siz[to[i]]]+md)%md); } int main() { n=read();k=read(); for(int i=2,u,v;i<=n;i++) { u=read();v=read(); add(u,v);add(v,u); } init(); for(int i=1;i<=n;i++) if(deg[i]>2) goto hell; return chain::mina(); hell:; dfs_pre(1); dfs_calc(1,ans=0); printf("%lld\n",ans*qpow(2,md-2)%md); return 0; }

[2018.6.22集訓]admirable-拆系數FFT-多項式相關