AOJ 0581 お土産購入計画 (Gifts)

2012年の情報オリンピック予選の6問目。当時解けなかったのは当然としても、今見ても結構つらい問題で、5年越しにようやくACしたのでなんとなくブログに書いておく。

問題文

まあ日本語なので読みましょう

お土産購入計画 | Aizu Online Judge

解法

流石にbitDPはしたくなって当然。

単純に考えると

dp[y][x][k][hash] := 「地点(x, y)におり後k回遡れる状態で、マンハッタン距離が3以下であるような周りの点の内辿ったことのある点の集合がhashで表される時の、それまでに買えたお土産の最大個数」

というdpを立てたくなる。

本質はその計算量で、この状態では、S = {マンハッタン距離が3以下であるような周りの点}としてO(HWK * 2^|S|)であり、|S|=24であることを考えると、概算としてこれは50*50*3 * 2^24 = 125829120000であって到底間に合わない。

ただこれには相当な無駄があって、例えば以下のような状態というのは考えなくて良い。

　　　　　　　　　　　

ここで、黒は今自分がいる点であり、オレンジは既に辿った点を表している。

このような状態になるためには、下と右のオレンジのラインの両方を辿り、黒に戻ってくるという遷移が必要で、これはKが6以上じゃないと不可能。

もっと変なケースを考えると一瞬怖い気持ちになるが、そもそもhashで持つのは既に辿った点の集合であるということに注意すると、更に考えなくていい状態が多いだろうと考えられる。

逆に考えなくてはならない状態の内、最もオレンジの個数が多いものは、こんな感じになる。

　　　　　　　　　　

これは例えば右下下左上上と動くことによって達成できる状態で、この時点でk=0となる。

本当に最大か？？もうちょっと増えない？？という気持ちに一瞬なるが、まず、遡る動きをすると、kの数が減り、記憶しておかなければならないオレンジの位置が減ってしまうので、出来るだけオレンジは右下に来るようにしてよい。で、この状態になるまでの動きには無駄がない(黒に戻ってくる動き以外の全てで、オレンジの数が増えている)。そして、仮にここから下や右に移動したとすると、オレンジの丸が増える時には、少なくとも1つ以上の元々あるオレンジがフェードアウトするので、これ以上増えることはない(はず(本当か？(一応__builtin_popcountでhash計測したら5が最大だったょ)))。

したがって、大きく見積もっても、2^24通りの内、考え無くてはならない状態の上界は

Σ[i=0..5]24Ci = 55454であって、55454 * 50 ** 2 * 3 * 4 = 1663620000 = 1.7*10^9な時点で、JOI予選本番中であれば、普通にGOサインを出す判断ができる。

また、この上界はかなり余分めに見積もっているのは明らかで、本当に発生しうる状態のみを計算するように真面目に書けば、普通に少なくとも10分の1ぐらいにはなりそうに思える(そもそもオレンジが5個あるケースというのは、さっきの2x3の長方形しかないはずで、24C5どころか8通りしかないように見える)。実際、適当に書いても1秒以内でAC出来た。

http://judge.u-aizu.ac.jp/onlinejudge/review.jsp?rid=2632622

ちなみにJOIの公式の解説には「考えなければならない近傍は高々12個で状態は2^12なので余裕です」と書いてあるが、正直これはまだ理解できていない。仮に考えなければならない近傍が12個だったとして、その近傍の位置が固定されていればそれは当然2^12だが、実際には固定されていないはずでコンビネーションが出てくるべきなわけで、2^12という見積もり方は嘘なんじゃないかと思うのだけど、誰か分かる人おしえてほしい。

符号でバグってそうだったら-Wsign-conversionも使おうねという話

TL;DR的にはタイトルまんまです。

CodeForcesでまあこんな感じで激冷えしてしまったんですが、

原因の一端に、自分がコンパイラオプションを正しく認識していなかった事もあった上、ほとんど日本語でも英語でも、似たような現象に悩まされてそうな人がいなかったのでブロガーしちゃう。

さて、本質的に何が原因でバグっていたかというと、

vector<int> a(3);
int i = -1;
a[(a.size() + i%a.size())%a.size()] = 1;

こういうコードがあったとして、これがa[2]=1になることを期待してしまっていたことに起因している。

