フラックスソートを使用する

フラックスソート (fluxsort) は、クイックソート型の上方向き分割と、マージソート型の補助メモリ上での仕上げを組み合わせた、安定かつ適応的な比較ソートである。

分岐予測ミスを抑える分割や、整列済み度合いを見るアナライザが特徴である。本記事の実装は簡略化しており、小区間の quadsort の代わりに挿入ソート/マージソートを用いる。

  1. アナライザ: 全体が昇順なら何もしない。厳密な降順なら反転して終了する。4 分割した各区間の隣接昇順ペアが半数超なら、その区間をマージソートで仕上げる。
  2. ピボット選択: 区間をほぼ等間隔に 9 点取り、3 組の三点中央値の中央値(準中央値)をピボットにする。
  3. 安定な二重書き込み分割: 要素を先頭から走査し、≤ ピボット を配列前方へ、> ピボット を後方へ、出現順を保ったまま書き分ける(本番の fluxsort では後者を swap 領域へ書く)。
  4. 等値の第二走査: 右側が空(すべて ≤ ピボット)なら、< ピボット= ピボット に分け直し、等値帯を再帰から外す。
  5. 不均衡フォールバック: 左右の長さ比が 1:16 を超えて偏ったら、両側をマージソートする。
  6. 小区間: 要素数が閾値未満なら挿入ソートで仕上げる。
procedure flux_stable_partition(A, pivot)
  L := empty; R := empty
  for each x in A in order
    if x ≤ pivot then append x to L else append x to R
  A := L concatenated with R
  return length(L)

procedure flux_sort_range(A)
  if length(A) < INSERTION_THRESHOLD then
    insertion_sort(A); return
  pivot := quasimedian_of_9(A)
  left := flux_stable_partition(A, pivot)
  right := length(A) - left
  if right = 0 then
    move keys < pivot to front, equals after them
    flux_sort_range(A[0 .. lt))
    return
  if left < length(A)/16 or right < length(A)/16 then
    merge_sort(A[0 .. left)); merge_sort(A[left .. end)); return
  flux_sort_range(A[0 .. left))
  flux_sort_range(A[left .. end))

procedure fluxsort(A)
  if A is sorted then return
  if A is reverse-sorted then reverse(A); return
  for each quarter Q of A
    if ordered_pairs(Q) > half then merge_sort(Q)
  if A is sorted then return
  flux_sort_range(A)

最良は整列済み検出により O(n)、平均・最悪は O(n log n) を狙う。補助配列に最大 O(n) を使う安定ソートである。

以下のデモでは視認性のため挿入閾値を 4、不均衡判定を 1/4 に緩めている。

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

ロムート分割型クイックソートはインプレースで不安定な分割を行う。フラックスソートは補助メモリ(または同等の二重書き込み)で安定分割し、分岐の少ない走査を前提にしている。

パターン撃退クイックソートも悪パターン対策のハイブリッドだが、不安定で補助メモリをほぼ使わない。フラックスソートは安定性を保ったまま分割とマージを行き来する。

マージソートティムソートは常にマージ中心である。フラックスソートはランダム寄りの入力では分割を主とし、整列度が高い区間や不均衡時だけマージ側へ寄せる。

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

Size Average time Maximum time Average memory Maximum memory
256 0.000008 0.000074 62 68
512 0.000018 0.000057 62 68
1024 0.000036 0.000096 70 76
2048 0.000074 0.000169 78 84
4096 0.000158 0.001303 94 100
8192 0.000337 0.001851 141 148
16384 0.000731 0.004769 58 64
32768 0.001557 0.003763 108 148
65536 0.003316 0.005650 364 404
131072 0.007068 0.011157 876 916
262144 0.015100 0.029708 1904 2052
計測に使用したコードを表示する

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



const FLUX_INSERTION_THRESHOLD: usize = 24;

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

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

fn flux_reverse(a: &mut [usize]) {
    let mut lo = 0;
    let mut hi = a.len();
    while lo + 1 < hi {
        hi -= 1;
        a.swap(lo, hi);
        lo += 1;
    }
}

/// Count ascending adjacent pairs (presortedness measure).
fn flux_ordered_pairs(a: &[usize]) -> usize {
    a.windows(2).filter(|w| w[0] <= w[1]).count()
}

fn flux_median3_idx(a: &[usize], i: usize, j: usize, k: usize) -> usize {
    let (x, y, z) = (a[i], a[j], a[k]);
    if x < y {
        if y < z {
            j
        } else if x < z {
            k
        } else {
            i
        }
    } else if x < z {
        i
    } else if y < z {
        k
    } else {
        j
    }
}

