スポンサーサイト

上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。

SRM 551 Div.1 Hard : SweetFruits

SRM の記事はいつもは書かないんだけど、公式の editorial がない ( 遅いだけ? ) のと、解法がきれいだったので書くことにした。

Tags : プログラミング TopCoder SRM

問題

n 個のフルーツがあり、それぞれ甘いか苦いか。
甘いフルーツについては、甘さが 0 以上の整数値で与えられる。
二つのフルーツを辺でつなぐことを n-1 回くり返して全域木を作る。
木の甘みとは、「少なくとも一つの甘いフルーツと隣接するフルーツ」の甘みの総和のこととする。
ただし、苦いフルーツは甘み 0 とみなす。
木の甘みが maxSweetness を超えないような全域木の作り方は何通りあるか? mod 10^9+7 で求めよ。

1 ≦ n ≦ 40

解答

木の甘みの定義に注意すると、全域木は
( 木の甘みに寄与するフルーツ )-( 苦いフルーツ )-( 木の甘みに寄与しないフルーツ )
という形になる。-はいくつかの辺の意味。
( 木の甘みに寄与するフルーツ ) を左の頂点、( 木の甘みに寄与しないフルーツ ) を右の頂点と呼ぶことにする。

とりあえず、左にどの頂点を割り振ったかを固定しよう。
L 個の頂点を割り振ったとする。
右の頂点数を R とする。これは R = n - L - ( 苦い頂点数 ) で求まる。
全域木の個数を求めるのに行列木定理を用いる。 ( 昔、紹介記事を書いた。 )
そのために、全域木を求めたいグラフを構成する。
・左の頂点どうし
・苦い頂点どうし
・左の頂点と苦い頂点
・右の頂点と苦い頂点
の間には辺を張る。
右の頂点どうしに辺を張らないのは、右の頂点は木の甘みに寄与してはいけないから。
このようにしてできたグラフの全域木の個数を求める。
でも、これだと余分なものまでカウントしてしまうことに気付く。
たとえば、n = 4, L = 2, R = 1 のとき、全域木は

という 3 種類があるけど、一番右のものは「左の頂点が甘みに寄与する」という条件を満たさない。

うまく条件に合うものだけを数え上げるためには一工夫いる。
num[L] を、左の頂点数が L のときの条件に合う全域木の個数とする。
num[L] = ( 行列木定理で得た値 ) - Σ[j=1 to L] LCj num[L-j]
が成立。この漸化式を用いると L=0 から順に計算できる。
行列木定理で出てきた全域木は、左の頂点が苦い頂点とだけつながっているケースがまずい。
なので、左の頂点のうち何個が苦い頂点とだけつながっているか ( この個数が j ) を各 j について調べる。
「左の頂点のうち j 個が苦い頂点とだけつながっている」全域木の個数は、j 個の左の頂点をないものと思ってカウントした値 num[L-j] に等しい。j 個の頂点の選び方は LCj 通りあるので、この漸化式が正しいことがわかる。

あとは、左に L 個の頂点を割り振る方法が何通りあるかがわかればよい。
n ≦ 40 なので反射的に半分全列挙したくなる。
半分全列挙したあと、素直に二分探索で組み合わせると、2^20・20・log(2^20) くらいかかりそうでやばい。
でも、定数倍改善すると最悪ケースで 1.9 sec くらいでぎりぎり間に合った。
しゃくとりとかすれば log が消えるような気もする。 ( よく考えてない )

ソースコード

#include<vector>
#include<algorithm>

#define rep(i,n) for(int i=0;i<(n);i++)

using namespace std;

typedef long long ll;

const ll M=1000000007;

const int N_BIN=100;
ll nCr[N_BIN+1][N_BIN+1];
void binom(){
rep(n,N_BIN+1) nCr[n][0]=1;
rep(n,N_BIN) rep(r,n+1) nCr[n+1][r+1]=(nCr[n][r+1]+nCr[n][r])%M;
}

ll xgcd(ll a,ll b,ll &x,ll &y){
if(b==0){ x=1; y=0; return a; }
ll g=xgcd(b,a%b,y,x); y-=a/b*x;
return g;
}

ll modinv(ll a,ll m){
ll x,y;
if(xgcd(a,m,x,y)==1) return (x+m)%m;
return -1;
}