バグった時点で、薄々signedからunsignedとかの暗黙な変換が起きてるんだろうなーとは思っていたので、-Wconversionオプションをつけてコンパイルとかしてみたものの、一切なんのwarningも出なかった。そのせいで、思い過ごしかと思って全然違うところの修正をしてしまったりした挙句9ペナぐらい食らってしまった。もちろん上のコードに対しても、少なくとも僕の環境(Ubuntu16.04デフォルトのg++ (GCC) 6.3.0)では-Wconversionオプションではwaningが出ない。

どうすればいいのかなあと思って検索してみた所、-Wsign-conversionなるものがあるらしいということに気がつき、試してみた所ちゃんと検出されてwarningが出た。

まあどうせ100行ぐらいのコードだから暗黙な変換が起こりそうな場所ぐらい覚えておけよという話なんだけど、変にコンパイラオプションに頼ろうとして痛い目を見た経験は中々ないし覚えておかないとまたやってしまいそう。

ちなみに、clangだと-Wconversionつければ勝手に-Wsign-conversionも有効になるっぽい。ただ-Wsign-conversionは-Wsign-conversionで、競プロだとsignedな値で配列参照することがあるだろうから、vector::atに対してsignedな添え字で参照することに対して毎回怒られるのはちょっと面倒かもなあという気もする。まあバグってそうだったらつけて見て、配列参照の部分は無視すればそこまで大変でもないかな？

MP法とKMP法の違い

大学生なのでデータ構造とアルゴリズムという授業を取っているんですが、今まで自分がKMP法だと思い込んでいたものがMP法だった上、ネット上の記事が無限にバラバラなのでメモしておく。

MP法

MP法はMorris-Pratt algorithmと呼ばれるもので、Knuthが登場しない。

MP法では、パターン文字列に対してborderと呼ばれる値を計算することにより、シフト量を定める。

ある文字列Sのborderとは、Sの部分文字列であって、Sの接頭辞かつ接尾辞であるような文字列Wのことである。

MP法では、i=0..|S|に対して、S[0:i-1]の最長border(の長さ)をO(|S|)で求め、それを元にシフト量を決める。具体的な最長borderの求め方については、snukeのこの記事が一番分かりやすいのでこれを見るといいと思う。

snuke.hatenablog.com

さて最長borderが求まれば、1-originとして、j文字目で不一致した場合のシフト量はj-border[j]とすればよい(配列インデックスと何文字目を表すかのoriginの違いによって、±1となる可能性があるので気をつけること)。

結局、MP法とはどちらかというと、borderを計算することの出来る線形アルゴリズムという側面が強く、別に文字列検索にだけ使うようなものではない。

例えば検索文字列自体を検索する用途以外でMP法を使う例としては以下の問題などがある。

codeforces.com

KMP法

KMP法はMP法にKnuthがくっつくやつで(は？)、上記のborderによるシフト量の決定には、まだ定数の改善余地があるため、それを改善する。

ある文字列Sの、接頭辞Pを考える。

すなわち、S = P + Xである(ここで+は文字列のconcatenate)。

ただし、 PやXが空文字列であるとちょっと困るので、P' = P + '_', X' = X + '_'というように、末尾に番兵があることとする。

Pのtagged border(あるいはstrict border、strong borderともいう) Wとは、Pのborder Wであって、かつP[|W|] != X[0]であるものを指す。

MP法におけるシフトの何に無駄があるかというと、S[i:i+j]が、pat[:j]までのj文字が一致していて、S[i+j]とpat[j]で初めて不一致が起きた時、当然シフトが発生して、次は、pat[border[j]:]とS[i+j:]が一致しているかを比較し始めるわけだが、仮にpat[j] == pat[border[j]]であったとすると、このシフトはまだ不十分であるということになる。

なぜなら、今、不一致が起きたからS[i+j] != pat[j]であって、かつpat[j] == pat[border[j]]であるから、明らかにpat[border[j]] != S[i+j]となって、pat[border[j]:]とS[i+j:]は先頭からして一致しない。

よって、borderの中でも、P[|W|] == X[0]を満たすようなborderは文字列検索においては考えなくて良いことが分かり、そうでないようなborderの中で最長のものをシフト量の決定時に使えば、少し検索の定数が減ることが分かる。

実装はMP法からほとんど変わらないため、コストのかかる改善ではない。

例えば上記のsnukeのコードを、以下のように変更するだけでいい。

A[0] = -1;
int j = -1;
for (int i = 0; i < S.size(); i++) {
  while (j >= 0 && S[i] != S[j]) j = A[j];

  j++;
  // 元のコードでは A[i+1] = j;
  if (S[i+1] == S[j]) A[i+1] = A[j]; 
  else A[i+1] = j;
}

