AtCoder LibraryをJavaに移植

特集

AtCoder Library (ACL)をJavaに移植したものです。
ただし、色々とコメントを付け足していたり、場合によっては勝手に機能追加していたりもします。
できたものから追記していきます。
（10/3完了目標に下方修正しているけどそれもちょっと厳しい・・・。）

既に他にやっている人たちの成果物がありますが、自分がやり始めようとした時には既に終わりかけているようだったのとか、自分の都合で改変したいことも出てくるかもしれないとか、最近あまり時間が取れず周りのペースに合わせられる自信がないとか、GitHubの使い方よくわかってないとかで、入っていくのは諦めました。

C++は完全に読めるわけではないので、自分が移植する上では、上記の人たちの成果物はかなり参考にさせてもらっています。
参考どころかだんだんただのパクリになっており、自力で行っている作業はコメントの追加だけなものも多々あります。

以下に記載していく自分の成果物も、CC0ライセンスで公開します。（って言うのが正しいのかよくわかってないけど、ご自由に使っていただいて構いませんが、責任は持ちません。）


※アサーションを有効にするには、実行時に「-ea」オプションを付ける必要あり。

Fenwick Tree

  • 冒頭に記載の先駆者たちが作っていた、配列を引数に初期化するものも拝借。
  • 引数がrだけのものも使う想定。今までそれでやってたので。
/**
 * Fenwick Tree（Binary Indexed Tree）
 */
class FenwickTree {
    private int n;
    private long[] data;

    /**
     * 長さnの配列（a[0]～a[n-1]）を作る。初期値は全て0。<br>
     * O(n)
     * 
     * @param n 配列の長さ
     */
    public FenwickTree(int n) {
        this.n = n;
        data = new long[n];
    }

    /**
     * 初期データを元にFenwick Treeを構成する。<br>
     * O(n)
     * 
     * @param data 初期データ
     */
    public FenwickTree(long[] data) {
        this(data.length);
        build(data);
    }

    /**
     * a[p] += x を行う。<br>
     * O(log n)
     * 
     * @param p 加算位置（0≦p＜n）
     * @param x 加算値
     */
    void add(int p, long x) {
        assert 0 <= p && p < n : "p=" + p;

        p++;
        while (p <= n) {
            data[p - 1] += x;
            p += p & -p;
        }
    }

    /**
     * a[l] + ... + a[r-1] を返す。<br>
     * O(log n)
     * 
     * @param l 開始位置（含む）    （0≦l≦r≦n）
     * @param r 終了位置（含まない）（0≦l≦r≦n）
     */
    long sum(int l, int r) {
        assert 0 <= l && l <= r && r <= n : "l=" + l + ", r=" + r;

        return sum(r) - sum(l);
    }

    /**
     * a[0] + ... + a[r-1] を返す。<br>
     * O(log n)
     * 
     * @param r 終了位置（含まない）（0≦r≦n）
     */
    long sum(int r) {
        assert 0 <= r && r <= n : "r=" + r;

        long s = 0;
        while (r > 0) {
            s += data[r - 1];
            r -= r & -r;
        }
        return s;
    }

    private void build(long[] dat) {
        System.arraycopy(dat, 0, data, 0, n);
        for (int i = 1; i <= n; i++) {
            int p = i + (i & -i);
            if (p <= n) {
                data[p - 1] += data[i - 1];
            }
        }
    }
}
提出

ACL Practice Contest　B - Fenwick Tree


Segment Tree（2020/9/30更新）

  • 冒頭に記載の先駆者たちの成果物ほぼそのまま。
  • 指定する型がInteger等の不変オブジェクトであればあまり気にする必要はないが、自前のオブジェクトの場合は、二項演算の戻り値を新しいオブジェクトにしないと壊れることに注意する（引数のオブジェクトをそのまま返すとかはNG）。
/**
 * セグメント木
 * <pre>
 * 【使用例】
 * ■最大
 * SegTree<Integer> st = new SegTree<>(n, 0, (v1, v2) -> Math.max(v1, v2));
 * ■最小
 * SegTree<Integer> st = new SegTree<>(n, Integer.MAX_VALUE, (v1, v2) -> Math.min(v1, v2));
 * ■和
 * SegTree<Integer> st = new SegTree<>(n, 0, (v1, v2) -> v1 + v2);
 * </pre>
 */
class SegTree<S> {
    private final int n; // 要素数
    private final int size; // SegTreeのデータサイズ
    private final BinaryOperator<S> op; // 二項演算
    private final S e; // 単位元
    private final S[] data;

    /**
     * 長さnの配列aを作る。初期値は全て単位元。<br>
     * O(n)
     * 
     * @param n 要素数
     * @param e 単位元
     * @param op 二項演算
     */
    @SuppressWarnings("unchecked")
    public SegTree(int n, S e, BinaryOperator<S> op) {
        this.n = n;
        int k = 1;
        while (k < n) k <<= 1;
        this.size = k;
        this.e = e;
        this.op = op;
        this.data = (S[]) new Object[size << 1];
        Arrays.fill(data, e);
    }

    /**
     * 長さdat.lengthの配列aをdatを元に作る。<br>
     * O(n)
     * 
     * @param dat 初期データ
     * @param e 単位元
     * @param op 二項演算
     */
    public SegTree(S[] dat, S e, BinaryOperator<S> op) {
        this(dat.length, e, op);
        build(dat);
    }

    private void build(S[] dat) {
        int l = dat.length;
        System.arraycopy(dat, 0, data, size, l);
        for (int i = size - 1; i > 0; i--) {
            data[i] = op.apply(data[i << 1 | 0], data[i << 1 | 1]);
        }
    }

    /**
     * a[p] = xとする。<br>
     * O(log n)
     * 
     * @param p 設定位置（0≦p＜n）
     * @param x 設定値
     */
    void set(int p, S x) {
        assert 0 <= p && p < n : "p=" + p;

        data[p += size] = x;
        p >>= 1;
        while (p > 0) {
            data[p] = op.apply(data[p << 1 | 0], data[p << 1 | 1]);
            p >>= 1;
        }
    }

    /**
     * a[p]を取得する。<br>
     * O(1)
     * 
     * @param p 取得位置（0≦p＜n）
     */
    S get(int p) {
        assert 0 <= p && p < n : "p=" + p;

        return data[p + size];
    }

    /**
     * 区間l～(r-1)の計算を行う。<br>
     * l＝rのときは単位元を返す。<br>
     * O(log n)
     * 
     * @param l 開始位置（含む）    （0≦l≦r≦n）
     * @param r 終了位置（含まない）（0≦l≦r≦n）
     */
    S prod(int l, int r) {
        assert 0 <= l && l <= r && r <= n : "l=" + l + ", r=" + r;

        S sumLeft = e;
        S sumRight = e;
        l += size;
        r += size;
        while (l < r) {
            if ((l & 1) == 1) sumLeft = op.apply(sumLeft, data[l++]);
            if ((r & 1) == 1) sumRight = op.apply(data[--r], sumRight);
            l >>= 1;
            r >>= 1;
        }
        return op.apply(sumLeft, sumRight);
    }

    /**
     * 全区間の計算を行う。<br>
     * O(1)
     */
    S allProd() {
        return data[1];
    }

