Mister雑記

競プロやります。

MisterのJAG夏合宿2019参加記

The atama 3人でJAG(ACM-ICPC OB/OGの会)が主催する夏合宿に参加してきました。

目次

Day 1

The atamaの集合は早い

初日は集合場所がオリンピックセンター(以降オリセン)だったので、万全を期してThe atamaは一旦合流してから行くことに。 結果全員が時間通りに集合して、11:00にオリセンに着いた。全員揃った時刻ではFAだった気がする。

Day1の貢献は少ない

その後人がぞろぞろと会場入りして、Day1のコンテストが始まった。 以降(謎の)守秘義務により、メンバー名を僕、S、Oで表すこととする。

コンテストリンク

  • いつも通りA問題をS、B問題をO、C問題を僕が読む
  • B問題をOが実装→AC
  • A問題をSが実装→AC
  • C問題の考察が終わる
    • セグ木が要るので慣れているOに実装を投げる→AC
  • D問題に考察を集中
    • サイクルに分解して中国剰余定理を使えば解けるね
    • 似た問題を見たことがあるということでSが実装→AC
  • Jが解かれているので考察する
    • 最大マッチングと連結成分じゃん
    • ライブラリ担当?のOが実装→AC
  • 実装中に僕がE問題を考察
    • 積を距離としたときBFSとかで矛盾を調べれば解ける
    • 素因数分解とか考えたが、計算量が削れない
  • SとOがH問題を考察
    • 「どの行にあるかは分かる」というところまでは来た
    • しかしそこから先が詰められず
  • 少し解かれているG問題とかを考えたが微塵も分からない
  • 再びE問題に考察を集中
    • とりあえず対数でOが実装→WA
    • S「ハッシュ的な方法でなんとかならない?」→MODを取るもWA
    • 非連結なケースで死ぬことに気づき修正→AC
      • ちなみにコンテスト後に対数のを修正して投げたが通らなかった
  • K問題が畳み込みで解けることに気づく
    • 僕「誰かFFT持ってますか?」S「いいえ」O「いいえ」 おしまい
    • 僕がうろ覚えで実装を試みるも、できるはずもなく終了

結果ABCDEJの15位(オンサイトは多分10位)。以下反省点。

  • 僕が実装をしなかった
  • 全員が全問題を読むべきだった
  • FFTくらい用意しておくべきだった
  • Eの考察に時間がかかりすぎた

振り返るとABCDEJを2:30くらいで片付けて、せめてK問題までは通したかったという印象。OがLi-Chao木を持っていたのでI問題もワンチャンあったかもしれない。

何にせよ実力不足、特に速度不足を痛感する結果となった。

314の夜は長い

コンテスト後は夕飯を除いてほとんど自由時間。談話室でボドゲに混ざらせていただいたりした。

なお僕が割り振られた314号室はTwitterでよく見る顔ぶれで、さらに何人か訪問者が来て結局2:00まで起きていた。ワイワイできてとても楽しかった。

Day2

夏合宿の朝食は遠い

7:00 - 9:00が朝食の時間なのだが、宿泊者が数百人いるのに食堂が1つしかないので非常に混む。というわけで早起きして6:55に食堂に着いたのだが、それでも結局10分待たされた。

ちなみにDay3はそれを踏まえて6:45に向かったところ、並んではいたが1回目の波で入ることができた。

Day2の精進は足りない

その後しばらく部屋で寝て、Day2のコンテストが始まった。

コンテストリンク

  • 難易度順がランダムということで、FA狙いで各位バラバラに問題を読むことに
  • 僕がE問題から読み始める
    • 面倒そうな式だな→いやギャグやんけこれ→AC
    • 5分ACで全体2位だったが、FAとは2分差だった
  • J問題がいけそうになったので僕が実装
    • サンプルで嘘であることに気がついて悲しくなる
  • SとOがA問題を考察
    • Sが実装→AC
  • 僕がD問題、SがJ問題、OがI問題を考察
  • D問題の解法が分かるも、計算量がヤバい
    • しばらく考えていたが、計算量の見積もりミスに気づく
    • 僕が実装→AC
  • I問題の考察が終わる
    • Oに実装をしてもらう
    • 大きいケースは良いが、小さいケースのシミュレートが大変
    • ジャッジがバグっているという情報が入るが、とりあえず投げておく
  • SがJ問題の解となる式を導出する
    • しかしこれの計算量がどうしても減らせない
    • 見積もり的にコンテスト中には間に合わないが、とりあえず埋め込みを走らせる
  • 他の問題も何も分からず、時間だけが過ぎていく
  • I問題のリジャッジが行われるもWA
    • コードを印刷するも何もミスは見当たらない
    • 大きいケースを走らせると、infという出力が出る
    • 累乗関数が名前衝突を起こしていることが発覚→AC
      • 3人のコーディングスタイルの違いが災いの元だった
  • 終了数分前にF問題が一般マッチングに帰着できることに気づいて終了