elseの部分は、元のコードと変わっていないので当たり前に受け取れる。

問題はif文の後のA[i+1] = A[j]という処理で、何故これでいいのかというと、

まず、帰納的に考えれば、j < i+1より、A[j]は既に正しく計算された値である。

したがって、S[A[j]] != S[j]が成立しており、S[i+1] == S[j]ならば、明らかにS[i+1] != S[A[j]]が成り立ち、MP法と同じ議論からこれは最長tagged borderでもある。

初め僕は「なんでwhileじゃなくてifなんだ？？」と思ってたけどまあこういう理由から、「A[i+1] = j; while (A[i+1] > -1 && S[i+1] == S[A[i+1]]) A[i+1] = A[A[i+1]];」とかしなくても良いということが分かる。

終わり。

Codeforces Round #415 (Div. 1) C. Find a car

dpの遷移が複雑すぎるだろ。

問題概要

縦10⁹, 横10⁹の合計10¹⁸個のマス目があり、上からi番目, 左からj番目のマス目を(i, j)と表す。マス目に以下のルールで数を書き込んでいく:

• マス(i, j)には、マス{(x, j) | 1 ≦ x < i}とマス{(i, y) | 1 ≦ y < j}に書かれている数全てを除いた上で、最小の自然数(0は含まれない)を書き込む。
• つまり、マス(i, j)には、そのマスの左と上を見て、書かれていない数の内、最小のものを書く。
• 当然(1, 1)は1とする

具体的には、以下のような形になる(問題文から引用)

以下の質問に対して Q回答えよ:

• x1, y1, x2, y2, kが与えられるので、長方形区間{(x, y) | x1 ≦ x ≦ x2 ∧ y1 ≦ y ≦ y2}について、k以下の数が書かれているマスについて、それらのマスに書かれている数の総和を答える。
• ただし総和はあまりにも大きな数になるので、10⁹+7で割った余りを答えれば良い。

制約

• 1 ≦ Q ≦ 10⁴
• 1 ≦ x1 ≦ x2 ≦ 10⁹
• 1 ≦ y1 ≦ y2 ≦ 10⁹
• 1 ≦ k ≦ 10⁹

解法

まず求めたい値をcalc(x1, y1, x2, y2)とおくと、calc(x1, y1, x2, y2) = calc(1, 1, x2, y2) - calc(1, 1, x1-1, y2) - calc(1, 1, x2, y1-1) + calc(1, 1, x1-1, y1-1)に分解できるので、全て(1, 1)から考えて良いことがわかる。 ここから少し考察すると、(1, 1)を含む任意の長方形状のマス目は簡単な操作で2倍に拡大できる事に気づく。

つまり、

こういうことで、 例えば

この性質を用いれば、まあ倍倍にしていけば、ぐらいで計算できそうだな、というのは分かる。ここからが問題。 基本的にdpをするんだけど、dpの遷移が大変なので本当に辛かった。 やはり遷移を直で200行ぐらい書くよりは、悩んでもいいからループでちゃんと遷移をまとめる方が良いなという気持ちになった。 具体的な方針としては、1x1の1個のマス目から初めて、目的の長方形へ伸ばしていって、その際に必要な値をdpとして伝搬させて行く。この際、最終的な最大値kから逆算して、各過程に置いて許される最大値というものを前計算しておく。 つまり例えばk = 7だったとすれば、1つ前の長方形ではk' = 4である必要がある、これは長方形を2倍した時、マス目の値も2倍か2倍-1されるから。 まずdpとして持つ状態は、

• どの部分のマス目を考えているか - 右端、下端、右下のマス、それ以外の4通りが必要になる
• 最大値かどうか - 桁dp的な考え方をする必要がある。つまり、上記のように、マス目のに書かれている数が各過程において可能な最大値である場合、場合によっては2倍すると不適となる場合がある。
• 例えば上記の例でマスの数が4であった場合、2倍すると8になるが、許される最大値が7である時、これは不適となる。

また、遷移の式は、dp[i+1][nj][nk] += dp[i][j][k]*2 - num[i][j][k]のようになるわけなので、実はマス目の総和だけでなくて、条件に適するマス目の個数というのも持っておく必要がある。ちなみに僕は誤読していて個数を求める問題だと思っていて200行直書き遷移を書き直す必要が出てきて詰んだ。

書いたコード