    /**
     * f(op(l～(r-1)))＝trueとなる最大のrを返す。<br>
     * （計算区間にrを追加するとfalseになる）<br>
     * O(log n)
     * 
     * @param l 左端（含む）（0≦l≦n）
     * @param f 判定関数（f(e)＝true）
     * @return 条件を満たす最大のr
     */
    int maxRight(int l, Predicate<S> f) {
        assert 0 <= l && l <= n : "l=" + l;
        assert f.test(e);

        if (l == n) return n;
        l += size;
        S sum = e;
        do {
            l >>= Integer.numberOfTrailingZeros(l);
            if (!f.test(op.apply(sum, data[l]))) {
                while (l < size) {
                    l = l << 1;
                    if (f.test(op.apply(sum, data[l]))) {
                        sum = op.apply(sum, data[l]);
                        l++;
                    }
                }
                return l - size;
            }
            sum = op.apply(sum, data[l]);
            l++;
        } while ((l & -l) != l);
        return n;
    }

    /**
     * f(op(l～(r-1)))＝trueとなる最小のlを返す。<br>
     * （計算区間に(l-1)を追加するとfalseになる）<br>
     * O(log n)
     * 
     * @param r 右端（含まない）（0≦r≦n）
     * @param f 判定関数（f(e)＝true）
     * @return 条件を満たす最小のl
     */
    int minLeft(int r, Predicate<S> f) {
        assert 0 <= r && r <= n : "r=" + r;
        assert f.test(e);

        if (r == 0) return 0;
        r += size;
        S sum = e;
        do {
            r--;
            while (r > 1 && (r & 1) == 1) r >>= 1;
            if (!f.test(op.apply(data[r], sum))) {
                while (r < size) {
                    r = r << 1 | 1;
                    if (f.test(op.apply(data[r], sum))) {
                        sum = op.apply(data[r], sum);
                        r--;
                    }
                }
                return r + 1 - size;
            }
            sum = op.apply(data[r], sum);
        } while ((r & -r) != r);
        return 0;
    }

    // **************** DEBUG **************** //

    private int indent = 6;

    void setIndent(int newIndent) {
        this.indent = newIndent;
    }

    @Override
    public String toString() {
        return toSimpleString();
    }

    /**
     * セグメント木全体の要素をツリー状に出力。
     */
    String toDetailedString() {
        return toDetailedString(1, 0);
    }

    private String toDetailedString(int k, int sp) {
        if (k >= size) return indent(sp) + data[k];
        StringBuilder sb = new StringBuilder();
        sb.append(toDetailedString(k << 1 | 1, sp + indent));
        sb.append("\n");
        sb.append(indent(sp) + data[k]);
        sb.append("\n");
        sb.append(toDetailedString(k << 1 | 0, sp + indent));
        return sb.toString();
    }

    private String indent(int n) {
        StringBuilder sb = new StringBuilder();
        while (n-- > 0) sb.append(' ');
        return sb.toString();
    }

    /**
     * n件の要素のみを、配列形式で出力。
     */
    String toSimpleString() {
        StringBuilder sb = new StringBuilder();
        sb.append('[');
        for (int i = 0; i < size; i++) {
            sb.append(data[i + size]);
            if (i < size - 1) sb.append(',').append(' ');
        }
        sb.append(']');
        return sb.toString();
    }
}
提出

ACL Practice Contest　J - Segment Tree


String（2020/9/30更新）

  • 冒頭に記載の先駆者たちの成果物ほぼそのまま。
  • 不要なものを消したりしたくらい。
/**
 * 文字列アルゴリズム詰め合わせ
 * 
 * <pre>
 * 【Suffix Array、LCP Arrayの例】
 * char[] s = "abcbcba".toCharArray();
 * int[] sa = StringAC.suffixArray(s);  // [6, 0, 5, 3, 1, 4, 2]
 * int[] la = StringAC.lcpArray(s, sa); // [1, 0, 1, 3, 0, 2]
 * 
 * 【saが表す文字列】
 * "a"
 * "abcbcba"
 * "ba"
 * "bcba"
 * "bcbcba"
 * "cba"
 * "cbcba"
 * 
 * 【Zアルゴリズムの例】
 * char[] s2 = "ababac".toCharArray();
 * int[] za = StringAC.zAlgorithm(s);  // [6, 0, 3, 0, 1, 0]
 * </pre>
 */
class StringAC {
    private static final int THRESHOLD_NAIVE = 50;

    /**
     * 長さnの文字配列sのSuffix Arrayとして、長さnの配列を返す。<br>
     * 戻り値の配列は、x＝0～(n-1)の順列で、「s[x]..s[n-1]」の昇順。<br>
     * O(n)
     */
    static int[] suffixArray(char[] s) {
        int n = s.length;
        int[] s2 = new int[n];
        for (int i = 0; i < n; i++) {
            s2[i] = s[i];
        }
        return sais(s2, 255);
    }

    /**
     * 長さnの数値配列sのSuffix Arrayとして、長さnの配列を返す。<br>
     * 戻り値の配列は、x＝0～(n-1)の順列で、「s[x]..s[n-1]」の昇順。<br>
     * 時間：O(n log n)、空間：O(n)
     */
    static int[] suffixArray(int[] s) {
        int n = s.length;
        int[] vals = Arrays.copyOf(s, n);
        Arrays.sort(vals);
        int p = 1;
        for (int i = 1; i < n; i++) {
            if (vals[i] != vals[i - 1]) {
                vals[p++] = vals[i];
            }
        }
        int[] s2 = new int[n];
        for (int i = 0; i < n; i++) {
            s2[i] = Arrays.binarySearch(vals, 0, p, s[i]);
        }
        return sais(s2, p);
    }

    /**
     * 長さnの数値配列sのSuffix Arrayとして、長さnの配列を返す。<br>
     * 戻り値の配列は、x＝0～(n-1)の順列で、「s[x]..s[n-1]」の昇順。<br>
     * 0≦sの全要素≦upper<br>
     * O(n + upper)
     */
    static int[] suffixArray(int[] s, int upper) {
        assert (0 <= upper);
        for (int d : s) {
            assert (0 <= d && d <= upper) : "d=" + d + ", upper=" + upper;
        }
        return sais(s, upper);
    }

    /**
     * 長さnの文字配列sのLCP Arrayとして、長さn-1の配列を返す。<br>
     * 戻り値の配列のi番目は、sa[i]と[i+1]の共通接頭辞の長さ。<br>
     * O(n)
     * 
     * @param s 元の文字配列
     * @param sa Suffix Array
     */
    static int[] lcpArray(char[] s, int[] sa) {
        int n = s.length;
        int[] s2 = new int[n];
        for (int i = 0; i < n; i++) {
            s2[i] = s[i];
        }
        return lcpArray(s2, sa);
    }

    /**
     * 長さnの数値配列sのLCP Arrayとして、長さn-1の配列を返す。<br>
     * 戻り値の配列のi番目は、sa[i]と[i+1]の共通接頭辞の長さ。<br>
     * O(n)
     * 
     * @param s 元の数値配列
     * @param sa Suffix Array
     */
    static int[] lcpArray(int[] s, int[] sa) {
        int n = s.length;
        assert (n >= 1);
        int[] rnk = new int[n];
        for (int i = 0; i < n; i++) {
            rnk[sa[i]] = i;
        }
        int[] lcp = new int[n - 1];
        int h = 0;
        for (int i = 0; i < n; i++) {
            if (h > 0) h--;
            if (rnk[i] == 0) {
                continue;
            }
            int j = sa[rnk[i] - 1];
            for (; j + h < n && i + h < n; h++) {
                if (s[j + h] != s[i + h]) break;
            }
            lcp[rnk[i] - 1] = h;
        }
        return lcp;
    }

