けんちょんの競プロ精進記録

競プロの精進記録や小ネタを書いていきます

プリューファーコード:全域木の数え上げ (頂点次数制約つき)

ARC 106 F に関連して、頂点次数制約のついた全域木の個数を求める問題がまさにあったので、その解説を。

問題概要 (New Year Contest 2015 E - ひも)

頂点数が  N であるような完全グラフの全域木であって、以下の条件を満たすものが何通りあるか、1000000007 で割ったあまりを求めよ。

  • 各頂点  v における次数が  a_{v} である

制約

  •  2 \le N \le 10^{5}
  •  1 \le a_{i} \le 3

 

解法

次数列の制約のない「全域木の数え上げ」は、Cayley の定理としてよく知られている!!


【Cayley の定理】
頂点数が  N であるような完全グラフの部分グラフであって、全域木をなすものの個数は

 N^{N-2}

で与えられる。


Cayley の定理の証明方法としては、たとえば以下の記事に 6 通りの方法が載っている。

joisino.hatenablog.com

 

次数制約つきの Cayley の定理

次数の指定があっても綺麗になる。


頂点数が  N であるような完全グラフの全域木であって、頂点  v の次数が  d_{v} ( d_{v} の総和が  2N-2) であるものの個数は

 \frac{(N-2)!}{(d_{1}-1)!(d_{2}-1)! \dots (d_{N}-1)!}

で与えられる。


むしろこれを補題として先に示すことで Cayley の定理を証明する方法も考えられる。補題がわかっていれば、Cayley の定理は次のように簡単に導出できる。

 \sum_{d_{i} \ge 1, d_{1} + \dots + d_{N} = 2N-2} \frac{(N-2)!}{(d_{1}-1)! \dots (d_{N}-1)!}
 = \sum_{d_{i} \ge 0, d_{1} + \dots + d_{N} = N-2} \frac{(N-2)!}{d_{1}! \dots d_{N}!}

の値が  N^{N-2} となることを示せばよい。そしてこの各項は、

 (x_{1} + \dots + x_{N})^{N-2}

の次数列が  d_{1}, \dots, d_{N} である項の係数と一致することに注意する。よって、 x_{1} = \dots = x_{N} = 1 とすることで計算できて、

 \sum_{d_{i} \ge 0, d_{1} + \dots + d_{N} = N-2} \frac{(N-2)!}{d_{1}! \dots d_{N}!} = N^{N-2}

となることが示された。

 

補題の証明

補題を示す。数学的帰納法を用いる。また、Cayley の定理の証明方法として有名な「プリューファーコード」による方法を活用する。

ブリューファーコードとは、全域木に対して、「各要素が  1, 2, \dots, N であるような数列  c_{1}, \dots, c_{N-2}」に一対一に対応付けるものである。全域木からプリューファーコードへの対応は次の通り。


  • 全域木の頂点数が 2 になるまで以下の操作を繰り返す (初期状態では数列は空とする)
    • 全域木の「葉」のうち、頂点番号が最小のものに隣接している頂点の番号を数列に末尾に追加する
    • その「葉」を削除する

下図は  N = 7 の木に対して、プリューファーコードを生成する一例を示している。

この「全域木」から「プリューファーコード」への写像は、「逆写像」が自然に作れる (後述) ので、全単射になっていることがわかる。このことから Cayley の定理 (全域木が  N^{N-2} 個) が直接導かれる。

さて、プリューファーコードを真似して、以下の場合に場合分けする

  • プリューファーコードの先頭の要素が頂点 1 のとき
  • プリューファーコードの先頭の要素が頂点 2 のとき
  • ...
  • プリューファーコードの先頭の要素が頂点 N のとき

数学的帰納法の仮定より、プリューファーコードの先頭要素が頂点  1 であるようなものの個数は

 \frac{(N-3)!}{(d_{1}-2)!(d_{2}-1)! \dots (d_{N}-1)!}

となる。一般に、先頭が頂点  v であるものの個数は、上の式の分母の積において、 v 番目のみが  (d_{v}-2)! の形をしていて、それ以外は  (d_{v} - 1)! の形をしているようなものになる。よって、分母を通分して総和をとると、

 \frac{(N-3)!( (d_{1}-1) + \dots + (d_{N}-1) )}{(d_{1}-1)! \dots (d_{N}-1)!}
 = \frac{(N-2)!}{(d_{1}-1)! \dots (d_{N}-1)!}

が導かれる。補題は示された。

 

プリューファーコードから木への対応

最後に、プリューファーコード  c_{1}, \dots, c_{N-2} から、頂点数  N の木を作るのは、次のようにしてできる。ちょうどいい感じの構築問題と言えそう。

  • まず、 c_{1}, c_{2}, \dots, c_{N-1}, c_{N-2} に含まれていない最小の値を  a_{1} として、2 頂点  a_{1}, c_{1} を結ぶ
  • 次に、 a_{1}, c_{2}, \dots, c_{N-1}, c_{N-2} に含まれていない最小の値を  a_{2} として、2 頂点  a_{2}, c_{2} を結ぶ
  • ...
  • 次に、 a_{1}, a_{2}, \dots, a_{N-1}, c_{N-2} に含まれていない最小の値を  a_{N-2} として、2 頂点  a_{N-2}, c_{N-2} を結ぶ
  • 最後に、 a_{1}, a_{2}, \dots, a_{N-2} に含まれていない 2 つの値を  d, e として、2 頂点  d, e を結ぶ

