AtCoder Regular Contest 018 D - 僕は友達が少ない

問題

arc018.contest.atcoder.jp

解法

解説を読みました。

AtCoder Regular Contest 018 解説 from AtCoder Inc.

www.slideshare.net

要するに「最小全域木のコストとそれを構成する場合の数を求めよ」という問題です。

解くための前提として, 行列木定理というのを知っていなければなりません。

各行，列がそれぞれ頂点に対応しており，対角成分にはその頂点の次数，非対角成分については枝がある部分に−1，ない部分に 0 を格納した |V|×|V| の正方行列をラプラシアン行列（グラフラプラシアン）と言います。で, 任意のグラフ G について, G の全域木の個数は, G に対応したラプラシアン行列 L の余因子の値に等しいという定理です(余因子は任意のものを選んで良い。下のコードでは (|V|-1, |V|-1) に関する余因子をとっている)。

これを使って問題を解きましょう。最小全域木は, 「コストが小さい順に辺を見ていき, その辺を付け加えても閉路が生じないなら貪欲に辺を加えていく」というアルゴリズムで構成できました。ということは, あるコストまでの辺を見た場合, 生成されているグラフの連結状況は一定であるはずです。

コストが小さい順に, 「辺のコストが c の辺の集合 G」というのを考えます。G の各連結成分についてラプラシアン行列を作ると, 全域木を作る場合の数を各連結成分で数えてくれます。
辺で連結になった頂点同士はひとつのまとまった頂点にして, この作業を連結成分が 1 つになるまで繰り返せば良いです。

typedef long long number;
typedef vector<number> vec;
typedef vector<vec> matrix;

const ll MOD = 1e9+7;

ll mod_pow(ll x, ll p, ll MOD) {
    ll a = 1;
    while (p) {
        if (p%2) a = a*x%MOD;
        x = x*x%MOD;
        p/=2;
    }
    return a;
}

// mod_inverse
ll mod_inverse(ll a, ll m) {
    return mod_pow(a, MOD-2, m);
}

// long long 専用 行列式を求める関数
number det(matrix A, const ll MOD) {
    const int n = A.size();
    assert(n == (int)A[0].size());
    ll ans = 1;
    for (int i = 0; i < n; i++) {
        int pivot = -1;
        for (int j = i; j < n; j++) if (A[j][i]) {
            pivot = j;
            break;
        }
        if (pivot == -1) return 0;
        if (i!=pivot) {
            swap(A[i], A[pivot]);
            ans *= -1;
        }
        ll inv = mod_inverse(A[i][i], MOD);
        for (int j = i+1; j < n; j++) {
            ll c = A[j][i]*inv % MOD;
            for (int k = i; k < n; k++) {
                A[j][k] = (A[j][k] - c*A[i][k])%MOD;
            }
        }
        (ans *= A[i][i]) %= MOD;
    }
    if (ans < 0) ans += MOD;
    return ans;
}

struct UnionFind {
    vector<int> par;
    int n, cnt;
    UnionFind(const int& x = 0) {init(x);}
    void init(const int& x) {par.assign(cnt=n=x, -1);}
    inline int find(const int& x) {return par[x] < 0 ? x : par[x] = find(par[x]);}
    inline bool same(const int& x, const int& y) {return find(x) == find(y);}
    inline bool unite(int x, int y) {
        if ((x = find(x)) == (y = find(y))) return false;
        --cnt;
        if (par[x] > par[y]) swap(x, y);
        par[x] += par[y];
        par[y] = x;
        return true;
    }
    inline int count() const {return cnt;}
    inline int count(int x) {return -par[find(x)];}
};

const int MAX = 100010;
vector<pii> Es[MAX];
vector<int> G[200];
int edge[200][200], deg[200], used[200];
vector<int> ap;

void dfs(int v) {
    used[v] = 1;
    ap.push_back(v);
    for (int u : G[v]) if (!used[u]) dfs(u);
}

int main() {
    int N, M;
    cin >> N >> M;
    UnionFind uf(N);
    for (int i = 0; i < M; i++) {
        int a, b, c;
        cin >> a >> b >> c;
        a--; b--;
        Es[c].emplace_back(a, b);
    }
    ll sum = 0, way = 1;
    for (int c = 1; c < MAX; c++) {
        vector<int> vs;
        for (pii& p : Es[c]) {
            p.first = uf.find(p.first);
            p.second = uf.find(p.second);
            if (p.first == p.second) continue;
            vs.push_back(p.first); vs.push_back(p.second);
        }
        if (vs.size() == 0) continue;
        sort(vs.begin(), vs.end());
        vs.erase(unique(vs.begin(), vs.end()), vs.end());
        int size = vs.size();
        for (int i = 0; i < size; i++) {
            G[i].clear();
            deg[i] = used[i] = 0;
            for (int j = 0; j < size; j++) edge[i][j] = 0;
        }
        for (pii p : Es[c]) {
            if (p.first == p.second) continue;
            if (!uf.same(p.first, p.second)) {
                uf.unite(p.first, p.second);
                sum += c;
            }
            int a = lower_bound(vs.begin(), vs.end(), p.first) - vs.begin();
            int b = lower_bound(vs.begin(), vs.end(), p.second) - vs.begin();
            G[a].push_back(b); G[b].push_back(a);
            deg[a]++; deg[b]++;
            edge[a][b]++; edge[b][a]++;
        }
        for (int i = 0; i < size; i++) if (!used[i]) {
            ap.clear();
            dfs(i);
            int L = ap.size();
            matrix mat(L-1, vec(L-1));
            for (int j = 0; j < L-1; j++) for (int k = 0; k < L-1; k++) {
                if (j==k) mat[j][k] = deg[ap[j]];
                else mat[j][k] -= edge[ap[j]][ap[k]];
            }
            (way *= det(mat, MOD)) %= MOD;
        }
    }
    cout << sum << " " << way << endl;
}