結果ADEIの20位(オンサイトで多分13位)。I問題のジャッジがバグってなければもっと上だったと思うが、それにしても悲しい結果だった。

反省点は「コーディングスタイルの違いがバラバラだった」くらいで、他の問題は現状どうしようもないという感覚が強い。シンプルに実力不足を痛感した。

懇親会の飯は多い

コンテスト後諸々を経て、19:00から懇親会が始まった。

とはいえ根本的にコミュ障なので、あまり人脈が広げられず悲しかった。それでも初めて話した人もいて楽しかった。

あと飯が競プロ関係のものが多かった。以下は哀愁漂うカツサンドくん。

懇親会後、風呂から部屋に戻るとこたつがめくんがABCのコードゴルフをしていた。僕はABCに参加する予定はなかったので観戦していたが、AtCoderのコーディングエリアにコードをベタ書きする様は見ていて圧巻だった。コンテスト全体でもちゃんと全完していて凄いなぁとなった。

ちなみにこの日の夜もなんやかんや1:30まで起きていた。人はなぜ。

Day3

Day3の嘘は怖い

部屋でダラダラしていたらチェックアウトがビリになっていた。全員とっくに出る準備はできていたんです。

それはさておきDay3のコンテストが始まった。

コンテストリンク

  • やや変則的だがSがテンプレを書きつつ僕がA問題を考察することに
    • 解けたので解法をSに伝える→AC
  • 実装中に僕がB問題を考察する
    • 順番は前後しないから二分探索だなとなる
    • greater or equal to」の「or equal to」を見逃して1WA→AC
  • 諸事情でPCが使えない時間があり、B問題の実装待ち中にSとC問題を考察する
    • 往復回数は区間グラフの彩色問題で、これは離散数学演習でやったため解ける
    • 最後の週が面倒だが、それを詰め切る前にPCが使えるようになって後は任せる
    • 詰め切れたらしくOが実装→AC
  • D問題が構築なのでThe atamaのatamaことSが取り組む
    • X = 0とA ^ O = Xが無理ということが分かる
    • Xのbitが全部1の場合に構築ができる
    • 他の場合も帰着できるかと思ったができなくて終了
  • E問題は見た感じフローっぽいのでグラフを考える
    • 流量下限を課すことでグラフが完成する(嘘)
    • 僕が最小費用流を実装をする間Oにグラフの構築を頼む
    • 実装が終わるもなんか上手くいかない
    • どうすれば直るのかよく分からないままコンテストが終わった
  • 後で聞いた話だが、終盤SはJ問題(これも構築)を考察していた
    • Nが奇数の場合が無理であることは分かっていた
    • Nが偶数の場合、N頂点の完全グラフを完全マッチングに分解できれば解けていた
    • コンテスト後に調べたら構築方法が出てきた。

結局ABCの20位(オンサイトで多分12位)。コンテスト後無力感に包まれていた。 なお他の問題もちゃんと読んでいて、合間合間で考察をしていたが全く分からなかった。

コンテスト後にhenoくんにE問題のフローが嘘であることを指摘してもらった。 そしてD問題はアドホック構築ではなく2-SATが想定解だった。僕たちの3時間は一体......。

万豚記の閉店は早い

そんなこんなで3日間に渡る合宿は終わりを告げた。 JAGの皆さん、writer陣の皆さん、本当にありがとうございました。 正直なところ、夏休みの頭にやってほしかった。

ちなみに合宿中に何度か「近くにある万豚記という店の炒飯が旨い」とOに熱弁されていたので、帰りにThe atama3人で寄ることにした。

ARC100-E 「Or Plus Max」 (700) 高速ゼータ変換解法

atcoder.jp

本記事は以下の記事の番外編です。

misteer.hatenablog.com

問題

略。

考察

misteer.hatenablog.com

