ABC331 勉強メモ

先日開催された、AtCoder Beginner Contest 331の勉強メモです。
問題への言及がありますので、今後解く予定の方はご注意ください。
使用言語はPython・提出はPyPy3です。

(追記)G問題を解いたので追記しました。

# 公式解説

https://atcoder.jp/contests/abc331/editorial



# A Tomorrow

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_a)

__解法__
`if d > D: d = 1` のように、月をまたぐときは日付を1日に変更しましょう。
`if d > D: d %= D` だとWAします。

<details><summary>WAの理由</summary>

`D = 1`のとき、`d %= D`だと翌日が0日になってしまいます。(終)
~~~:sample input (sample output: 2023 6 1)
12 1
2023 5 1
~~~
~~~python: WA例(入力は上記のサンプル)
M,D = map(int,input().split())
y,m,d = map(int,input().split())
d += 1
if d > D: d %= D; m += 1
if m > M: m %= M; y += 1
print(y,m,d)  #2023 6 0
~~~

</details>


# B Buy One Carton of Milk

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_b)

__解法__
各パック数の全探索を行います。計算量は$O(N^3)$です。
2種類のパックの購入数を固定して、残りのパックは必要最小限の量だけ買うことにしてもよいです。計算量は$O(N^2)$です。
DPの場合、状態の持ち方はいろいろと考えられます。
`DP1[x]: 卵をx個買うために必要な最小金額`
`DP2[i][x]: i種類目のパックまで見たとき、卵をx個買うために必要な最小金額`
`DP3[c]: c円払ったとき、卵を購入できる最大数`
のような方法がありますが、最初のDPがもっとも書きやすいと思います。
計算量はパック種類数を$p(=3)$・1パックの最高金額を$X(<10^4)$として、順に$O(pN)$, $O(pN^2)$, $O(pNX)$です。
提出例: [全探索 O(N^2)](https://atcoder.jp/contests/abc331/submissions/48177306), [DP1](https://atcoder.jp/contests/abc331/submissions/48177452), [DP2](https://atcoder.jp/contests/abc331/submissions/48177554), [DP3](https://atcoder.jp/contests/abc331/submissions/48177670)

# C Sum of Numbers Greater Than Me

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_c)

__解法__
毎回`sum(a for a in A if a > A[i])`のように数えていると$O(N^2)$でTLEします。
なので見方を変えて、index順ではなくA[i]の数字順に見ることにします。
そのためにはAに含まれる数字の昇順(または降順)に整頓できればよく、これはdefaultdictが便利です。

~~~python: defaultdictでまとめる
from collections import defaultdict
D = defaultdict(list)  #D[i] のようにアクセスすると [] が用意される

#同じ数字をまとめる
for i in range(N):
    num = A[i]
    D[num].append(A[i])

#降順探索
for num in sorted(D.keys(), reverse = True):  #大きい→小さい
    for i in D[num]:
        print(A[i] == num)  #True
~~~
計算量はdefaultdictのソート部分がもっとも重く、$O(NlogN)$です。
なお本番はセグメント木を使いましたが、計算量が悪化するので非推奨です。
提出例: [変数で保持](https://atcoder.jp/contests/abc331/submissions/48177925), [SegTree](https://atcoder.jp/contests/abc331/submissions/48177940)


# D Tile Pattern

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_d)

__解法__
二次元累積和を取って、黒く塗られたマスを高速に数え上げる用意をします。
二次元累積和の原理はimosさんの解説が勉強になります。

https://imoz.jp/algorithms/imos_method.html

二次元累積和を知っていればあとは実装を頑張るだけです。
区間の上下左右斜め中央の9方向を見てもよいですが、公式解説の通り4方向を見た方が実装が楽･･･だと思います。
9方向を見る場合、「上下左右」がそれぞれ`A`, `C+1`, `B`, `D+1`に対応します。
見比べるとデバッグしやすいです。
提出例: [9way](https://atcoder.jp/contests/abc331/submissions/48152836), [4way](https://atcoder.jp/contests/abc331/submissions/48152827)

# E Set Meal

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_e)

__解法__
Bを降順にソートしましょう。
`C = sorted([(B[i],i) for i in range(M)], reverse = True)`
とすると、`C[j][0]`には副菜の値段, `C[j][1]`には元々の添字が格納されます。