まあ本当に遷移は頑張るという感じで、こういう問題は嫌いですが、本番中に1発で通せるべきだよなあと思い反省の意を込めてこの記事を書いておきます

Codeforces Round #415 (Div. 1) C. Find a car

みんなのプロコン 2017 予選 D

みんなになれなかった。

問題概要

読んで: D: 工場 - 「みんなのプロコン」 | AtCoder

解法

正攻法はそのままこのクエリを実現しようとすれば結合則成り立つので値を2,3個持ってくださいというやつみたいです(たしかにね。)。

クエリの時系列と日数の時系列の2次元なので、こういう時はとりあえず日数の時系列でクエリをソートして、segment treeのインデックスをクエリのインデックスにすることを考える。

D = D1, D2, ... , Dqを順に代入しつつ以下の値を考える(ただしD1 ~ Dqは昇順にソートしたクエリの日付)。

num[qidx] := qidx番目までのクエリを考慮した時、D日目までに販売した賞品の数

これをセグ木によって持つことを考える。

Dqidx = Dとなるようなqidx番目のクエリを処理することを考える。

経理質問(回答クエリ)については、num[qidx]を答えられれば良いので、この値をセグ木で実現できるなら特に考える必要がない。

注文追加(更新クエリ)について考える。

qidx番目のクエリによる更新が反映されるべきなのは、num[qidx] ~ num[Q]であって、それぞれのi ∈ [qidx, Q]について、num[i] = max(num[i]+Aqidx, K*D)という更新を行う。

K*DとAqidxは定数でnum[i]はiについて流石に単調増加なので、num[t] > K*D-Aqidxなる最小のtを二分探索で見つければ、t未満のインデックスについては区間に対してAを加算、t以上のインデックスについては区間に対してK*Dを代入をすればよく、これは、seg木にnum[i]の代わりにnum[i]-num[i-1]をもたせれば、Lazy Propagationによる区間代入をサポートした、区間和を求めるseg木によって実現できるので解ける。めちゃくちゃバグる。

Codeforces Round #176 (Div. 1) D. Tourists

JOIerお家芸感。なんかこれはJOIerは結構見た瞬間簡単だなあと思う気がする。おそらくDではない。

原点にある点がn個あり、i個目の点は時刻q[i]に原点を出発し、x軸の正方向に沿って速さ1で進む。また、それとは別にm個の線分が存在して、i番目の線分は時刻t[i]に区間[l[i], r[i]]を覆うように一瞬の内に出現する。各点について、点が線分に覆われながら移動している時間の長さを求めよ。

なんか本当に簡単で、div1Dにしては簡単だと思うんだけど、1つの線分(l, r, t)に注目した時、q[i] >= t-lであるような点iはその線分によって、r-l秒覆われていることになる。そして、これはどう見てもqに関して単調性がある、すなわち、q[j] >= q[i] >= t-lであれば、点jもr-l秒覆われることになるから、点をq[i]で昇順にsortしておけば、二分探索により、r-l秒覆われる点が全て分かる上、imos法を用いれば、全ての線分について考えたときの具体的な秒数の合計も一瞬で分かる。ここで少し怖くなるのは、t-r+1 <= q[i] < t-lであるような点iについては、点が移動している途中で線分が出現するので、具体的には、q[i]+r-t秒だけ微妙に覆われることになる。しかしこれも同様に二分探索とimos法によって対応可能で、どうするかというと、二分探索でt-r+1 <= q[i] < t-lを満たすiの範囲を求めた上で、定数r-tを普通のimos法で計算(A[i]にi番目の点のこの定数に関する累積和が入るとする)、定数1をimos法で計算する(P[i]にこの累積和が入るとする)と、微妙に覆われた時間の累計はP[i]*q[i] + A[i]で求まる。1次関数の累積を求めているイメージ。基本的にはこれだけなんだけど、後は線分に重複があるところだけ解決しないといけない。これは、各座標において、その座標を覆う線分の内、tが最小であるものを求めればよい。これは結構明らかだと思う(tが最小であるような線分以外を無視すれば重複は解決され、かつ最も覆われる点が多い)。実際に各座標についてやると当然計算量がやばいが、これは簡単なイベント走査で解決できる。線分(l, r, t)についてイベント(l, t)とイベント(r, t)を作成し、イベントの集合を1つ目の要素に関してソート。multisetを1つ使って(l, t)が来たらtをinsert、(r, t)が来たらtをeraseしていけば、連続する2つのイベント(x1, t1), (x2, t2)について、区間[x1, x2]におけるtの最小値は、イベント(x1, t1)を処理し終わった後のmultisetの最小値となる。これも多分典型。なので各イベントの処理後に新しい線分の集合に、それらの値を追加してやれば、重複を取り除いた線分が全て求められて、解ける。このやり方ではどうやっても、重複を取り除いた線分は流石に自明にO(m)個しか存在しないので、かなり適当に追加しても問題ない。