以降 X = \{ 0, 1, \cdots, n - 1 \}とする。 また S \subseteq Xに対して Sを二進表記したものを B_S \in {\mathbb N}としたとき、 A_S := A_{B_S}という略記を使うことにする。

上の記事で述べたように、各 S \subseteq Xについて「 T \subseteq Sの中で A_Tの大きい方から2つ」が欲しい。そこで、以下のような加法の定義された集合 Gを考える。

  •  g \in Gについて g \subseteq {\cal P}(X), |g| \leq 2
  •  g, g' \in Gについて、 g + g'を「 g \cup g'から A_Sが大きい方2つのみを残した集合」として定義する。

これは結合則・交換則を満たすため高速ゼータ変換を施すことができる。初期値は f(S) = \{ S \}

実装例

atcoder.jp

struct G {
    vector<int> idx;
    static vector<int> a;
 
    G& operator+=(const G& opp) {
        std::copy(opp.idx.begin(), opp.idx.end(),
                  std::back_inserter(idx));
        std::sort(idx.begin(), idx.end(), [&](Int i, Int j) {
            return a[i] == a[j] ? i < j : a[i] > a[j];
        });
        idx.erase(std::unique(idx.begin(), idx.end()), idx.end());
        idx.resize(2);
        return *this;
    }
};
vector<int> G::a;
 
int main() {
    int n;
    std::cin >> n;
    vector<int> a(1 << n);
    for (auto& x : a) std::cin >> x;
    G::a = a;
 
    vector<G> g(1 << n);
    for (int s = 0; s < (1 << n); ++s) {
        g[s].idx.push_back(s);
    }
 
    // 高速ゼータ変換 下位集合版
    for (int k = 0; k < n; ++k) {
        for (int b = 0; b < (1 << n); ++b) {
            if ((b >> k) & 1) {
                g[b] += g[b ^ (1 << k)];
            }
        }
    }
 
    Int max = 0;
    for (int s = 1; s < (1 << n); ++s) {
        max = std::max(max, a[g[s].idx[0]] + a[g[s].idx[1]]);
        std::cout << max << std::endl;
    }
    return 0;
}

高速ゼータ・メビウス変換

ググれば大量の記事が出てくるが、知見の整理ということで自分でも書いてみる。

ゼータ・メビウス変換とは

「高速」ゼータ・メビウス変換とは名前の通りゼータ・メビウス変換を高速化したものなので、まずはそちらを理解する必要がある。

定義: ゼータ変換、メビウス変換

 Xを有限集合、 Gを加法が定義された集合1とする。

このとき写像 f: {\cal P}(X) \to Gから、以下を満たす写像 g: {\cal P}(X) \to Gを求めることをゼータ変換という2

 \displaystyle g(S) = \sum_{T \supseteq S} f(T)

なお、この式中の T \supseteq S T \subseteq Sにしたものも同様にゼータ変換と呼ぶことにする3

そしてこの逆変換、すなわち g: {\cal P}(X) \to Gから以下を満たす f: {\cal P}(X) \to Gを求めることをメビウス変換という4

 \displaystyle f(S) = \sum_{T \supseteq S} (-1)^{|T| - |S|}g(T)

この式中の T \supseteq S T \subseteq Sにすると、上で包含関係の向きを逆にした場合のゼータ変換の逆変換になる。

この変換がゼータ変換の逆変換であることは全く自明ではないが、この記事ではその証明は省略する5


以降、集合 Xとして \{ 1, 2, \cdots, n \}のみを考え、関数 fを長さ 2^nの配列として表現する。これと対応して、以降 f(\{ 2, 4, 5 \}) f(11010)のように表現する。

ゼータ変換の高速化

ここからゼータ変換、特に実装しやすい g(S) = \sum_{T \subseteq S} f(T)の方を高速に行うことを考える。メビウス変換も符号を少し変えるだけでほぼ同様に実現できる。

計算量の評価は、加算が行われる回数 C (Complexity)によって行う。

愚直解

まず愚直に求めようとすると、以下のようなコードになるだろう。

for (int s = 0; s < (1 << n); ++s) {
    for (int t = 0; t < (1 << n); ++t) {
        if ((s & t) == t) g[s] += f[t];
    }
}

 S T Xの部分集合全体を回り、 S \cap T = T、すなわち T \subseteq Sのときだけ g(S) f(T)を加算する。定義通りである。

 Xの部分集合は計 2^n個あることから、このコードの計算量は C = 2^n \cdot 2^n = 4^nとなる。

