ACL Contest 1 A,B,C問題メモ

A - Reachable Towns

問題

$xy$ 平面上に $N$ 個の点 $(x_1,y_1),...,(x_N,y_N)$
- $x$ 座標が $1,2,...,N$ の箇所に点が1個ずつある
- $y$ 座標が $1,2,...,N$ の箇所に点が1個ずつある
$k=1～N$ について、以下の問いに答えよ
- 今いる点から「$x$ 座標も $y$ 座標も小さい点」と「$x$ 座標も $y$ 座標も大きい点」に移動できる
- 点 $k$ から出発して、好きなだけ移動を繰り返して到達できる点の個数を（点 $k$ を含めて）求めよ
$1 \le N \le 2 \times 10^5$

解法

ある点から到達できる点のグループは、グループ内のどの点から出発しても同じグループに到達できる。
なので、Union-Find（Disjoint Set Union）を使えばよさそう。

ただ問題は、行き来できる点をどのように見つければよいか。
各点から、他の点が移動できる点かどうかを全て確かめていては $O(N^2)$ で無理。

$x,y$ 座標が大きい点の方から、小さい方向に移動できる点を探す。

  ↙●
  ●  
      ●
          ↙●
        ●  ↙
          ●
              ●
x １２３４５６７

まず点を $x$ 座標でソートし、その順に、該当する点を探していく。

自分より後に処理される点は、当然 $x$ 座標が大きいので、該当する点ではあり得ない。
自分より前に処理された中で、$y$ 座標が小さい点に対して、Uniteすればよい。

これは、std::set やFenwick Tree などで管理できる。

ただ、この方法も、以下のようなケースではUniteする対象がどんどん増え、結局 $O(N^2)$ かかってしまう。

            ...
          ●
        ●
      ●
    ●
  ●  
x １２３４５６７

$x=5$ まで処理された時点では、$x=1～4$ の点は既にUniteされているので、5から辺を張りに行くのはどれか1点でよい。
そう考えると、一番左下の点のみ残しておけばいい具合に無駄を減らせそう。

「$x$ 座標も $y$ 座標も自分より小さい点が無い点」のみをFenwickTreeに登録し、探索対象とする。
その点がハブとなってくれて、Uniteされるグループは正しいまま、探す対象が少なくなる。

（Pythonでは std::set が無いから代用とは言え）いきなりA問題から FenwickTreeとUnionFindの合わせ技が出てくるのは飛ばしてるなあ、と思ったけど、実は使わなくても解ける方法があるらしい。先入観よ。

Python3

import os
import sys

import numpy as np


def solve(inp):
    UNIONFIND_TABLE = []

    def unionfind_init(n):
        UNIONFIND_TABLE.append(np.full(n, -1, dtype=np.int64))
        return len(UNIONFIND_TABLE) - 1

    def unionfind_getroot(ins, x):
        table = UNIONFIND_TABLE[ins]
        stack = []
        while table[x] >= 0:
            stack.append(x)
            x = table[x]
        for y in stack:
            table[y] = x
        return x

    def unionfind_unite(ins, x, y):
        table = UNIONFIND_TABLE[ins]
        r1 = unionfind_getroot(ins, x)
        r2 = unionfind_getroot(ins, y)
        if r1 == r2:
            return
        d1 = table[r1]
        d2 = table[r2]
        if d1 <= d2:
            table[r2] = r1
            table[r1] += d2
        else:
            table[r1] = r2
            table[r2] += d1

    def unionfind_find(ins, x, y):
        return unionfind_getroot(ins, x) == unionfind_getroot(ins, y)

    def unionfind_getsize(ins, x):
        table = UNIONFIND_TABLE[ins]
        return -table[unionfind_getroot(ins, x)]

    FENWICK_TREE = []
    FENWICK_LOGN = []

    def fenwick_init(n):
        log_n = 0
        m = n
        while m:
            log_n += 1
            m >>= 1
        FENWICK_TREE.append(np.zeros(n + 1, dtype=np.int64))
        FENWICK_LOGN.append(log_n)
        return len(FENWICK_TREE) - 1

    def fenwick_add(ins, i, x):
        arr = FENWICK_TREE[ins]
        n = arr.size - 1
        while i <= n:
            arr[i] += x
            i += i & -i

    def fenwick_sum(ins, i):
        arr = FENWICK_TREE[ins]
        result = 0
        while i > 0:
            result += arr[i]
            i ^= i & -i
        return result

    def fenwick_lower_bound(ins, x):
        arr = FENWICK_TREE[ins]
        log_n = FENWICK_LOGN[ins]
        n = arr.size - 1
        sum_ = 0
        pos = 0
        for i in range(log_n, -1, -1):
            k = pos + (1 << i)
            if k < n and sum_ + arr[k] < x:
                sum_ += arr[k]
                pos += 1 << i
        return pos + 1

    n = inp[0]

    bit = fenwick_init(n + 2)
    uft = unionfind_init(n)
    xxx = inp[1::2]
    yyy = inp[2::2]
    idx = np.argsort(xxx)

    y_dict = np.zeros(n + 1, np.int64)
    for i in range(n):
        y_dict[yyy[i]] = i

    for i in idx:
        y = yyy[i]
        k = fenwick_sum(bit, y + 1)

        if k == 0:
            fenwick_add(bit, y + 1, 1)
            continue

        while k > 0:
            u = fenwick_lower_bound(bit, k) - 1
            j = y_dict[u]
            unionfind_unite(uft, i, j)
            k -= 1

    ans = np.zeros(n, np.int64)
    for i in range(n):
        ans[i] = unionfind_getsize(uft, i)

    return ans


