競技プログラミングで頻出の「ダブリング」を解説する

競技プログラミングでは頻繁に出てくる「ダブリング」という手法について説明しようと思います。競プロをはじめて間もない人や、競プロ外の人に向けて書きたいと思います。

最初に予防線を張っておきますが、内容が正しいかどうかは保証しません。

繰り返し二乗法

繰り返し二乗法という有名なアルゴリズムがあります。

例えば、3の100億(10^10)乗を計算せよと言われた時に、 1回1回計算していたのでは時間がいくらあっても足りません。しかし繰り返し二乗法を使えば、log(100億)くらいの計算量で計算出来るようになります。

具体的にどういう仕組みかを説明するためにより小さな場合として3の11乗を計算するとした時に、

3^11 = (3^8) x (3^2) x (3^1)

と3^(2^k)の積に分解出来るならば、 11を1011と2進数で表した時の1の数分だけで計算が終わることになります。（a^bはaのb乗の意味です）

100億というのは、2進数で表すと 1001010100000010111110010000000000 だそうです。従って、もともと100億という時間かかっていたものが、ビットの数の11という時間だけで済むことになります。単純にいえば10億倍高速化したということです。 これがアルゴリズムの力です。

この手法は、行列の累乗にも応用出来ます。

ポイントは、「何乗」の部分を2進数表記して各々に分解するという点です。

ダブリング

繰り返し二乗法というのは、ダブリングの特殊な場合です。掛け算のような単純な演算だけでなく、ダブリングを使うことが出来る性質を持つ演算であれば適用することが出来ます。木構造で共通の祖先を求めるLCA(Lowest Common Ancestor)というアルゴリズムもダブリングの考えを用います。

上で挙げた繰り返し二乗法の例を一般化してみましょう。

3をかけるというのを関数f(x)=3xとして、xの初期値idを1とします。すると、3のk乗というのは、f(f(f(f(……f(id)))))と、fをk回適用するのと同じことになります。

そして、ダブリングというのはこれを、

ffffffff . ff . f id = f3 . f1 . f0 id = g id

のように分解して、先に関数gを計算するということです。

従って、

関数fkが計算出来る
関数fkが合成可能である

という点が要件となります。

関数が合成可能とは

このうち、「関数が合成可能である」というのはどういうことかというと、 f(x)の出力が環境に依らないつまり「純粋関数」である場合に、合成可能であることがいえます。これは、f(x)=3xなどの例では自明ですが、そのような場合だけでなく、例えば、xに1を入れたら必ず4を返し、xに4を入れたら3を返すようなわけのわからない関数であっても、常にその出力が決まっているのであればこの要件を満たします。

この性質を利用した問題にyukicode No.1013: ○マス進むがあります。

競プロでは、値が大きくなる問題では、大きな素数のmodをとるということが要求されることがあります。これは、繰り返し二乗法の問題では決まってそうなります（だって3の100億乗とか64bitを軽く超えますからね）。この時にもダブリングが使えるというのは、modをとったとしても純粋関数の性質を失っていないからだと説明出来ます。

フレームワーク化

以上の概念を実装したのが以下になります。rust-comp-snippets

リンクした問題も、このフレームワークを使ってAC済です。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
trait Doublable {
    type T: std::fmt::Debug;
    fn x0(&self) -> Self::T;
    fn f(&self) -> Self::T;
    fn ap(&self, f: &Self::T, x: &Self::T) -> Self::T;
    fn inv(&self, x: &Self::T) -> Self::T;
}
struct Doubling<D: Doublable> {
    d: D,
    f_table: Vec<D::T>,
}
impl <D: Doublable> Doubling<D> {
    pub fn new(d: D, maxbit: usize) -> Self {
        let mut f = vec![d.f()];
        for i in 1..=maxbit {
            let x = d.x0();
            let fx = d.ap(&f[i-1], &x);
            let ffx = d.ap(&f[i-1], &fx);
            f.push(d.inv(&ffx));
        }
        Doubling {
            d: d,
            f_table: f,
        }
    }
    pub fn pow(&self, k: i64) -> D::T {
        let mut k = k;
        let mut res = self.d.x0();
        let mut i = 0;
        while k > 0 {
            if k & 1 == 1 {
                res = self.d.ap(&self.f_table[i], &res);
            }
            k >>= 1;
            i += 1;
        }
        res
    }
}