部分集合の列挙

しかしこのコードでは Tの探索に無駄が多すぎる。できることなら Xではなく Sの部分集合全体だけを回したいが、実はこれは以下のコードで実現できる。

for (int s = 0; s < (1 << n); ++s) {
    for (int t = s; t >= 0; t = ((t - 1) & s)) {
        g[s] += f[t];
    }
}

この手法はビット演算 (bit 演算) の使い方を総特集! 〜 マスクビットから bit DP まで 〜に載っている。

 Sの部分集合は 2^{|S|}個あり、 Xのサイズ kの部分集合は {}_n C_k個あることから、このコードの計算量は C = \sum_{k = 0}^{n} {}_n C_k \cdot 2^k = 3^nとなる6

高速ゼータ変換

さらに高速化して、 C = n \cdot 2^{n - 1}でゼータ変換を行うのが本記事の主題である高速ゼータ変換と呼ばれる手法である。

方針は「下位要素から順に拡張していく」というもの。 まず dp_{S, k}を「 1 \sim k Sの部分集合、それより上は Sと一致しているもの全部の総和」と定義する。言葉では分かりにくいので例を挙げると、

 \displaystyle dp_{\color{red}{110}\color{blue}{101}, 3} =
f(\color{red}{110}\color{blue}{000}) +
f(\color{red}{110}\color{blue}{001}) +
f(\color{red}{110}\color{blue}{100}) +
f(\color{red}{110}\color{blue}{101})

といった具合である。定義より、 dp_{S, 0} = f(S) dp_{S, n} = g(S)

これを k = 1から k = nまで順に更新していく。このとき、 k \in Sか否かで場合分けをする。具体例で考えてみると、

 \displaystyle dp_{\color{red}{10}\color{green}{0}\color{blue}{10}, 3} =
f(\color{red}{10}\color{green}{0}\color{blue}{00}) +
f(\color{red}{10}\color{green}{0}\color{blue}{10}) =
dp_{\color{red}{10}\color{green}{0}\color{blue}{10}, 2} \\
dp_{\color{red}{10}\color{green}{1}\color{blue}{10}, 3} =
f(\color{red}{10}\color{green}{0}\color{blue}{00}) +
f(\color{red}{10}\color{green}{0}\color{blue}{10}) +
f(\color{red}{10}\color{green}{1}\color{blue}{00}) +
f(\color{red}{10}\color{green}{1}\color{blue}{10}) =
dp_{\color{red}{10}\color{green}{0}\color{blue}{10}, 2} +
dp_{\color{red}{10}\color{green}{1}\color{blue}{10}, 2}

のように、 k \in Sの場合は dp_{S \backslash \{ k \}, k - 1}も足す必要がある。

これを実装に落とし込むと以下の通り。

for (int s = 0; s < (1 << n); ++s) {
    dp[s][0] = f[s];
}

for (int k = 1; k <= n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        dp[s][k] = dp[s][k - 1];
        if ((s >> (k - 1)) & 1) {
            dp[s][k] += dp[s ^ (1 << (k - 1))][k - 1];
        }
    }
}

for (int s = 0; s < (1 << n); ++s) {
    g[s] = dp[s][n];
}

集合及びDPの添字は1-indexedなのに対してbit演算は0-indexedなので、シフト数を1減らしてやる必要がある。

配列の使い回し

さらに更新について考えてみると、DPテーブルは1次元のものを使い回せることに気づく。すなわち、以下のような実装ができる。

for (int s = 0; s < (1 << n); ++s) {
    dp[s] = f[s];
}

for (int k = 1; k <= n; ++k) {
    // この時点でdp[s]は上の実装のdp[s][k-1]と一致
    for (int s = 0; s < (1 << n); ++s) {
        if ((s >> (k - 1)) & 1) {
            dp[s] += dp[s ^ (1 << (k - 1))];
        }
    }
}

for (int s = 0; s < (1 << n); ++s) {
    g[s] = dp[s];
}

ここで「 dp_{S \backslash \{ k \}} dp_{S}より先に更新されているが大丈夫なのか」となるが、 k \not \in S \backslash \{ k \}なので、この週ではそもそも dp_{S \backslash \{ k \}}は更新されていない。よって問題ない。

さらにkを1ずらしたりDP配列として直にgを使ったりして整理することで、他の記事でもよく見られるような実装になる。