if sys.argv[-1] == 'ONLINE_JUDGE':
    from numba.pycc import CC

    cc = CC('my_module')
    cc.export('solve', '(i8[:],)')(solve)
    cc.compile()
    exit()

if os.name == 'posix':
    # noinspection PyUnresolvedReferences
    from my_module import solve
else:
    from numba import njit

    solve = njit('(i8[:],)', cache=True)(solve)
    print('compiled', file=sys.stderr)

inp = np.fromstring(sys.stdin.read(), dtype=np.int64, sep=' ')
ans = solve(inp)
print('\n'.join(map(str, ans)))

B - Sum is Multiple

問題

正整数 $N$ が与えられる
$(1+2+...+k)$ が $N$ の倍数となるような、最小の $k$ を求めよ
$1 \le N \le 10^{15}$

解法

シンプルですごい。

まずパッと見、$(1+2+...+k)=\dfrac{k(k+1)}{2}$ と式変形できる。
よって、「$k(k+1)$ が $2N$ の倍数となる」と言い換えられるのだが、ここで重要なのは $k$ と $k+1$ は互いに素であるということ。

つまり、$2N$ が $2^a3^b5^c$ などと素因数分解されたとして、$k$ と $k+1$ の両方に2や3が素因数として入ることは無い。

$2^a$ も $3^b$ も $5^c$ も、それぞれ $k$ と $k+1$ のどちらかに全入れしなければならない。

素因数のユニークな個数なんて、素数を小さい方からかけていってもまぁ15個もしない内に $10^{15}$ を越えるので、たかだか15個の数についてどっちに入れるか、$2^{15}$ 通り試すのは十分可能である。

さて、$k$ が $A=2^a5^c$ の倍数、$k+1$ が $B=3^b$ の倍数と割り振ったとして、これについての答えを求める。
（ここでの割り振りは一例であって、重要ではない）

中国剰余定理の考え方を使って、以下のように求められる。

$k+1=Bm$（$m$:整数）と表すとして、これを $k \equiv 0 \mod{A}$ に代入すると、$Bm \equiv 1 \mod{A}$ となる。

ここで、$A$ と $B$ は互いに素なので、$B$ は逆元 $B^{-1}$ を持つ。これを使えば $m \equiv B^{-1} \mod{A}$ となり、 $m$ は、$B^{-1}$ をベースとして、それに $A$ の倍数を足し引きした数に絞られる。

正整数の中での最小値を求めたいので、$0 \le m \lt A$ である $m$ を使って $k+1$ を求めればよい。

ただし、$A=1$ や $B=1$ の場合は、うまく逆元を計算できない（というか、何でもよい）ので、個別に考える。

$A=1$ の時、$k \equiv 0 \mod{1}$ かつ $k+1 \equiv 0 \mod{2N}$ である。よって、$k=2N-1$ が最小となる。
$B=1$ も同様に考えると、$k=2N$ となる。

Python3

def mod_inv(a, MOD):
    """ 前提: aとMODが互いに素であること """
    b = MOD
    u = 1
    v = 0
    while b:
        t = a // b
        a -= t * b
        u -= t * v
        a, b = b, a
        u, v = v, u
    u %= MOD
    return u


def prime_factorize(n):
    i = 2
    res = {}
    while i * i <= n:
        if n % i != 0:
            i += 1
            continue
        res[i] = 0
        while n % i == 0:
            res[i] += 1
            n //= i
        i += 1
    if n != 1:
        res[n] = 1
    return res


