左傾ヒープソートを使用する

左傾ヒープソート (leftist heap sort) は、要素を左傾ヒープへ挿入したあと、最小値を繰り返し取り出して昇順にする整列である。

左傾ヒープは、根が最小(または最大)となる二分木である。ヒープ条件に加え、各節点の ヌルパス長null path length, NPL)について左子の方が右子以上になるよう子を並べ替える(左傾性)。 NPL は「その節点から右の子だけを辿って最初の欠損に至るまでの辺数」で、葉は 0、空木は -1 とおく。右背骨の長さは O(log n) に抑えられ、合併が右背骨に沿って進むため速い。

  1. 合併: 2 本のヒープの根を比較し、キーの大きい方を小さい方の右部分木と再帰的に合併する。 終わったら左右の NPL を見て、左傾性が崩れていれば左右を入れ替え、根の NPL を「右子の NPL + 1」に更新する。最悪 O(log n)
  2. 挿入: 単一節点のヒープを既存ヒープと合併する。
  3. 抽出: 根を外し、左右の子ヒープを合併して新しい根を得る。最悪 O(log n)
  4. 書き戻し: 取り出したキーを配列の先頭から順に書けば昇順になる。
procedure npl(H)
  if H is empty
    return -1
  return H.npl

procedure merge(H1, H2)
  if H1 is empty
    return H2
  if H2 is empty
    return H1
  if H1.key > H2.key
    swap H1, H2
  H1.right = merge(H1.right, H2)
  if npl(H1.left) < npl(H1.right)
    swap H1.left, H1.right
  H1.npl = npl(H1.right) + 1
  return H1

procedure leftist_heap_sort(A)
  H = empty leftist heap
  for x in A
    H = merge(H, singleton(x))
  for i from 0 to length(A) - 1
    A[i] = H.key
    H = merge(H.left, H.right)

最悪時間計算量は O(n log n) であり、節点用に O(n) の追加記憶域が要る(インプレースではない)。等値キーの相対順序は合併時の規約に依存し、一般に不安定である。実装は二分木の合併として素直だが、ポインタ経由の木はキャッシュ効率では配列上のヒープソートに劣りやすい。

優先度付きキューとしての左傾ヒープは、合併を右背骨に沿って書ける点が二項ヒープより単純で、フィボナッチヒープほど複雑な遅延操作も要らない。整列用途ではその操作を「すべて挿入してからすべて取り出す」形に固定したものが左傾ヒープソートである。

類似アルゴリズムとの相違点

ヒープソートは配列上の二分ヒープをインプレースで縮める。左傾ヒープはポインタの二分木で、合併と左傾性の修復が中心になる。

二項ヒープソートは次数の異なる二項木を二進加算のように結合する。左傾ヒープは単一の二分木を保ち、NPL で右背骨の高さを抑える。

ペアリングヒープソートは多分岐木の合併と子の二パス・ペアリングを使う。左傾ヒープは常に二分木で、合併は右背骨の再帰である。

弱ヒープソートは配列上の不完全木と逆ビットで比較回数を抑える。ヒープ同士の合併を第一級には扱わない。

計算時間量および空間計算量を計測する

Size Average time Maximum time Average memory Maximum memory
256 0.000018 0.000194 8 8
512 0.000040 0.000261 16 16
1024 0.000086 0.000288 32 32
2048 0.000189 0.000407 64 64
4096 0.000424 0.001954 128 128
8192 0.000960 0.001728 256 256
16384 0.002188 0.004864 512 512
32768 0.005174 0.121067 1024 1024
65536 0.012176 0.046988 2048 2048
131072 0.032271 0.094475 4096 4096
262144 0.087888 0.504213 8192 8192
計測に使用したコードを表示する

set -euo pipefail

WORKDIR="$(mktemp -d)"
trap 'rm -rf "$WORKDIR"' EXIT

cat > "$WORKDIR/Dockerfile" <<'EOF'
FROM rust:1.95.0

WORKDIR /app

RUN mkdir -p src

RUN cat > Cargo.toml <<'CARGO'
[package]
name = "rust-benchmark"
version = "0.1.0"
edition = "2021"

[profile.release]
lto = true
codegen-units = 1
panic = "abort"
CARGO

RUN cat > src/main.rs <<'RUST'
use std::{
    alloc::{GlobalAlloc, Layout, System},
    env,
    process::Command,
    sync::atomic::{AtomicUsize, Ordering},
    time::{Duration, Instant},
};

/// Counts live heap bytes and the high-water mark so auxiliary sort buffers
/// (swap Vecs, etc.) are measured as explicit heap growth during the sort.
struct TrackingAllocator;

static LIVE_BYTES: AtomicUsize = AtomicUsize::new(0);
static PEAK_BYTES: AtomicUsize = AtomicUsize::new(0);

fn record_alloc(size: usize) {
    let live = LIVE_BYTES.fetch_add(size, Ordering::Relaxed) + size;
    PEAK_BYTES.fetch_max(live, Ordering::Relaxed);
}