for (int s = 0; s < (1 << n); ++s) {
    g[s] = f[s];
}

for (int k = 0; k < n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        if ((s >> k) & 1) {
            g[s] += g[s ^ (1 << k)];
        }
    }
}

長い道のりだったが、ようやく高速ゼータ変換のコードに辿り着いた。

このコードの計算量は、各 S \subseteq Xについて |S|回加算更新が行われるため、 C = \sum_{k = 0}^{n} {}_n C_k \cdot k = n \cdot 2^{n - 1}となる7。加算が O(1)なら、 n = 18くらいなら安定して間に合うだろう( C \simeq 2.4 \times 10^6)。

高速ゼータ・メビウス変換の実装例まとめ

最後にまとめとして、各種変換の実装例を載せておく。 なお配列 f, gはコード外で宣言済みとする。

(注※) verifyはしていないので、コピペは自己責任で。実装ミスがあれば報告していただけると助かります。

高速ゼータ変換 上位集合版

 g(S) = \sum_{T \supseteq S} f(T)

for (int s = 0; s < (1 << n); ++s) {
    g[s] = f[s];
}

for (int k = 0; k < n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        if (!((s >> k) & 1)) {
            g[s] += g[s ^ (1 << k)];
        }
    }
}

先程とは逆に、 k \not \in Sのときに dp_{S \cup \{ k \}}を加える。この遷移は Sが大きい方から小さい方に加えられるので、更新順序的に問題ない。

ゼータ変換 下位集合版

 g(S) = \sum_{T \subseteq S} f(T)

for (int s = 0; s < (1 << n); ++s) {
    g[s] = f[s];
}

for (int k = 0; k < n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        if ((s >> k) & 1) {
            g[s] += g[s ^ (1 << k)];
        }
    }
}

先程長々と解説した通り。

メビウス変換 上位集合版

 f(S) = \sum_{T \supseteq S} (-1)^{|T| - |S|} g(T)

for (int s = 0; s < (1 << n); ++s) {
    f[s] = g[s];
}

for (int k = 0; k < n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        if (!((s >> k) & 1)) {
            f[s] -= f[s ^ (1 << k)];
        }
    }
}

非自明ではあるが、実はゼータ変換の更新の符号を反転させるだけでいい。  kが抜けると dp_{S \cup \{k\}}の各項について |S|が1減るので符号が反転する、みたいなお気持ち。

メビウス変換 下位集合版

 f(S) = \sum_{T \subseteq S} (-1)^{|S| - |T|} g(T)

for (int s = 0; s < (1 << n); ++s) {
    f[s] = g[s];
}

for (int k = 0; k < n; ++k) {
    for (int s = 0; s < (1 << n); ++s) {
        if ((s >> k) & 1) {
            f[s] -= f[s ^ (1 << k)];
        }
    }
}

先と同じく符号を反転させるだけでいい。


  1. ゼータ変換では、 \sumを使うため Gは可換半群であるべきだと思われる。しかし代数学にはあまり自信がないので注釈に留めておく。

  2. 別な表現として、 g = \zeta fを満たすような写像 \zeta: G^{{\cal P}(X)} \to G^{{\cal P}(X)}をゼータ変換、(存在すれば)その逆写像 \zeta^{-1}メビウス変換と呼ぶこともできる。が、却って分かりにくくなりそうなので割愛。

  3. この辺は記事によってはメビウス変換になっていたりして曖昧な部分だが、ここでは高速ゼータ変換/高速メビウス変換に従うことにする。

  4. メビウス変換では逆元を必要とするため、 Gはより強い性質をもつアーベル群である必要があると思われる。

  5. 指数時間アルゴリズム入門の53ページに式変形が載っている。 f gが逆になっているので注意。

  6. 2つ目の等号は、二項定理で (1 + 2)^nを展開すると成り立つことが分かる。

  7. この2つ目の等号を示す方法として、例えば (1 + x)^nを二項定理を使って展開し、辺々を微分してから x = 1を代入するなどがある。

ARC100-E 「Or Plus Max」 (700)

atcoder.jp

問題

長さ 2^Nの数列 (A_i) (0-indexed)が与えられる。

 1 \leq k \lt 2^Nなる各整数 kに対し、 \max \{ A_i + A_j \mid i \neq j, \; i \lor j \leq k \}を求めよ。ここで \lorはbitwise orを意味する。

