三点中央値型クイックソートを使用する

三点中央値型クイックソート (quick sort with median-of-three pivot) は、部分配列の左端・中央・右端の 3 要素から中央値となるものをピボットに選び、分割の偏りを抑える。

ロムート分割型クイックソートが右端固定だと昇順・降順に近い入力で最悪計算量 O(n²) になりやすいのに対し、三点中央値は整列済みやほぼ整列済みでも端の極端な値を避けやすい。分割そのものはロムート方式を使い、違いは「どの位置をピボットにするか」にある。

ランダムピボット型クイックソートが乱択で期待値を上げるのに対し、こちらは決定的な 3 点比較で同じ目的に近づける。

  1. 三点中央値: 部分配列 lo … hiA[lo]A[mid]A[hi]mid = lo + ⌊(hi - lo) / 2⌋)を比べ、値の中央値にあたる添字を選ぶ。
  2. ロムート分割: 選んだ要素を右端へ移し、片方向走査で「未満」と「以上」に分ける。返り値 p がピボットの最終位置になる。
  3. 再帰: lo … p - 1p + 1 … hi に同じ処理を繰り返す。十分短い区間は挿入ソートで仕上げる。
procedure median_of_three_quick_sort(A, lo, hi)
  if lo >= hi then
    return
  if hi - lo < INSERTION_THRESHOLD then
    insertion_sort(A, lo, hi)
    return
  m = median_of_three_index(A, lo, hi)
  swap(A[m], A[hi])
  p = lomuto_partition(A, lo, hi)
  median_of_three_quick_sort(A, lo, p - 1)
  median_of_three_quick_sort(A, p + 1, hi)

procedure median_of_three_index(A, lo, hi)
  mid = lo + floor((hi - lo) / 2)
  // A[lo], A[mid], A[hi] のうち値の中央値の添字を返す
  ...

procedure lomuto_partition(A, lo, hi)
  pivot = A[hi]
  i = lo
  for j from lo to hi - 1
    if A[j] < pivot then
      swap(A[i], A[j])
      i = i + 1
  swap(A[i], A[hi])
  return i

整列済み配列では中央が中央値になりやすく、右端固定より分割がバランスしやすい。それでも巧妙に構成した入力では最悪計算量は O(n²) のままである。一般に不安定である。

平均計算量は O(n log n)、追加空間は再帰スタックの O(log n) 程度を見込む。

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

ロムート分割型クイックソートは分割手順が同じだが、ピボットを右端固定にする。本アルゴリズムはその右端へ移す前に三点中央値を選ぶ点が違う。

ランダムピボット型クイックソートは一様乱択で期待計算量を上げる。三点中央値は乱数を使わず、左端・中央・右端の比較だけで偏りを抑える。

ホーア分割型クイックソートは両端ポインタの分割で、返り値の意味と再帰区間の切り方が異なる。三点中央値はホーア分割とも組み合わせられるが、ここではロムート分割との組み合わせを示す。

三分割クイックソートは等値帯をその場で確定する。デュアルピボットクイックソートはピボットを 2 つ使い 3 区間に分ける。いずれもピボット選びの三点中央値とは直交する改良である。

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

Size Average time Maximum time Average memory Maximum memory
256 0.000005 0.000032 57 64
512 0.000011 0.000054 62 68
1024 0.000024 0.000122 58 64
2048 0.000052 0.000271 58 64
4096 0.000112 0.000543 62 68
8192 0.000241 0.000898 65 72
16384 0.000520 0.000878 62 68
32768 0.001119 0.001794 58 64
65536 0.002382 0.003276 58 64
131072 0.005057 0.006681 73 80
262144 0.010675 0.013158 62 68
計測に使用したコードを表示する

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::{
    env,
    process::Command,
    time::{Duration, Instant},
};
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;
        }
    }
}

fn partition_at(a: &mut [usize], lo: usize, hi: usize, pivot_idx: usize) -> usize {
    a.swap(pivot_idx, hi);
    let pivot = a[hi];
    let mut i = lo;
    for j in lo..hi {
        if a[j] < pivot {
            a.swap(i, j);
            i += 1;
        }
    }
    a.swap(i, hi);
    i
}



fn quick_median_of_three_idx(a: &[usize], lo: usize, hi: usize) -> usize {
    let mid = lo + (hi - lo) / 2;
    let (x, y, z) = (a[lo], a[mid], a[hi]);
    if (x <= y && y <= z) || (z <= y && y <= x) {
        mid
    } else if (y <= x && x <= z) || (z <= x && x <= y) {
        lo
    } else {
        hi
    }
}

fn quick_median_of_three_sort_range(a: &mut [usize], lo: usize, hi: usize) {
    if hi <= lo {
        return;
    }
    if hi - lo < 16 {
        insertion_sort(&mut a[lo..=hi]);
        return;
    }
    let pivot_idx = quick_median_of_three_idx(a, lo, hi);
    let p = partition_at(a, lo, hi, pivot_idx);
    if p > 0 {
        quick_median_of_three_sort_range(a, lo, p - 1);
    }
    quick_median_of_three_sort_range(a, p + 1, hi);
}

fn quick_median_of_three_sort(a: &mut [usize]) {
    if let Some(hi) = a.len().checked_sub(1) {
        quick_median_of_three_sort_range(a, 0, hi);
    }
}


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

    quick_median_of_three_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 memory_usage_kb() -> usize {
    // VmHWM (peak RSS, KiB). Reported memory subtracts a per-size baseline that only
    // holds the input array, so the table reflects auxiliary space during sorting.
    let contents = std::fs::read_to_string("/proc/self/status")
        .unwrap_or_default();

    for line in contents.lines() {
        if let Some(rest) = line.strip_prefix("VmHWM:") {
            let kb = rest
                .split_whitespace()
                .next()
                .unwrap_or("0")
                .parse::<usize>()
                .unwrap_or(0);

            return kb;
        }
    }

    0
}

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

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

fn run_baseline(size: usize) -> usize {
    let _hold = input_array(size, 1);
    memory_usage_kb()
}

fn run_once(size: usize, seed: usize) -> (u128, usize) {
    let mut array = input_array(size, seed as u64);

    let start = Instant::now();

    benchmark_sort(&mut array);

    let elapsed = start.elapsed();
    let mem = memory_usage_kb();

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

    (micros(elapsed), mem)
}

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

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 == "--baseline-once") {
        run_baseline_child(&args);
        return;
    }
    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 baseline_output = Command::new(env::current_exe().expect("failed to find current executable"))
            .arg("--baseline-once")
            .arg(size.to_string())
            .output()
            .expect("failed to run benchmark baseline process");

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

        let baseline_stdout = String::from_utf8(baseline_output.stdout)
            .expect("baseline process returned non-UTF-8 output");
        let baseline_mem = baseline_stdout
            .split_whitespace()
            .next()
            .expect("missing baseline memory usage")
            .parse::<usize>()
            .expect("invalid baseline memory usage");

        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 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;
            }

            let aux_mem = mem.saturating_sub(baseline_mem);

            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