振動ソートを使用する

振動ソート (oscillating merge sort / oscillating sort) は、後退読み取り可能なテープ向けに提案された外部マージ整列である。入力の分配(ラン生成)とマージを交互に進め、バランスマージのように全ランを先に配り切ってからマージパスへ入るのではなく、途中でマージを挟みながら大きなランを組み立てる。

テープ本数を \(n\)(入力 1 本 + 作業 \(n - 1\) 本)とすると、各マージはおおむね \(n - 2\) 方向になる。本稿のデモとベンチマークでは \(n = 5\)、すなわち 3 方向マージ(WAY = 3)と、作業テープ 4 本分に相当するレベル管理を採用する。

物理テープの代わりに「マージ段(レベル)」ごとのラン列をベクタで持ち、同レベルに WAY 本たまったら次レベルへ 1 本まとめる。バッチ長は \(3^0, 3^1, 3^2, \ldots\) のように累乗で伸びる。

  1. 初期ラン生成: 配列を固定長(デモは 3、計測は 32)の区間に区切り、各区間を内部整列してレベル 0 のランとする。
  2. 分配と振動: 入力ランを最大 WAY 本ずつレベル 0 へ置き、そのたびに「同レベルに WAY 本あるか」を見てマージする。分配とマージを交互に繰り返すのが振動の由来である。
  3. k 方向マージ: WAY 本の昇順ランの先頭を比較し、最小(同値ならより左のラン)を出力へ確定する。結果は 1 段上のレベルへ 1 ランとして積む。
  4. 仕上げ: 入力が尽きたあと、残った各レベルのランを同じ k 方向マージで 1 本にまとめて配列へ書き戻す。

後退読み取りそのものはメモリ上のシミュレーションでは省略する。テープ実装では、直前に書き出した昇順ランを後ろから読むことで巻き戻しを避けつつ次のマージに渡す、という点が歴史的な利点だった。

procedure merge_k_way(runs[0..k))
  heads[i] = 0 for each run i
  while some run still has unread elements
    pick run i with smallest heads[i] value
      (ties: smallest i, for stability)
    append runs[i][heads[i]] to output
    heads[i] = heads[i] + 1
  return output

procedure collapse(by_level, lv)
  while length(by_level[lv]) >= WAY
    batch = take WAY runs from by_level[lv]
    append merge_k_way(batch) to by_level[lv + 1]
    collapse(by_level, lv + 1)

procedure oscillating_merge_sort(A)
  pending = create_runs(A, run_size)
  by_level[0] = empty list
  while pending is not empty
    for up to WAY runs from pending
      append run to by_level[0]
      collapse(by_level, 0)
  leftover = concatenate all by_level[*]
  while length(leftover) > 1
    k = min(WAY, length(leftover))
    append merge_k_way(leftover[0 .. k)) to leftover; drop those k
  copy leftover[0] back into A

パス数はおおよそ \(\log_{n-2}(N / r)\)(\(N\) は要素数、\(r\) は初期ラン長)なので、全体の時間は \(O(N \log N)\) である。作業領域はランとマージバッファに依存し、テープ実装では追加メモリはほぼ定数、本稿のメモリ実装では \(O(N)\) のバッファを使う。同値を左ラン優先で取れば安定ソートになる。

デモでは上段をレベル0、下段を常設のレベル1作業列とする。マージで選んだ値は下段へ1本ずつ移し、マージが終わったらその結果を上段へ1本ずつ戻す。後退読み取り可能なテープでは、分配とマージを交互に進めることでドライブの空き時間と巻き戻しを抑えられる、というのが振動ソートの歴史的な狙いである。

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

マージソート多方向マージソートは、通常すべての分割(または初期ラン)を用意してから段階的にマージする。振動ソートは「少し分配してはマージする」点で位相が違う。

ポリフェーズマージソートもテープ本数が少ない外部整列向けだが、フィボナッチ分布で全ランを先に配置してからパスを回す。振動ソートは分布を完了させず、\(n - 2\) 方向のバッチを累乗で積み上げる。

本稿のレベル管理は、テープ上の「後ろから読む」詳細を省略した教育用モデルである。実テープでは読み書き方向の反転そのものが巻き戻し削減の本体になる。

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