const int N_MAX=40;
int moddet(int n,int A[N_MAX][N_MAX],int M){
rep(i,n) rep(j,n) {
A[i][j]%=M;
if(A[i][j]<0) A[i][j]+=M;
}

int res=1;
rep(i,n){
int pivot=-1;
for(int k=i;k<n;k++) if(A[k][i]!=0) pivot=k;
if(pivot==-1) return 0;

if(pivot!=i){
rep(j,n) swap(A[i][j],A[pivot][j]);
res=(ll)res*(M-1)%M;
}

for(int k=i+1;k<n;k++) {
ll c=A[k][i]*modinv(A[i][i],M)%M;
for(int j=i;j<n;j++){
A[k][j]-=c*A[i][j]%M;
if(A[k][j]<0) A[k][j]+=M;
}
}
}
rep(i,n) res=(ll)res*A[i][i]%M;

return res;
}

const int V_MAX=N_MAX;
int matrix_tree(int n,const bool G[V_MAX][V_MAX],int M){
int deg[V_MAX]={};
rep(u,n) rep(v,n) if(G[u][v]) deg[v]++;

static int H[V_MAX][V_MAX];
rep(u,n-1) rep(v,n-1) H[u][v]=(u==v?deg[u]:(G[u][v]?-1:0));

return moddet(n-1,H,M);
}

class SweetFruits{
public:
int countTrees(vector<int> sweet,ll maxsweet){
binom();

int n=sweet.size(),n_bitter=0;
rep(i,n) if(sweet[i]==-1) {
n_bitter++;
sweet.erase(sweet.begin()+i);
n--;
i--;
}

// 苦い頂点がないときは例外処理 ( アルゴリズム的には以下のものでまとめて扱えるけど、ぎりぎり TLE する )
if(n_bitter==0){
ll sum=0;
rep(i,n) sum+=sweet[i];
ll res=0;
if(sum<=maxsweet){
res=1;
rep(i,n-2) res=res*n%M;
}
return res;
}

int n1=(n+1)/2,n2=n-n1;
vector<ll> harf[21];
rep(S,1<<n1){
int pc=0;
ll sum=0;
rep(i,n1) if(S&1<<i) pc++, sum+=sweet[i];
if(sum<=maxsweet) harf[pc].push_back(sum);
}
rep(i,n1+1) sort(harf[i].begin(),harf[i].end());

ll cnt[41]={}; // cnt[i] := ( 苦くない頂点たちのサイズ i の部分集合で、甘みの和が maxsweet 以下になるようなものの個数 )
rep(S,1<<n2){
int pc=0;
ll sum=0;
rep(i,n2) if(S&1<<i) pc++, sum+=sweet[n1+i];
rep(i,n1+1){
cnt[pc+i]+=upper_bound(harf[i].begin(),harf[i].end(),maxsweet-sum)-harf[i].begin();
cnt[pc+i]%=M;
}
}

// 左の頂点は甘みの総和に寄与する ( 別の甘い頂点とくっついている )
// 右の頂点は甘みの総和に寄与しない ( 甘くない頂点とだけくっついている )
ll num[41]; // num[i] := ( 左の頂点に i 個, 右の頂点に n-i 個を割り振ったときの全域木の個数 )
rep(i,n+1){ // i : 左の頂点数, n-i : 右の頂点数
// グラフを構築
bool G[40][40]={};
// 右の頂点-苦い頂点の間に辺を張る
for(int u=i;u<n;u++) for(int v=n;v<n+n_bitter;v++) G[u][v]=G[v][u]=true;
// 左の頂点-左の頂点間に辺を張る
rep(u,i) rep(v,i) G[u][v]=true;
// 左の頂点-苦い頂点間に辺を張る
rep(u,i) for(int v=n;v<n+n_bitter;v++) G[u][v]=G[v][u]=true;
// 苦い頂点-苦い頂点間に辺を張る
for(int u=n;u<n+n_bitter;u++) for(int v=n;v<n+n_bitter;v++) G[u][v]=true;
// 自己ループを消す
rep(u,n+n_bitter) G[u][u]=false;

num[i]=matrix_tree(n+n_bitter,G,M);
rep(j,i){
num[i]-=nCr[i][i-j]*num[j]%M; // 左の頂点のうち i-j 個がソロで苦い頂点にくっついているケースを除く
if(num[i]<0) num[i]+=M;
}
}

ll ans=0;
rep(i,n+1) ans=(ans+cnt[i]*num[i])%M; // 左の頂点に i 個割り当てる方法は cnt[i] 通りある ( 残りは自動的に左の頂点 )
return ans;
}
};

スポンサーサイト

コメントの投稿

非公開コメント

プロフィール

fura2

Author : fura2
数学・コンピュータを中心に、考えたこと・やったことを書いていきます。

誤植等を含め、間違いはご指摘いただければ幸いです。

FC2カウンター
検索フォーム
最新記事
最新コメント
最新トラックバック
月別アーカイブ
リンク
RSSリンクの表示
上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。