やおよろず

ゆるふわな内容を書いても良い。と言うか何でも良い。そんな自由度の高いブログをめざします。

Educational Codeforces Roundをそんなにやっていけなかった話

この記事は カレーのち ぴょこりんクラスタ Advent Calendar 2018 - Adventar のために書いたものです。

本当は今日までにスマブラでＶＩＰ部屋に行きたかったんですが無理でした。精進します。

はじめに

ぴょこりんクラスタアドベントカレンダーってあるじゃないですか。
今もまさにそれ用の記事を書いているわけですけれども。
これまで自分は参加していなかったし、去年や一昨年は、２５日分の記事を３人で埋めていたと言うから恐ろしいですよね。

自分は今年になってはじめての参戦となりますが、
一昨日の記事でも触れた夏の大合宿(仮)に参加したあたりで、アドベントカレンダーの記事を書くことが確定していました。

２５日を４人で割ってもおよそ６日*1。
そんなに沢山ネタを作るのは大変だなぁ、と思ってましたが、
「途中経過でもええんやで」という誰かしらのアドバイスを受けて思いついたのが続き物ネタ。
そのための布石として実行したのが、12/5に投稿した、Educational Codeforces Roundをやっていきたい話 - やおよろずでした。

本題

と言うわけで今回は、実際に12/5の記事を執筆して以降、 自分がどのくらいEducational Codeforces Roundをやっていけたかの話をしたいと思います。

成果

以下が成果となります。ご査収ください。

f:id:yaoshimax:20181222044132p:plain
10問解きました

12/5時点*2ではRound1 A,B*3 しか解いていませんでしたが、
その後にRound1のC、Round2のABCD、Round3のABCDEを解き、いい感じに階段状の図を作る事に成功しました*4。
「なんだ１０問か…」と思われるかもしれませんが、submission数で言えば29回。
ならせば１ヶ月毎日１submissionしたと言えるので、かなり頑張った方だと思います*5。

以下では、Educational Codeforces Roundのいくつかの問題について、解説とか苦労話などを書いておきます。

８回誤答したRound 2 D「Area of Two Circles' Intersection」

問題ページはこちら

問題概要

２つの円があります。円の中心座標と半径はx_1, y_1, r_1と x_2, y_2, r_2として入力で与えられます。
この円の共有部分の面積を求めて出力してください。

解き方とか

数学ゲーですね。

  • 円が交わらないときは当然０
  • 片方の円が完全にもう片方を包含する時は、小さい円の面積
  • 上記以外の時、交点１つと中心座標(x1,y1), (x2,y2)をがなす角度 \theta_1, \theta_2 を用いて、 {r_1}^2 (\theta_1 - \sin \theta_1 \cos \theta_1) + {r_2}^2(\theta_2- \sin \theta_2 \cos \theta_2) が解

で求まります。 \theta_1 とかは余弦定理などを使えば計算できます*6。

が、誤差が1e-6以内である必要があるため、適当な計算をして居ると誤差で不正解扱いになります。

f:id:yaoshimax:20181222051309p:plain
誤差と格闘している様子

  • 多倍長整数とか勝手にやってくれるしpythonの方が良いかな？と思ったが、pythonだと誤差が解決できなかった
    • 多分pythonのfloatが53bitなのが悪い*7。この問題、c++でlong doubleを使う必要があるらしい*8。
  • {r_1}^2 \theta_1 + {r_2}^2\theta_2 - {r_1}^2 \sin \theta_1 \cos \theta_1 - {r_2}^2 \sin \theta_2 \cos \theta_2と計算すると r_1, r_2の差が大きい時に誤差で解がずれる

あたりが特にはまりどころでした。
ネット上から他の方の解答を引っ張ってきて途中計算dumpしながらデバッグするくらいにはしんどかったです。

個人的に教育的っぽかったRound 3 D「Gadgets for dollars and pounds」

問題ページはこちら

問題概要（意訳）

m個の商品がある。これらの価格はドルまたはポンドで与えられている。
あなたはs円までの予算で、m個のうちから k個の商品を購入したい。
ただし、あなたが所持できるのは日本円だけ。
商品購入時に、必要額のドル・ポンドを円と交換する必要がある（予め円をドルに替えて翌日以降持ち越し、とかは不可）。
ちなみにあなたは、ここから先、 i 日後には1ドル a_i円、１ポンド b_i円で変換できる事は知っている（ 0 \le i \le n-1）。
最速で何日までに k個の商品を購入する事ができるのかを調べ、何日にどの商品を購入すべきかと合わせて出力せよ

解き方とか
動的計画法で頑張る方針（これだと解けない）

１日ずつ変換レートを見ながら商品の買う・買わないを判断しつつ、
i 日目までに商品集合 M を買う値段の最小値を dp[ i ][ M ]円として保持する事を考えます
コードに起こそうとすると以下のような感じでしょうか。

