−目次

DISCO presents ディスカバリーチャンネルコードコンテスト2020 本戦 A～C問題メモ

DISCO presents ディスカバリーチャンネルコードコンテスト2020 本戦 A～C問題メモ

DISCO presents ディスカバリーチャンネルコードコンテスト2020 本戦

ちょっと難しめの典型手法を一ひねりしたような問題が並んでいる。

多分これを使うんだろうなということは想像できても、それでいい証明はなかなか難しい。

A - Div/de

問題

紙に $N$ 個の正整数 $A_1,A_2,...,A_N$ が書かれている
先攻と後攻に分かれてゲームを行う
自分の手番の時、1以外の整数を1つ選び、その約数（自身を除く）に書き換える
書き換えられなくなった方の負け
互いが最善を尽くすとき、どちらが必勝か答えよ
$1 \le N \le 100$
$1 \le A_i \le 100$

解法

石取りゲームNimのアレンジ。

石取りゲーム(Nim)

各整数についてGrundy数を求める。Grundy数とは、石取りゲームの1つの山（今回は1つの整数）について、そこから遷移できない状態の中で最小の非負整数。

「1」のGrundy数は0とする
2以上の数のGrundy数は、約数のGrundy数の最大値+1となる
- 自身の約数から遷移できる状態へは当然、自身からも遷移できるため、約数のGrundy数以下になることは無い

$O(N \sqrt{A_i})$ なので、全て約数列挙して確かめても間に合う。

後はGrundy数を全てXORして、0なら後手必勝、それ以外なら先手必勝。

適当に $A_i$ を決めると先手必勝になるケースが多いのに、サンプルでは後手必勝ケースを多めに入れているのがいやらしい?

        
              def make_divisors(n):
    i = 1
    while i * i <= n:
        if n % i == 0:
            yield i
            yield n // i  # 留意: nが平方数の時重複して返される
        i += 1
 
grundy = [-1] * 101
grundy[1] = 0
 
for i in range(2, 101):
    can = set(range(101))
    for d in make_divisors(i):
        if d != i:
            can.discard(grundy[d])
    grundy[i] = min(can)
 
n = int(input())
aaa = list(map(int, input().split()))
 
p = 0
for a in aaa:
    p ^= grundy[a]
 
print('Yes' if p != 0 else 'No')

B - Hawker on Graph

問題

$N$ 頂点 $M$ 辺の有向重み付きグラフがある
$i$ 番目の辺は頂点 $u_i→v_i$ を重み $w_i$ で繋ぐ
はじめ $W$ 円を持ち、頂点 $S$ を始点として、グラフ上を $K$ 回ちょうど移動する
- 同じ辺・頂点は何回でも通ってよい
- 終了時はどの頂点にいてもよい
- 移動するたびに、辺の重みだけ所持金が加算される
- 負の重みなら所持金は減るが、結果が0未満になる場合は0円となる
$K$ 回ちょうどの移動を行えるか判定し、可能なら終了時に可能な最大所持金を求めよ
$2 \le N \le 200$
$0 \le M \le N(N-1)$
$0 \le W \le 10^{16}$
$1 \le K \le 10^9$
$-10^7 \le w_i \le 10^7$
実行制限時間: 6sec

解法

解説AC。演算子カスタマイズ行列累乗。

行列積による最短距離の計算

所持金とかは一旦無視して、普通の有向グラフ上の最短距離を考える。

グラフ上を $K$ 回移動した時の全頂点間最短距離は、トロピカル半環上の演算で、隣接行列の $K$ 乗によって求めることができるらしい。

グラフの隣接行列を $A$ （移動不可は∞）、トロピカル上の乗算を $\odot$ で表すとして $K$ 回の移動後の最短距離は以下で求められる。

初期状態
単位行列         隣接行列のK乗
 0 ∞ ∞
∞  0 ∞    ⊙    A^⊙K
∞ ∞  0

半環
トロピカル幾何学
動的計画法を実現する代数〜トロピカル演算でグラフの最短経路を計算する〜 - Qiita
Abc009 - D 解説半環のざっくりした意味・意義