unsafe impl GlobalAlloc for TrackingAllocator {
    unsafe fn alloc(&self, layout: Layout) -> *mut u8 {
        let ptr = System.alloc(layout);
        if !ptr.is_null() {
            record_alloc(layout.size());
        }
        ptr
    }

    unsafe fn dealloc(&self, ptr: *mut u8, layout: Layout) {
        LIVE_BYTES.fetch_sub(layout.size(), Ordering::Relaxed);
        System.dealloc(ptr, layout);
    }

    unsafe fn alloc_zeroed(&self, layout: Layout) -> *mut u8 {
        let ptr = System.alloc_zeroed(layout);
        if !ptr.is_null() {
            record_alloc(layout.size());
        }
        ptr
    }

    unsafe fn realloc(&self, ptr: *mut u8, layout: Layout, new_size: usize) -> *mut u8 {
        let new_ptr = System.realloc(ptr, layout, new_size);
        if !new_ptr.is_null() {
            LIVE_BYTES.fetch_sub(layout.size(), Ordering::Relaxed);
            record_alloc(new_size);
        }
        new_ptr
    }
}

#[global_allocator]
static GLOBAL: TrackingAllocator = TrackingAllocator;
const MIN_POWER: u32 = 8;
const MAX_POWER: u32 = 18;
const RUNS: usize = 8192;


struct LeftistNode {
    key: usize,
    npl: i32,
    left: Option<Box<LeftistNode>>,
    right: Option<Box<LeftistNode>>,
}

fn npl(node: &Option<Box<LeftistNode>>) -> i32 {
    match node {
        None => -1,
        Some(n) => n.npl,
    }
}

fn merge(
    a: Option<Box<LeftistNode>>,
    b: Option<Box<LeftistNode>>,
) -> Option<Box<LeftistNode>> {
    match (a, b) {
        (None, x) | (x, None) => x,
        (Some(mut x), Some(mut y)) => {
            if x.key > y.key {
                std::mem::swap(&mut x, &mut y);
            }
            x.right = merge(x.right.take(), Some(y));
            if npl(&x.left) < npl(&x.right) {
                std::mem::swap(&mut x.left, &mut x.right);
            }
            x.npl = npl(&x.right) + 1;
            Some(x)
        }
    }
}

fn insert_key(heap: Option<Box<LeftistNode>>, key: usize) -> Option<Box<LeftistNode>> {
    let node = Box::new(LeftistNode {
        key,
        npl: 0,
        left: None,
        right: None,
    });
    merge(heap, Some(node))
}

fn extract_min(heap: Option<Box<LeftistNode>>) -> (Option<usize>, Option<Box<LeftistNode>>) {
    let Some(root) = heap else {
        return (None, None);
    };
    let key = root.key;
    let left = root.left;
    let right = root.right;
    (Some(key), merge(left, right))
}

fn leftist_heap_sort(a: &mut [usize]) {
    let mut heap = None;
    for &key in a.iter() {
        heap = insert_key(heap, key);
    }
    for slot in a.iter_mut() {
        let (key, next) = extract_min(heap);
        heap = next;
        *slot = key.expect("leftist heap exhausted early");
    }
}


fn benchmark_sort(array: &mut [usize]) {

    leftist_heap_sort(array);

}

fn is_non_decreasing(a: &[usize]) -> bool {
    a.windows(2).all(|w| w[0] <= w[1])
}

fn same_multiset(a: &[usize], b: &[usize]) -> bool {
    if a.len() != b.len() {
        return false;
    }

    let mut left = a.to_vec();
    let mut right = b.to_vec();
    left.sort_unstable();
    right.sort_unstable();
    left == right
}

fn check_correctness_case(label: &str, mut input: Vec<usize>) {
    let original = input.clone();

    benchmark_sort(&mut input);

    if !is_non_decreasing(&input) {
        panic!("correctness case {}: output is not sorted", label);
    }

    if !same_multiset(&input, &original) {
        panic!("correctness case {}: elements were lost or added", label);
    }
}

fn few_unique_values(size: usize, unique: usize, seed: u64) -> Vec<usize> {
    let mut state = seed;

    (0..size)
        .map(|_| {
            state ^= state << 13;
            state ^= state >> 7;
            state ^= state << 17;
            (state as usize % unique) + 1
        })
        .collect()
}