def solve(n):
    n *= 2
    pf = prime_factorize(n)
    l = len(pf)
    pfk = []
    for i, p in enumerate(pf):
        pfk.append(p ** pf[p])

    ans = 10 ** 18

    for bit in range(1 << l):
        a = 1
        for d in range(l):
            if bit & (1 << d):
                a *= pfk[d]
        b = n // a

        if a == 1:
            k = n - 1
        elif b == 1:
            k = n
        else:
            m = mod_inv(b, a)
            k = b * m - 1
        ans = min(ans, k)

    return ans


n = int(input())
ans = solve(n)
print(ans)

C - Moving Pieces

問題

$N \times M$ の盤面
- . が空きマス
- o がコマ（複数個あり得る）
- # が壁
コマを1つ選んで、右または下に動かすことを繰り返す
- 壁や、他のコマがいるマスには移動できない
どの駒も動かせなくなったら終わり
全ての駒の移動回数の総和の最大値を求めよ
$1 \le N,M \le 50$
$1 \le コマの個数 \le 100$

解法

どのコマをどのマスで終わらせるか、の組み合わせで最小費用流。

各コマは右も下も外周か # のマス（以降、「行き止まりマス」）のいずれかを目指すんだけど、同じマスを目指すコマが複数あると、1つ以外は手前で止まって移動回数を減らさざるを得ない。

最終状態は、行き止まりマスの多くにはコマが置かれ、その周辺にも溢れたコマが密集している形となる。

その「周辺」がたとえば1本道だと、ちょっとのコマが集まっただけで詰まってしまい、より手前で止まらねばならず、移動回数の損失は大きくなる。
一方、周辺が開けていると、多くのコマが集まっても比較的、損失の増え方は小さい。

はじめ「コマと目指す行き止まりマス」の対応を考えたが、それだとこの違いを上手く吸収できない。

「コマと最終的なマス」の対応を考えることで、1対1のマッチングを求めればよくなり、最小費用流に持ち込める。

仮想的な始点・終点（各1つずつ）、コマ（oの数だけ）、各マス（$N \times M$）、の分だけ頂点を作る
始点→コマ→そのコマが移動できるマス→終点と通るように辺を張る
- 始点からコマに、容量1、コスト0の辺を張る
- 各コマから経路探索し、到達できるマスに、容量1、コスト「-移動距離」の辺を張る
- 各マスから終点に、容量1、コスト0の辺を張る

最大スコアを最小費用流によって求めたいときは、正負を逆転させれば使えるようになる。
これで始点から終点に、コマの個数分の流量を流したときの最小費用を求めれば、それが答えとなる。

負辺ができてしまうが、始点から終点への全ての経路で負辺（コマ→マスに張られる辺）を通るのは1回固定なので、適当に下駄を履かせておけばよい。

ところで、この方針ではコマ同士の、途中経路での衝突が考慮されていない。
「あるコマから到達できるマスであっても、他のコマを考慮すると邪魔されて動けなくなることはないか？」という点が気になる。

初期状態      最終状態
                            このような実際は不可能な
①②....  →  ....②①      マッチングが発生してしまう

が、衝突する場合はマッチング結果の対応を入れ替えると解消でき、総移動回数は変わらないので問題ない。
適宜入れ替えれば問題ない、ということが言えたので、実際に矛盾しないマッチングを求める必要は無い。

Python3

import os
import sys
from heapq import heappop, heappush

import numpy as np