ある集合とそれに対する2種類の演算がいくつかの条件を満たせば、何でも半環と呼ぶことができる。半環なら○○してもよい、と証明済みの何らかの処理があれば、半環であることさえ調べればその処理を安心して適用できる（いちいち個別に証明しなくてよい）のが、半環という概念を持ってくる利点である。

今回の場合、グラフの移動を行列の積で表現するが、半環では結合法則が成り立ちどこから先に計算してもよいので、累乗の計算に安心してダブリングが使える（ということなんだと思う。この辺の関係性はよく理解できてない）。

ダブリング - sataniC++

詳細は解説pdfとか、Webの記事とかに頼るとして、とりあえず計算に必要な情報は以下。

	一般的な代数	トロピカル
零元	0	$\infty$
単位元	1	0
加法	$A+B$	$\min(A,B)$
乗法	$A \times B$	$A+B$

ただし今回の場合、最大値を求めたいので、minでなくmaxを使い、 $\infty$ も $-\infty$ となる。

破産の考慮

今回は、途中で所持金が0になるかも知れないことを考慮する必要がある。（以降、「破産」と呼ぶ）

悩ましいのは、途中で破産してでも、それ以降高い金額の獲得が可能となるグラフに移った方がよくなる、なども考えられること。

Ｓ--1→○-- -1億 -→○--100→○
↑1   1↓           ↑100 100↓
○←1--○           ○←100--○

そのため単純に最も稼げる（重みの和が大きくなる）ルートを探すだけではダメ。かといって始点から順番に計算するとなると、-1億の辺を渡る段階ではその後に本当にそれ以上稼げるのか、 $K$ 回移動できるのかは見通せない。

通った経路が同じ場合、1回でも破産をしていれば、最初の所持金がどうあれ最終結果は同じになる。（破産を考慮せずに描いた所持金の推移グラフの、最低値を記録してからの増分が、最終的な所持金となる）

￥| ／＼                    | ／＼
  |     ＼                  |     ＼          ／
  |_______＼／＼____／  →  |_______＼／＼__／___
  |             ＼／        |

なので、1回でも破産をしたかどうかでわけるのがよさそう。

「破産していない世界線」「破産を経験した世界線」2つの状態を同時に考えつつ、全点間距離を求める必要がある。

これが、解説pdfの $\max(w+a,b)$ の $a,b$ に相当する。実装上は別個の行列として持ち、更新時に情報をマージする。

$A_{k,i,j}$ …初期資金が潤沢にあると仮定し一度も破産しない場合の、 $k$ 回かけて $i$ から $j$ に移動した場合の最高獲得金額
$B_{k,i,j}$ … $k$ 回かけて $i$ から $j$ への移動中に1回以上破産している場合に得られる最高獲得金額（ $i$ への移動で破産する場合も含む）

実際は $k$ はダブリングで更新していき、 $i,j$ の2次元の $N \times N$ 行列で持つ。

$A$ について、途中で $-W$ を下回ったらそれ以降意味ないじゃないか、という気が一瞬するが、これは別に始点からの移動のみに使うのではなく、「ある程度移動して稼いだ後、さらに $i$ から $j$ に行く」際の情報にもなる。そこで破産するかどうかは未知数なため、 $A$ に下限はない。

$A$ の更新は簡単で、ただ行列積を行えばよい。破産条件が無い場合の最良経路を求める場合と変わらない。

$A_{k+l} = A_{k} \odot A_{l}$

$B$ の更新は、以下のようになる。

$B_{k+l} = \max(B_{k} \odot A_{l}, Z_{k} \odot B_{l})$

ただし $Z_{k}$ は、 $k$ 回で到達可能な頂点組のコストを0にリセット、他は $-\infty$ にしたものとなる。 $\max$ は、2つの行列の各要素毎に大きい方を取る処理とする。

言葉で表すと、以下のようになる。

$B_{k} \odot A_{l}$ … 前半の $k$ 回中に破産してから、後半の $l$ 回では破産せず移動した
$Z_{k} \odot B_{l}$ … 後半の $l$ 回中に破産した

Python

計算量 $O(N^3 \log{K})$ はそのまま代入して約 $2.4 \times 10^8$ 、しかも行列アクセス多用となると、さすがに6秒制限ではPythonでは無理、PyPyでも厳しい。（上手い書き方があるかも知れないが）