上図のように、プリューファーコードに登場する値が、ちょうど木の「内点」になることがわかる。そして上の手続きは

「指定された内点に対して、まだ使っていない最小番号の頂点を葉にしてくっつける」

という意味を持っていることがわかる。

 

コード

補題が示されれば、それを実装するだけ。ただし

 a_{1} + \dots + a_{N} = 2N-2

とならない場合は 0 通りとすることに注意。

#include <bits/stdc++.h>
using namespace std;

// modint
template<int MOD> struct Fp {
    long long val;
    constexpr Fp(long long v = 0) noexcept : val(v % MOD) {
        if (val < 0) val += MOD;
    }
    constexpr int getmod() const { return MOD; }
    constexpr Fp operator - () const noexcept {
        return val ? MOD - val : 0;
    }
    constexpr Fp operator + (const Fp& r) const noexcept { return Fp(*this) += r; }
    constexpr Fp operator - (const Fp& r) const noexcept { return Fp(*this) -= r; }
    constexpr Fp operator * (const Fp& r) const noexcept { return Fp(*this) *= r; }
    constexpr Fp operator / (const Fp& r) const noexcept { return Fp(*this) /= r; }
    constexpr Fp& operator += (const Fp& r) noexcept {
        val += r.val;
        if (val >= MOD) val -= MOD;
        return *this;
    }
    constexpr Fp& operator -= (const Fp& r) noexcept {
        val -= r.val;
        if (val < 0) val += MOD;
        return *this;
    }
    constexpr Fp& operator *= (const Fp& r) noexcept {
        val = val * r.val % MOD;
        return *this;
    }
    constexpr Fp& operator /= (const Fp& r) noexcept {
        long long a = r.val, b = MOD, u = 1, v = 0;
        while (b) {
            long long t = a / b;
            a -= t * b, swap(a, b);
            u -= t * v, swap(u, v);
        }
        val = val * u % MOD;
        if (val < 0) val += MOD;
        return *this;
    }
    constexpr bool operator == (const Fp& r) const noexcept {
        return this->val == r.val;
    }
    constexpr bool operator != (const Fp& r) const noexcept {
        return this->val != r.val;
    }
    friend constexpr istream& operator >> (istream& is, Fp<MOD>& x) noexcept {
        is >> x.val;
        x.val %= MOD;
        if (x.val < 0) x.val += MOD;
        return is;
    }
    friend constexpr ostream& operator << (ostream& os, const Fp<MOD>& x) noexcept {
        return os << x.val;
    }
    friend constexpr Fp<MOD> modpow(const Fp<MOD>& r, long long n) noexcept {
        if (n == 0) return 1;
        if (n < 0) return modpow(modinv(r), -n);
        auto t = modpow(r, n / 2);
        t = t * t;
        if (n & 1) t = t * r;
        return t;
    }
    friend constexpr Fp<MOD> modinv(const Fp<MOD>& r) noexcept {
        long long a = r.val, b = MOD, u = 1, v = 0;
        while (b) {
            long long t = a / b;
            a -= t * b, swap(a, b);
            u -= t * v, swap(u, v);
        }
        return Fp<MOD>(u);
    }
};

const int MOD = 1000000007;
using mint = Fp<MOD>;

// Binomial coefficient
template<class T> struct BiCoef {
    vector<T> fact_, inv_, finv_;
    constexpr BiCoef() {}
    constexpr BiCoef(int n) noexcept : fact_(n, 1), inv_(n, 1), finv_(n, 1) {
        init(n);
    }
    constexpr void init(int n) noexcept {
        fact_.assign(n, 1), inv_.assign(n, 1), finv_.assign(n, 1);
        int MOD = fact_[0].getmod();
        for(int i = 2; i < n; i++){
            fact_[i] = fact_[i-1] * i;
            inv_[i] = -inv_[MOD%i] * (MOD/i);
            finv_[i] = finv_[i-1] * inv_[i];
        }
    }
    constexpr T com(int n, int k) const noexcept {
        if (n < k || n < 0 || k < 0) return 0;
        return fact_[n] * finv_[k] * finv_[n-k];
    }
    constexpr T fact(int n) const noexcept {
        if (n < 0) return 0;
        return fact_[n];
    }
    constexpr T inv(int n) const noexcept {
        if (n < 0) return 0;
        return inv_[n];
    }
    constexpr T finv(int n) const noexcept {
        if (n < 0) return 0;
        return finv_[n];
    }
};

BiCoef<mint> bc;

int main() {
    int N;
    cin >> N;
    bc.init(N);
    long long sum = 0;
    mint res = bc.fact(N-2);
    for (int i = 0; i < N; ++i) {
        int a; cin >> a;
        res *= bc.finv(a-1);
        sum += a-1;
    }
    cout << (sum == N-2 ? res : 0) << endl;
}