スムースソートを使う

スムースソート (smooth sort) は、1981年にエドガー・ダイクストラが提案した 比較ソート である。

ヒープソート (heap sort) と同様に 最悪時間計算量は O(n log n) で、追加の補助記憶域は O(1) に抑えられる インプレース なアルゴリズムだが、入力がすでに昇順に近いほど作業量が少なくなる適応的 (adaptive) 性質 を持つ。

通常の二分ヒープでは、根が配列の一端にあり、最大値を取り出すたびに 根と末尾を交換して ヒープを底側から縮めていく。この配置は「ヒープとしての整理」と「昇順に書き出す向き」が噛み合わず、すでにソート済みの入力でも一度ヒープ状に混ぜ直すため、最良でも必ず Ω(n log n) 程度の仕事が残る。

スムースソートでは、二分木ではなく レオナルド数 に基づく レオナルド木 を積み重ねた 森 を維持する。 各木はヒープ条件を満たし、さらに 各木の根の値は左から右へ弱い順（非減少） として並べる。 そうすると 常に右端が全体の最大 になり、末尾から確定した最大値を「取り出す」処理が すでに整っている入力ではほとんど余計な比較を要しない。

レオナルド数 L(k) は

  • L(0) = L(1) = 1
  • L(k) = L(k-2) + L(k-1) + 1（k ≥ 2）

で定められ、1, 1, 3, 5, 9, 15, … と続く。サイズ L(k) のレオナルド木は、大きい方の子が左、小さい方が右になるよう 二つのより小さいレオナルド木と根 で構成される。

任意の長さ n は、高々 O(log n) 個 の互いに異なるレオナルド数の和として表せる。スムースソートはこの性質を利用し、森に含まれる木の本数を常に対数オーダーに抑える。

すでに昇順に近い場合は O(n) に近い時間で済み、定数倍は入力の乱れに依存する。最悪の漸近性能を維持しつつ、順序の良い入力ではヒープソートより反応が柔らかい のが特徴である。名前は、乱れ具合に応じて実行時間が 滑らかに O(n) に近づくことから付けられている。

スムースソートの処理の概略は、単一レオナルド木での根の沈下 sift_in と、複数の木があるときに根が隣へ繰られないよう整える interheap_sift、および ビットパターン＋末尾木の順 でレオナルド森を表す符号化に沿って二段フェーズを進める、という構成で整理できる。

procedure sift_in(A, rootIdx, size)
  // インデックス size のレオナルド木 1 本の中でヒープ条件を満たすよう根から下げる。
  // L[k]=L[k-2]+L[k-1]+1、子方向は二分木のサイズにより左右どちらかへ潜る。
  if size < 2 then return
  tmp = A[rootIdx]
  r = rootIdx
  sz = size
  loop
    right = r - 1
    left = right - L[sz - 2]
    if A[right] < A[left]
      candidate = left
      nsz = sz - 1
    else
      candidate = right
      nsz = sz - 2
    if A[candidate] <= tmp then break
    A[r] = A[candidate]
    r = candidate
    sz = nsz
    if sz <= 1 then break
  A[r] = tmp

procedure interheap_sift(A, rootIdx, heap_state)
  // heap_state は「どの桁に木の根があるか」（マスク）と「現在見ている木の順 order」（オフセット）などを束ねたもの。
  tmp = A[rootIdx]
  r = rootIdx
  state = heap_state anchored at rootIdx
  loop while mask(state) <> 1
    maxValue = tmp
    if order(state) > 1 then
      right = r - 1
      left = right - L[order(state) - 2]
      maxValue = max(maxValue, A[left], A[right])
    next = r - L[order(state)]
    if A[next] <= maxValue then break
    A[r] = A[next]
    r = next
    // マスクを調整して左隣の木へカーソル移動
    march_one_tree_left(state)
  A[r] = tmp
  sift_in(A, r, order(state))

procedure smooth_sort(A)
  n = length(A)
  if n <= 1 then return
  heap_state ← initial_singleton_forest_encoding()
  // 第 1 段階（ヒープ化）
  for i = 1 to n - 1
    advance_leonardo_forest(heap_state, i)
    if wide_bottom_condition(i, heap_state, n) then sift_in(A, i, order_at_tail(heap_state)) else interheap_sift(A, i, heap_state)
  // 第 2 段階（右端から順に確定・森の更新）
  for i = n - 1 down to 2
    if order_at_tail(heap_state) < 2 then shrink_trailing_trees(heap_state)
    else
      split_rightmost_into_two_children(heap_state, i)
      for child_root c in newborn_pair_roots()
        interheap_sift(A, c, heap_state)

実装はマーティン・ノブラウチ・レブエルタ氏の smoothsort.c と同系のビット列＋オフセット表現に沿っている（レオナルド数表と sift_in / interheap_sift）。デモでは 浮動穴の連鎖代入 などの微最適化は省略し、配列の内容が変わるたびにバーを更新して追いやすくしている。

実装の複雑さと定数倍の大きさから、汎用ライブラリの sort として採用されることは稀だが、アルゴリズム史と適応的整列の理解には価値がある。