NumPyでは、通常の行列積は np.dot() で一発だが、演算をカスタマイズする場合もNumPyに関数があれば一括計算させる記法がある。

vectorization - Numpy matrix multiplication with custom dot product - Stack Overflow

トロピカル演算
c = np.max(a[..., np.newaxis] + b[np.newaxis, ...], axis=1)

乗算を「等しければ1, 異なれば0」という演算に置きかえた場合
c = np.sum(a[..., np.newaxis] == b[np.newaxis, ...], axis=1)

これを使えば一番重たいパートを高速化でき、TLEは防げる。

次の問題は、何故か1ケースだけWAになってしまうこと。想像だが、制約が大きすぎてfloat64では精度が足りなくなってしまうのが原因かも知れない。

今回 $\infty$ を表現するため手軽なfloat64型を選んだが、 float64の仮数部は符号を含めず52bitで、 $2^{53} \simeq 9 \times 10^{15}$ までしか区別できない。制約では $10^{16}$ が最大なので、ギリ足りない。 float128を使えば解決すると思われるが、すると今度は演算が重くなりすぎてTLE。

提出時はエスパーをしてそのケースだけ-1すると無理矢理通った（よくない）。

int64型にすると精度は足りる。 -INFが何回も足されるとオーバーフローしてしまうので、1回足されただけではオーバーフローしない数（ $-2^{62}$ ）にしておき、1回毎に-INFを下回る値は-INFにリセットする、という処理を入れればこちらも通った。……が、ぱっと見、何やってるかわかりにくい。

Numpyは「値がINFの場合は処理を飛ばす」など柔軟な措置が執りにくいのはもどかしい。他には1回毎にMODをとらないとオーバーフローする累積和とか。（とはいえPythonで演算が高速に出来るのは大変助かるというのは大前提）

        
          float型での実装（アドホックな処理は略）
          
              
              import sys
 
import numpy as np
 
 
def tropical_dot(a, b, c, d):
    e = np.max(a[..., np.newaxis] + c[np.newaxis, ...], axis=1)
    f = np.max(b[..., np.newaxis] + c[np.newaxis, ...], axis=1)
    a[a != -np.inf] = 0  # a,cはもう使わないので, 経路情報だけ利用
    g = np.max(a[..., np.newaxis] + d[np.newaxis, ...], axis=1)
    return e, np.maximum(f, g)
 
 
n, m, y, s, k = map(int, input().split())
graph = np.full((n, n), -np.inf, dtype=np.float64)
restart = np.full((n, n), -np.inf, dtype=np.float64)
 
for line in sys.stdin:
    u, v, w = map(int, line.split())
    u -= 1
    v -= 1
    graph[u, v] = w
    restart[u, v] = 0
 
da = np.full((n, n), -np.inf, dtype=np.float64)
db = np.full((n, n), -np.inf, dtype=np.float64)
np.fill_diagonal(da, 0)
np.fill_diagonal(db, 0)
tc = graph
td = restart
 
while k:
    if k & 1:
        da, db = tropical_dot(da, db, tc, td)
    tc, td = tropical_dot(tc, td, tc, td)
    k >>= 1
 
s -= 1
if da[s].max() == -np.inf:
    print(-1)
else:
    das = int(np.around(da[s].max()))
    dbs = int(np.around(db[s].max()))
    print(max(das + y, dbs))

            

        
      

C - Smaller-Suffix-Free Sequences

問題

数列 $T=(T_1,T_2,...,T_k)$ が「Smaller-Suffix-Free」であるとは、以下の条件を満たすことを指す
- $T$ の全てのSuffix $(T_1,...,T_k), (T_2,...,T_k), ..., (T_k)$ を考える
- この中で、 $T$ 自体が辞書順で最小である
長さ $N$ の数列 $A_1,A_2,...,A_N$ が与えられる
各 $i=1～N$ につき、以下を求めよ
- $(A_i,A_{i+1},...,A_j)$ がSmaller-Suffix-Freeとなる最大の $j$
$1 \le N \le 2 \times 10^5$
$1 \le A_i \le 2 \times 10^5$

解法

$A$ の $i～j$ 要素目から成る数列を $A_{i,j}$ と表す。