    /**
     * 入力の長さをnとして、長さnの配列を返す。<br>
     * 戻り値の配列のi番目は、「s[0]..s[n-1]」と「s[i]..s[n-1]」の共通接頭辞の長さ。<br>
     * O(n)
     */
    static int[] zAlgorithm(char[] s) {
        int n = s.length;
        if (n == 0) return new int[0];
        int[] z = new int[n];
        for (int i = 1, j = 0; i < n; i++) {
            int k = j + z[j] <= i ? 0 : Math.min(j + z[j] - i, z[i - j]);
            while (i + k < n && s[k] == s[i + k]) k++;
            z[i] = k;
            if (j + z[j] < i + z[i]) j = i;
        }
        z[0] = n;
        return z;
    }

    /**
     * 入力の長さをnとして、長さnの配列を返す。<br>
     * 戻り値の配列のi番目は、「s[0]..s[n-1]」と「s[i]..s[n-1]」の共通接頭辞の長さ。<br>
     * O(n)
     */
    int[] zAlgorithm(int[] s) {
        int n = s.length;
        if (n == 0) return new int[0];
        int[] z = new int[n];
        for (int i = 1, j = 0; i < n; i++) {
            int k = j + z[j] <= i ? 0 : Math.min(j + z[j] - i, z[i - j]);
            while (i + k < n && s[k] == s[i + k]) k++;
            z[i] = k;
            if (j + z[j] < i + z[i]) j = i;
        }
        z[0] = n;
        return z;
    }

    private static int[] sais(int[] s, int upper) {
        int n = s.length;
        if (n == 0) return new int[0];
        if (n == 1) return new int[] { 0 };
        if (n == 2) {
            if (s[0] < s[1]) {
                return new int[] { 0, 1 };
            } else {
                return new int[] { 1, 0 };
            }
        }
        if (n < THRESHOLD_NAIVE) {
            return saNaive(s);
        }

        int[] sa = new int[n];
        boolean[] ls = new boolean[n];
        for (int i = n - 2; i >= 0; i--) {
            ls[i] = s[i] == s[i + 1] ? ls[i + 1] : s[i] < s[i + 1];
        }

        int[] sumL = new int[upper + 1];
        int[] sumS = new int[upper + 1];

        for (int i = 0; i < n; i++) {
            if (ls[i]) {
                sumL[s[i] + 1]++;
            } else {
                sumS[s[i]]++;
            }
        }

        for (int i = 0; i <= upper; i++) {
            sumS[i] += sumL[i];
            if (i < upper) sumL[i + 1] += sumS[i];
        }

        Consumer<int[]> induce = lms -> {
            Arrays.fill(sa, -1);
            int[] buf = new int[upper + 1];
            System.arraycopy(sumS, 0, buf, 0, upper + 1);
            for (int d : lms) {
                if (d == n) continue;
                sa[buf[s[d]]++] = d;
            }
            System.arraycopy(sumL, 0, buf, 0, upper + 1);
            sa[buf[s[n - 1]]++] = n - 1;
            for (int i = 0; i < n; i++) {
                int v = sa[i];
                if (v >= 1 && !ls[v - 1]) {
                    sa[buf[s[v - 1]]++] = v - 1;
                }
            }
            System.arraycopy(sumL, 0, buf, 0, upper + 1);
            for (int i = n - 1; i >= 0; i--) {
                int v = sa[i];
                if (v >= 1 && ls[v - 1]) {
                    sa[--buf[s[v - 1] + 1]] = v - 1;
                }
            }
        };

        int[] lmsMap = new int[n + 1];
        Arrays.fill(lmsMap, -1);
        int m = 0;
        for (int i = 1; i < n; i++) {
            if (!ls[i - 1] && ls[i]) {
                lmsMap[i] = m++;
            }
        }

        int[] lms = new int[m];
        {
            int p = 0;
            for (int i = 1; i < n; i++) {
                if (!ls[i - 1] && ls[i]) {
                    lms[p++] = i;
                }
            }
        }

        induce.accept(lms);

        if (m > 0) {
            int[] sortedLms = new int[m];
            {
                int p = 0;
                for (int v : sa) {
                    if (lmsMap[v] != -1) {
                        sortedLms[p++] = v;
                    }
                }
            }
            int[] recS = new int[m];
            int recUpper = 0;
            recS[lmsMap[sortedLms[0]]] = 0;
            for (int i = 1; i < m; i++) {
                int l = sortedLms[i - 1], r = sortedLms[i];
                int endL = (lmsMap[l] + 1 < m) ? lms[lmsMap[l] + 1] : n;
                int endR = (lmsMap[r] + 1 < m) ? lms[lmsMap[r] + 1] : n;
                boolean same = true;
                if (endL - l != endR - r) {
                    same = false;
                } else {
                    while (l < endL && s[l] == s[r]) {
                        l++;
                        r++;
                    }
                    if (l == n || s[l] != s[r]) same = false;
                }
                if (!same) {
                    recUpper++;
                }
                recS[lmsMap[sortedLms[i]]] = recUpper;
            }

            int[] recSA = sais(recS, recUpper);

            for (int i = 0; i < m; i++) {
                sortedLms[i] = lms[recSA[i]];
            }
            induce.accept(sortedLms);
        }
        return sa;
    }

    private static int[] saNaive(int[] s) {
        int n = s.length;
        int[] sa = new int[n];
        for (int i = 0; i < n; i++) {
            sa[i] = i;
        }
        insertionsort(sa, (l, r) -> {
            while (l < n && r < n) {
                if (s[l] != s[r]) return s[l] - s[r];
                l++;
                r++;
            }
            return -(l - r);
        });
        return sa;
    }

    private static void insertionsort(int[] a, IntBinaryOperator comparator) {
        final int n = a.length;
        for (int i = 1; i < n; i++) {
            final int tmp = a[i];
            if (comparator.applyAsInt(a[i - 1], tmp) > 0) {
                int j = i;
                do {
                    a[j] = a[j - 1];
                    j--;
                } while (j > 0 && comparator.applyAsInt(a[j - 1], tmp) > 0);
                a[j] = tmp;
            }
        }
    }
}
提出

ACL Practice Contest　I - Number of Substrings


Math（2020/9/22更新）

  • 使うときは、MathACクラスごと貼り付けてもよいが、ほぼメソッド単体で貼り付けても使えるようにしたい考え。（invModやcrt使用時は、invGcdも一緒に貼り付ける必要があるが）
  • そのため、safe_mod程度のものはメソッド化せず、ベタ実装してしまっている。
  • （2020/9/22）CRTも作成。
/**
 * 数論的アルゴリズム詰め合わせ
 */
class MathAC {

    /**
     * xのn乗 mod m<br>
     * 制約：0≦n, 1≦m<br>
     * O(log n)
     */
    static long powMod(long x, long n, int m) {
        assert 0 <= n : "n=" + n;
        assert 1 <= m : "m=" + m;

        if (m == 1) {
            return 0;
        }
        long r = 1;
        long y = x % m;
        if (y < 0) {
            y += m;
        }
        while (n > 0) {
            if ((n & 1) == 1) {
                r = r * y % m;
            }
            y = y * y % m;
            n >>= 1;
        }
        return r;
    }