制約

  •  1 \leq N \leq 18
  •  1 \leq A_i \leq 10^9

考察

以降、 j \subseteq i \iff i \land j = j ( \landはbitwise and)という記法を用いる。これは整数を「立っているbitからなる集合」と見たときに、集合における意味と一致する。また同様の理屈で、 |i|を立っているbitの数とする(絶対値ではない)。


 \max \{ A_i + A_j \mid i \lor j \lt k \}は既に求まっているので、代わりに \max \{ A_i + A_j \mid i \lor j = k \}を求めればよい。  i \lor j = kのとき i, j \subseteq kなので、「 i \subseteq kなる iを全列挙して上から2つを取る」という解法が考えられる。しかしこの計算量は O(2^{|k|})となるので遅すぎる。そこでこれを以下のような手法で高速化する1


例えばもっと簡単な問題として、各 kについて X_k = \max \{ A_i \mid i \subseteq k \}を求めたいとする。このとき、 X_kは以下の式で効率的に求められる。

 \displaystyle X_k = \max(A_k, \max \{ X_l \mid l \subset k, \; |l| = |k| - 1 \})

つまりそれまでの結果を使うことで、 kからbitを1つ下ろしたものだけを列挙すれば十分となる。この工夫によって、計算量は愚直な場合の O(2^{|k|})から O(|k|)に落ちる。


この手法を今回の問題に適用する。今回は大きい方から2つの要素が欲しいので、各 kについて「 i \subseteq kなる iで、 A_iが大きい方から2つ」という集合 Y_kを求めることにする。

ここで X_kでは最大値自体を持っていたのに対して、 Y_kは最大値を取るindexを持っていることに注意。これは A_i = A_jかつ Y_k = \{ i, j \}となるような場合があるためである。

このとき Y_kは、集合

 \displaystyle k \cup \bigcup_{l \subseteq k, |l| = |k| - 1} Y_l

から2つを取ることで求められる。この集合のサイズは O(N)なので、全体の計算量は O(N 2^N)となる。実装時は要素の重複を消し忘れないように注意。

実装例

atcoder.jp

int main() {
    Int n;
    std::cin >> n;
    Int m = (1 << n);
 
    Vector<Int> a(m);
    for (auto& x : a) std::cin >> x;
 
    Vector<Vector<Int>> max(m);
    max[0] = {0};
    // kの部分bit全体で、aが大きい方から2つ(0だけ1つ)
 
    Int ans = 0;
    for (Int k = 1; k < m; ++k) {
        auto& idx = max[k];
        idx.push_back(k);
 
        // kからbitを1つだけ下ろしたlについて、max[l]をマージ
        for (Int i = 0; i < n; ++i) {
            if ((k >> i) & 1) {
                Int l = k - (1 << i);
                for (auto x : max[l]) {
                    idx.push_back(x);
                }
            }
        }
 
        // a[i]の大きさによってソート
        // 後でuniqueで重複を消すため、a[i]が同じ場合はiについてソート
        std::sort(idx.begin(), idx.end(), [&](Int i, Int j) {
            return (a[i] == a[j] ? i < j : a[i] > a[j]);
        });
 
        // 重複を消して大きい方から2つを残す
        idx.erase(std::unique(idx.begin(), idx.end()), idx.end());
        idx.resize(2);
 
        ans = std::max(ans, a[idx[0]] + a[idx[1]]);
        std::cout << ans << std::endl;
    }
    return 0;
}

補足

misteer.hatenablog.com


  1. この手法は高速ゼータ変換と似て非なるものである。高速ゼータ変換を使った解法は本記事末尾のリンクを参照。

ARC095-E 「Symmetric Grid」 (700)

atcoder.jp

問題

 H \times Wのアルファベットからなるグリッドが与えられる。

この行と列を自由に入れ替えることで、点対称にできるか判定せよ。

制約

  •  1 \leq H, W \leq 12

考察

脳死 O(HW \times H! W!)をやりたくなるが、当然それは許されない。

そこで行の並び方だけ全列挙することにする。すると、「互いに逆の並びになっている列のペアで両端から埋めていく」という貪欲によって、列の並び方を全列挙することなく可能性判定ができる。

ただし Wが奇数の場合は注意が必要。このとき真ん中の列に来るものは他とペアになっている必要はないが、代わりに自分自身とペアになる、すなわち回文になっている必要がある。