def solve(n, m, field):
    MINCOSTFLOW_LINKS = []
    INF = 10 ** 10

    def mincostflow_init(n):
        """ n: 頂点数 """
        lst = [[0]]
        lst.clear()
        MINCOSTFLOW_LINKS.append([lst.copy() for _ in range(n)])
        return len(MINCOSTFLOW_LINKS) - 1

    def mincostflow_add_link(ins, frm, to, capacity, cost):
        """ インスタンスID, 辺始点頂点番号, 辺終点頂点番号, 容量, コスト """
        links = MINCOSTFLOW_LINKS[ins]
        links[frm].append([to, capacity, cost, len(links[to])])
        links[to].append([frm, 0, -cost, len(links[frm]) - 1])

    def mincostflow_flow(ins, s, t, quantity):
        """ インスタンスID, フロー始点頂点番号, フロー終点頂点番号, 要求流量 """
        links = MINCOSTFLOW_LINKS[ins]
        n = len(links)
        res = 0
        potentials = np.zeros(n, dtype=np.int64)
        dist = np.full(n, INF, dtype=np.int64)
        prev_v = np.full(n, -1, dtype=np.int64)
        prev_e = np.full(n, -1, dtype=np.int64)

        while quantity:
            dist.fill(INF)
            dist[s] = 0
            que = [(0, s)]

            while que:
                total_cost, v = heappop(que)
                if dist[v] < total_cost:
                    continue
                for i, (u, cap, cost, _) in enumerate(links[v]):
                    new_cost = dist[v] + potentials[v] - potentials[u] + cost
                    if cap > 0 and new_cost < dist[u]:
                        dist[u] = new_cost
                        prev_v[u] = v
                        prev_e[u] = i
                        heappush(que, (new_cost, u))

            # Cannot flow quantity
            if dist[t] == INF:
                return -1

            potentials += dist

            cur_flow = quantity
            v = t
            while v != s:
                cur_flow = min(cur_flow, links[prev_v[v]][prev_e[v]][1])
                v = prev_v[v]
            quantity -= cur_flow
            res += cur_flow * potentials[t]

            v = t
            while v != s:
                link = links[prev_v[v]][prev_e[v]]
                link[1] -= cur_flow
                links[v][link[3]][1] += cur_flow
                v = prev_v[v]

        return res

    nm = (n + 2) * (m + 2)
    m2 = m + 2

    starts = np.where(field == 1)[0]
    s_size = starts.size

    ins = mincostflow_init(nm + s_size + 2)
    s = nm + s_size
    t = s + 1

    stack = np.zeros(10 ** 7, np.int64)
    SENTINEL = INF - 1

    for i in range(s_size):
        mincostflow_add_link(ins, s, nm + i, 1, 0)

        stack[0] = starts[i]
        stack[1] = 0
        sl = 0
        sr = 2
        stacked = np.zeros(nm, np.int8)
        stacked[starts[i]] = 1
        while sl < sr:
            v = stack[sl]
            c = stack[sl + 1]
            sl += 2

            mincostflow_add_link(ins, nm + i, v, 1, SENTINEL - c)

            if stacked[v + 1] == 0 and field[v + 1] != 2:
                stack[sr] = v + 1
                stack[sr + 1] = c + 1
                stacked[v + 1] = 1
                sr += 2
            if stacked[v + m2] == 0 and field[v + m2] != 2:
                stack[sr] = v + m2
                stack[sr + 1] = c + 1
                stacked[v + m2] = 1
                sr += 2

    for i in range(nm):
        if field[i] != 2:
            mincostflow_add_link(ins, i, t, 1, 0)

    f = mincostflow_flow(ins, s, t, s_size)

    ans = SENTINEL * s_size - f
    return ans


SIGNATURE = '(i8,i8,i1[:],)'
if sys.argv[-1] == 'ONLINE_JUDGE':
    from numba.pycc import CC

    cc = CC('my_module')
    cc.export('solve', SIGNATURE)(solve)
    cc.compile()
    exit()

if os.name == 'posix':
    # noinspection PyUnresolvedReferences
    from my_module import solve
else:
    from numba import njit

    solve = njit(SIGNATURE, cache=True)(solve)
    print('compiled', file=sys.stderr)

n, m = map(int, sys.stdin.readline().split())
field = '#' * (m + 3) + '##'.join(sys.stdin.read().split()) + '#' * (m + 3)
field = np.fromiter(map('.o#'.index, field), np.int8)
ans = solve(n, m, field)
print(ans)

よりよい解法

コマ1つ1つについて各マスへのコストを求めなくても、通常のグリッド問題のように、隣接するマスにコスト $-1$ の辺を張る方法でも、最小費用流が可能である。

そうすると各コマで辺を共有できるので、高速になる。

先ほどの方法では、コマの数 $p$、壁でないマスの数が $NM$ に比例するとして、だいたい $p+NM+pNM$ 本の辺が必要
- $p \le 100, NM \le 2500$ より、ざっくり25万本
この方法では、だいたい $p+3NM$ 程度でいける
- ざっくり7500本

ただ、その場合、負辺を防ぐために下駄を履かせたくても、各コマによって負辺を通る回数が違う＝足される下駄の回数が異なるので、正しく求められない。

以下のように考えると下駄の回数を統一した上で辺を非負にできる。

全てのコマは、$(1,1)$ から $(N,M)$ への $N+M-2$ 回の移動分のコストがデフォルトでかかる（下駄）
その中で、実際におこなった移動分だけを0にできる
$(N+M-2) \times コマの個数$ 分の下駄が履かされた状態なので、ここから最小費用を引くと、答えとなる

これは、以下のように実装する。