    /**
     * xのmod mでの逆元<br>
     * 制約：1≦m, gcd(x, m)＝1<br>
     * O(log m)
     */
    static long invMod(long x, int m) {
        assert 1 <= m : "m=" + m;
        long[] z = invGcd(x, m);
        assert z[0] == 1 : "gcd(x, m) != 1";
        return z[1];
    }

    private static long[] invGcd(long a, long b) {
        a = a % b;
        if (a < 0) {
            a += b;
        }
        if (a == 0) {
            return new long[] { b, 0 };
        }

        long s = b, t = a;
        long m0 = 0, m1 = 1;
        while (t > 0) {
            long u = s / t;
            s -= t * u;
            m0 -= m1 * u;
            long tmp = s;
            s = t;
            t = tmp;
            tmp = m0;
            m0 = m1;
            m1 = tmp;
        }
        if (m0 < 0) {
            m0 += b / s;
        }
        return new long[] { s, m0 };
    }

    /**
     * <pre>
     * 中国剰余定理
     * x≡r[i] (mod m[i])を解く。
     * O(n log lcm(m[i]))
     * 
     * ■制約
     * |r|＝|m|
     * 1≦m[i]
     * lcm(m[i])がlongに収まる。
     * 
     * ■戻り値
     * x≡y (mod z) となる{y, z} (0≦y＜z＝lcm(m[i]))
     * 答えがない場合は{0, 0}
     * n＝0の時は{0, 1}
     * </pre>
     */
    static long[] crt(long[] r, long[] m) {
        assert r.length == m.length : "|r|=" + r.length + ", |m|=" + m.length;

        int n = r.length;
        long r0 = 0, m0 = 1;
        for (int i = 0; i < n; i++) {
            assert 1 <= m[i] : "m[" + i + "]=" + m;

            long r1 = r[i] % m[i];
            if (r1 < 0) {
                r1 += m[i];
            }
            long m1 = m[i];
            if (m0 < m1) {
                long tmp = r0;
                r0 = r1;
                r1 = tmp;
                tmp = m0;
                m0 = m1;
                m1 = tmp;
            }
            if (m0 % m1 == 0) {
                if (r0 % m1 != r1) {
                    return new long[] { 0, 0 };
                }
                continue;
            }

            long[] ig = invGcd(m0, m1);
            long g = ig[0], im = ig[1];
            long u1 = m1 / g;
            if ((r1 - r0) % g != 0) {
                return new long[] { 0, 0 };
            }

            long x = (r1 - r0) / g % u1 * im % u1;
            r0 += x * m0;
            m0 *= u1;
            if (r0 < 0) {
                r0 += m0;
            }
        }
        return new long[] { r0, m0 };
    }

    /**
     * floor((a * i + b) / m) の i＝0～n-1 の総和<br>
     * y = (a/m)x + b/m、x軸、y軸、x = nで囲まれた領域の格子点の数（x軸上とx = n上は含まない）<br>
     * 制約：0≦n≦10^9, 1≦m≦10^9, 0≦a,b＜m<br>
     * O(log(n + m + a + b))
     */
    static long floorSum(long n, long m, long a, long b) {
        long ans = 0;
        if (a >= m) {
            ans += (n - 1) * n * (a / m) / 2;
            a %= m;
        }
        if (b >= m) {
            ans += n * (n / m);
            b %= m;
        }

        long ymax = (a * n + b) / m;
        long xmax = (ymax * m - b);
        if (ymax == 0) {
            return ans;
        }
        ans += (n - (xmax + a - 1) / a) * ymax;
        ans += floorSum(ymax, a, m, (a - xmax % a) % a);
        return ans;
    }
}
提出（CRT）

ACL Contest 1　B - Sum is Multiple

提出（Floor Sum）

ACL Practice Contest　C - Floor Sum


Convolution（2021/3/23更新）

  • 冒頭に記載の先駆者たちのものをほぼそのまま拝借。
  • 任意modで使えるconvolutionLLの方は、練習用問題でTLEが出たまま放置中。
  • （2021/3/23）ABCで使ったらTLEが出てしまったため、定数倍改善を試みたところ、引数や戻り値を全体的にlong→intに変更したら2～3割ほど早くなった。ソースは提出のリンク先参照。
/**
 * 畳み込み
 */
class Convolution {

    /**
     * 畳み込み<br>
     * 制約：n = |a|+|b|として、n - 1 ≦ 2^c かつ (mod - 1)が2^cで割り切れる<br>
     * O(n log n + log mod)
     *
     * @param a 数列
     * @param b 数列
     * @param mod NTT用素数（998244353, 1053818881, 1004535809, ...）
     * @return ret[i] = sum(j=0～i) {a[j] * b[i-j]}
     */
    static long[] convolution(long[] a, long[] b, int mod) {
        int n = a.length;
        int m = b.length;
        if (n == 0 || m == 0) {
            return new long[0];
        }

        int z = 1 << ceilPow2(n + m - 1);
        {
            long[] na = new long[z];
            long[] nb = new long[z];
            System.arraycopy(a, 0, na, 0, n);
            System.arraycopy(b, 0, nb, 0, m);
            a = na;
            b = nb;
        }

        int g = primitiveRoot(mod);
        long[] sume = sumE(mod, g);
        long[] sumie = sumIE(mod, g);

        butterfly(a, sume, mod);
        butterfly(b, sume, mod);
        for (int i = 0; i < z; i++) {
            a[i] = a[i] * b[i] % mod;
        }
        butterflyInv(a, sumie, mod);
        a = Arrays.copyOf(a, n + m - 1);

        long iz = powMod(z, mod - 2, mod);
        for (int i = 0; i < n + m - 1; i++)
            a[i] = a[i] * iz % mod;
        return a;
    }

    /**
     * 畳み込み<br>
     * 制約：n = |a|+|b|として、n - 1 ≦ 2^24<br>
     * O(n log n)
     *
     * @param a 数列
     * @param b 数列
     * @param mod 任意のmod
     * @return ret[i] = sum(j=0～i) {a[j] * b[i-j]}
     */
    static long[] convolutionLL(long[] a, long[] b, int mod) {
        int n = a.length;
        int m = b.length;
        if (n == 0 || m == 0)
            return new long[0];

        int mod1 = 754974721;
        int mod2 = 167772161;
        int mod3 = 469762049;

        long[] c1 = convolution(a, b, mod1);
        long[] c2 = convolution(a, b, mod2);
        long[] c3 = convolution(a, b, mod3);

        int retSize = c1.length;
        long[] ret = new long[retSize];
        int[] mods = { mod1, mod2, mod3, mod };
        for (int i = 0; i < retSize; ++i) {
            ret[i] = garner(new long[] { c1[i], c2[i], c3[i] }, mods);
        }
        return ret;
    }

    private static int primitiveRoot(int m) {
        if (m == 2) return 1;
        if (m == 167772161) return 3;
        if (m == 469762049) return 3;
        if (m == 754974721) return 11;
        if (m == 998244353) return 3;

        int[] divs = new int[20];
        divs[0] = 2;
        int cnt = 1;
        int x = (m - 1) / 2;
        while (x % 2 == 0) {
            x /= 2;
        }
        for (int i = 3; (long) (i) * i <= x; i += 2) {
            if (x % i == 0) {
                divs[cnt++] = i;
                while (x % i == 0) {
                    x /= i;
                }
            }
        }
        if (x > 1) {
            divs[cnt++] = x;
        }
        for (int g = 2;; g++) {
            boolean ok = true;
            for (int i = 0; i < cnt; i++) {
                if (powMod(g, (m - 1) / divs[i], m) == 1) {
                    ok = false;
                    break;
                }
            }
            if (ok) return g;
        }
    }

