Maximum closure problem

Maximum closure problem の学習まとめ。以下のサイトを参考にした。

Closure problem - Wikipedia
http%3A%2F%2Friot.ieor.berkeley.edu%2F~dorit%2Fpub%2Fscribe%2Flec11%2FLec11posted.pdf&usg=AFQjCNEhubnfNuotTV0_wrbkGGdEF9UHMA&sig2=ewamZFzAoKgnMP5FaFV5cA&bvm=bv.151325232,d.amc *1

Maximum closure problemとは、ノードに重み(正負あり)の付いた有向グラフがあたえられたとき、重みの合計値が最大になるグラフのClosedなサブセットを求める問題である。
Closedのグラフとは、グラフの全ノードについて、その下流ノードがすべてそのグラフに属するもの。たとえば下の図だと、(a)はClosedだが(b)はClosedではない。
f:id:yambe2002:20170331101834p:plain
(*1より抜粋)


求め方
次のようなS-Tグラフを作成する。
・Sから正の重みがついたノードへ、その重み付きのエッジを追加
・Tから負の重みがついたノードへ、その重み付きのエッジを追加
・それ以外のエッジは、重み∞にする
そしてS-Tの最小カットを行い、S側の集合が求めるサブセットになる。

例えば次のような元グラフがあったとする。
f:id:yambe2002:20170331101642p:plain

これをS-Tグラフ化したものが以下になる。(点線が最小カット)
f:id:yambe2002:20170331101655p:plain
S側の集合{b, c, d, e, g}が最大のClosedセットであり、重みの合計値は6になっている。

有向グラフ {G} について、{ (s \cup S, t \cup T) }{ G_{s, t} }の最小カットであれば、{S}{G}の重みの合計が最大になっているclosedセットである。({ G_{s, t} }{ G }から作成した{s-t}グラフ)

証明
{ W^+ \equiv { j \in V | w_i > 0 }}{ W^- \equiv { j \in V | w_i \leq 0  }}とする。
f:id:yambe2002:20170331103755p:plain
(*1より抜粋)

{ W^{+}}は定数なので、 { C( s \cup S, t \cup T ) }を最小化するということは、{ \sum_{ i \in S } W_i }の最大化と同じことがわかる。

最小s-tカットの求め方は、以下のスライドが分かりやすかった。
https://www.slideshare.net/KuoE0/acmicpc-minimum-cut

また、Dinicという方法で最小カット(最大フロー)を求める方法はこちらが参考になる。
https://www.slideshare.net/KuoE0/acmicpc-dinics-algorithm


DinicのC#実装

public class Dinic
{
    public class Edge
    {
        public int V;
        public double C;
        public Edge N;
        public Edge B;
    }

    Edge[] _pool = new Edge[0];
    int _topIdx;

    Edge[] _adj = new Edge[0];

    public int MaxIdx;
    int _minV;

    public Dinic()
    {
        Expand(ref _pool, 128, true);
        Expand(ref _adj, 128, false);
        Init();
    }

    public void Init()
    {
        for (int i = 0; i < _adj.Length; i++) _adj[i] = null;
        _topIdx = 0;
        MaxIdx = 0;
        _minV = Int32.MaxValue;
    }

    void Expand(ref Edge[] edges, int size, bool fillWithEmptyEdge = false)
    {
        if (edges.Length <= size)
        {
            var newEdges = new Edge[Math.Max(size, edges.Length * 2)];
            for (int i = 0; i < edges.Length; i++) newEdges[i] = edges[i];
            if (fillWithEmptyEdge)
            {
                for (int i = edges.Length; i < newEdges.Length; i++)
                {
                    newEdges[i] = new Edge();
                }
            }
            edges = newEdges;
        }
    }

    public void AddEdge(int u, int v, double c, double bc = 0)
    {
        MaxIdx = Math.Max(MaxIdx, Math.Max(u + 1, v + 1));
        _minV = Math.Min(_minV, Math.Min(u, v));

        Expand(ref _pool, _topIdx + 1, true);
        Expand(ref _adj, MaxIdx);

        _pool[_topIdx].V = v;
        _pool[_topIdx].C = c;
        _pool[_topIdx].N = _adj[u];
        _adj[u] = _pool[_topIdx];
        _topIdx++;

        _pool[_topIdx].V = u;
        _pool[_topIdx].C = bc;
        _pool[_topIdx].N = _adj[v];
        _adj[v] = _pool[_topIdx];
        _topIdx++;

        _adj[u].B = _adj[v];
        _adj[v].B = _adj[u];

        if (u == v)
        {
            _adj[u].N.B = _adj[u];
            _adj[v].B = _adj[v].N;
        }
    }

    int[] _d;
    int[] _q;
    int _qh;
    int _qt;
    Edge[] _cur;

    public double GetMaxFlow(int s, int t)
    {
        MaxIdx = Math.Max(MaxIdx, Math.Max(s + 1, t + 1));
        _minV = Math.Min(_minV, Math.Min(s, t));

        _d = new int[MaxIdx];
        _q = new int[MaxIdx];
        _cur = new Edge[MaxIdx];

        double f = 0;

        while (ReLabel(s, t))
        {
            for (int i = 0; i < MaxIdx; i++) _cur[i] = _adj[i];
            f += Augument(s, t, s, double.MaxValue);
        }

        return f;
    }

    bool ReLabel(int s, int t)
    {
        for (int i = 0; i < MaxIdx; i++) _d[i] = Int32.MaxValue;
        _d[_q[_qh = _qt = 0] = t] = 0;
        while (_qh <= _qt)
        {
            var u = _q[_qh++];

            for (var i = _adj[u]; i != null; i = i.N)
            {
                if (i.B.C != 0 && _d[i.V] > _d[u] + 1)
                {
                    _d[i.V] = _d[u] + 1;
                    if ((_q[++_qt] = i.V) == s)
                    {
                        return true;
                    }
                }
            }
        }
        return false;
    }

    double Augument(int s, int t, int u, double e)
    {
        if (u == t) return e;
        double f = 0;

        for (var i = _cur[u]; i != null; i = i.N)
        {
            if (i.C > 0 && _d[u] == _d[i.V] + 1)
            {
                var df = Augument(s, t, i.V, Math.Min(e, i.C));
                if (df > 0)
                {
                    i.C -= df;
                    i.B.C += df;
                    e -= df;
                    f += df;
                }
            }
            if (!(e > 0)) break;
        }

        return f;
    }

    // (S, T)
    public Tuple<List<int>, List<int>> GetStCut(int s)
    {
        var used_flg = new bool[MaxIdx];

        Dfs(s, used_flg);

        var ret_s = new List<int>();
        var ret_t = new List<int>();
        for (int i = _minV; i < MaxIdx; i++)
        {
            if (!used_flg[i])
            {
                ret_t.Add(i);
            }
            else
            {
                ret_s.Add(i);
            }
        }

        return new Tuple<List<int>, List<int>>(ret_s, ret_t);
    }

    void Dfs(int u, bool[] used_flg)
    {
        if (used_flg[u]) return;
        used_flg[u] = true;

        for (var i = _adj[u]; i != null; i = i.N)
        {
            if (i.C > 1e-6) Dfs(i.V, used_flg);
        }
    }
}