for day in [0..n-1] do
  for all M, which is subset of [1, ..., m] do
    for all U, which is subset of M do
       cost = i日目に商品集合Ｕを買うコスト(円);
       dp[i][M] = min(dp[i][M], dp[i-1][M\U] + cost)
    end
    if Mのサイズがk以上でdp[i][M]がs円以下 do
       //i日目までに商品はk個以上買える。
       print Mの買い方
       return
    end
  end 
end

この場合、 m個の商品のうちどれを買ったか・買わないかの情報（商品集合 M ）を保持しなければなりません。 この状態数は 2mとなりますが、 m は最大で105なので、そんな情報持てません。

ちゃんとした解き方

と言うわけでどうするかと言うと、注目すべき性質としては以下があります。

  • ドル額で d 円のとある商品を x 日目までに買おうとするなら、最適な購入費用は d \times \min_{i\le x} a_i である。
    • 要するに、 x 日目までの一番ドルが安い日に買っておけば良い
  • ポンド額で p 円のとある商品を x 日目までに買おうとするなら、最適な購入費用は ｐ \times \min_{i\le x} b_i である。
  • 「 x 日目までに s 円の予算で商品が k 個以上買える」が真ならば、任意の x 以上整数 yにおいて「 y までに s 円の予算で商品が k 個以上買える」 も真

以上の条件が成り立つ事により、以下のような二分探索ができます。

left = -1
right = n
while(left+1 < n) do
 mid = (left+right)/2
 price = []
 for all 商品 do
  price.append(商品を mid日目までに購入するとした場合の最安値)
 end
 sort(price)
 sum = 0  // mid日目までにk個の商品を買うとしたら最安値でいくらになるか
 for i in [0,...k-1] do
   sum += price[i]
 end
 if( sum <=  s ) do
  right = mid
 else
  left = mid
 end
end

// right日が答えなので、その日迄に買うk個の商品の購入の仕方を出力する

この手によって、計算量を O (m* \log(m) \log (n)) におとす事ができるので、 m の最大値が 105とかでも許されます。

実装がダルかった、Round 3 Ｅ「Minimum spanning tree for each edge」

問題ページはこちら

問題概要

n 頂点 m 枝の重み付き無向グラフが与えられます。
このグラフは連結であること、多重辺やセルフループがないことは保証されています。
枝 e_0, e_1, ... , e_{m-1} それぞれについて、 c_i = \min (e_i を含むグラフの全域木のコスト) を求めてください。

解き方
基本方針

最小全域木を計算する代表的なアルゴリズムに、プリムのアルゴリズムがあります*9。 プリム法 - Wikipedia

平たく言ってしまえば、重みの軽い枝は、ループを作らない限り、片っ端から最小全域木を構成するものと考えて良い。
ってアルゴリズムです。

f:id:yaoshimax:20181222130705p:plain
プリムのアルゴリズムを説明しようとした雑な図

今回は「ある枝 e をつかった全域木の最小コスト」を枝ごとに求めなければなりませんが、
そう言う全域木は、
枝 e の端点を u,v とすると「グラフの最小全域木上で u から vまでのパスを考えたときに、一番コストの高い枝を、枝 eと置き換える」
事によって作る事ができます。

f:id:yaoshimax:20181222135138p:plain
雑に説明を試みる図

全域木上の２頂点間のパス中の最大辺を求める方法

というわけで、最小全域上の２頂点間のパスの最大辺を求める必要が出るんですが。
頂点数が最大で 2* 105 なので、毎回愚直に深さ優先探索とかして求めようとすると時間がかかるし、
全頂点対間のデータをテーブルに保持するのは頂点対が2 * 1010くらいあるので無理です。

これは、木と２頂点 u, v が与えられた時に、一番近しい共通親（Lowest common ancestor）を求めるアルゴリズムを応用することで効率よく計算できます。

まずLowest common ancestorを効率良く計算するにはどうするかと言うと、 前処理で各頂点 n について１番目の親、２番目の親、４番目の親、…2k番目の親をparents[n][k]に保持しておいて*10、 以下のような擬似コードにより、不必要な所をすっ飛ばしながら親を見つける、というものです。

function LCA(int u, int v){
  if(depth[v] < depth[u]) swap(u,v);
  // ここでdepth[v]がdepth[u]に等しくなるまでvの親をたどる（省略）

  while(u!=v && parents[u][0]!=parents[v][0]){
    // 親が一致しないギリギリである2^k番目の祖先までu,vを更新する
    int k =-1;
    while(parents[u][k+1]!=parents[u][k+1]) k++;
    u = parents[u][k];
    v = parents[v][k];
  }
  if(u!=v) return parents[u][0];
  else return u;
}

これがどうして全域木上の２頂点パス中の最大辺を求める事に繋がるか？というと、 各頂点 nについて、１番目の親、２番目の親、…2k番目の親をたどる迄の最大辺をmaxEdge[n][k]と保存しておけば 以下のような擬似コードが使えるのです。