    private static long powMod(long x, long n, int m) {
        if (m == 1) return 0;
        long r = 1;
        long y = x % m;
        while (n > 0) {
            if ((n & 1) != 0) r = (r * y) % m;
            y = (y * y) % m;
            n >>= 1;
        }
        return r;
    }

    private static int ceilPow2(int n) {
        int x = 0;
        while ((1L << x) < n) {
            x++;
        }
        return x;
    }

    private static long[] sumE(int mod, int g) {
        long[] sum_e = new long[30];
        long[] es = new long[30];
        long[] ies = new long[30];
        int cnt2 = Integer.numberOfTrailingZeros(mod - 1);
        long e = powMod(g, (mod - 1) >> cnt2, mod);
        long ie = powMod(e, mod - 2, mod);
        for (int i = cnt2; i >= 2; i--) {
            es[i - 2] = e;
            ies[i - 2] = ie;
            e = e * e % mod;
            ie = ie * ie % mod;
        }
        long now = 1;
        for (int i = 0; i < cnt2 - 2; i++) {
            sum_e[i] = es[i] * now % mod;
            now = now * ies[i] % mod;
        }
        return sum_e;
    }

    private static long[] sumIE(int mod, int g) {
        long[] sum_ie = new long[30];
        long[] es = new long[30];
        long[] ies = new long[30];

        int cnt2 = Integer.numberOfTrailingZeros(mod - 1);
        long e = powMod(g, (mod - 1) >> cnt2, mod);
        long ie = powMod(e, mod - 2, mod);
        for (int i = cnt2; i >= 2; i--) {
            es[i - 2] = e;
            ies[i - 2] = ie;
            e = e * e % mod;
            ie = ie * ie % mod;
        }
        long now = 1;
        for (int i = 0; i < cnt2 - 2; i++) {
            sum_ie[i] = ies[i] * now % mod;
            now = now * es[i] % mod;
        }
        return sum_ie;
    }

    private static void butterfly(long[] a, long[] sumE, int mod) {
        int n = a.length;
        int h = ceilPow2(n);

        for (int ph = 1; ph <= h; ph++) {
            int w = 1 << (ph - 1), p = 1 << (h - ph);
            long now = 1;
            for (int s = 0; s < w; s++) {
                int offset = s << (h - ph + 1);
                for (int i = 0; i < p; i++) {
                    long l = a[i + offset];
                    long r = a[i + offset + p] * now % mod;
                    a[i + offset] = (l + r) % mod;
                    a[i + offset + p] = (l - r + mod) % mod;
                }
                int x = Integer.numberOfTrailingZeros(~s);
                now = now * sumE[x] % mod;
            }
        }
    }

    private static void butterflyInv(long[] a, long[] sumIE, int mod) {
        int n = a.length;
        int h = ceilPow2(n);

        for (int ph = h; ph >= 1; ph--) {
            int w = 1 << (ph - 1), p = 1 << (h - ph);
            long inow = 1;
            for (int s = 0; s < w; s++) {
                int offset = s << (h - ph + 1);
                for (int i = 0; i < p; i++) {
                    long l = a[i + offset];
                    long r = a[i + offset + p];
                    a[i + offset] = (l + r) % mod;
                    a[i + offset + p] = (mod + l - r) * inow % mod;
                }
                int x = Integer.numberOfTrailingZeros(~s);
                inow = inow * sumIE[x] % mod;
            }
        }
    }

    private static long garner(long[] c, int[] mods) {
        int n = c.length + 1;
        long[] cnst = new long[n];
        long[] coef = new long[n];
        Arrays.fill(coef, 1);
        for (int i = 0; i < n - 1; i++) {
            int m1 = mods[i];
            long v = (c[i] - cnst[i] + m1) % m1;
            v = v * powMod(coef[i], m1 - 2, m1) % m1;

            for (int j = i + 1; j < n; j++) {
                long m2 = mods[j];
                cnst[j] = (cnst[j] + coef[j] * v) % m2;
                coef[j] = (coef[j] * m1) % m2;
            }
        }
        return cnst[n - 1];
    }
}
提出（long版）

ACL Practice Contest　F - Convolution

提出（int版）

ACL Practice Contest　F - Convolution
ABC196　F - Substring 2　※convolution2（任意mod）だとTLE


Disjoint Set Union

  • 自分が元々持っていたUnion Findでは連結成分の数を取得できるようにしていたので、その機能（numメソッド）を追加。
  • groupsメソッドは本当にO(n)？　removeIfがどのように動くのかよくわからないが、毎回後ろの要素を詰める動きをしてしまうのであれば、O(n^2)になりそうな気がする。　でもとりあえずA問題の提出に付けてみても30msほどしか増えなかったので、大丈夫なのかな・・・。　まああまり使わなそうなので、実際使ってTLE喰らってから直すでもいいかな。
/**
 * Disjoint Set Union（Union Find）
 */
class DSU {
    private int n;
    private int[] parentOrSize;
    private int num;

    /**
     * n頂点0辺の無向グラフを作る。<br>
     * O(n)
     * 
     * @param n 頂点数
     */
    public DSU(int n) {
        this.n = n;
        this.parentOrSize = new int[n];
        Arrays.fill(parentOrSize, -1);
        num = n;
    }

    /**
     * 辺(a, b)を足す。<br>
     * ならしO(α(n))
     * 
     * @param a 頂点番号（0≦a＜n）
     * @param b 頂点番号（0≦b＜n）
     * @return 代表元
     */
    int merge(int a, int b) {
        assert 0 <= a && a < n;
        assert 0 <= b && b < n;

        int x = leader(a);
        int y = leader(b);
        if (x == y) {
            return x;
        }
        if (-parentOrSize[x] < -parentOrSize[y]) {
            int tmp = x;
            x = y;
            y = tmp;
        }
        parentOrSize[x] += parentOrSize[y];
        parentOrSize[y] = x;
        num--;
        return x;
    }

    /**
     * 頂点a, bが連結かどうか。<br>
     * ならしO(α(n))
     * 
     * @param a 頂点番号（0≦a＜n）
     * @param b 頂点番号（0≦b＜n）
     * @return true：連結　false：非連結
     */
    boolean same(int a, int b) {
        assert 0 <= a && a < n;
        assert 0 <= b && b < n;

        return leader(a) == leader(b);
    }

    /**
     * 頂点aの属する連結成分の代表元を返す。<br>
     * ならしO(α(n))
     * 
     * @param a 頂点番号（0≦a＜n）
     * @return 代表元
     */
    int leader(int a) {
        assert 0 <= a && a < n;

        if (parentOrSize[a] < 0) {
            return a;
        } else {
            return parentOrSize[a] = leader(parentOrSize[a]);
        }
    }

    /**
     * 頂点aの属する連結成分の要素数を返す。<br>
     * ならしO(α(n))
     * 
     * @param a 頂点番号（0≦a＜n）
     * @return 要素数
     */
    int size(int a) {
        assert 0 <= a && a < n;

        return -parentOrSize[leader(a)];
    }