__heapq解__
(主菜, 最高額の副菜) のペアからはじめて、副菜の降順に探索します。
heapqには(合計金額, 主菜の添字, 副菜の添字)を格納し、合計金額が __大きい__ 順に(主菜, 副菜)の組を調べます。
ここで、Pythonのheapqは値の __小さい__ 順にソートされる点に注意してください。
なので実装上、heapqには ( __-1 *__ 合計金額, 主菜, 副菜) と正負反転した値を入れる必要があります。
計算量は副菜ソートに$O(MlogM)$、ヒープ化に$O(N)$、ヒープの処理が$O(LlogN)$、悪い食べ合わせの判定が$O(L)$です。

なお、処理開始前のヒープ化を忘れるとサンプル3で落ちます。
~~~python: ヒープ化忘れの例
import heapq
Q = [9,8,7,6,5]
heapq.heappop(Q)  #9

R = [9,8,7,6,5]
heapq.heapify(R)  #O(N)
heapq.heappop(R)  #5
~~~
提出例: [heapq(620ms)](https://atcoder.jp/contests/abc331/submissions/48153603), [SortedListで代用(753ms)](https://atcoder.jp/contests/abc331/submissions/48178905), [SortedSetで代用(778ms)](https://atcoder.jp/contests/abc331/submissions/48178907)

__二分探索・主菜の全探索__
heapqを使わず、すべての主菜に対して「選べる最高額の副菜は？」の判定を行います。副菜が選べた時点で次の主菜に切り替えることで、この判定は$O(N+L)$と線形時間で可能です。

別解として、解の二分探索に載せることもできます。
具体的には「○円以上の定食が存在するか？」の判定を行います。この判定は上記の判定法で可能です。
計算量は最大金額をXとして、$O((N+L)logX)$になりますが十分高速です。実装例では練習として解の二分探索で実装しています。
提出例: [二分探索・主菜全探索(651ms)](https://atcoder.jp/contests/abc331/submissions/48154361)

__セグメント木・sortedcontainers__
ユーザ解説の通りです。練習で実装しました。
解法の詳細は公式解説のリンクをご確認ください。

https://atcoder.jp/contests/abc331/editorial

なお、本問はSortedSet・SortedListどちらでもACできます。
提出例: [SegTree(241ms)](https://atcoder.jp/contests/abc331/submissions/48153600), [SortedList(573ms)](https://atcoder.jp/contests/abc331/submissions/48178930), [SortedSet(578ms)](https://atcoder.jp/contests/abc331/submissions/48178932)


# F Palindrome Query

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_f)

__解法__
前から読んだ文字列と後ろから読んだ文字列が一致すれば回文です。
この文字列の一致判定には、Rolling Hashを用いることにします。
法や基数の選択法、高速な実装法はkeymoonさんの解説が詳細です。

https://qiita.com/keymoon/items/11fac5627672a6d6a9f6

<details><summary>Rolling Hash 個人用メモ</summary>

(追記)誤字を修正しました。

ローリングハッシュとは、文字列や数列を○進数とみなすことで、10進数の整数値で文字列の一致判定を行えるようにするアルゴリズムです。
変換元の進数を基数(base)、10進数型整数の丸め範囲を法(mod)と呼んで管理します。
~~~:Rolling Hash  base = 100
数列: [2, 3, 4, 5, 7]を100進数表記とみなし、ハッシュ値に変換します。
[2, 3]          → 203
      [4, 5, 7] →     40507
   [3, 4, 5]    →   30405
[2, 3, 4, 5, 7] → 203040507

前計算すれば累積和の要領で、狙った範囲のハッシュ値をO(1)で取得できます。
      [4, 5, 7] →     40507 のハッシュ値が知りたい場合
[2, 3, 4, 5, 7] → 203040507 : 全幅のハッシュ値から
[2, 3](0, 0, 0) → 203000000 : 203 * 10^6 を引くことで
      [4, 5, 7] →     40507 : 取得できます。
~~~
ハッシュ値の定義はいろいろあるようですが、ここでは
```math:
hash(S) = \sum_{i=0}^{N-1}S[i] × base^{N-i-1} (\% mod)
```
の計算式を用います(%は余りを表します)。

さて、ローリングハッシュではハッシュの衝突が起きえます。
これは「本来は別の数列なのに、同じ数列と誤評価してしまう」現象で、法が小さいほど起きやすいです。
たとえば、(base,mod = 10, 331)を考えます(331は素数です)。
この場合、長さ3・ハッシュ値100の数列を取ってみても[1, 0, 0] と [4, 3, 1] と [7, 6, 2]の3種類があり、これらは同じ数列だと誤判定されてしまいます。
対策としては法を64bit型素数にする・32bit型素数で2個以上の基数を使う などが挙げられます。32bit型素数1個だと不十分です(衝突確率の評価は[wikipediaの誕生日攻撃の記事](https://ja.wikipedia.org/wiki/%E8%AA%95%E7%94%9F%E6%97%A5%E6%94%BB%E6%92%83)にあります)。
なお基数は乱択が推奨されますが、AtCoderではhackがないので定数でもかまいません。(終)

</details>

__セグメント木__
ローリングハッシュの計算式を見ると、2つの数列を連結したときのハッシュ値は$O(1)$で計算できることが分かります。
~~~:hash値の計算法  base = 100
[2, 3]          : 203 *10^6 +
      [4, 5, 7] :     40507 =
[2, 3, 4, 5, 7] : 203040507
~~~
具体的には (現在のハッシュ値, 数列長) を保持していれば区間ハッシュ値が計算できます。
さらに隣接する数列同士の演算順序は問わず、かつ長さ0の数列を(0,0)と定義することで単位元として扱えることから、ローリングハッシュはセグメント木に乗せることができます。
実装上は「前から読んだ文字列」と「後ろから読んだ文字列」の2種類を別管理して、クエリごとにハッシュ値を取得すればよいです。
注意点として、法に32bit型素数(998244353, 2^31-1 など)を用いる場合、ハッシュの衝突率が高いので最低2周はローリングハッシュを行いましょう。
計算量はセグメント木の準備に$O(N)$、クエリは1周あたり$O(QlogN)$です。
提出例: [2^31-1 2周(1218ms)](https://atcoder.jp/contests/abc331/submissions/48176319), [2^61-1 1周(1168ms)](https://atcoder.jp/contests/abc331/submissions/48176664)

__平方分割__
一応平方分割でも解けますが、非想定解です。Pythonだと定数倍を相当詰める必要があります。
計算量は前準備に$O(N)$、クエリは1周あたり$O(Q√N)$です。[提出例(2547ms)](https://atcoder.jp/contests/abc331/submissions/48175374)

<details><summary>実装の諸注意</summary>

区間を$√N$個のバケットに分割し、各バケットごとに全幅のハッシュ値を計算しておきます。
1. 更新クエリでは、更新後に該当バケットの全幅ハッシュ値を再計算します。
計算量は一点更新が$O(1)$、ハッシュ値の再計算が$O(√N)$です。
2. 取得クエリでは、左端 + 全幅ハッシュ地帯 + 右端 と分割してハッシュの計算を行います。
計算量は左端・全幅地帯・右端どれも$O(√N)$です。

ただし、Pythonでは定数倍高速化が必須です。高速化案を示します。
1. メルセンヌ素数を用いた高速演算。
詳細はkeymoonさんの解説記事をご確認ください。
1. 回文判定幅の削減。
前から半分と後ろから半分が一致していれば回文なので、ハッシュ計算幅を半分に減らせます。
1. 端のバケットの計算高速化。
ほとんど全域にまたがる区間ハッシュ値は、愚直計算するよりも、全幅ハッシュ値から不要部分を削除する計算法のほうが高速です。
1. 基数の冪乗の前計算。
前処理$O(√N)$、冪乗計算ごと$O(1)$にできます。
1. `input = sys.stdin.readline`で入力高速化。
ただし、文字列の末尾に改行文字`\n`が入るので注意してください。

(追記)上記改善案をすべて行うとより高速になりました。 [提出例(2157ms)](https://atcoder.jp/contests/abc331/submissions/48185862)

その他、基数が小さいほど演算コストが下がる現象も確認しています。
基数が10^3以下だと特に高速で、10^9以上だと一気に遅くなりました。(終)


</details>


# G Collect Them All

[問題文はこちら](https://atcoder.jp/contests/abc331/tasks/abc331_g)

__解説__
コンプガチャの期待値を求める問題です。難しすぎます。
[公式解説](https://atcoder.jp/contests/abc331/editorial/7763)を読んで通しましたが、難しすぎて解説ができません。がんばってください。
注意点としてnumpyの高速フーリエ変換だと実行時間に間に合わないので、数論変換(NTT)のライブラリを自作して畳み込みを行ってください。
具体的な実装はizu_noriさんの解説が勉強になります。

https://qiita.com/izu_nori/items/1c5cdef0500ffa0276f5

提出例: [NTT(1700ms)](https://atcoder.jp/contests/abc331/submissions/48199434)

<details><summary>個人用 解法メモ</summary>

[kyopro_friendsさんの公式解説](https://atcoder.jp/contests/abc331/editorial/7763)の1番目、包除原理による解法を追います。
なお数学は高校数学レベルなので、稚拙な用語が頻出します。ご容赦ください。



__用語の再確認__
基本的に公式解説に準じますが、「ノートの数字」を「賞品」と読み替えます。
$f(S) = \sum_{i∈S} \frac{C_i}{N}$ ･･･ 賞品集合Sに含まれる賞品が当たる確率
$U = {1,2, ... , M}$ ･･･ 全賞品の集合(全体集合)

__使用する極限計算__
1. $\sum_{n=0}^∞ a^n = \frac{1}{1-a}$ $(0\leq a<1)$
証明: [無限等比級数の収束・発散(高校数学の美しい物語)](https://manabitimes.jp/math/1055)
$S_n = \sum_{i=0}^n a^i$として、$(a-1)S_n = a^{n+1} - 1$から導けます。
1. $\lim_{n\to \infty} na^n = 0$ $(0\leq a<1)$
証明: $\frac {1}{a} = 1+t$ $(t>0)$とすると、$(1+t)^x>1+xt+{}_xC_2 t^2$ から
$\frac{x}{(1+t)^x} < \frac{x}{1 + xt +{}_xC_2 t^2} < \frac{1}{t+\frac{x-1}{2}t} → 0(x→∞)  $ となります。
xが整数でない場合も、うまく整数に丸めることで証明できます。


__$p(n)$の計算__
高々n回で全商品が揃う確率は、$p(n) = \sum_{S⊂U} (-1)^{M-|S|}f(S)^n$です。

<details><summary>包除原理の説明</summary>

簡単のため、賞品が等確率で出る場合の場合の数を考えます。
3種類の賞品{a,b,c}を、4回以内に全種揃える場合の数を数えると
{a,b,c}から4個選び、一列に並べる場合の数 : $3^4$通り
\- {a,b}から4個選んだ場合(cがない場合)を除外 : 各$2^4$通り
\- {b,c}から4個選んだ場合を除外
\- {c,a}から4個
\+ {a}から4個選んだ場合が引きすぎなので戻す : 各$1^4$通り
\+ {b}から4個
\+ {c}から4個
\- {}から4個 (0通り) = 36通り
となります。賞品ごとの確率の偏りがある場合も、場合の数に確率傾斜をつければよいです。
これを一般化すると上記の式になります。(終)

</details>

さて、ちょうどn回目にすべての賞品が揃った確率を$g(n)$とすると、
$g(n) = p(n) - p(n-1)$ ･･･ n-1回目には揃っていないが、n回目には揃っている確率
と定義できます。
さらに、n回以内にすべての賞品が揃う期待値を$h(n)$とすると、
$h(n) = \sum_{x=1}^n x × g(n)$ ･･･ x回目にはじめて揃った確率と回数の積和
となります。本問で求める期待値は$\lim_{n\to \infty} h(n)$です。

$h(n)$の式を整理します。
$h(n)$ 
$= 1×(p(1)-p(0)) + 2×(p(2)-p(1)) + ･･･ + n×(p(n)-p(n-1))$
$= -p(0) - p(1) - ･･･ - p(n-1) + np(n)$
$= np(n) - \sum_{x=0}^{n-1}p(x)$ ･･･(1)

ここで、 $p(n) = \sum_{S⊂U} (-1)^{M-|S|}f(S)^n$ のうち $S=U$を外に出すと
$p(n) = f(U) + \sum_{S\subsetneqq U} (-1)^{M-|S|}f(S)^n $
となります。特にΣの内部をXとして変形すると、$f(U)=1$から
$np(n) = nf(U) + nX = n + nX$ ･･･(2)
$nX = \sum_{S\subsetneqq U} (-1)^{M-|S|} × n×f(S)^n$ ･･･(3)
となります。ここで、$S\subsetneqq U$に対して$0 \leq f(S)<1$であり、
$\lim_{n\to \infty} na^n = 0$ $(0\leq a<1)$ ･･･(4)
から、$n→∞$のとき(3)は0に収束するため、最終的には無視できます。

(1)に(2)を代入して変形を続けると
$h(n) = np(n) - \sum_{x=0}^{n-1}p(x)$ ･･･(1)
$= -nX + n - \sum_{x=0}^{n-1}p(x)$
$= -nX + \sum_{x=0}^{n-1}(1-p(x))$ ･･･ 公式解説 等号1個目
$= -nX + \sum_{x=0}^{n-1} \sum_{S\subsetneqq U}(-1)^{M-1-|S|}f(S)^x$
$= -nX + \sum_{S\subsetneqq U}(-1)^{M-1-|S|}\sum_{x=0}^{n-1} f(S)^x$ ･･･ 部分和なので変形可能

求める期待値は$h(n)$の極限なので、(4)と$\sum_{n=0}^∞ a^n$の極限公式を用いて
$\lim_{n\to \infty} h(n) =\sum_{S\subsetneqq U}(-1)^{M-1-|S|}\frac{1}{1-f(S)}$ ･･･ (5)
と変形できました。さらにf(S)を$\frac{1}{N}$ごとの値で刻むと、
```math:
\sum_{k=0}^{N-1}\frac{N}{N-k} \sum_{S|f(S)=\frac{k}{N}} (-1)^{M-1-|S|}
```
とできます。Σの内部はDPで計算できます。
`DP[m][k]: m番目の賞品までの集合のうち、f(S)=k/Nを満たすSに対するΣの値`
と定義すると、これは$O(NM)$で計算できます。
ここまででようやく指数時間から脱却できました。

ここからさらに計算量を落とす方法が形式的べき級数･･･でしょうか？
DP[m][k] = $[X^k] (-1)^{M-1} \prod_{i=1}^{m} (1-X^{C_i})$
で確かに計算できるのですが、もはや着想のきっかけすらつかめません。

ともあれ、畳み込みで計算を高速化できるようなのでやります。
高速フーリエ変換を行うと1回あたり$O(NlogN)$で畳み込みが可能です。
さらにマージテクの要領で、数列長が短いもの同士で優先的に畳み込みを行うことで計算量を落とせるようです。
これで$O(N(logN)^2)$に落ちたので、畳み込みが高速であればACできます。(終)

<details><summary>雰囲気で計算量評価</summary>

マージテク併用の畳み込みで、$O(N(logN)^2)$か$O(NlogNlogM)$になるようです。
しかしこの部分の解説資料が存在しなかったので、個人的に納得できる範囲で解釈を行いました。
証明ではないのでご理解ください。

部分問題として、以下の命題を考えます。
> 長さN, 総和Mの自然数数列A $[A_1, A_2, ･･･ , A_N]$がある。
> 現在の最小要素a,bを取り出し、a+bを戻す操作を繰り返す。コストが$(a+b)log(a+b)$かかる。
> 最終的に数列長が1となったときのコスト総和は、$O(NlogNlogM)$である。

$aloga + blogb \leq (a+b)log(a+b) \leq (a+b)logN$から、マージ結果のa+b部分の総和が$O(NlogM)$であればよいです。
ここでマージ過程を二分木に見立てます。黄色の部分がa+bの部分、すなわちコストが発生する部分です。

![ABC331G 二分木参考図.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3624214/e9fa4d51-e377-0f1f-00b0-676fa9b0660a.png)

黄色の部分がすべて図の太線部分(深さが$logM$になる部分)以上に収まるように移動できればよく、直感的には収まりそうです。
「子の値は、親の値未満・親のペアの値以下である」ことを利用して、子を(必要なら部分木ごと)親方向か、親のペア方向に移動させるスライドパズルを行えば収まる気がします。

気がするだけで間違っていた場合、ご容赦ください。(終)

</details>


</details>








# おわりに
おわりです。
いろいろなときかたがあるとおもしろいなとおもいました。