線形RMQ 分割幅34のSparseテーブル

参考

おそらく、以下の Codeforces 解説が最も簡単な線形な静的RMQ実装である。これを自分で使いやすいようにインデックスを逆順にして実装したい。
codeforces.com

分割の幅はO(log n)なら良いが、34の定数にできるのでそれで実装する。

幅33以下の場合

左閉じ右開きの区間(L, R) が R-L ≦ 33 の場合を答えるために小さいRMQを作る。
インデックス i から、連続33要素を見たときに、接頭最小値になっているインデックスに印をつける。

この i より後ろの32要素を32-bitの整数で保存する。このビット表現を全てのiについて求めたいが、スライド最小値のアルゴリズムによって線形時間で求めることができる（ただし、各ビット演算は定数時間でできるとする）。
33よりも小さい範囲でRMQを求めたい場合は、ビットマスクで外側の範囲を0で上書きする。

L, Rが共に34の倍数の場合

34要素ごとに最小値を取り出した数列を作るとこれの長さはn/34である。この数列で単純なSparseテーブルを構築するとL,Rが34の倍数の場合にRMQを答えることができる。
実装は定数34を用いているが、log nで分割されたSparseテーブルの構築の時間計算量・空間計算量はO(n)である。

一般のL, Rの場合

R-L ≦ 33の場合は前述の通り計算できる。
そうでない場合は

// 34x は L を34の倍数に切り上げた値。
int x = (L + 33) / 34;
// 34y は R を34の倍数に切り下げた値。
int y = R / 34;

最小値の候補は3つ

• 区間(L, L+33) の最小値を小さいRMQで求める
• 区間(34x, 34y) の最小値を34で分割したSparseテーブルで求める
• 区間(R-33, R)の最小値を小さいRMQで求める

この3つのうち、最小値が区間(L, R)の最小値である。

コード

judge.yosupo.jp

template<class T, class Comp=less<T> > struct LinearRmq {
    static const int WIDTH = 34;
    static const int BIT_LEN = sizeof(unsigned) * CHAR_BIT;
    static_assert(BIT_LEN == 32);

    vector<T> a;
    vector<unsigned> mask;
    vector<int> table;
    Comp comp;

    LinearRmq(const vector<T> &a_): LinearRmq(a_.begin(), a_.end()) {}

    template<class Iter> LinearRmq(Iter begin, Iter end): a(begin, end), mask(a.size()) {
	int n = a.size();
	unsigned cur = 0;
	for (int i=n, b; i--;) {
	    while (cur && select(i, i + 1 + (b = __builtin_ctz(cur))) == i) { cur ^= 1u << b; }
	    mask[i] = cur;
	    cur = cur << 1 | 1;
	}
	if (int len = n / WIDTH) {
	    int levels = log2plus1(len);
	    table.resize(len * levels, -1);
	    for (int i=0; i<len; i++) { table[i] = select(i*WIDTH, small(i*WIDTH+1)); }
	    for (int s=0; (2<<s)<=len; s++) for (int i=0; i+(2<<s)<=len; i++) {
		table[len*(s+1)+i] = select(table[len*s+i], table[len*s+i+(1<<s)]);
	    }
	}
    }

    int find_index(int l, int r) const {
	assert(0 <= l && l < r && r <= (int)a.size());
	if (r - l < WIDTH) {
	    return l + log2plus1(mask[l] & bit_mask(r - l - 1));
	}
	int ans = small(l);
	int lx = (l + WIDTH - 1) / WIDTH, rx = r / WIDTH;
	if (lx != rx) {
	    int len = (int)a.size() / WIDTH;
	    int s = BIT_LEN - 1 - __builtin_clz(rx - lx);
	    ans = select(ans, select(table[len*s+lx], table[len*s+rx-(1<<s)]));
	}
	ans = select(ans, small(r - (WIDTH - 1)));
	return ans;
    }

    T find_value(int l, int r) const {
	return a[find_index(l, r)];
    }

    int select(int i, int j) const {
	return comp(a[i], a[j])? i: j;
    }

private:

    static inline constexpr unsigned bit_mask(int n) {
	return n? ~0u >> (BIT_LEN - n): 0;
    }

    // x == 0 => 0;
    // o.w.   => floor(log2(x)) + 1;
    static inline constexpr int log2plus1(unsigned x) {
	return x? BIT_LEN - __builtin_clz(x): 0;
    }

    // argmin_{i| l <= i < l+WIDTH} a[i];
    int small(int l) const {
	return l + log2plus1(mask[l]);
    }
};

注意

線形RMQは構成の計算量がスパーステーブルより優れているが、出力命令の中での比較をする回数が異なる。

• 単純なスパーステーブル: 1回
• この記事の線形RMQ: 3回

RMQを問題解決のパーツとして使う時に、出力が高速であるスパーステーブルの方が有利になりやすい。