    /**
     * 連結成分の数を返す。<br>
     * O(1)
     * 
     * @return 連結成分数
     */
    int num() {
        return num;
    }

    /**
     * グラフを連結成分に分けた情報を返す。<br>
     * O(n)
     * 
     * @return 「1つの連結成分の頂点番号のリスト」のリスト
     */
    List<List<Integer>> groups() {
        int[] leaderBuf = new int[n];
        int[] groupSize = new int[n];
        for (int i = 0; i < n; i++) {
            leaderBuf[i] = leader(i);
            groupSize[leaderBuf[i]]++;
        }
        List<List<Integer>> result = new ArrayList<>(n);
        for (int i = 0; i < n; i++) {
            result.add(new ArrayList<>(groupSize[i]));
        }
        for (int i = 0; i < n; i++) {
            result.get(leaderBuf[i]).add(i);
        }
        result.removeIf(List::isEmpty);
        return result;
    }
}
提出

ACL Practice Contest　A - Disjoint Set Union
ACL Contest 1　A - Reachable Towns


Max Flow（2020/9/18追記）

  • ACLをほぼそのまま移植したつもり。
  • ACLでは容量をintかlong longか指定できるようになっているようだが、めんどくさいのでlong決め打ち。
/**
 * 最大流
 */
class MaxFlow {
    private final int n;
    private List<int[]> pos;
    private List<List<Edge2>> g;

    class Edge {
        /** 有向辺の始点 */
        final int from;
        /** 有向辺の終点 */
        final int to;
        /** 最大容量 */
        long cap;
        /** 流量 */
        long flow;

        public Edge(int from, int to, long cap, long flow) {
            this.from = from;
            this.to = to;
            this.cap = cap;
            this.flow = flow;
        }
    }

    private class Edge2 {
        final int to, rev;
        long cap;

        public Edge2(int to, int rev, long cap) {
            this.to = to;
            this.rev = rev;
            this.cap = cap;
        }
    }

    /**
     * n頂点0辺のグラフを作る。<br>
     * O(n)
     * 
     * @param n 頂点数
     */
    public MaxFlow(int n) {
        this.n = n;
        pos = new ArrayList<>();
        g = new ArrayList<>(n);
        for (int i = 0; i < n; i++) {
            g.add(new ArrayList<>());
        }
    }

    /**
     * fromからtoへ最大容量cap、流量0の辺を追加する。<br>
     * ならしO(1)
     * 
     * @param from 有向辺の始点（0≦from＜n）
     * @param to 有向辺の終点（0≦to＜n）
     * @param cap 最大容量（0≦cap）
     * @return 何番目に追加された辺か
     */
    int addEdge(int from, int to, long cap) {
        assert 0 <= from && from < n : "from=" + from;
        assert 0 <= to && to < n : "to=" + to;
        assert 0 <= cap : "cap=" + cap;

        int m = pos.size();
        pos.add(new int[] { from, g.get(from).size() });
        g.get(from).add(new Edge2(to, g.get(to).size(), cap));
        g.get(to).add(new Edge2(from, g.get(from).size() - 1, 0));
        return m;
    }

    /**
     * i番目に追加された辺を取得する。<br>
     * O(1)
     * 
     * @param i 辺のインデックス（0≦i＜辺数）
     * @return 辺情報
     */
    Edge getEdge(int i) {
        assert 0 <= i && i < pos.size() : "i=" + i + ", size=" + pos.size();

        Edge2 e = g.get(pos.get(i)[0]).get(pos.get(i)[1]);
        Edge2 re = g.get(e.to).get(e.rev);
        return new Edge(pos.get(i)[0], e.to, e.cap + re.cap, re.cap);
    }

    /**
     * 全ての辺を取得する。<br>
     * O(辺数)
     * 
     * @return 辺情報リスト
     */
    List<Edge> edges() {
        int m = pos.size();
        List<Edge> result = new ArrayList<>(m);
        for (int i = 0; i < m; i++) {
            result.add(getEdge(i));
        }
        return result;
    }

    /**
     * i番目に追加された辺の容量、流量を変更する。他の辺は変更しない。<br>
     * O(1)
     * 
     * @param i 辺のインデックス（0≦i＜辺数）
     * @param newCap 変更後の容量（0≦newFlow≦newCap）
     * @param newFlow 変更後の流量（0≦newFlow≦newCap）
     */
    void changeEdge(int i, long newCap, long newFlow) {
        assert 0 <= i && i < pos.size() : "i=" + i + ", size=" + pos.size();
        assert 0 <= newFlow && newFlow <= newCap : "newCap=" + newCap + ", newFlow=" + newFlow;

        Edge2 e = g.get(pos.get(i)[0]).get(pos.get(i)[1]);
        Edge2 re = g.get(e.to).get(e.rev);
        e.cap = newCap - newFlow;
        re.cap = newFlow;
    }

    /**
     * 最大流。頂点sからtへ流せるだけ流す。<br>
     * O(n^2・m)
     * 
     * @param s 始点（0≦s＜n）
     * @param t 終点（0≦t＜n）
     * @return 流せた量
     */
    long flow(int s, int t) {
        return flow(s, t, Long.MAX_VALUE);
    }

    /**
     * 最大流。頂点sからtへflowLimitに達するまで流せるだけ流す。<br>
     * O(n^2・m)
     * 
     * @param s 始点（0≦s＜n）
     * @param t 終点（0≦t＜n）
     * @param flowLimit 流量制限
     * @return 流せた量
     */
    long flow(int s, int t, long flowLimit) {
        assert 0 <= s && s < n : "s=" + s;
        assert 0 <= t && t < n : "t=" + t;

        int[] level = new int[n];
        int[] iter = new int[n];
        Queue<Integer> que = new ArrayDeque<>();

        long flow = 0;
        while (flow < flowLimit) {
            bfs(s, t, level, que);
            if (level[t] == -1) {
                break;
            }
            Arrays.fill(iter, 0);
            while (flow < flowLimit) {
                long f = dfs(t, flowLimit - flow, s, iter, level);
                if (f <= 0) {
                    break;
                }
                flow += f;
            }
        }
        return flow;
    }

    private void bfs(int s, int t, int[] level, Queue<Integer> que) {
        Arrays.fill(level, -1);
        level[s] = 0;
        que.clear();
        que.add(s);
        while (!que.isEmpty()) {
            int v = que.poll();
            for (Edge2 e : g.get(v)) {
                if (e.cap == 0 || level[e.to] >= 0) {
                    continue;
                }
                level[e.to] = level[v] + 1;
                if (e.to == t) {
                    return;
                }
                que.add(e.to);
            }
        }
    }

    private long dfs(int v, long up, int s, int[] iter, int[] level) {
        if (v == s) {
            return up;
        }
        long res = 0;
        int levelv = level[v];
        while (iter[v] < g.get(v).size()) {
            Edge2 e = g.get(v).get(iter[v]++);
            Edge2 re = g.get(e.to).get(e.rev);
            if (levelv <= level[e.to] || re.cap == 0) {
                continue;
            }
            long d = dfs(e.to, Math.min(up - res, re.cap), s, iter, level);
            if (d <= 0) {
                continue;
            }
            e.cap += d;
            re.cap -= d;
            res += d;
            if (res == up) {
                break;
            }
        }
        return res;
    }