fn run_correctness_checks() {
    check_correctness_case("empty", vec![]);
    check_correctness_case("single", vec![42]);
    check_correctness_case("duplicates", vec![3, 1, 3, 2, 1, 2]);
    check_correctness_case("sorted", vec![1, 2, 3, 4, 5]);
    check_correctness_case("reverse", vec![5, 4, 3, 2, 1]);
    check_correctness_case("all_equal", vec![7, 7, 7, 7]);
    check_correctness_case("skewed_range", vec![1_000_000, 2, 1_000_001, 1, 999_999]);
    // Static-buffer Grail skips the in-buffer build when key collection is sparse
    // (ideal_buffer = false). Exercising that path catches regressions in buffer gating.
    check_correctness_case(
        "few_keys_len16",
        vec![2, 2, 2, 2, 2, 2, 2, 2, 4, 3, 1, 2, 3, 4, 1, 4],
    );
    // Seed 0 is a fixed point of the xorshift below, so it would degenerate into
    // yet another all-equal case instead of a 4-value mix. Start at 1.
    for seed in 1..=32 {
        check_correctness_case(
            &format!("few_keys_len32_seed_{seed}"),
            few_unique_values(32, 4, seed),
        );
    }
}


fn shuffled(size: usize, seed: u64) -> Vec<usize> {
    let mut v: Vec<usize> = (1..=size).collect();

    let mut state = seed;

    for i in (1..size).rev() {
        state ^= state << 13;
        state ^= state >> 7;
        state ^= state << 17;

        let j = (state as usize) % (i + 1);

        v.swap(i, j);
    }

    v
}

fn micros(d: Duration) -> u128 {
    d.as_micros()
}

fn input_array(size: usize, seed: u64) -> Vec<usize> {
    shuffled(size, seed)
}

/// Peak heap growth during `benchmark_sort`, in KiB (explicit buffers such as swap).
fn run_once(size: usize, seed: usize) -> (u128, usize) {
    let mut array = input_array(size, seed as u64);

    let base_bytes = LIVE_BYTES.load(Ordering::Relaxed);
    PEAK_BYTES.store(base_bytes, Ordering::Relaxed);

    let start = Instant::now();

    benchmark_sort(&mut array);

    let elapsed = start.elapsed();
    let peak_bytes = PEAK_BYTES.load(Ordering::Relaxed);
    let aux_kb = peak_bytes.saturating_sub(base_bytes) / 1024;

    let expected: Vec<usize> = (1..=size).collect();
    if array != expected {
        panic!(
            "sort failed with seed {} for size {}",
            seed,
            size
        );
    }

    (micros(elapsed), aux_kb)
}

fn run_child(args: &[String]) {
    let size = args[2].parse::<usize>().expect("invalid size");
    let seed = args[3].parse::<usize>().expect("invalid seed");
    let (elapsed_us, mem) = run_once(size, seed);
    println!("{} {}", elapsed_us, mem);
}

fn main() {
    let args: Vec<String> = env::args().collect();
    if args.get(1).is_some_and(|arg| arg == "--run-once") {
        run_child(&args);
        return;
    }

    run_correctness_checks();

    println!(
        "| {:>10} | {:>15} | {:>15} | {:>15} | {:>15} |",
        "Size",
        "Average time",
        "Maximum time",
        "Average memory",
        "Maximum memory"
    );

    println!(
        "|{:-<11}:|{:-<16}:|{:-<16}:|{:-<16}:|{:-<16}:|",
        "",
        "",
        "",
        "",
        ""
    );

    for power in MIN_POWER..=MAX_POWER {
        let size = 1usize << power;

        let mut total_time: u128 = 0;
        let mut max_time: u128 = 0;

        let mut total_mem: usize = 0;
        let mut max_mem: usize = 0;

        for seed in 1..=RUNS {
            let output = Command::new(env::current_exe().expect("failed to find current executable"))
                .arg("--run-once")
                .arg(size.to_string())
                .arg(seed.to_string())
                .output()
                .expect("failed to run benchmark child process");

            if !output.status.success() {
                panic!(
                    "benchmark child process failed: {}",
                    String::from_utf8_lossy(&output.stderr)
                );
            }

            let stdout = String::from_utf8(output.stdout)
                .expect("child process returned non-UTF-8 output");
            let mut fields = stdout.split_whitespace();
            let elapsed_us = fields
                .next()
                .expect("missing elapsed time")
                .parse::<u128>()
                .expect("invalid elapsed time");
            let aux_mem = fields
                .next()
                .expect("missing memory usage")
                .parse::<usize>()
                .expect("invalid memory usage");

            total_time += elapsed_us;

            if elapsed_us > max_time {
                max_time = elapsed_us;
            }

            total_mem += aux_mem;

            if aux_mem > max_mem {
                max_mem = aux_mem;
            }
        }

        let avg_time = total_time / RUNS as u128;
        let avg_mem = total_mem / RUNS;

        println!(
            "| {:>10} | {:>15} | {:>15} | {:>15} | {:>15} |",
            size,
            format!("{}.{:06}", avg_time / 1_000_000, avg_time % 1_000_000),
            format!("{}.{:06}", max_time / 1_000_000, max_time % 1_000_000),
            avg_mem,
            max_mem
        );
    }
}
RUST

RUN cargo build --release

CMD ["./target/release/rust-benchmark"]
EOF

docker build -t rust-benchmark "$WORKDIR"
docker run --rm --init rust-benchmark