各バケットに区間のindexを持たせる平方分割
結論
平方分割の実装はかなりきれいになりました
— りあん (@rian_tkb) 2020年9月30日
本質部分は backet 内の private method 6 個を書き換えれば良いですhttps://t.co/v2Rx6coYay
まぁ本質以外もちょこちょこ書きかえなきゃいけなくて結局そこでバグらせるんですが……
この記事について
この記事は、平方分割の簡単かもしれない実装方法を紹介するものです。
平方分割って何? という人には、以下の 3 つの記事がわかりやすいのでおすすめです。
www.slideshare.net caddi.tech kujira16.hateblo.jp
(本当はこの記事でも平方分割についてわかりやすい包括的な説明をしようと思っていたのですが力尽きました…)
実装方法
多くの記事では、主に以下のような実装がなされることが多いです。
[TODO: ここによくある実装を書く]
ただ、どうしてもインデックスをいじりつつ更新等をしていくことになるため、実装が煩雑になったりバグの原因になったりすることが多いです。
例えば、区間更新、区間 Max を平方分割で解くときに、区間 Max のクエリに対し、区間の一部が被っているバケットについてはそのバケット内の全ての値を更新した上で Max を求め、区間に完全に包含されているバケットについては値の更新はせずに区間に対する更新を適用した上での Max を求める、と言った処理をすることになり、それをインデックスをいじっている中で行うのはかなりつらいです(個人の感想です)。
実装方針
ここで、以下のような実装をすることを考えます。
- 問題を解く関数側では、とりあえずそのバケットに区間に包含されているかは全く気にせず全てのバケットに対し更新、計算を行わせる。
- バケット側で区間との包含関係(全く被っていない、完全に包含されている、一部が被っている)を求め、それぞれに応じた処理を行う。
これはかなりセグメント木に似てますね。 実際に平方分割は、セグメント木が 2 分木なのに対する √N 分木と考えることができるので、同じような実装ができるのも妥当とも言えそうです。
具体的な実装
ここでは、区間 Add 区間 Min の以下の問題を例にコードを示していきます。
なお、この問題は遅延セグメント木などを用いることでより良い計算量で解くことが可能であることに注意してください。
骨組み
今回は、平方分割 class sqrt_decomp
の中に、1 つのバケットを表すような class bucket
を用意します。
class sqrt_decomp { class bucket { int l, r; vector<int> v; public: bucket(const vector<int>& a, int l, int r) : l(l), r(r) { v = vector<int>(r - l); for (int i = 0; i < r - l; ++i) { v[i] = a[i + l]; } } }; vector<bucket> v; public: sqrt_decomp(const vector<int>& a, int bucket_size) { v = vector<bucket>(); for (int i = 0; i < (int)a.size(); i += bucket_size) { v.emplace_back(a, i, min(i + bucket_size, (int)a.size())); } } };
更新クエリ
まずは更新を書いてみます。
a_l, a_{l+1}, ... , a_{r-1}
に x
を加算するクエリを考えます。
バケット全体が更新区間に含まれる場合
まず、そのバケット全体が更新区間に含まれている場合を考えます。そのとき、当然その区間の値には全て x
が加算されます。
ただ、そのバケット内の値を一つずつ変更していくと結局全体の計算量 になってしまうので、「そのバケット全体にいくつ加算されたか」という情報を持っておくことにします。
int add_val; void update_all(int x) { add_val += x; }
バケットの一部が更新区間に含まれる場合
次に、そのバケットの一部が更新区間に含まれる場合を考えます。 この場合は、以下の 3 つの処理を順番に行う必要があります(演算の性質によっては一部の更新をサボれることもありますが、ここでは割愛します)。
- バケット全体に対する更新を適用する (
reflect_update
) - 今回のクエリでの更新を適用する
- バケット全体に対する回答クエリで答えるべき値(今回は最小値)を求める (
reflect_calc
)
1 番でバケット全体に対する更新を適用したら、きちんと「そのバケット全体にいくつ加算されたか」を初期化することを忘れないようにしましょう。
vector<int> v; int min_val; int add_val; void reflect_update() { for (auto& i : v) { i += add_val; } add_val = 0; } void reflect_calc() { min_val = INT_MAX; for (auto& i : v) { min_val = min(min_val, i); } } void update_partial(int l, int r, int x) { reflect_update(); for (int i = l; i < r; ++i) { v[i] += x; } reflect_calc(); }
呼び出し方法
上で作った関数を呼び出せば更新クエリは完了です。
呼び出し元からはそのバケットとクエリ区間との位置関係は気にせず全てのバケットに対して更新を走らせて、実際にどういった更新を行うかは中で判定するようにすると、問題に応じて書きかえる部分がシンプルになるので楽だと思います。
int l, r; void update(int s, int t, int x) { if (r <= s || t <= l) { return; } else if (s <= l && r <= t) { update_all(x); } else { update_partial(max(s, l) - l, min(t, r) - l, x); } }
ここまで書いたコードを全て bucket
class の中に放り込み、sqrt_bucket
class の方に全ての bucket に対して update
を走らせる関数を書けば変更クエリは完成です。
vector<bucket> v; void update(int l, int r, int x) { for (auto& i : v) { i.update(l, r, x); } }
回答クエリ
次に回答クエリ(今回は最小値)を書いていきます。
a_l, a_{l+1}, ... , a_{r-1}
の最小値を求めるクエリを考えます。
バケット全体が回答区間に含まれる場合
バケット全体がクエリの区間に含まれている場合、予め計算した最小値に対し、バケット全体への加算を足し合わせればその部分の最小値が得られます。
int min_val; int add_val; int calc_all() { return min_val + add_val; }
バケットの一部が回答区間に含まれる場合
バケットの一部がクエリの区間に含まれている場合、更新の時と同じように以下のような 3 つの処理を順番に行います。
- バケット全体に対する更新を適用する (
reflect_update
) - バケット全体に対する回答クエリで答えるべき値(今回は最小値)を求める (
reflect_calc
) - 今回のクエリの区間に対する答えるべき値(最小値)を求める
vector<int> v; int calc_partial(int l, int r) { reflect_update(); reflect_calc(); int res = INT_MAX; for (int i = l; i < r; ++i) { res = min(res, v[i]); } return res; }
実は…?
実は回答クエリにおいて必ずしもバケット全体に更新を適用する必要はなく、以下のようにクエリと被っている区間のみ一時的に更新を評価するだけでも良いです。
この場合、更新クエリが一部分に被ったときにのみしか更新が伝播しないため、上記の実装よりも若干高速になることが多いです。
vector<int> v; int calc_partial(int l, int r) { int res = INT_MAX; for (int i = l; i < r; ++i) { res = min(res, v[i] + add_val); } return res; }
呼び出し方法
更新クエリと同様に、クエリ区間との位置関係は気にせず全てのバケットに対してクエリを走らせる方針を取ります。
int l, r; int calc(int s, int t) { if (r <= s || t <= l) { return INT_MAX; } else if (s <= l && r <= t) { return calc_all(); } else { return calc_partial(max(s, l) - l, min(t, r) - l); } }
ここまで書いたコードを全て bucket
class の中に放り込み、sqrt_bucket
class の方に全ての bucket に対して calc
を走らせる関数を書けば回答クエリも完成です。
vector<bucket> v; int calc(int l, int r) { int res = INT_MAX; for (auto& i : v) { res = min(res, i.calc(l, r)); } return res; }
提出例
以上のコードをまとめたものが以下になります。 コード自体は長いけど部分部分は簡単……?
class sqrt_decomp { class bucket { int l, r; vector<int> v; int min_val; int add_val; void reflect_update() { for (auto& i : v) { i += add_val; } add_val = 0; } void reflect_calc() { min_val = INT_MAX; for (auto& i : v) { min_val = min(min_val, i); } } void update_partial(int l, int r, int x) { reflect_update(); for (int i = l; i < r; ++i) { v[i] += x; } reflect_calc(); } void update_all(int x) { add_val += x; } int calc_partial(int l, int r) { reflect_update(); reflect_calc(); int res = INT_MAX; for (int i = l; i < r; ++i) { res = min(res, v[i]); } return res; } int calc_all() { return min_val + add_val; } public: bucket(const vector<int>& a, int l, int r) : l(l), r(r), add_val(0) { v = vector<int>(r - l); for (int i = 0; i < r - l; ++i) { v[i] = a[i + l]; } reflect_calc(); } void update(int s, int t, int x) { if (r <= s || t <= l) { return; } else if (s <= l && r <= t) { update_all(x); } else { update_partial(max(s, l) - l, min(t, r) - l, x); } } int calc(int s, int t) { if (r <= s || t <= l) { return INT_MAX; } else if (s <= l && r <= t) { return calc_all(); } else { return calc_partial(max(s, l) - l, min(t, r) - l); } } }; vector<bucket> v; public: sqrt_decomp(const vector<int>& a, int bucket_size) { v = vector<bucket>(); for (int i = 0; i < (int)a.size(); i += bucket_size) { v.emplace_back(a, i, min(i + bucket_size, (int)a.size())); } } void update(int l, int r, int x) { for (auto& i : v) { i.update(l, r, x); } } int calc(int l, int r) { int res = INT_MAX; for (auto& i : v) { res = min(res, i.calc(l, r)); } return res; } };
提出例です https://onlinejudge.u-aizu.ac.jp/solutions/problem/DSL_2_H/review/4879712/rian_tkb/C++17
抽象化について
せっかくなので抽象化して中身をいじらずに様々な問題に適用できるようにしたいのですが、平方分割は遅延セグメント木などに比べて自由度が高く、上手い抽象化方法が考えれらていない状態です。 上手い抽象化が出来た方はぜひ教えていただけると助かります。
類題
- 区間代入区間 Min (AOJ)
- 区間 Add 区間 Sum (AOJ)
- 区間代入区間 Sum (AOJ)
- 上のコードとほぼ同じようなコードで AC できると思います。
- こちらも遅延セグメント木で解くことができます。
- Range K-th Smallest (Library Checker)
- 蟻本の平方分割の項に載っている問題です。が、制約がかなり大きいので通りません(かなり頑張って高速化をすれば通るかもしれませんが…)
- 記事中では
reflect_calc
において「答えそのもの」を求めていますが、実は「答えを高速に求めるための情報」を求めるようにすることでより多様な問題に対応することができ、この問題もそのような問題の一つです。 - この問題では「各バケット内の値をソートしたもの」を求めておくことで、「そのバケット内に
x
以下の数が何個存在するか」という、クエリに依存する値を高速に求めることができます。 - TLE ですがコードです https://judge.yosupo.jp/submission/25166
- 素直に領域木を書いたりウェーブレット行列を使ったりしましょう。
- HUPC2020 day1-G Freqs (AOJ)
- Range K-th Smallest で用いた「
x
以下の数が何個存在するか」というクエリに chmin, chmax, add で更新しつつ答える問題です。 - そのバケット全体に chmin, chmax, add が飛んできても元の要素の順序関係が崩れないので、各要素に更新を伝播しなくとも二分探索で「
x
以下の数が何個存在するか」を求めることができます。 - そのバケットの一部分に chmin, chmax, add が飛んできたら、そのバケット全体を更新してからもう一度ソートしたものを再計算すれば良いです。その操作にかかる計算量は で、一部分に更新クエリが被るようなバケットの個数は 1 つの更新クエリに対し高々 2 個なので、クエリあたりの計算量は となり、間に合います。
- バケット 1 つあたりのサイズに対して がかかるため、バケットサイズは 100 ~ 200 程度にすると速いです(計算量的には変わらないはず…?)。
- コードです https://onlinejudge.u-aizu.ac.jp/solutions/problem/3170/review/4880051/rian_tkb/C++17
- Range K-th Smallest で用いた「
- Range Chmin Chmax Add Range Sum (Library Checker)
- segment tree beats で解けることで有名な問題です。
- 上の問題のコードを少し変えるとこちらも解くことができます。
- 「ソートしたもの」の他に「ソートしたものにおける累積和」を持っておくことで、「chmin, chmax により値が変わらなかったもの」の総和を で求めることができ、上の問題と同じ計算量で解くことができます。
- コードです https://judge.yosupo.jp/submission/25173