    /**
     * 最小カット<br>
     * limitなしのflow(s, t)をちょうど1回読んだ後に呼ぶ。<br>
     * O(n + m)
     * 
     * @param s 始点（0≦s＜n）
     * @return true：始点から到達可能な頂点
     */
    boolean[] minCut(int s) {
        assert 0 <= s && s < n : "s=" + s;

        boolean[] visited = new boolean[n];
        Queue<Integer> que = new ArrayDeque<>();
        que.add(s);
        visited[s] = true;
        while (!que.isEmpty()) {
            int p = que.poll();
            for (Edge2 e : g.get(p)) {
                if (e.cap > 0 && !visited[e.to]) {
                    visited[e.to] = true;
                    que.add(e.to);
                }
            }
        }
        return visited;
    }
}
提出

ACL Practice Contest　D - Maxflow


Min Cost Flow（2020/9/18追記）

  • ACLをほぼそのまま移植したつもり。
  • ACLでは容量やコストをintかlong longか指定できるようになっているようだが、めんどくさいのでlong決め打ち。
  • 戻り値はpairの代わりに長さ2の配列。
/**
 * 最小費用流
 */
class MinCostFlow {
    private final int n;
    private List<int[]> pos;
    private List<List<Edge2>> g;

    class Edge {
        /** 有向辺の始点 */
        final int from;
        /** 有向辺の終点 */
        final int to;
        /** 最大容量 */
        long cap;
        /** 流量 */
        long flow;
        /** コスト */
        long cost;

        public Edge(int from, int to, long cap, long flow, long cost) {
            this.from = from;
            this.to = to;
            this.cap = cap;
            this.flow = flow;
            this.cost = cost;
        }
    }

    private class Edge2 {
        final int to, rev;
        long cap, cost;

        public Edge2(int to, int rev, long cap, long cost) {
            this.to = to;
            this.rev = rev;
            this.cap = cap;
            this.cost = cost;
        }
    }

    private class Q {
        final long key;
        final int to;

        public Q(long key, int to) {
            this.key = key;
            this.to = to;
        }
    }

    /**
     * n頂点0辺のグラフを作る。<br>
     * O(n)
     * 
     * @param n 頂点数
     */
    public MinCostFlow(int n) {
        this.n = n;
        pos = new ArrayList<>();
        g = new ArrayList<>(n);
        for (int i = 0; i < n; i++) {
            g.add(new ArrayList<>());
        }
    }

    /**
     * fromからtoへ最大容量cap、コストcostの辺を追加する。<br>
     * ならしO(1)
     * 
     * @param from 有向辺の始点（0≦from＜n）
     * @param to 有向辺の終点（0≦to＜n）
     * @param cap 最大容量（0≦cap）
     * @param cost コスト（0≦cost）
     * @return 何番目に追加された辺か
     */
    int addEdge(int from, int to, long cap, long cost) {
        assert 0 <= from && from < n : "from=" + from;
        assert 0 <= to && to < n : "to=" + to;
        assert 0 <= cap : "cap=" + cap;
        assert 0 <= cost : "cost=" + cost;

        int m = pos.size();
        pos.add(new int[] { from, g.get(from).size() });
        g.get(from).add(new Edge2(to, g.get(to).size(), cap, cost));
        g.get(to).add(new Edge2(from, g.get(from).size() - 1, 0, -cost));
        return m;
    }

    /**
     * i番目に追加された辺を取得する。<br>
     * O(1)
     * 
     * @param i 辺のインデックス（0≦i＜辺数）
     * @return 辺情報
     */
    Edge getEdge(int i) {
        assert 0 <= i && i < pos.size() : "i=" + i + ", size=" + pos.size();

        Edge2 e = g.get(pos.get(i)[0]).get(pos.get(i)[1]);
        Edge2 re = g.get(e.to).get(e.rev);
        return new Edge(pos.get(i)[0], e.to, e.cap + re.cap, re.cap, e.cost);
    }

    /**
     * 全ての辺を取得する。<br>
     * O(辺数)
     * 
     * @return 辺情報リスト
     */
    List<Edge> edges() {
        int m = pos.size();
        List<Edge> result = new ArrayList<>(m);
        for (int i = 0; i < m; i++) {
            result.add(getEdge(i));
        }
        return result;
    }

    /**
     * 最小費用流。頂点sからtへ流せるだけ流す。<br>
     * Fを流量、mを辺数として、O(F(n + m) log n)
     * 
     * @param s 始点（0≦s＜n）
     * @param t 終点（0≦t＜n）
     * @return [0]：流せた量、[1]：その時のコスト
     */
    long[] flow(int s, int t) {
        return flow(s, t, Long.MAX_VALUE);
    }

    /**
     * 最小費用流。頂点sからtへflowLimitに達するまで流せるだけ流す。<br>
     * Fを流量、mを辺数として、O(F(n + m) log n)
     * 
     * @param s 始点（0≦s＜n）
     * @param t 終点（0≦t＜n）
     * @param flowLimit 流量制限
     * @return [0]：流せた量、[1]：その時のコスト
     */
    long[] flow(int s, int t, long flowLimit) {
        List<long[]> result = slope(s, t, flowLimit);
        return result.get(result.size() - 1);
    }

    /**
     * <pre>
     * 流量とコストの関係の折れ線を取得する。
     * 戻り値の最初の要素は(0, 0)
     * 戻り値の[0]、[1]は共に狭義単調増加
     * 辺のコストの最大をxとして、最後の要素の[0]はx
     * 
     * ■制約
     * sからtへ流したフローの流量がlongに収まる。
     * 流したフローのコストの総和がlongに収まる。
     * 0≦nx≦8×10^18＋1000
     * </pre>
     * 
     * @param s 始点（0≦s＜n、s≠t）
     * @param t 終点（0≦t＜n、s≠t）
     * @return ＜[0]：流量、[1]：その時のコスト＞のリスト
     */
    List<long[]> slope(int s, int t) {
        return slope(s, t, Long.MAX_VALUE);
    }

    /**
     * <pre>
     * 流量とコストの関係の折れ線を取得する。
     * 戻り値の最初の要素は(0, 0)
     * 戻り値の[0]、[1]は共に狭義単調増加
     * 辺のコストの最大をxとして、最後の要素の[0]はmin(x, flowLimit)
     * 
     * ■制約
     * sからtへ流したフローの流量がlongに収まる。
     * 流したフローのコストの総和がlongに収まる。
     * 0≦nx≦8×10^18＋1000
     * </pre>
     * 
     * @param s 始点（0≦s＜n、s≠t）
     * @param t 終点（0≦t＜n、s≠t）
     * @param flowLimit 流量制限
     * @return ＜[0]：流量、[1]：その時のコスト＞のリスト
     */
    List<long[]> slope(int s, int t, long flowLimit) {
        assert 0 <= s && s < n : "s=" + s;
        assert 0 <= t && t < n : "t=" + t;
        assert s != t : "s=t=" + s;

        long[] dual = new long[n];
        long[] dist = new long[n];
        int[] pv = new int[n];
        int[] pe = new int[n];
        boolean[] vis = new boolean[n];

        long flow = 0;
        long cost = 0;
        long prevCost = -1;
        List<long[]> result = new ArrayList<>();
        result.add(new long[] { flow, cost });
        while (flow < flowLimit) {
            if (!dualRef(s, t, dual, dist, pv, pe, vis)) {
                break;
            }
            long c = flowLimit - flow;
            for (int v = t; v != s; v = pv[v]) {
                c = Math.min(c, g.get(pv[v]).get(pe[v]).cap);
            }
            for (int v = t; v != s; v = pv[v]) {
                Edge2 e = g.get(pv[v]).get(pe[v]);
                e.cap -= c;
                g.get(v).get(e.rev).cap += c;
            }
            long d = -dual[s];
            flow += c;
            cost += c * d;
            if (prevCost == d) {
                result.remove(result.size() - 1);
            }
            result.add(new long[] { flow, cost });
            prevCost = cost;
        }
        return result;
    }