いかにもSuffix Arrayを使いそうな問題だが、開始位置 $i$ だけでなく終端位置 $j$ も変化するので、単純に元の数列 $A$ に対してSuffix Arrayを構築しただけではどうにもならなそうに見える。

しかし、例えば $i=3$ について考えるとして、Suffix Arrayが辞書順に

a7 a8
...
a3 a4 a5 a6 a7 a8
...
a5 a6 a7 a8
a4 a5 a6 a7 a8
a6 a7 a8

などとなっているとする。

$i=3$ 以降ではじめて $A_{3,N}$ より辞書順で前に来るのが $A_{7,N}$ である時、これだけから $A_{3,j}$ がsmaller-suffix-freeを満たせる最大の $j$ は、6であると言えてしまう。

以下の2つを言えれば、それが証明できる。

$j = 6$ ではsmaller-suffix-freeを満たす
$j \ge 7$ の全てで満たさなくなる

1.

$i=3, j=6$ の場合のSuffixはこの4つである。この辞書順を考える。

a3 a4 a5 a6
a4 a5 a6
a5 a6
a6

元のSuffix Arrayでは $A_{3,N}$ が一番先に来る。

ここで $A_{3,6}$ と各 $A_{k,6}$ で、長さを $A_{k,6}$ の方に揃えた同士での比較をする。（つまり、 $A_{3,9-k}$ と $A_{k,6}$ で比較する）

元のSuffix Arrayから当然、 $(a_3) \le (a_6)$ 、 $(a_3,a_4) \le (a_5,a_6)$ など、 $A_{3,9-k}$ の方が小さいか同じ。

小さい場合は、疑いなく $A_{3,6} \lt A_{k,6}$ と言える。

唯一、等しくなる場合に、より長さが短い $A_{k,6}$ の方が逆転して先になってしまうが……

どういう場合に等しくなってしまうかというと、 $A_{3,4}=A_{5,6}$ となったとして、数列を末尾まで復元してみる。

a3 a4 a5 a6 a7 a8
a5 a6 a7 a8

ここで、 $A_{3,4}=A_{5,6}$ なので、差が生じているとしたらそれ以降である。それ以降にのみ着目すると、

a5 a6 a7 a8
a7 a8

$A_{3,N} \lt A_{5,N}$ を満たすには、 $A_{5,N} \lt A_{7,N}$ でなければならない。しかし、元のSuffix Arrayでは、 $A_{7,N} \lt A_{3,N} \lt A_{5,N}$ なので矛盾する。

よって等しくなることはなく、 $A_{k,6} \ (k=4～6)$ は全て $A_{3,6}$ より辞書順で後、つまり $A_{3,6}$ はsmaller-suffix-freeを満たすと言える。

2.

どんな $j \ge 7$ を取ろうと、 $A_{7,j} \lt A_{3,j}$ となることを示す。

1.と同様に、長さを $A_{7,j}$ の方に揃えて比較をすると、

a7 a8 ... aj
a3 a4 ... a[j-4]

Suffix Arrayでの順序より $A_{7,j} \le A_{3,j-4}$ である。今回の場合、長さが $A_{7,j}$ の方が短いので、たとえ同じになったとしてもこの順序が覆ることはない。

以上、具体的な例だけ取り上げたので一般化できてないが、こんな感じで証明でき、最初に述べたとおり、 $i$ 以降にはじめて $A_{k,N} < A_{i,N}$ となる $k$ に対し、 $k-1$ が $i$ に対する答えとなる。

Suffix Array

SA-IS法というのを使うと速いらしい。

SA-IS 法のメモ - まめめも

感想

Suffix Arrayっぽさはわかっても、そこから解答となる特徴にどうたどり着けるかちょっと見えない。さらにその証明(1)も「prefixが一致した場合は証明が無理そうに見えるが、実はそのままいくと矛盾」という信じて突き進まなければいけない壁がある。

コード（長いので折りたたみ）

        
          
          
              
              from collections import Counter