始点、終点、空きマスの数だけの頂点を用意
始点から、各コマが初期で置かれているマス $(x_i,y_i)$ へ辺を張る
- 容量は1、コストは $x_i+y_i-2$　（$(1,1)$ からの移動分）
壁以外のマス $(x_j,y_j)$ から、終点へ辺を張る
- 容量は1、コストは $N+M-x_j-y_j$　（$(N,M)$ への移動分）
マスから右か下に隣接する移動できるマスに辺を張る
- 容量は∞、コストは $0$

Python3

import os
import sys
from heapq import heappop, heappush

import numpy as np


def solve(n, m, field):
    MINCOSTFLOW_LINKS = []
    INF = 10 ** 18

    def mincostflow_init(n):
        """ n: 頂点数 """
        lst = [[0]]
        lst.clear()
        MINCOSTFLOW_LINKS.append([lst.copy() for _ in range(n)])
        return len(MINCOSTFLOW_LINKS) - 1

    def mincostflow_add_link(ins, frm, to, capacity, cost):
        """ インスタンスID, 辺始点頂点番号, 辺終点頂点番号, 容量, コスト """
        links = MINCOSTFLOW_LINKS[ins]
        links[frm].append([to, capacity, cost, len(links[to])])
        links[to].append([frm, 0, -cost, len(links[frm]) - 1])

    def mincostflow_flow(ins, s, t, quantity):
        """ インスタンスID, フロー始点頂点番号, フロー終点頂点番号, 要求流量 """
        links = MINCOSTFLOW_LINKS[ins]
        n = len(links)
        res = 0
        potentials = np.zeros(n, dtype=np.int64)
        dist = np.full(n, INF, dtype=np.int64)
        prev_v = np.full(n, -1, dtype=np.int64)
        prev_e = np.full(n, -1, dtype=np.int64)

        while quantity:
            dist.fill(INF)
            dist[s] = 0
            que = [(0, s)]

            while que:
                total_cost, v = heappop(que)
                if dist[v] < total_cost:
                    continue
                for i, (u, cap, cost, _) in enumerate(links[v]):
                    new_cost = dist[v] + potentials[v] - potentials[u] + cost
                    if cap > 0 and new_cost < dist[u]:
                        dist[u] = new_cost
                        prev_v[u] = v
                        prev_e[u] = i
                        heappush(que, (new_cost, u))

            # Cannot flow quantity
            if dist[t] == INF:
                return -1

            potentials += dist

            cur_flow = quantity
            v = t
            while v != s:
                cur_flow = min(cur_flow, links[prev_v[v]][prev_e[v]][1])
                v = prev_v[v]
            quantity -= cur_flow
            res += cur_flow * potentials[t]

            v = t
            while v != s:
                link = links[prev_v[v]][prev_e[v]]
                link[1] -= cur_flow
                links[v][link[3]][1] += cur_flow
                v = prev_v[v]

        return res

    nm = (n + 2) * (m + 2)
    m2 = m + 2

    starts = np.where(field == 1)[0]
    s_size = starts.size

    ins = mincostflow_init(nm + 2)
    s = nm
    t = s + 1

    for i in starts:
        x, y = divmod(i, m2)
        mincostflow_add_link(ins, s, i, 1, x + y - 2)

    for i in range(nm):
        if field[i] == 2:
            continue
        x, y = divmod(i, m2)
        mincostflow_add_link(ins, i, t, 1, n + m - x - y)
        if field[i + 1] != 2:
            mincostflow_add_link(ins, i, i + 1, INF, 0)
        if field[i + m2] != 2:
            mincostflow_add_link(ins, i, i + m2, INF, 0)

    f = mincostflow_flow(ins, s, t, s_size)

    ans = (n + m - 2) * s_size - f
    return ans


SIGNATURE = '(i8,i8,i1[:],)'
if sys.argv[-1] == 'ONLINE_JUDGE':
    from numba.pycc import CC

    cc = CC('my_module')
    cc.export('solve', SIGNATURE)(solve)
    cc.compile()
    exit()

if os.name == 'posix':
    # noinspection PyUnresolvedReferences
    from my_module import solve
else:
    from numba import njit

    solve = njit(SIGNATURE, cache=True)(solve)
    print('compiled', file=sys.stderr)

n, m = map(int, sys.stdin.readline().split())
field = '#' * (m + 3) + '##'.join(sys.stdin.read().split()) + '#' * (m + 3)
field = np.fromiter(map('.o#'.index, field), np.int8)
ans = solve(n, m, field)
print(ans)

目次

ACL Contest 1 A,B,C問題メモ

A - Reachable Towns

問題

解法

B - Sum is Multiple

問題

解法

C - Moving Pieces

問題

解法

よりよい解法