Size Average time (s) Maximum time (s) Average memory (KiB) Maximum memory (KiB)
256 0.000052 0.000473 7 7
512 0.000072 0.012272 15 15
1024 0.000114 0.000665 27 27
2048 0.000221 0.001313 61 61
4096 0.000421 0.007946 135 135
8192 0.000852 0.029075 209 209
16384 0.001841 0.015629 544 544
32768 0.003724 0.017902 1089 1089
65536 0.006871 0.029669 2172 2172
131072 0.014895 0.125831 4330 4330
262144 0.023055 0.332180 8655 8655
計測に使用したコードを表示する

#!/usr/bin/env swift
import Foundation

// This standalone Swift driver creates the same temporary Docker build
// context as the former shell wrapper.  The benchmark program itself remains
// embedded below so readers can copy one complete, reproducible file.
struct BenchmarkError: Error, CustomStringConvertible {
    let message: String

    var description: String { message }

    init(_ message: String) {
        self.message = message
    }
}

func runCommand(_ executable: String, _ arguments: [String]) throws {
    let process = Process()
    process.executableURL = URL(fileURLWithPath: "/usr/bin/env")
    process.arguments = [executable] + arguments
    process.standardInput = FileHandle.standardInput
    process.standardOutput = FileHandle.standardOutput
    process.standardError = FileHandle.standardError

    do {
        try process.run()
    } catch {
        throw BenchmarkError("Could not start \(executable): \(error)")
    }
    process.waitUntilExit()
    guard process.terminationStatus == 0 else {
        throw BenchmarkError(
            "Command failed (\(process.terminationStatus)): " +
            "\(executable) \(arguments.joined(separator: " "))"
        )
    }
}