from itertools import accumulate
from operator import itemgetter
 
 
def induced_sort(n, aaa, lsa, lmss, a_acc, a_dic, ak):
    """
    :param lsa: L-S Array
    :param lmss: Leftmost S の配列（順番により結果は変化する）
    :param a_acc: 頭文字がa_dic[a]の要素がSA上で並び始めるindex
    :param a_dic: aaaの各要素に昇順に割り当てる連番
    :param ak: aaaに出現する要素の種類数（len(a_dic)）
    :return:
    """
 
    sa = [-1] * n
 
    counter = [0] * ak
    for i in reversed(lmss):
        di = a_dic[aaa[i]]
        sa[a_acc[di + 1] - 1 - counter[di]] = i
        counter[di] += 1
 
    counter = [0] * ak
    for i in sa:
        if i <= 0:
            continue
        if lsa[i - 1] == 0:
            continue
        di = a_dic[aaa[i - 1]]
        sa[a_acc[di] + counter[di]] = i - 1
        counter[di] += 1
 
    counter = [0] * ak
    for i in reversed(sa):
        if i <= 0:
            continue
        if lsa[i - 1] == 1:
            continue
        di = a_dic[aaa[i - 1]]
        sa[a_acc[di + 1] - 1 - counter[di]] = i - 1
        counter[di] += 1
 
    return sa
 
 
def sa_is(n: int, aaa: list):
    """
    :param n: aaaの長さ
    :param aaa: Suffix Arrayを求めたい数列
               （初回呼び出し時は末尾に最小ダミー要素を1つ追加すること）
    :return:
    """
 
    a_cnt = Counter(aaa)
    a_kind = sorted(a_cnt.keys())
    a_dic = {a: i for i, a in enumerate(a_kind)}
    a_acc = list(accumulate([0] + list(map(a_cnt.get, a_kind))))
    ak = len(a_kind)
 
    # LS array
    lsa = [0] * n  # 0:S, 1:L
    for i in range(n - 2, -1, -1):
        if aaa[i] == aaa[i + 1]:
            lsa[i] = lsa[i + 1]
        elif aaa[i] > aaa[i + 1]:
            lsa[i] = 1
 
    # LeftMost S
    lmss = [i for i, (lsl, lsr) in enumerate(zip(lsa, lsa[1:]), start=1)
            if lsl == 1 and lsr == 0]
 
    sa = induced_sort(n, aaa, lsa, lmss, a_acc, a_dic, ak)
 
    lms_set = set(lmss)
    sa = [i for i in sa if i in lms_set]
 
    nums = [-1] * n
    nums[sa[0]] = num = 0
 
    for i, j in zip(sa, sa[1:]):
 
        for d in range(n):
            is_id_lms = i + d in lms_set
            is_jd_lms = j + d in lms_set
            if aaa[i + d] != aaa[j + d] or is_id_lms != is_jd_lms:
                num += 1
                break
            elif d > 0 and (is_id_lms or is_jd_lms):
                break
 
        nums[j] = num
 
    nums = [i for i in nums if i != -1]
 
    if num < len(nums) - 1:
        sa = sa_is(len(nums), nums)
        lmss2 = list(map(lmss.__getitem__, sa))
    else:
        lmss2 = [lmss[i] for i, j in sorted(enumerate(nums), key=itemgetter(1))]
 
    sa = induced_sort(n, aaa, lsa, lmss2, a_acc, a_dic, ak)
 
    return sa
 
 
class MinBit:
    def __init__(self, n, SENTINEL):
        self.size = n
        self.tree = [SENTINEL] * (n + 1)
        self.SENTINEL = SENTINEL
 
    def get_min(self, i):
        s = self.SENTINEL
        while i > 0:
            s = min(s, self.tree[i])
            i ^= i & -i
        return s
 
    def add(self, i, x):
        while i <= self.size:
            self.tree[i] = min(x, self.tree[i])
            i += i & -i
 
 
def solve(n, aaa):
    sa = sa_is(n + 1, aaa + [-1])
 
    bit = MinBit(n, n)
    ans = [0] * n
    for s in sa[1:]:
        k = n - s
        ans[s] = bit.get_min(k)
        bit.add(k, s)
    return ans
 
 
n = int(input())
aaa = list(map(int, input().split()))
print('\n'.join(map(str, solve(n, aaa))))

            

        
      

−目次

DISCO presents ディスカバリーチャンネル コードコンテスト2020 本戦 A～C問題メモ