これで O(HW^2 H!)に落ちるが、まだ際どい。そこで行の並べ方を枝刈りする。点対称なグリッドにおいて、その対称な位置にある行を交換してもこれは点対称である。よって、両端から「上の行の方が小さくなるように」行番号を2つずつ決めていけば、それで可能性を全列挙できたことになる。

これで O(HW^2 H!!)に落ち、 12!! = 12 \cdot 10 \cdot \cdots \cdot 2 = 46080なので割と際どいが間に合う。


以上の解法はEditorialそのままだが、難しいのは実装の方である。まず「両端から行番号を決めていく」というのが、いざ実装しようとすると結構面倒なことが分かる。そして何よりも「 H, Wが奇数」のケースが本当に面倒くさい。

通すまでに相当な数のバグを踏んだので、以下の実装例に自分が踏んだバグを一通り書いておく。

実装例

atcoder.jp

以下では無名構造体内部で探索するような実装をしているが、ただ「グローバル変数を増やしたくない」というだけの理由なので、気にならない人は普通にいくつか関数を宣言すればいいと思われる。

// 各関数は上から下に流れるようになっている
// ! で始まるコメントは自分が踏んだバグ
struct {
    Int h, w;
    Vector<Int> col;
    Vector<String> s, t;
    // sは入力 tは行を並び替えた後、転置したもの
    Bool possible;

    void operator()() {
        possible = false;
        std::cin >> h >> w;

        s.resize(h), t.resize(w);
        for (auto& x : s) std::cin >> x;
        for (auto& y : t) y.resize(h);  // ! tのresizeし忘れ

        col.resize(h);
        enumcol_center();
    }

    // 真ん中の行の列挙 
    // ! これを後回しにして、中心の全列挙をミスる
    void enumcol_center() {
        if (h % 2 == 0) {
            enumcol(0, 0);
        } else {
            for (Int x = 0; x < h; ++x) {
                col[h / 2] = x;
                enumcol(1 << x, 0);
            }
        }
    }

    // bは固定済みの行 iはいくつのペアを固定したか
    void enumcol(Int b, Int i) {
        Int x;
        for (x = 0; x < h; ++x) {
            if (!((b >> x) & 1)) break;
        }
        if (x == h) {
            align();
            return;  // ! ここでreturnしないとhを固定し始める
        }

        for (Int y = x + 1; y < h; ++y) {
            if (!((b >> y) & 1)) {
                col[i] = x;
                col[h - i - 1] = y;
                enumcol(b | (1 << x) | (1 << y), i + 1);
            }
        }
    }

    // colに従ってtを埋める
    void align() {
        for (Int y = 0; y < w; ++y) {
            for (Int x = 0; x < h; ++x) {
                t[y][x] = s[col[x]][y];
            }
        }
        check_center();
    }

    // 真ん中の行の列挙 
    // ! enumcol_center()と同じく
    void check_center() {
        if (w % 2 == 0) {
            check(0);
        } else {
            for (Int x = 0; x < w; ++x) {
                // 回文判定
                String z = t[x];
                std::reverse(z.begin(), z.end());
                if (z == t[x]) check(1 << x);
            }
        }
    }

    // 両端から、互いに逆なペアを決めていく
    // bは決定済みの行番号
    void check(Int b) {
        Int x;
        for (x = 0; x < w; ++x) {
            if (!((b >> x) & 1)) break;
        }

        if (x == w) {
            possible = true;
            return;  // ! ここで枝刈りしないと、全文字が同じケースで非常に時間がかかる
                     //   なんならYESを出力してexit()を呼んでもいい
        }

        String z = t[x];
        std::reverse(z.begin(), z.end());

        for (Int y = x + 1; y < w; ++y) {
            if (!((b >> y) & 1)) {
                if (t[y] == z) {
                    check(b | (1 << x) | (1 << y));
                    return;
                }
            }
        }
    }
} solve;

int main() {
    solve();
    std::cout << (solve.possible ? "YES" : "NO") << std::endl;
    return 0;
}

ちなみに、行列挙を「next_permutation()で枝刈りする」として手を抜くとTLEになる。

JSC2019予選 - D 「Classified」 (600)

atcoder.jp

問題

 N頂点完全グラフの各辺に対して、以下を満たすように正整数を割り当てる。

  • 頂点 vと正整数 Lを任意に取る。
  •  vから Lが割り当てられた辺のみを辿って vに戻るとき、そのパスの長さは常に偶数。