    private boolean dualRef(int s, int t, long[] dual, long[] dist, int[] pv, int[] pe, boolean[] vis) {
        Arrays.fill(dist, Long.MAX_VALUE);
        Arrays.fill(pv, -1);
        Arrays.fill(pe, -1);
        Arrays.fill(vis, false);
        PriorityQueue<Q> que = new PriorityQueue<>((o1, o2) -> Long.compare(o1.key, o2.key));
        dist[s] = 0;
        que.add(new Q(0, s));
        while (!que.isEmpty()) {
            int v = que.poll().to;
            if (vis[v]) {
                continue;
            }
            vis[v] = true;
            if (v == t) {
                break;
            }
            for (int i = 0; i < g.get(v).size(); i++) {
                Edge2 e = g.get(v).get(i);
                if (vis[e.to] || e.cap == 0) {
                    continue;
                }
                long cost = e.cost - dual[e.to] + dual[v];
                if (dist[e.to] - dist[v] > cost) {
                    dist[e.to] = dist[v] + cost;
                    pv[e.to] = v;
                    pe[e.to] = i;
                    que.add(new Q(dist[e.to], e.to));
                }
            }
        }
        if (!vis[t]) {
            return false;
        }

        for (int v = 0; v < n; v++) {
            if (vis[v]) {
                dual[v] -= dist[t] - dist[v];
            }
        }
        return true;
    }
}
提出

ACL Practice Contest　E - MinCostFlow
ACL Contest 1　C - Moving Pieces　※O(K+NM)頂点、O(KNM)辺解法
ACL Contest 1　C - Moving Pieces　※O(NM)頂点、O(K+NM)辺解法


Strongly Connected Component（2020/9/23追記）

  • ACLをだいたいそのまま移植したつもり。
  • 頂点が何番目の強連結成分に属しているかを取得できる配列gidはprivateにせず。
/**
 * 強連結成分分解
 */
class SCC {
    /** gid[i]：頂点iが何番目の強連結成分に属しているか */
    int[] gid;

    private final int n;
    private List<Edge> edges;

    private int nowOrd, groupNum;
    private int[] start, low, ord;
    private Edge[] elist;
    private Deque<Integer> visited;

    private class Edge {
        final int from, to;

        public Edge(int from, int to) {
            this.from = from;
            this.to = to;
        }
    }

    /**
     * n頂点0辺のグラフを作る。<br>
     * O(n)
     * 
     * @param n 頂点数
     */
    public SCC(int n) {
        this.n = n;
        edges = new ArrayList<>();
    }

    /**
     * fromからtoへ有向辺を追加する。<br>
     * ならしO(1)
     * 
     * @param from 有向辺の始点（0≦from＜n）
     * @param to 有向辺の終点（0≦to＜n）
     */
    void addEdge(int from, int to) {
        assert 0 <= from && from < n : "from=" + from;
        assert 0 <= to && to < n : "to=" + to;

        edges.add(new Edge(from, to));
    }

    /**
     * <pre>
     * 強連結成分分解
     * 辺の本数をmとして、O(n + m)
     * 
     * ・全ての頂点がちょうど1つずつ、内側の配列のどれかに含まれる。
     * ・内側の配列（1つの強連結成分）内での頂点の順序は未定義。
     * ・外側の配列はトポロジカルソートされている。
     * </pre>
     * 
     * @return 強連結成分ごとの頂点の配列
     */
    int[][] scc() {
        start = new int[n + 1];
        for (Edge e : edges) {
            start[e.from + 1]++;
        }
        for (int i = 1; i <= n; i++) {
            start[i] += start[i - 1];
        }
        elist = new Edge[edges.size()];
        int[] counter = Arrays.copyOf(start, n + 1);
        for (Edge e : edges) {
            elist[counter[e.from]++] = e;
        }

        nowOrd = 0;
        groupNum = 0;
        visited = new ArrayDeque<>();
        low = new int[n];
        ord = new int[n];
        gid = new int[n];
        Arrays.fill(ord, -1);

        for (int i = 0; i < n; i++) {
            if (ord[i] == -1) {
                dfs(i);
            }
        }
        for (int i = 0; i < n; i++) {
            gid[i] = groupNum - 1 - gid[i];
        }

        int[] counts = new int[groupNum];
        for (int x : gid) {
            counts[x]++;
        }
        int[][] groups = new int[groupNum][];
        for (int i = 0; i < groupNum; i++) {
            groups[i] = new int[counts[i]];
        }
        for (int i = 0; i < n; i++) {
            int idx = gid[i];
            groups[idx][--counts[idx]] = i;
        }
        return groups;
    }

    private void dfs(int v) {
        low[v] = ord[v] = nowOrd++;
        visited.add(v);
        for (int i = start[v]; i < start[v + 1]; i++) {
            int to = elist[i].to;
            if (ord[to] == -1) {
                dfs(to);
                low[v] = Math.min(low[v], low[to]);
            } else {
                low[v] = Math.min(low[v], ord[to]);
            }
        }
        if (low[v] == ord[v]) {
            while (true) {
                int u = visited.pollLast();
                ord[u] = n;
                gid[u] = groupNum;
                if (u == v) {
                    break;
                }
            }
            groupNum++;
        }
    }
}
提出

ACL Practice Contest　G - SCC


2-SAT（2020/9/27追記）

  • 使う際は、SCCも一緒に貼り付ける。
/**
 * 2-Satisfiability
 */
class TwoSat {
    /** 最後に呼んだsatisfiableの条件を満たす割当 */
    boolean[] ans;

    private final int n;
    private final SCC scc;

    /**
     * n変数の2-SATを作る。<br>
     * O(n)
     * 
     * @param n 変数の数
     */
    public TwoSat(int n) {
        this.n = n;
        this.scc = new SCC(2 * n);
        ans = new boolean[n];
    }

    /**
     * (x[i] = f) or (x[j] = g) という条件を足す。<br>
     * 制約：0≦i,j＜n<br>
     * ならしO(1)
     */
    void addClause(int i, boolean f, int j, boolean g) {
        assert 0 <= i && i < n : "i=" + i;
        assert 0 <= j && j < n : "j=" + j;

        scc.addEdge(2 * i + (f ? 0 : 1), 2 * j + (g ? 1 : 0));
        scc.addEdge(2 * j + (g ? 0 : 1), 2 * i + (f ? 1 : 0));
    }

    /**
     * 全ての条件を満たす割当が存在するかどうかを判定する。<br>
     * 複数回呼ぶことも可能。<br>
     * 足した条件の数をmとして、O(n + m)
     * 
     * @return true：割当が存在する
     */
    boolean satisfiable() {
        scc.scc();
        int[] id = scc.gid;
        for (int i = 0; i < n; i++) {
            if (id[2 * i] == id[2 * i + 1]) {
                return false;
            }
            ans[i] = id[2 * i] < id[2 * i + 1];
        }
        return true;
    }
}
提出

ACL Practice Contest　H - Two SAT