パワーソートを使用する

パワーソート (powersort) は、すでに整列している連続区間(ラン)を見つけ、隣接ランの中点から求めたパワーに従ってマージ順を決めながらマージソートの枠組みで併合する。

ティムソートと同様のラン検出を使うが、どのラン同士をいつマージするかの方針だけを差し替えた改良版と捉えられる。

Python 3.11 以降の list.sort() では、ラン検出や短いランの拡張などティムソート由来の仕組みを保ったまま、マージ方針がパワーソートに置き換わっている。

  1. ランの検出: 左から昇順または厳密な降順の連続区間を見つける。降順ランは反転して昇順にそろえる。
  2. ランの拡張: 長さが最小ラン長 min_run 未満なら挿入ソートで伸ばす(デモでは見やすさのため 4 に固定)。
  3. パワーの計算: 隣接する 2 ランの「中点位置」から、理想のマージ木上でのノードの深さに相当する整数パワーを求める。
  4. スタックに従ったマージ: 未マージのランをスタックに積み、新しいパワーがスタック先端より小さくなるまで左側のランと現在のランを併合する。
  5. 仕上げ: 入力末尾まで進んだあと、スタックに残ったランをすべてマージして全体を整列する。

ティムソートはスタック上端 3 本の長さ関係を見る経験則でマージ順を決めていた。

パワーソートは 2 ランの中点だけからパワーを 1 つ計算し、ほぼ最適な二分マージ木に沿う順序で併合する。ラン検出や短いランの拡張はティムソートと同じ思想だが、マージ順の決め方だけがパワーという 1 つの整数に集約される点が特徴である。

理論上、既存ラン長 (L₁, …, Lᵣ) に対する適応性は、加法項 O(n) を除き、最適に近い。

procedure node_power(n, b1, e1, b2, e2)
  n1 := e1 - b1
  n2 := e2 - b2
  a := (b1 + n1/2) / n
  b := (b2 + n2/2) / n
  p := 0
  while floor(a · 2^p) = floor(b · 2^p) do
    p := p + 1
  return p

procedure powersort(A)
  S := empty stack of (run, power)
  b1 := 0; e1 := first_run_end(A, 0)
  while e1 < length(A)
    b2 := e1; e2 := first_run_end(A, b2)
    P := node_power(length(A), b1, e1, b2, e2)
    while S is not empty and S.top().power > P
      (b1, e1) := merge(S.pop().run, A[b1..e1))
    S.push((A[b1..e1), P))
    b1 := b2; e1 := e2
  while S is not empty
    (b1, e1) := merge(S.pop().run, A[b1..e1])

ティムソートと同様の適応型だが、マージ順の設計で非効率な併合を避けやすく、安定ソートである。

実装の複雑さはティムソートの 3 本ルールより読み取りやすく、理論的な保証も強い一方、パワー計算やマージ用バッファなどオーバーヘッドは残る。

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

マージソートの下界に近い適応性を目指すしていて既存ラン長の分布に応じて併合順を選ぶ。ポリフェーズマージソートはテープ本数が少ない外部整列向きである。

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