/// Quasimedian of 9: median of three medians-of-three sampled across the range.
fn flux_quasimedian9(a: &[usize]) -> usize {
    let n = a.len();
    if n < 9 {
        return a[n / 2];
    }
    let step = n / 8;
    let i0 = 0;
    let i1 = step;
    let i2 = step * 2;
    let i3 = step * 3;
    let i4 = step * 4;
    let i5 = step * 5;
    let i6 = step * 6;
    let i7 = step * 7;
    let i8 = n - 1;
    let m0 = flux_median3_idx(a, i0, i1, i2);
    let m1 = flux_median3_idx(a, i3, i4, i5);
    let m2 = flux_median3_idx(a, i6, i7, i8);
    a[flux_median3_idx(a, m0, m1, m2)]
}

fn flux_merge(a: &mut [usize], swap: &mut [usize]) {
    let n = a.len();
    if n <= 1 {
        return;
    }
    let mid = n / 2;
    flux_merge(&mut a[..mid], swap);
    flux_merge(&mut a[mid..], swap);
    let (left, right) = a.split_at(mid);
    let mut i = 0;
    let mut j = 0;
    let mut k = 0;
    while i < left.len() && j < right.len() {
        if left[i] <= right[j] {
            swap[k] = left[i];
            i += 1;
        } else {
            swap[k] = right[j];
            j += 1;
        }
        k += 1;
    }
    while i < left.len() {
        swap[k] = left[i];
        i += 1;
        k += 1;
    }
    while j < right.len() {
        swap[k] = right[j];
        j += 1;
        k += 1;
    }
    a.copy_from_slice(&swap[..n]);
}

/// Stable dual-destination partition: `≤ pivot` stay toward the front of `a`,
/// `> pivot` are collected after them (educational stand-in for fluxsort’s
/// main/swap split). Returns the length of the left partition.
fn flux_stable_partition(a: &mut [usize], swap: &mut [usize], pivot: usize) -> usize {
    let n = a.len();
    swap[..n].copy_from_slice(a);
    let mut left = 0usize;
    for i in 0..n {
        if swap[i] <= pivot {
            left += 1;
        }
    }
    let mut l = 0usize;
    let mut r = left;
    for i in 0..n {
        let x = swap[i];
        if x <= pivot {
            a[l] = x;
            l += 1;
        } else {
            a[r] = x;
            r += 1;
        }
    }
    left
}

fn flux_partition_sort(a: &mut [usize], swap: &mut [usize]) {
    let n = a.len();
    if n <= 1 {
        return;
    }
    if n < FLUX_INSERTION_THRESHOLD {
        insertion_sort(a);
        return;
    }

    let pivot = flux_quasimedian9(a);
    let left = flux_stable_partition(a, swap, pivot);
    let right = n - left;

    // All keys ≤ pivot: filter equals out so recursion makes progress
    // (fluxsort’s “second sweep” for generic / low-cardinality data).
    if right == 0 {
        swap[..n].copy_from_slice(a);
        let mut lt = 0usize;
        for i in 0..n {
            if swap[i] < pivot {
                a[lt] = swap[i];
                lt += 1;
            }
        }
        let mut eq = lt;
        for i in 0..n {
            if swap[i] == pivot {
                a[eq] = swap[i];
                eq += 1;
            }
        }
        if lt > 1 {
            flux_partition_sort(&mut a[..lt], swap);
        }
        return;
    }

    // Worst-case guard: one side < 1/16 of the other → mergesort both sides.
    let unbalanced = left > 0 && (left < n / 16 || right < n / 16);

    if unbalanced {
        flux_merge(&mut a[..left], swap);
        flux_merge(&mut a[left..], swap);
        return;
    }

    if left > 1 {
        flux_partition_sort(&mut a[..left], swap);
    }
    if right > 1 {
        flux_partition_sort(&mut a[left..], swap);
    }
}

fn flux_analyze(a: &mut [usize], swap: &mut [usize]) -> bool {
    let n = a.len();
    if n <= 1 {
        return true;
    }
    if flux_is_sorted(a) {
        return true;
    }
    if flux_is_reverse_sorted(a) {
        flux_reverse(a);
        return true;
    }

    // Four-segment presortedness: if more than half the adjacent pairs in a
    // segment are ordered, finish that segment with mergesort (stand-in for
    // quadsort). Remaining disorder is handled by partitioning afterward.
    let q = n / 4;
    if q >= 2 {
        let bounds = [0, q, q * 2, q * 3, n];
        for s in 0..4 {
            let lo = bounds[s];
            let hi = bounds[s + 1];
            if hi - lo < 2 {
                continue;
            }
            let pairs = hi - lo - 1;
            if flux_ordered_pairs(&a[lo..hi]) * 2 > pairs {
                flux_merge(&mut a[lo..hi], swap);
            }
        }
        if flux_is_sorted(a) {
            return true;
        }
    }
    false
}

fn flux_sort(a: &mut [usize]) {
    let n = a.len();
    if n <= 1 {
        return;
    }
    let mut swap = vec![0usize; n];
    if flux_analyze(a, &mut swap) {
        return;
    }
    flux_partition_sort(a, &mut swap);
}


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

    flux_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