ちなみに初め誤読して、「点は線分に覆われなくなったらその場で止まる。」という条件があると思っていたんだけど、それでも解けるので考えてみてください。

Codeforces Round #176 (Div. 1)

codeforces.com

はいサイアクの1完。

Bが分からない時は本当に辛いンゴねえ

A. nが与えられる。1~nの順列P[1]~P[n]が「∀i(1 <= i <= n) P[P[i]] = n+1-i」を満たす。このようなPを1つ出力せよ。そのようなPがなければ-1と出力すること。

順列なのでグラフではdirected "cycle"(cyclicではない) graphの集合と対応する。グラフにおいては、上の条件は、頂点iから有向辺を2本辿った先の頂点が常にn+1-iであることを表している。頂点n+1-iから2本辿った先が頂点iにならないといけないから、そのような条件を満たす順列は、大きさ1または4のdirected cycle graphのみで構成されている(ただしnが奇数でかつ頂点(n+1)/2についてのみ大きさ1は許される)。

よって与えられたnをまず4で割って、剰余が2, 3であればそのようなPはない。

後は適当にやればよくて、nが奇数であれば、ans[(n+1)/2] = (n+1)/2にだけ気をつけて、(k)→(k+1)→(n+1-k)→(n-k)→(k)みたいな閉路を考えて答えを構成してやればよい。

なんかめっちゃバグった。アホ。

B. 1~nの順列P[1]..P[n]と自然数kをとって、P[2]P[3]...P[k]P[1], P[k+2]P[K+3]...P[2k]P[k+1], ...を新しい順列として返す写像f(P, k)を定義する(つまり順列をk個ごとで区切って、区切った中で先頭の数を最後に移動させる)。L[1] = [1, 2, 3, ..., n]として、L[k] = f(L[k-1], k)と置く。L[n]を求めよ。

まず最後に移動させる操作の回数は、Σ[i=2..N]ceil(N/i) = O(NlogN)であるから、L[1][1] ~ L[n][n]のn^2個の数の内、いい感じに最後に移動させる数だけを処理することが出来れば、つまり1つの区間内のシフトをO(1)で出来ればO(NlogN)が達成できることはわかる。

ここからが問題で、僕はL[n]におけるi番目の数がL[1]においてどこにあったかというのを独立に考えようとしたのでドツボにハマった。つまり、L[1] ~ L[n]を正方形に並んだn^2個の数だというイメージで見た時、僕は縦に区切って処理しようとしてしまったんだけど、横に区切らないといけなかった。横に区切る、つまりL[k]からL[k+1]を一度に求めようという気持ちでいくとちゃんと出来て、これは、最後に移動させられる数を無視してみると、左に1個ずつズレていくだけだから、このズレを実際にはズラさずに固定してやると良い。どういうことかというと、配列L[k]から配列L[k+1]を作る時、本来ならL[k+1][i-1] = L[k][i]とやるわけだけど、定数1分ずつ全てがズレるのだから、i' = i-1と置いてやれば、L[k+1][i'] = L[k][i]と出来て、L[k+1]を新しく生成するのではなく、L[k]をinplaceに置き換えることにすると(dpのメモリ削減みたいなもん)、1ズレる数に関しては値を参照したりコピーする必要がなくなる。後は、最後に移動させる必要のある数だけ、L[k+1][i+k-1] = L[k][i]みたいな感じにする必要があるから、i'とinplaceを考えると、後ろからL[i'+k] = L[i]で更新していくだけでよく、上の評価式からこの更新の回数はO(NlogN)であるから解ける。

ちなみにL[n][i]がL[1]においてどこにあるかという考え方だと、O(NlogN(logN+loglogN))は達成できるんだけど、N=10^6でしかも旧版の問題に対する処置としてTLEが1/2の1秒だから普通に間に合わない。

自分用メモ: 流石にどっちで対処するかぐらいちょっと考えましょう、というか結構dpの変数置き換えはよくやるし好きなんだからこれぐらい思いつくべき