Size Average time Maximum time Average memory Maximum memory
256 0.000009 0.000129 2 2
512 0.000020 0.000093 4 4
1024 0.000043 0.000167 8 8
2048 0.000105 0.000344 16 16
4096 0.000234 0.000535 32 32
8192 0.000452 0.001962 64 64
16384 0.000926 0.001786 128 128
32768 0.001931 0.003279 256 256
65536 0.004336 0.015979 512 512
131072 0.008890 0.015049 1024 1024
262144 0.017345 0.026970 2048 2048
計測に使用したコードを表示する

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 as AtomicOrdering},
    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, AtomicOrdering::Relaxed) + size;
    PEAK_BYTES.fetch_max(live, AtomicOrdering::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(), AtomicOrdering::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(), AtomicOrdering::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;
fn insertion_sort(a: &mut [usize]) {
    for i in 1..a.len() {
        let mut j = i;
        while j > 0 && a[j - 1] > a[j] {
            a.swap(j - 1, j);
            j -= 1;
        }
    }
}



#[derive(Clone, Copy)]
struct PowerRun {
    lo: usize,
    hi: usize,
    power: u32,
}

fn node_power(n: usize, b1: usize, e1: usize, b2: usize, e2: usize) -> u32 {
    let a = (b1 as f64 + (e1 - b1) as f64 / 2.0) / n as f64;
    let b = (b2 as f64 + (e2 - b2) as f64 / 2.0) / n as f64;
    let mut p = 0u32;
    while (a * 2f64.powi(p as i32)).floor() == (b * 2f64.powi(p as i32)).floor() {
        p += 1;
    }
    p
}

fn merge_power_runs(a: &mut [usize], left: PowerRun, right: PowerRun) -> PowerRun {
    let lo = left.lo;
    let hi = right.hi;
    let mid = left.hi + 1;
    let mut merged = Vec::with_capacity(hi - lo + 1);
    let (mut l, mut r) = (left.lo, mid);
    while l <= left.hi && r <= right.hi {
        if a[l] <= a[r] {
            merged.push(a[l]);
            l += 1;
        } else {
            merged.push(a[r]);
            r += 1;
        }
    }
    merged.extend_from_slice(&a[l..=left.hi]);
    merged.extend_from_slice(&a[r..=right.hi]);
    a[lo..=hi].copy_from_slice(&merged);
    PowerRun { lo, hi, power: 0 }
}

fn prepare_power_run(a: &mut [usize], start: usize, min_run: usize) -> usize {
    let n = a.len();
    let mut i = start + 1;
    if i < n && a[i - 1] > a[i] {
        while i < n && a[i - 1] > a[i] {
            i += 1;
        }
        a[start..i].reverse();
    } else {
        while i < n && a[i - 1] <= a[i] {
            i += 1;
        }
    }
    let end = (start + min_run).min(n).max(i);
    insertion_sort(&mut a[start..end]);
    end
}

fn power_sort(a: &mut [usize]) {
    const MIN_RUN: usize = 32;
    let n = a.len();
    if n <= 1 {
        return;
    }
    let mut stack: Vec<PowerRun> = Vec::new();
    let mut b1 = 0usize;
    let mut e1 = prepare_power_run(a, 0, MIN_RUN);
    while e1 < n {
        let b2 = e1;
        let e2 = prepare_power_run(a, b2, MIN_RUN);
        let p = node_power(n, b1, e1, b2, e2);
        while stack.last().is_some_and(|top| top.power > p) {
            let top = stack.pop().unwrap();
            let cur = PowerRun {
                lo: b1,
                hi: e1 - 1,
                power: 0,
            };
            let merged = merge_power_runs(a, top, cur);
            b1 = merged.lo;
            e1 = merged.hi + 1;
        }
        stack.push(PowerRun {
            lo: b1,
            hi: e1 - 1,
            power: p,
        });
        b1 = b2;
        e1 = e2;
    }
    while let Some(top) = stack.pop() {
        let cur = PowerRun {
            lo: b1,
            hi: e1 - 1,
            power: 0,
        };
        let merged = merge_power_runs(a, top, cur);
        b1 = merged.lo;
        e1 = merged.hi + 1;
    }
}


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

    power_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),
        );
    }
    // Small-input cutoffs (insertion sort below 32 elements, etc.) hide duplicate-key
    // bugs in the recursive path, so repeat the duplicate cases at the smallest
    // benchmark size, which every algorithm must handle within reasonable time.
    check_correctness_case("all_equal_len256", vec![7; 256]);
    for seed in 1..=4 {
        check_correctness_case(
            &format!("few_keys_len256_seed_{seed}"),
            few_unique_values(256, 4, seed),
        );
    }
    // Blit's equal-key second sweep used to copy the whole range into a fixed
    // 512-element swap; lengths above that must still sort without panicking.
    check_correctness_case("all_equal_len600", vec![7; 600]);
    for seed in 1..=4 {
        check_correctness_case(
            &format!("few_keys_len2048_seed_{seed}"),
            few_unique_values(2048, 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 bytes (explicit buffers such as swap).
/// Kept in bytes so the parent can average before rounding; converting to KiB here
/// would truncate sub-KiB buffers to 0 in every run and hide them from the average.
fn run_once(size: usize, seed: usize) -> (u128, usize) {
    let mut array = input_array(size, seed as u64);

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

    let start = Instant::now();

    benchmark_sort(&mut array);

    let elapsed = start.elapsed();
    let peak_bytes = PEAK_BYTES.load(AtomicOrdering::Relaxed);
    let aux_bytes = peak_bytes.saturating_sub(base_bytes);

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

    (micros(elapsed), aux_bytes)
}

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;
        // Memory is summed in bytes and converted to KiB once, after averaging.
        let avg_mem_kb = total_mem / RUNS / 1024;
        let max_mem_kb = max_mem / 1024;

        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_kb,
            max_mem_kb
        );
    }
}
RUST

RUN cargo build --release

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

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