このような割り当てで、割り当てられる正整数の最大値が最小のものを出力せよ。

制約

  •  2 \leq N \leq 500

考察

まず条件を分かりやすいように言い換えると、

  • 任意の正整数 Lについて、 Lが割り当てられた辺による誘導部分グラフが二部グラフをなす。

となる。「奇数長の閉路を持たない」と「二部グラフである」が同値であることが肝。

まず1を割り当てる辺を決める。 極力辺を使った方が後々楽なので、辺数が最多となる割り当てがしたい。 そのような割り当てが何かと考えると、これは N頂点を半々に分割した完全二部グラフとなる。

次に2を割り当てるが、残った辺に注目するとこれは2つのクリークをなしている。よってこれらに再帰的に上と同じ割り当てをしてやればいい。

1つの正整数を割り当てる度に頂点数が半減するので、最大値は \lceil \log_2 N \rceilになる。 最小性の証明は思いつかなかったが、まぁ最適だろうと直感を信じて投げたら通った。

実装例

atcoder.jp

using namespace std;

int main() {
    int N;
    cin >> N;

    vector<vector<int>> ans(N, vector<int>(N, -1));

    // 頂点集合[l, r)からなるクリークにd以上を割り当てる
    function<void(int, int, int)> dfs =
        [&](int l, int r, int d) {
            if (r - l <= 1) return;

            int m = (l + r) / 2;
            // [l, m)と[m, r)に分割して完全に部グラフを構築
            for (int i = l; i < m; ++i) {
                for (int j = m; j < r; ++j) {
                    ans[i][j] = d;
                }
            }

            dfs(l, m, d + 1);
            dfs(m, r, d + 1);
        };
    dfs(0, N, 1);

    for (int v = 0; v < N; ++v) {
        for (int u = v + 1; u < N; ++u) {
            cout << ans[v][u] << " ";
        }
        cout << endl;
    }
    return 0;
}

JSC2019予選 - C 「Cell Inversion」 (500)

atcoder.jp

問題

長さ 2NのBとWからなる文字列 Sが与えられる。この文字列に対して、以下の操作を N回行う。

  • まだ選んでいない2つのindexを選ぶ( l, rとする)。
  •  S[l, r]を全て反転させる。

最終的に全部Wになるような操作列の個数を求めよ。

制約

  •  1 \leq N \leq 10^5

考察

 N個のマッチングを決めてしまえば、最終的な文字列は選ぶ順番に依存しない。 そこでマッチングを左から決めていくことを考える。

 S_1に着目すると、これは絶対に一度だけ反転される。よって S_1がWである場合、答えは 0になる。  S_1がBである場合は1とペアになるindexがあるが、これは後で決めることにする。

次に S_2に着目すると、これの反転する回数は以下の2パターンがある。

  • 1のペアを2にする場合は1回。
  • 2と3以降のどれかをペアにする場合は2回。

このうち S_2がWになるのはいずれか一方のみである。

これと同じことが任意のindexで言える。 つまり、 S_iを見ている時点で右端を保留しているペアが K個ある場合、

  • 保留中のペアのうちいずれか1つの右端を iにする場合は K回。
    • どれと iをペアにするかが K通りある。
    •  Kは1減る。
  •  iをペアの左端にする場合は K + 1回。
    •  Kは1増える。

このうち S_iがWになるのはいずれか一方のみである。 ただし K = 0の場合前者は選べない、つまり条件を満たすマッチングはない。

これを先頭から繰り返せばマッチングの数が求まる。 ただし最終的に K \neq 0になったときは、これはinvalidなマッチングなので答えは0になることに注意。

最後にマッチングをどの順に適用するか、つまり N!を掛ければ答えとなる。

実装例

atcoder.jp

using namespace std;
using mint = ModInt<1000000007>;

int main() {
    int N;
    string S;
    cin >> N >> S;

    int K = 0;
    mint ans = 1;
    for (auto c : S) {
        if ((c == 'W') == (K % 2 == 0)) {
            // 保留中のどれかの右端にする
            ans *= K;
            if ((--K) < 0) break;
        } else {
            // iを左端にする 右端は保留
            ++K;
        }
    }

    if (K != 0) {
        // invalidなマッチング
        cout << 0 << endl;
        return 0;
    }

    // N!を掛けて出力
    for (int i = 1; i <= N; ++i) ans *= i;
    cout << ans << endl;
    return 0;
}