これは最後の一押しがしんどかった。
https://atcoder.jp/contests/diverta2019/tasks/diverta2019_f
問題
N頂点M辺の連結無向グラフが与えられる。
M辺に1~Mの重さを1本ずつ割り当てるM!通りのうち、最初の(N-1)本が最小全域木を成すようなケースの重みの総和を求めよ。
解法
最初に、最初の(N-1)本のうちいくつかの辺がすでに最小全域木に含まれることがわかっている場合、それ以上追加しても連結状態が変わらないような辺の数を列挙しておこう。
f(mask) := 最小全域木に含まれることが確定した辺の集合がmaskで表されるとして、それ以上追加しても連結状態が変化しないような辺の数の和
これは高速ゼータ変換の要領でO(N*2^N)で求められる。
最小全域木を成す(N-1)本内における重さ順が定まっているとき、その組み合わせと条件を満たす重みの総和を求めることを考えよう。
これを応用すればbitdpの要領で、重さ順を(N-1)!通り探索しなくてもO(N*2^N)で全並び順を列挙できる。
元のグラフの辺を、最小全域木を成すべき(N-1)本を黒頂点、残りを白頂点とみなした無向グラフを考える。
(N-1)本の重さが定まっている前提なので、黒頂点は軽い辺から重い辺に順次有向辺を張ったとする。
仮にここでは頂点番号を1~(N-1)としておく。
1~mの黒頂点に対応する辺が最小全域木に含まれることが確定して、初めて追加しても連結状態に変化を与えない(=m個の黒頂点群以上の重さであれば最小全域木に含まれない)辺に対応する白頂点は、mから辺を張ることにしよう。
そのような頂点数は、f(1~m)-f(1~(m-1))となる。
こうしてM頂点のDAGができる。求めるものが「並び順の組み合わせ」であれば、これはDAGのトポロジカルソート順である。
しかしここでは求めるのは黒頂点の重さの総和である。
重い順に頂点を確定されることを考える。
m番の黒頂点の重さを確定するには、(m+1)番の黒頂点と、m番につながる白頂点の重さを先に確定させなければならない。
そこで、
(N-1)番につながる白頂点群→(N-1)番の黒頂点→(N-2)番につながる白頂点群→(N-2)番の黒頂点…
の順で挿入DPの要領で重さを確定させていくことを考えよう。
状態を以下のtupleで表すとする。
(n,b,c,s) := n個目までの頂点の重さ順が確定したとき、黒頂点はb個含まれていて、その組み合わせはc通りであり、黒頂点の重みの総和はsである
それぞれの頂点を追加するとき、状態は以下の通り遷移する。
- 白頂点を追加するとき、どこに追加してもよいのでその位置は(n+1)通り。よって(n.b.c.s)→(n+1,b,(n+1)c,(n+2)s)と状態遷移する。
- sが(n+2)sになるのは直感的じゃないので説明を加える。まず白頂点の追加の組み合わせが(n+1)通りなので、総和は(n+1)sになるのはわかる。残りsはどこから来るのか。
- 白頂点を手前に挿入することで、既存の黒頂点の重みが後ろにずれる可能性があるので、その分を考慮しなければならない。仮に今黒頂点がi番目にあるとき、白頂点の挿入位置のうち黒頂点が後ろにずれるのはi通りである。よって全体ではこのようなiの総和分ずれるわけだが、その値は定義を考えるとそもそもsと一致するので、s追加される。
- 黒頂点は、確定済みのどの頂点よりも軽くないといけないので、(n.b.c.s)→(n+1,b+1,c,s+c*(b+1))と状態遷移する。
あとは上記処理で最終的なsを答えればよい。
int N,M; int A[400],B[400]; vector<int> E[20]; vector<int> path[20][20]; int rev[20][20]; int ne[1<<20]; ll C[1<<20],S[1<<20]; ll mo=1000000007; const int NUM_=400001; static ll fact[NUM_+1],factr[NUM_+1],inv[NUM_+1]; ll modpow(ll a, ll n = mo-2) { ll r=1;a%=mo; while(n) r=r*((n%2)?a:1)%mo,a=a*a%mo,n>>=1; return r; } ll comb(ll N_, ll C_) { if(C_<0 || C_>N_) return 0; return factr[C_]*fact[N_]%mo*factr[N_-C_]%mo; } void solve() { int i,j,k,l,r,x,y; string s; inv[1]=fact[0]=factr[0]=1; for (int i=2;i<=NUM_;++i) inv[i] = inv[mo % i] * (mo - mo / i) % mo; for (int i=1;i<=NUM_;++i) fact[i]=fact[i-1]*i%mo, factr[i]=factr[i-1]*inv[i]%mo; MINUS(rev); cin>>N>>M; FOR(i,M) { cin>>A[i]>>B[i]; A[i]--; B[i]--; if(i<N-1) { E[A[i]].push_back(B[i]); E[B[i]].push_back(A[i]); } rev[A[i]][B[i]]=rev[B[i]][A[i]]=i; } FOR(i,N) { path[i][i].push_back(i); queue<int> Q; Q.push(i); while(Q.size()) { x=Q.front(); Q.pop(); FORR(e,E[x]) if(path[i][e].empty()) { path[i][e]=path[i][x]; path[i][e].push_back(e); Q.push(e); } } } FOR(i,M) { int mask=0; FOR(j,path[A[i]][B[i]].size()-1) mask |= 1<<(rev[path[A[i]][B[i]][j]][path[A[i]][B[i]][j+1]]); ne[mask]++; } FOR(i,N-1) FOR(j,1<<(N-1)) if(j&(1<<i)) ne[j]+=ne[j^(1<<i)]; C[(1<<(N-1))-1]=1; int mask; for(int mask=(1<<(N-1))-1;mask>0;mask--) { int num=(N-1)-__builtin_popcount(mask); FOR(j,N-1) if((mask&(1<<j))) { int CN=M-ne[mask]; int W=ne[mask]-ne[mask^(1<<j)]-1; ll wc=C[mask]*fact[CN+W]%mo*factr[CN]%mo; ll ws=S[mask]*fact[CN+W+1]%mo*factr[CN+1]%mo; (C[mask^(1<<j)]+=wc)%=mo; (S[mask^(1<<j)]+=ws+(num+1)*wc)%=mo; } } cout<<S[0]<<endl; }
まとめ
本番、高速ゼータ変換とかは思いついたのに、sの状態遷移が思いつかなかった。