function maxEdgeBetween(int u, int v){
  if(depth[v] < depth[u]) swap(u,v);
  ans = -1;
  // ここでdepth[v]がdepth[u]に等しくなるまでvの親をたどる（省略）

  while(u!=v && parents[u][0]!=parents[v][0]){
    // 親が一致しないギリギリである2^k番目の祖先までu,vを更新する
    int k =-1;
    while(parents[u][k+1]!=parents[u][k+1]) k++;
    // 更新する際に最大辺を更新する
    ans = max(ans, maxEdge[u][k]);
    ans = max(ans, maxEdge[v][k]);
    u = parents[u][k];
    v = parents[v][k];
  }
  if(u!=v){
    ans = max(ans, maxEdge[u][0]);
    ans = max(ans, maxEdge[v][0]);
  }
  return ans
}
解き方まとめ

と言うわけでこの問題については、

  • グラフの最小全域木とそのコスト S 求める
  • 最小全域木をたどりながらparentsとmaxEdgeを作る
  • 各枝 e = (u,v)について、S+weight(e) + maxEdgeBetween(u,v)を出力する

みたいな事をすれば、 O(m \log(n)) の計算量で解く事ができます。

プリムのアルゴリズムの実装に加えてLowest common ancestorを解く実装もする必要があるので糞だるかったです。

ＥＦ解けなかったけどＧが解けた。Round 55 G「Petya and Graph」

問題ページはこちら

問題概要

n 頂点 m 枝の無向グラフが与えられます。
各枝には w_i、各頂点には a_iの重みが与えられています。

頂点の部分集合Uにより定義できる誘導部分グラフについて、
得点を「Ｕに含まれる枝の重みの総和ーＵに含まれる頂点の重みの総和」で定義します。
この最大値を求めてください。

解き方

何となく燃やす埋める問題っぽさがあるなぁと思いながら考えてたらグラフの最小カットに帰着できました。

始点 s,終点 t、 元のグラフの枝 e_1,e_2,e_3... に対応する頂点集合 Ｗ、
元のグラフの頂点 v_1,v_2,v_3... に対応する頂点集合 Ｖ
を持つようなグラフを考えて枝を以下のように張ります。

  • 枝 e_i の端点 u,vから枝に対して、重み無限大の辺
  • 始点 sから各頂点に対して、 a_iの重みの辺
  • 各枝から終点 tに対して、 w_iの重みの辺

f:id:yaoshimax:20181222144946p:plain
雑な説明（作成都合で始点sが右側なのは許して）

こうすると求める値は、このグラフの最小カットＣを用いて、 \sum w_i - C と表せます。

なんで？を雰囲気で説明すると、
始点～頂点の位置で辺をカットした頂点の集合を V_1、辺をカットしなかった頂点の集合を V_2とすると、
枝～終点の位置では、 V_2に属する点を端点に持つ全ての枝～終点の経路をカットしなければいけないですよね。
これは、 V_1= Uとみなした際に、

  • 前者のカットによる値が、「Ｕに含まれる頂点の重みの総和が得点から減点される分」
  • 後者のカットによる値が、「Ｕに含まれない頂点を端点にもつ枝は、そもそもの得点を得られない分」

にそれぞれ対応していると捉える事ができます。
よって、この減点分を最小にするのが得点を最大にする事と等価といえるのです。

ここまで考察できてしまえば、このグラフの最小カットは sから tまでの最大フローに等しい*11 ので、
後は最大フローのライブラリで殴るだけですね。

終わりに

競プロ勢じゃない方にこの記事ってどのくらい伝わるんですかね…分かりづらかったらすみません（ぷよぷよぶり２回目）
ともあれ、それなりに歯ごたえのある問題も解けたし、
このあたりをサクサク解けるようなライブラリ化したいなーと言うモチベにもなったので、個人的には良い機会でした。
これからも競プロやっていきます*12。

*1:合宿当時まだ他の方の参戦は決まって居なかったので、４人で計算してます

*2:当時の記事のdemo画面などでご確認ください

*3:本質でないので脚注に書きますが、Round55のABCDGも解いてあります

*4:別に階段が好きな訳ではないです

*5:詭弁です。本当はもっと頑張りたかったと思っています

*6:実際は咄嗟に余弦定理を思い出すことができず、数式ごにょごにょして余弦定理求めるところから始めた

*7:テキトーに調べて53bitってなんで？と思いながらこう書きましたが、倍精度の仮数部が52bitなので妥当な値な気がしてきた

*8:pythonでもnumpy使えばlong doubleできそうですが競プロで使えるか分からないので未確認

*9:あとクラスカルもありますが。今回は使いません

*10:2011-12-05 - (iwi) ｛ 反省します - TopCoder部 の Doublingの図とか見るとわかりやすいかも？

*11:最大フロー最小カット定理 - Wikipedia

*12:でもそろそろ競プロ以外の趣味コーディングもしたいんですよね。UnityかReactあたり使うヤツももうちょっとやりたい