do {
    // The UUID avoids collisions when two benchmark copies are run at once.
    let workdir = FileManager.default.temporaryDirectory
        .appendingPathComponent("swift-sort-benchmark-\(UUID().uuidString)")
    try FileManager.default.createDirectory(at: workdir, withIntermediateDirectories: true)
    defer { try? FileManager.default.removeItem(at: workdir) }

    // A raw Swift string is used so the nested main.swift keeps its own
    // interpolation expressions such as \(seed) until Docker compiles it.
    let dockerfile = #"""
FROM swift:6.0

WORKDIR /app

RUN cat > alloc_track.c <<'ALLOC'
#define _GNU_SOURCE
#include <dlfcn.h>
#include <malloc.h>
#include <stdatomic.h>
#include <stddef.h>
#include <stdint.h>
#include <stdlib.h>
#include <string.h>

static atomic_size_t live_bytes = 0;
static atomic_size_t peak_bytes = 0;

static void *(*real_malloc)(size_t) = NULL;
static void *(*real_calloc)(size_t, size_t) = NULL;
static void *(*real_realloc)(void *, size_t) = NULL;
static void (*real_free)(void *) = NULL;

static void init_reals(void) {
    if (real_malloc) {
        return;
    }
    real_malloc = (void *(*)(size_t))dlsym(RTLD_NEXT, "malloc");
    real_calloc = (void *(*)(size_t, size_t))dlsym(RTLD_NEXT, "calloc");
    real_realloc = (void *(*)(void *, size_t))dlsym(RTLD_NEXT, "realloc");
    real_free = (void (*)(void *))dlsym(RTLD_NEXT, "free");
}

static void record_alloc(size_t size) {
    size_t live = atomic_fetch_add(&live_bytes, size) + size;
    size_t peak = atomic_load(&peak_bytes);
    while (live > peak) {
        if (atomic_compare_exchange_weak(&peak_bytes, &peak, live)) {
            break;
        }
    }
}

void alloc_track_reset_peak(void) {
    atomic_store(&peak_bytes, atomic_load(&live_bytes));
}

size_t alloc_track_live(void) { return atomic_load(&live_bytes); }
size_t alloc_track_peak(void) { return atomic_load(&peak_bytes); }

void *malloc(size_t size) {
    init_reals();
    void *p = real_malloc(size);
    if (p) {
        record_alloc(malloc_usable_size(p));
    }
    return p;
}

void *calloc(size_t nmemb, size_t size) {
    init_reals();
    void *p = real_calloc(nmemb, size);
    if (p) {
        record_alloc(malloc_usable_size(p));
    }
    return p;
}

void *realloc(void *ptr, size_t size) {
    init_reals();
    size_t old_size = 0;
    if (ptr) {
        old_size = malloc_usable_size(ptr);
    }
    void *p = real_realloc(ptr, size);
    if (p) {
        atomic_fetch_sub(&live_bytes, old_size);
        record_alloc(malloc_usable_size(p));
    } else if (size == 0) {
        atomic_fetch_sub(&live_bytes, old_size);
    }
    return p;
}

void free(void *ptr) {
    init_reals();
    if (ptr) {
        atomic_fetch_sub(&live_bytes, malloc_usable_size(ptr));
        real_free(ptr);
    }
}

ALLOC

RUN cat > main.swift <<'SWIFT'
import Foundation
#if canImport(Glibc)
import Glibc
#elseif canImport(Darwin)
import Darwin
#endif

@_silgen_name("alloc_track_live") func alloc_track_live() -> Int
@_silgen_name("alloc_track_peak") func alloc_track_peak() -> Int
@_silgen_name("alloc_track_reset_peak") func alloc_track_reset_peak()

extension UnsafeMutableBufferPointer where Element == Int {
    func swapAt(_ i: Int, _ j: Int) {
        let t = self[i]; self[i] = self[j]; self[j] = t
    }
}

let MIN_POWER: Int = 8
let MAX_POWER: Int = 18
let RUNS: Int = 8192


/// Oscillating merge sort (Sobel 1962), in-memory pedagogical simulation.
///
/// Classic tape form uses n drives (1 input + n−1 work) and an (n−2)-way merge.
/// Here n = 5 ⇒ WAY = 3 work inputs per merge and NUM_TAPES = WAY + 1 = 4 work
/// tapes in the narrative; we store runs by merge level instead of physical tapes.
/// Distribution of initial runs and merging are interleaved: after every WAY new
/// level-0 runs, collapse any level that has WAY pending runs into the next level
/// (batches grow as powers of WAY). Remaining runs are flushed with the same
/// k-way merge at the end.
fileprivate let RUN_SIZE = 32
fileprivate let WAY = 3

fileprivate func merge_k_way(_ runs: [[Int]]) -> [Int] {
    let k = runs.count
    var heads = [Int](repeating: 0, count: k)
    let total = runs.reduce(0) { $0 + $1.count }
    var out = [Int]()
    out.reserveCapacity(total)

    while true {
        var best: (Int, Int)? = nil
        for i in 0..<k {
            if heads[i] < runs[i].count {
                let v = runs[i][heads[i]]
                if let (bi, bv) = best {
                    if v < bv || (v == bv && i < bi) {
                        best = (i, v)
                    }
                } else {
                    best = (i, v)
                }
            }
        }
        guard let (i, v) = best else {
            break
        }
        out.append(v)
        heads[i] += 1
    }
    return out
}

fileprivate func create_runs(_ a: UnsafeMutableBufferPointer<Int>, _ run_size: Int) -> [[Int]] {
    var runs: [[Int]] = []
    var i = 0
    while i < a.count {
        let end = min(i + run_size, a.count)
        var run = Array(a[i..<end])
        run.sort()
        runs.append(run)
        i = end
    }
    return runs
}

fileprivate func collapse_from(_ by_level: inout [[ [Int] ]], _ level: Int) {
    while level < by_level.count && by_level[level].count >= WAY {
        var batch: [[Int]] = []
        batch.reserveCapacity(WAY)
        for _ in 0..<WAY {
            batch.append(by_level[level].removeFirst())
        }
        let merged = merge_k_way(batch)
        let next = level + 1
        while by_level.count <= next {
            by_level.append([])
        }
        by_level[next].append(merged)
        collapse_from(&by_level, next)
    }
}

func oscillating_merge_sort(_ a: inout [Int]) {
    a.withUnsafeMutableBufferPointer { oscillating_merge_sort($0) }
}

func oscillating_merge_sort(_ a: UnsafeMutableBufferPointer<Int>) {
    if a.count <= 1 {
        return
    }

    let pending = create_runs(a, RUN_SIZE)
    if pending.count <= 1 {
        if let r = pending.first {
            for i in 0..<a.count {
                a[i] = r[i]
            }
        }
        return
    }

    var by_level: [[[Int]]] = [[]]
    var index = 0
    while index < pending.count {
        var placed = 0
        while placed < WAY && index < pending.count {
            by_level[0].append(pending[index])
            index += 1
            placed += 1
            collapse_from(&by_level, 0)
        }
    }

    var leftover: [[Int]] = []
    for level in 0..<by_level.count {
        leftover.append(contentsOf: by_level[level])
    }
    while leftover.count > 1 {
        let k = min(WAY, leftover.count)
        let batch = Array(leftover.prefix(k))
        leftover.removeFirst(k)
        leftover.append(merge_k_way(batch))
    }

    let result = leftover.first ?? []
    for i in 0..<a.count {
        a[i] = result[i]
    }
}


func benchmark_sort(_ array: inout [Int]) {

    oscillating_merge_sort(&array)

}

func is_non_decreasing(_ a: [Int]) -> Bool {
    guard a.count >= 2 else { return true }
    for i in 1..<a.count {
        if a[i - 1] > a[i] { return false }
    }
    return true
}

func same_multiset(_ a: [Int], _ b: [Int]) -> Bool {
    if a.count != b.count {
        return false
    }

    var left = a
    var right = b
    left.sort()
    right.sort()
    return left == right
}

func check_correctness_case(_ label: String, _ input: [Int]) {
    var input = input
    let original = input

    benchmark_sort(&input)

    if !is_non_decreasing(input) {
        fatalError("correctness case \(label): output is not sorted")
    }

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

// Skip cases larger than the algorithm's measured size cap (MAX_POWER). That
// cap exists because larger inputs are impractically slow; forcing them here
// would stall the published measurement script before any table rows print.
func check_correctness_case_within_limit(_ label: String, _ input: [Int]) {
    if input.count > (1 << MAX_POWER) {
        return
    }
    check_correctness_case(label, input)
}

func few_unique_values(_ size: Int, _ unique: Int, _ seed: UInt64) -> [Int] {
    var state = seed
    var result = [Int]()
    result.reserveCapacity(size)
    for _ in 0..<size {
        state ^= state << 13
        state ^= state >> 7
        state ^= state << 17
        result.append(Int(state % UInt64(unique)) + 1)
    }
    return result
}

func run_correctness_checks() {
    check_correctness_case("empty", [])
    check_correctness_case("single", [42])
    check_correctness_case("duplicates", [3, 1, 3, 2, 1, 2])
    check_correctness_case("sorted", [1, 2, 3, 4, 5])
    check_correctness_case("reverse", [5, 4, 3, 2, 1])
    check_correctness_case("all_equal", [7, 7, 7, 7])
    check_correctness_case("skewed_range", [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",
        [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(
            "few_keys_len32_seed_\(seed)",
            few_unique_values(32, 4, UInt64(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", [Int](repeating: 7, count: 256))
    for seed in 1...4 {
        check_correctness_case(
            "few_keys_len256_seed_\(seed)",
            few_unique_values(256, 4, UInt64(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.
    // Respect MAX_POWER so algorithms with a low measured-size cap (slow,
    // sleep) do not hang here for minutes or months.
    check_correctness_case_within_limit("all_equal_len600", [Int](repeating: 7, count: 600))
    for seed in 1...4 {
        check_correctness_case_within_limit(
            "few_keys_len2048_seed_\(seed)",
            few_unique_values(2048, 4, UInt64(seed))
        )
    }
}


func shuffled(_ size: Int, seed: UInt64) -> [Int] {
    guard size > 0 else { return [] }

    var v = Array(1...size)
    var state = seed

    if size > 1 {
        for i in stride(from: size - 1, through: 1, by: -1) {
            state ^= state << 13
            state ^= state >> 7
            state ^= state << 17

            let j = Int(state % UInt64(i + 1))
            v.swapAt(i, j)
        }
    }

    return v
}

func micros(_ d: Duration) -> UInt64 {
    let c = d.components
    let fromSeconds = UInt64(c.seconds) * 1_000_000
    let fromAttos = UInt64(max(0, c.attoseconds / 1_000_000_000_000))
    return fromSeconds + fromAttos
}

func padLeft(_ value: String, _ width: Int) -> String {
    if value.count >= width {
        return value
    }
    return String(repeating: " ", count: width - value.count) + value
}

func formatSeconds(_ micros: UInt64) -> String {
    let whole = micros / 1_000_000
    let frac = micros % 1_000_000
    let fracStr = padLeft(String(frac), 6).replacingOccurrences(of: " ", with: "0")
    return "\(whole).\(fracStr)"
}

func input_array(_ size: Int, seed: UInt64) -> [Int] {
    shuffled(size, seed: 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.
func run_once(size: Int, seed: Int) -> (UInt64, Int) {
    var array = input_array(size, seed: UInt64(seed))

    let baseBytes = alloc_track_live()
    alloc_track_reset_peak()

    let start = ContinuousClock.now

    benchmark_sort(&array)

    let elapsed = ContinuousClock.now - start
    let peakBytes = alloc_track_peak()
    let auxBytes = max(0, peakBytes - baseBytes)

    let expected: [Int] = size > 0 ? Array(1...size) : []
    if array != expected {
        fatalError("sort failed with seed \(seed) for size \(size)")
    }

    return (micros(elapsed), auxBytes)
}

func run_child(_ args: [String]) {
    let size = Int(args[2])!
    let seed = Int(args[3])!
    let (elapsedUs, mem) = run_once(size: size, seed: seed)
    print("\(elapsedUs) \(mem)")
}

let args = CommandLine.arguments
if args.count > 1 && args[1] == "--run-once" {
    run_child(args)
} else {
    run_correctness_checks()

    let tableHeader =
        "| \(padLeft("Size", 10)) | " +
        "\(padLeft("Average time (s)", 16)) | " +
        "\(padLeft("Maximum time (s)", 16)) | " +
        "\(padLeft("Average memory (KiB)", 20)) | " +
        "\(padLeft("Maximum memory (KiB)", 20)) |"
    print(tableHeader)
    print("|-----------:|-----------------:|-----------------:|---------------------:|---------------------:|")

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

        var totalTime: UInt64 = 0
        var maxTime: UInt64 = 0

        var totalMem = 0
        var maxMem = 0

        for seed in 1...RUNS {
            let process = Process()
            process.executableURL = URL(fileURLWithPath: args[0])
            process.arguments = ["--run-once", "\(size)", "\(seed)"]
            let stdout = Pipe()
            let stderr = Pipe()
            process.standardOutput = stdout
            process.standardError = stderr

            do {
                try process.run()
            } catch {
                fatalError("failed to run benchmark child process: \(error)")
            }
            process.waitUntilExit()

            if process.terminationStatus != 0 {
                let err = String(data: stderr.fileHandleForReading.readDataToEndOfFile(), encoding: .utf8) ?? ""
                fatalError("benchmark child process failed: \(err)")
            }

            let data = stdout.fileHandleForReading.readDataToEndOfFile()
            let stdoutText = String(data: data, encoding: .utf8) ?? ""
            let fields = stdoutText.split(whereSeparator: \.isWhitespace)
            guard fields.count >= 2,
                  let elapsedUs = UInt64(fields[0]),
                  let auxMem = Int(fields[1]) else {
                fatalError("invalid child process output: \(stdoutText)")
            }

            totalTime += elapsedUs
            if elapsedUs > maxTime {
                maxTime = elapsedUs
            }

            totalMem += auxMem
            if auxMem > maxMem {
                maxMem = auxMem
            }
        }

        let avgTime = totalTime / UInt64(RUNS)
        // Memory is summed in bytes and converted to KiB once, after averaging.
        let avgMemKb = totalMem / RUNS / 1024
        let maxMemKb = maxMem / 1024

        let tableRow =
            "| \(padLeft(String(size), 10)) | " +
            "\(padLeft(formatSeconds(avgTime), 16)) | " +
            "\(padLeft(formatSeconds(maxTime), 16)) | " +
            "\(padLeft(String(avgMemKb), 20)) | " +
            "\(padLeft(String(maxMemKb), 20)) |"
        print(tableRow)
    }
}
SWIFT

RUN clang -O2 -fPIC -shared alloc_track.c -o liballoc_track.so -ldl

RUN swiftc -Ounchecked -whole-module-optimization \
    main.swift \
    -o swift-benchmark \
    -L. -lalloc_track \
    -Xlinker -rpath -Xlinker /app

ENV LD_PRELOAD=/app/liballoc_track.so
CMD ["./swift-benchmark"]
"""#
    try dockerfile.write(
        to: workdir.appendingPathComponent("Dockerfile"),
        atomically: true,
        encoding: .utf8
    )

    // Keeping build and run as separate child processes preserves Docker's
    // normal output and the original image tag used by the benchmark skill.
    try runCommand("docker", ["build", "-t", "swift-benchmark", workdir.path])
    try runCommand("docker", ["run", "--rm", "--init", "swift-benchmark"])
} catch {
    fputs("\(error)\n", stderr)
    exit(1)
}