Skip to content

Performance: high-memory mode slower and larger than C on long pairs; port-quirk optimizations #31

Open

Description

@cjfields

Tracking issue for Rust-vs-C performance. Output must stay byte-identical to C by default; validate every change with tests/differential/run.sh, and log results in .claude/PERFORMANCE_LOG.md. Items can be split off into their own issues as they're picked up.

Status (2026-10-15): a full Rust-vs-C sweep showed Rust faster than C, or at parity, on almost every config, with output identical everywhere: 5 metrics × 5 memory modes × 5 datasets, plus ends-free/extension and all five heuristics. Data in the comment below. The remaining gaps are mostly in the paths added or rewritten for C parity after the initial optimization pass: med/low, BiWFA, heuristics and corrected mode. Update 2026-10-17: a three-machine sweep at fc51088 (M2, EPYC, Xeon) gave a median Rust/C of 0.87 / 0.69 / 0.80, with output identical everywhere. Since then, #52 fixed the Xeon affine2p CIGAR gap and #58 the slab cost after long alignments.

Checklist

Done

Open: measured gaps

  • High-memory CIGAR peak memory is up to 2× C on long pairs. 10 kb affine2p: 418 vs 214 MB; affine: 77 vs 45 MB; 2 kb affine with match < 0: 66 vs 45 MB. Time is at parity. Find where the extra comes from (slab, CIGAR buffer, I1/D1 below).
  • BiWFA sequence handling. Each sub-problem re-copies the sequences with init_ascii, and the reverse copy is a scalar byte loop. C copies them once per alignment and then only moves the bounds and end-of-sequence sentinels (wavefront_sequences_set_bounds). BiWFA is at parity with C (0.92–1.05×) on long or divergent pairs, rather than ahead.
  • wfa-adaptive on divergent pairs: 1.02–1.12× C's time (edit and affine2p, 20% divergence). Not yet examined.
  • Corrected mode with distant first ends: when C stops at a wrong, Known issue (C parity): ends-free with whole-length free ends and match < 0 returns the no-overlap alignment #38-style first end (random flanks, whole-length free ends), corrected mode still costs +15–55% (short pairs up to +160%). Most of that is needed for optimality. A tighter early upper bound (e.g. from a cheap seed alignment) could make the pruning from Speed up Compat::Corrected ends-free (fast kernel, exact pruning); fix #41 #42 effective sooner.
  • dada2-rs workload (Optimize the WFA backend on long (PacBio) reads dada2-rs#51): ends-free with a match reward, corrected, banded, bounded free ends, CIGAR, many short-to-5 kb pairs. Being benchmarked on the dada2-rs side as length × divergence, NW vs WFA, with affine gap_opening = 0 (current) vs gap-linear. Banded WFA's cost grows with divergence while banded NW's doesn't, so the crossover divergence is the key number. Library-side findings come back here.

Open: measurement

Checked, no action needed

  • Per-component victim outputs for I1/D1 (and I2/D2), as C's wavefront_compute_allocate_output does: tried (patch in tmp/bench/victim-outputs.patch), with identical output but no measurable change in time, instructions or peak memory on the M2. I/D components with all-null inputs are rare after the first steps. Not merged. (It was suspected for the Xeon CIGAR gap, which x86: claim affine2p output lines with PREFETCHW; drop input prefetches #52 fixed.)
  • Independent I/D trim_ends (Banded alignments don't terminate as END_UNREACHABLE like C (hang in score-only mode) #4): per component, as in C (wavefront_compute_process_ends). Nothing to gain.
  • Per-step med/low checks (Port C's piggyback backtrace: --wfa-memory=med|low, buffer compaction #28): per step, not per cell; the plain kernels and extend loops have none. Under 1%.
  • Ends-free begin-gap null steps (endsfree_null_m): uses the slab free list, not a fresh allocation, and is rare. Same as C.
  • C-sized components (1000 × 1000) and the CIGAR buffer: allocated once, grown only when needed, with no per-alignment memset.
  • BiWFA sub-aligners sized 10 × 10: created once and kept. Only a few early grow steps. (The per-sub-problem sequence copy above is the real cost.)

Original description (2026-10-13)

Symptom

Rust is faster than C on the 1 kb benchmark (PERFORMANCE_LOG.md), but in high-memory (CIGAR) mode it falls behind on longer pairs. Measurements on wfa.utest.seq (incl. its 10 kb pairs), M2 Pro, single run, Rust main vs C 3282f83:

Configuration Rust C Rust max RSS C max RSS
gap-affine2p, -1,4,6,2,24,1, ends-free 0.1,0,0.1,0 12.7 s 4.3 s 5.3 GB 3.7 GB
gap-affine2p (default penalties) 2.5 s 1.7 s 3.7 GB 1.8 GB
gap-affine 0.41 s 0.35 s 972 MB 285 MB

Under parallel load, the first case timed out in the differential harness (#30).

Main cause (known target)

In non-modular (CIGAR) mode, every wavefront is allocated over the full historic diagonal range, [-(plen + scope), tlen + scope]. C allocates each wavefront at its own lo..hi plus margins (wavefront_compute_allocate_output, slab tight mode). Memory therefore grows with sequence length × score rather than with the wavefront widths. It's already noted as the next memory target in PERFORMANCE_LOG.md (#8).

Other candidates from the port's C-parity work

Each was kept for output identity; look for an equivalent that is cheaper:

Process

Activity

  1. cjfields commented on Oct 14, 2026

    @cjfields
    MemberAuthor

    PR #36 fixes the main target: wavefronts are now allocated over C's growing range, with a tight slab in CIGAR mode. Long pairs now match C on memory and time. For example, utest affine 927 → 275 MB (C 272) and affine2p 3.5 → 1.8 GB (C 1.75). The ends-free affine2p case with match < 0 goes from 18.8 s to 3.4 s, against C's 4.2 s.

    Another candidate from that work: Rust always allocates the I1/D1 (and I2/D2) outputs, while C uses the victim wavefront when both of their inputs are null (wavefront_compute_allocate_output).

  2. cjfields commented on Oct 15, 2026

    @cjfields
    MemberAuthor

    Rust-vs-C sweep, 2026-10-15 (main 0268377 vs C 3282f83, Apple Silicon M-series; C built by tests/differential/build_c.sh, i.e. without -march=native, Rust with target-cpu=native).

    • Measures. Instructions retired (whole process, so input parsing and output are included; that skews the 150 bp set) and alignment time (Time.Alignment), each the best of 3 runs at load ~1.5–2.5, plus max RSS. Ratios are Rust/C; below 1 means Rust is better.
    • Data (C generate_dataset):
      • s150_2: 50k × 150 bp, 2%
      • m1k_5: 20k × 1 kb, 5%
      • m1k_20: 5k × 1 kb, 20%
      • m1k_indel: 2k × 1 kb, 5% plus 5 indels of up to 50 bp
      • l10k_5: 300 × 10 kb, 5%
      • lng2k: 300 × 2 kb
    • Penalties. Defaults; "m-1" means a match reward (-1,4,6,2, and -1,4,6,2,24,1 for affine2p). Heuristic parameters: wfa-adaptive 10,50,1, xdrop/zdrop 100,1, banded -50,50.
    • Output was identical to C in all 301 valid configs. The 16 edit + X/Z-drop configs are rejected by both tools, so they're excluded. A — in the time column is a run over a minute, which the parser didn't read.

    Summary

    • High-memory CIGAR and score-only, all metrics: 0.13–1.02× time; Rust is faster almost everywhere.
    • BiWFA (ultralow): 0.35–1.05× time; at parity on long or divergent pairs.
    • med/low:
      • affine/affine2p use 1.17–1.31× C's instructions, though time is still 0.63–0.95×;
      • indel is 1.12–1.34× C's time;
      • edit/linear: 1.0–1.3× instructions, 0.6–1.03× time.
    • Ends-free/extension, including with match < 0: 0.58–0.93× time (med 0.63–0.93×).
    • Heuristics: 0.64–1.0× time, except wfa-adaptive at 20% divergence (1.02–1.12×).
    • Memory: close to C, except high-memory CIGAR on long pairs, up to 2× (10 kb affine2p 418 vs 214 MB).
    Core: metrics × memory modes × datasets (125 configs)
    Case R/C instr. R/C time C ms Rust ms C MB Rust MB
    s150_2 edit cigar 0.14 0.17 74 12 1 2
    s150_2 edit score 0.06 0.13 70 9 1 2
    s150_2 edit med 0.18 0.27 76 20 1 2
    s150_2 edit low 0.18 0.26 76 20 1 2
    s150_2 edit ultralow 0.19 0.35 95 33 1 2
    s150_2 indel cigar 0.16 0.19 75 14 1 2
    s150_2 indel score 0.07 0.14 71 10 1 2
    s150_2 indel med 0.19 0.28 76 21 1 2
    s150_2 indel low 0.19 0.28 77 22 1 2
    s150_2 indel ultralow 0.21 0.37 96 36 1 2
    s150_2 linear cigar 0.16 0.22 78 17 1 3
    s150_2 linear score 0.08 0.19 76 14 1 2
    s150_2 linear med 0.21 0.32 85 27 1 2
    s150_2 linear low 0.21 0.33 85 28 1 2
    s150_2 linear ultralow 0.23 0.43 113 48 1 3
    s150_2 affine cigar 0.24 0.41 109 45 1 4
    s150_2 affine score 0.19 0.39 108 42 1 2
    s150_2 affine med 0.37 0.72 127 91 1 2
    s150_2 affine low 0.37 0.70 128 89 1 2
    s150_2 affine ultralow 0.38 0.63 189 120 1 4
    s150_2 affine2p cigar 0.24 0.43 114 48 1 5
    s150_2 affine2p score 0.21 0.45 113 51 1 2
    s150_2 affine2p med 0.39 0.75 132 98 2 3
    s150_2 affine2p low 0.39 0.74 131 97 2 3
    s150_2 affine2p ultralow 0.47 0.74 239 178 2 6
    m1k_5 edit cigar 0.61 0.78 132 103 1 2
    m1k_5 edit score 0.51 0.73 114 83 1 2
    m1k_5 edit med 0.86 0.76 286 216 1 2
    m1k_5 edit low 0.86 0.76 286 216 1 2
    m1k_5 edit ultralow 0.65 0.84 190 160 1 2
    m1k_5 indel cigar 0.65 0.82 157 129 1 2
    m1k_5 indel score 0.57 0.77 137 106 1 2
    m1k_5 indel med 0.88 1.14 204 232 1 2
    m1k_5 indel low 0.88 1.12 206 230 1 2
    m1k_5 indel ultralow 0.69 0.87 218 190 2 3
    m1k_5 linear cigar 0.65 0.79 196 155 2 3
    m1k_5 linear score 0.60 0.76 177 135 1 2
    m1k_5 linear med 0.84 0.89 396 352 1 2
    m1k_5 linear low 0.87 0.92 442 408 1 2
    m1k_5 linear ultralow 0.69 0.81 297 242 2 3
    m1k_5 affine cigar 0.84 0.93 898 834 4 6
    m1k_5 affine score 0.84 0.90 877 791 1 2
    m1k_5 affine med 1.23 0.80 4170 3347 1 2
    m1k_5 affine low 1.22 0.87 4610 4028 1 2
    m1k_5 affine ultralow 0.85 0.89 1510 1343 3 5
    m1k_5 affine2p cigar 0.83 0.94 3440 3238 9 13
    m1k_5 affine2p score 0.81 0.93 3320 3090 2 3
    m1k_5 affine2p med 1.24 0.76 27350 20687 3 5
    m1k_5 affine2p low 1.24 0.83 29380 24374 3 4
    m1k_5 affine2p ultralow 0.93 1.01 16000 16159 6 10
    m1k_20 edit cigar 0.83 0.95 192 183 2 3
    m1k_20 edit score 0.86 0.95 174 166 1 2
    m1k_20 edit med 1.15 0.62 757 473 1 2
    m1k_20 edit low 1.15 0.63 753 475 1 2
    m1k_20 edit ultralow 0.84 0.98 234 229 2 3
    m1k_20 indel cigar 0.86 0.98 249 244 2 3
    m1k_20 indel score 0.91 0.98 230 226 1 2
    m1k_20 indel med 1.19 1.31 387 505 1 2
    m1k_20 indel low 1.19 1.31 388 507 1 2
    m1k_20 indel ultralow 0.87 1.00 300 300 2 3
    m1k_20 linear cigar 0.86 0.97 284 276 2 4
    m1k_20 linear score 0.90 0.94 270 254 1 2
    m1k_20 linear med 1.01 0.81 1080 873 1 2
    m1k_20 linear low 1.05 0.89 1200 1071 1 2
    m1k_20 linear ultralow 0.86 0.95 367 348 2 3
    m1k_20 affine cigar 0.93 1.02 1540 1571 8 15
    m1k_20 affine score 0.93 0.96 1500 1433 1 2
    m1k_20 affine med 1.29 0.78 9790 7619 2 3
    m1k_20 affine low 1.26 0.85 10820 9187 2 3
    m1k_20 affine ultralow 0.91 0.93 2570 2380 3 6
    m1k_20 affine2p cigar 0.84 0.99 8230 8115 42 57
    m1k_20 affine2p score 0.84 0.95 7610 7228 4 5
    m1k_20 affine2p med 1.24 — — 53078 13 14
    m1k_20 affine2p low 1.23 — — 64797 12 13
    m1k_20 affine2p ultralow 0.94 1.01 31350 31528 9 15
    m1k_indel edit cigar 0.86 0.96 81 78 2 3
    m1k_indel edit score 0.87 0.96 77 74 1 2
    m1k_indel edit med 1.19 0.61 480 293 1 2
    m1k_indel edit low 1.19 0.61 478 292 1 2
    m1k_indel edit ultralow 0.87 0.98 99 97 2 3
    m1k_indel indel cigar 0.87 0.97 72 70 2 3
    m1k_indel indel score 0.89 0.95 69 66 1 2
    m1k_indel indel med 1.22 1.28 113 144 1 2
    m1k_indel indel low 1.22 1.28 113 144 1 2
    m1k_indel indel ultralow 0.88 0.98 89 87 2 3
    m1k_indel linear cigar 0.86 0.96 83 80 2 4
    m1k_indel linear score 0.88 0.93 80 74 1 2
    m1k_indel linear med 1.00 0.68 418 282 1 2
    m1k_indel linear low 1.04 0.73 472 342 1 2
    m1k_indel linear ultralow 0.86 0.96 119 114 2 3
    m1k_indel affine cigar 0.90 0.99 231 229 6 9
    m1k_indel affine score 0.90 0.98 217 214 1 2
    m1k_indel affine med 1.27 0.69 1550 1070 2 3
    m1k_indel affine low 1.25 0.76 1730 1310 2 3
    m1k_indel affine ultralow 0.89 0.91 413 377 3 6
    m1k_indel affine2p cigar 0.84 0.97 722 700 17 19
    m1k_indel affine2p score 0.82 0.92 721 666 2 4
    m1k_indel affine2p med 1.24 0.73 6910 5035 6 8
    m1k_indel affine2p low 1.23 0.79 7380 5814 6 8
    m1k_indel affine2p ultralow 0.92 1.00 2960 2968 7 12
    l10k_5 edit cigar 0.91 0.97 81 78 3 6
    l10k_5 edit score 0.94 0.98 73 71 1 2
    l10k_5 edit med 1.28 0.74 252 186 1 3
    l10k_5 edit low 1.28 0.74 254 187 1 3
    l10k_5 edit ultralow 0.92 1.00 89 89 3 5
    l10k_5 indel cigar 0.94 1.02 103 105 7 9
    l10k_5 indel score 0.97 1.01 96 97 1 2
    l10k_5 indel med 1.32 1.33 164 218 1 3
    l10k_5 indel low 1.32 1.34 162 218 1 3
    l10k_5 indel ultralow 0.95 1.05 123 129 2 4
    l10k_5 linear cigar 0.94 1.01 110 111 7 11
    l10k_5 linear score 0.98 0.98 103 101 1 2
    l10k_5 linear med 1.08 0.97 352 341 2 3
    l10k_5 linear low 1.12 1.03 406 418 1 3
    l10k_5 linear ultralow 0.93 1.00 150 150 2 4
    l10k_5 affine cigar 0.94 0.97 1040 1014 45 77
    l10k_5 affine score 0.95 0.96 852 816 2 3
    l10k_5 affine med 1.31 0.89 4860 4334 10 11
    l10k_5 affine low 1.28 0.95 5720 5438 8 9
    l10k_5 affine ultralow 0.94 1.03 1240 1277 3 7
    l10k_5 affine2p cigar 0.85 0.97 5460 5284 214 418
    l10k_5 affine2p score 0.85 0.99 4570 4542 12 13
    l10k_5 affine2p med 1.24 0.73 41530 30289 67 74
    l10k_5 affine2p low 1.23 0.78 47660 37135 57 64
    l10k_5 affine2p ultralow 0.93 0.99 13060 12913 14 20
    Ends-free / extension (72 configs)
    Case R/C instr. R/C time C ms Rust ms C MB Rust MB
    m1k_5 affine ends-free cigar 0.82 0.88 2290 2016 4 6
    m1k_5 affine ends-free score 0.82 0.85 2250 1921 1 2
    m1k_5 affine ends-free med 1.17 0.63 12300 7699 2 3
    m1k_5 affine extension cigar 0.82 0.88 1240 1096 4 6
    m1k_5 affine extension score 0.79 0.83 1180 985 1 2
    m1k_5 affine extension med 1.17 0.80 4440 3549 1 2
    m1k_5 affine2p ends-free cigar 0.79 0.86 7180 6200 9 13
    m1k_5 affine2p ends-free score 0.78 0.85 7060 5987 2 3
    m1k_5 affine2p ends-free med 1.18 0.73 50930 37312 3 5
    m1k_5 affine2p extension cigar 0.79 0.86 4570 3949 9 13
    m1k_5 affine2p extension score 0.77 0.82 4540 3715 2 3
    m1k_5 affine2p extension med 1.19 0.75 28420 21275 3 5
    m1k_5 affine2p m-1 ends-free cigar 0.83 0.81 10960 8909 17 27
    m1k_5 affine2p m-1 ends-free score 0.80 0.77 11010 8513 2 4
    m1k_5 affine2p m-1 ends-free med 1.26 0.87 58010 50535 5 7
    m1k_5 affine2p m-1 extension cigar 0.82 0.79 7280 5744 17 25
    m1k_5 affine2p m-1 extension score 0.79 0.74 7240 5328 2 4
    m1k_5 affine2p m-1 extension med 1.26 0.93 33300 30901 5 6
    m1k_5 edit ends-free cigar 0.74 0.81 499 403 1 2
    m1k_5 edit ends-free score 0.74 0.80 474 380 1 2
    m1k_5 edit ends-free med 1.01 0.61 2140 1303 1 2
    m1k_5 edit extension cigar 0.62 0.76 165 126 1 2
    m1k_5 edit extension score 0.54 0.77 142 109 1 2
    m1k_5 edit extension med 0.84 0.75 320 238 1 2
    m1k_5 linear ends-free cigar 0.76 0.81 743 603 2 3
    m1k_5 linear ends-free score 0.76 0.81 725 588 1 2
    m1k_5 linear ends-free med 0.93 0.68 2580 1754 1 2
    m1k_5 linear extension cigar 0.71 0.83 280 233 2 3
    m1k_5 linear extension score 0.62 0.69 227 155 1 2
    m1k_5 linear extension med 0.87 0.90 475 426 1 2
    lng2k affine ends-free cigar 0.84 0.85 282 241 21 17
    lng2k affine ends-free score 0.84 0.89 262 234 1 2
    lng2k affine ends-free med 1.19 0.68 1470 999 4 5
    lng2k affine extension cigar 0.85 0.89 164 147 11 17
    lng2k affine extension score 0.84 0.84 158 133 1 2
    lng2k affine extension med 1.20 0.78 692 538 3 4
    lng2k affine2p ends-free cigar 0.79 0.93 1040 964 54 63
    lng2k affine2p ends-free score 0.79 0.89 925 822 5 7
    lng2k affine2p ends-free med 1.19 0.72 6900 4981 20 22
    lng2k affine2p extension cigar 0.79 0.91 733 665 49 60
    lng2k affine2p extension score 0.79 0.88 679 595 5 5
    lng2k affine2p extension med 1.19 0.73 4600 3363 18 15
    lng2k affine2p m-1 ends-free cigar 0.82 0.88 1730 1517 81 128
    lng2k affine2p m-1 ends-free score 0.82 0.78 1400 1090 6 8
    lng2k affine2p m-1 ends-free med 1.26 0.86 7930 6790 22 24
    lng2k affine2p m-1 extension cigar 0.82 0.85 1320 1125 69 124
    lng2k affine2p m-1 extension score 0.81 0.77 1060 817 6 8
    lng2k affine2p m-1 extension med 1.27 0.87 5640 4911 20 23
    lng2k edit ends-free cigar 0.82 0.83 65 54 2 3
    lng2k edit ends-free score 0.83 0.83 63 52 1 2
    lng2k edit ends-free med 1.10 0.58 292 171 1 2
    lng2k edit extension cigar 0.81 0.84 25 21 2 3
    lng2k edit extension score 0.82 0.82 23 19 1 2
    lng2k edit extension med 1.08 0.64 78 49 1 2
    lng2k linear ends-free cigar 0.83 0.84 88 75 3 6
    lng2k linear ends-free score 0.83 0.83 86 72 1 2
    lng2k linear ends-free med 0.99 0.68 333 227 2 2
    lng2k linear extension cigar 0.87 0.75 34 25 3 5
    lng2k linear extension score 0.87 0.71 31 22 1 2
    lng2k linear extension med 1.00 0.73 91 66 1 2
    m1k_5 affine m-1 ends-free cigar 0.90 0.75 7180 5370 11 16
    m1k_5 affine m-1 ends-free score 0.87 0.71 7140 5086 1 2
    m1k_5 affine m-1 ends-free med 1.27 0.79 30170 23931 2 3
    m1k_5 affine m-1 extension cigar 0.88 0.77 3670 2824 10 15
    m1k_5 affine m-1 extension score 0.83 0.73 3600 2640 1 2
    m1k_5 affine m-1 extension med 1.25 0.89 11030 9813 2 3
    lng2k affine m-1 ends-free cigar 0.91 0.79 1100 873 45 66
    lng2k affine m-1 ends-free score 0.90 0.73 914 665 2 3
    lng2k affine m-1 ends-free med 1.28 0.83 3790 3153 6 7
    lng2k affine m-1 extension cigar 0.91 0.79 666 523 40 65
    lng2k affine m-1 extension score 0.89 0.72 554 397 2 3
    lng2k affine m-1 extension med 1.27 0.89 1830 1637 4 5
    Heuristics (104 configs)
    Case R/C instr. R/C time C ms Rust ms C MB Rust MB
    m1k_5 affine wfa-adaptive cigar 0.78 0.92 583 535 4 6
    m1k_5 affine wfa-adaptive score 0.77 0.88 582 514 1 2
    m1k_5 affine wfa-adaptive med 1.09 0.96 1120 1077 1 2
    m1k_5 affine wfa-adaptive ultralow 0.85 0.97 1220 1181 3 5
    m1k_5 affine xdrop cigar 0.85 0.93 389 363 2 5
    m1k_5 affine xdrop score 0.83 0.90 406 366 1 2
    m1k_5 affine xdrop med 1.16 0.83 1120 930 1 2
    m1k_5 affine xdrop ultralow 0.92 0.93 832 772 2 5
    m1k_5 affine zdrop cigar 0.91 0.98 519 506 2 5
    m1k_5 affine zdrop score 0.81 0.91 392 356 1 2
    m1k_5 affine zdrop med 1.16 0.92 1170 1080 1 2
    m1k_5 affine zdrop ultralow 0.90 0.92 777 713 2 5
    m1k_5 affine banded-static cigar 0.79 0.91 618 559 4 6
    m1k_5 affine banded-static score 0.78 0.84 624 521 1 2
    m1k_5 affine banded-static med 1.19 0.77 2640 2037 1 2
    m1k_5 affine banded-static ultralow 0.85 0.88 1340 1173 3 5
    m1k_5 affine banded-adaptive cigar 0.79 0.89 608 541 4 6
    m1k_5 affine banded-adaptive score 0.78 0.85 596 506 1 2
    m1k_5 affine banded-adaptive med 1.19 0.76 2580 1957 1 2
    m1k_5 affine banded-adaptive ultralow 0.85 0.88 1350 1184 3 5
    m1k_5 edit wfa-adaptive cigar 0.57 0.77 108 84 1 2
    m1k_5 edit wfa-adaptive score 0.41 0.73 88 65 1 2
    m1k_5 edit wfa-adaptive med 0.70 0.86 146 125 1 2
    m1k_5 edit wfa-adaptive ultralow 0.65 0.90 188 169 1 2
    m1k_5 edit banded-static cigar 0.62 0.79 132 104 1 2
    m1k_5 edit banded-static score 0.52 0.76 115 87 1 2
    m1k_5 edit banded-static med 0.86 0.75 291 217 1 2
    m1k_5 edit banded-static ultralow 0.67 0.84 195 164 1 2
    m1k_5 edit banded-adaptive cigar 0.62 0.78 136 106 1 2
    m1k_5 edit banded-adaptive score 0.52 0.75 115 86 1 2
    m1k_5 edit banded-adaptive med 0.86 0.75 285 214 1 2
    m1k_5 edit banded-adaptive ultralow 0.67 0.84 193 162 1 2
    m1k_5 affine2p wfa-adaptive cigar 0.85 1.00 1810 1801 9 13
    m1k_5 affine2p wfa-adaptive score 0.79 0.91 1810 1642 2 3
    m1k_5 affine2p wfa-adaptive med 1.24 0.98 4440 4331 3 4
    m1k_5 affine2p wfa-adaptive ultralow 0.86 0.95 4900 4633 6 10
    m1k_5 affine2p xdrop cigar 0.86 0.91 1270 1157 4 9
    m1k_5 affine2p xdrop score 0.83 0.87 1310 1142 2 3
    m1k_5 affine2p xdrop med 1.26 0.87 4910 4282 3 4
    m1k_5 affine2p xdrop ultralow 0.86 0.88 2710 2391 2 3
    m1k_5 affine2p zdrop cigar 0.86 0.92 950 875 4 8
    m1k_5 affine2p zdrop score 0.84 0.94 850 800 2 3
    m1k_5 affine2p zdrop med 1.22 0.84 3710 3113 3 4
    m1k_5 affine2p zdrop ultralow 0.88 0.96 1750 1677 2 3
    m1k_5 affine2p banded-static cigar 0.79 0.89 1550 1374 9 13
    m1k_5 affine2p banded-static score 0.73 0.80 1480 1182 2 3
    m1k_5 affine2p banded-static med 1.22 0.77 10040 7696 3 4
    m1k_5 affine2p banded-static ultralow 0.90 1.03 8530 8787 6 10
    m1k_5 affine2p banded-adaptive cigar 0.79 0.88 1540 1357 9 13
    m1k_5 affine2p banded-adaptive score 0.73 0.78 1500 1163 2 3
    m1k_5 affine2p banded-adaptive med 1.22 0.76 10170 7726 3 4
    m1k_5 affine2p banded-adaptive ultralow 0.90 1.02 8680 8852 6 10
    m1k_20 affine wfa-adaptive cigar 0.96 1.05 604 636 8 10
    m1k_20 affine wfa-adaptive score 0.97 1.02 603 615 1 2
    m1k_20 affine wfa-adaptive med 1.23 0.84 1720 1452 2 2
    m1k_20 affine wfa-adaptive ultralow 0.99 1.04 1960 2035 3 6
    m1k_20 affine xdrop cigar 0.85 0.94 96 90 2 5
    m1k_20 affine xdrop score 0.84 0.91 99 90 1 2
    m1k_20 affine xdrop med 1.16 0.80 297 238 1 2
    m1k_20 affine xdrop ultralow 0.92 0.93 204 189 1 2
    m1k_20 affine zdrop cigar 0.92 0.97 128 124 2 5
    m1k_20 affine zdrop score 0.82 0.90 96 87 1 2
    m1k_20 affine zdrop med 1.17 0.87 304 264 2 2
    m1k_20 affine zdrop ultralow 0.90 0.93 193 180 1 2
    m1k_20 affine banded-static cigar 0.85 0.93 478 447 8 10
    m1k_20 affine banded-static score 0.84 0.88 474 417 1 2
    m1k_20 affine banded-static med 1.22 0.70 2440 1703 2 3
    m1k_20 affine banded-static ultralow 0.88 0.92 1640 1507 3 6
    m1k_20 affine banded-adaptive cigar 0.85 0.92 505 463 8 10
    m1k_20 affine banded-adaptive score 0.84 0.86 493 424 1 2
    m1k_20 affine banded-adaptive med 1.22 0.70 2420 1688 2 3
    m1k_20 affine banded-adaptive ultralow 0.88 0.89 1640 1464 3 6
    m1k_20 edit wfa-adaptive cigar 0.79 1.06 92 98 2 3
    m1k_20 edit wfa-adaptive score 0.84 1.07 76 82 1 2
    m1k_20 edit wfa-adaptive med 0.91 0.84 183 154 1 2
    m1k_20 edit wfa-adaptive ultralow 0.88 1.12 169 189 2 3
    m1k_20 edit banded-static cigar 0.78 0.94 118 110 2 3
    m1k_20 edit banded-static score 0.81 0.91 105 96 1 2
    m1k_20 edit banded-static med 1.05 0.64 422 272 1 2
    m1k_20 edit banded-static ultralow 0.84 0.96 207 199 2 3
    m1k_20 edit banded-adaptive cigar 0.78 0.95 119 113 2 3
    m1k_20 edit banded-adaptive score 0.81 0.93 105 98 1 2
    m1k_20 edit banded-adaptive med 1.05 0.64 422 272 1 2
    m1k_20 edit banded-adaptive ultralow 0.84 0.97 209 203 2 3
    m1k_20 affine2p wfa-adaptive cigar 0.98 1.08 2120 2296 23 30
    m1k_20 affine2p wfa-adaptive score 0.94 1.03 2070 2136 2 3
    m1k_20 affine2p wfa-adaptive med 1.27 0.84 9170 7708 3 5
    m1k_20 affine2p wfa-adaptive ultralow 0.96 1.01 10130 10243 7 12
    m1k_20 affine2p xdrop cigar 0.86 0.92 315 289 4 8
    m1k_20 affine2p xdrop score 0.83 0.87 322 279 2 3
    m1k_20 affine2p xdrop med 1.26 0.86 1250 1080 3 4
    m1k_20 affine2p xdrop ultralow 0.86 0.87 668 580 2 3
    m1k_20 affine2p zdrop cigar 0.87 0.92 238 218 4 8
    m1k_20 affine2p zdrop score 0.84 0.95 209 199 2 3
    m1k_20 affine2p zdrop med 1.22 0.82 958 788 3 4
    m1k_20 affine2p zdrop ultralow 0.88 0.94 440 413 2 3
    m1k_20 affine2p banded-static cigar 0.80 0.90 1240 1111 23 30
    m1k_20 affine2p banded-static score 0.74 0.80 1200 956 2 3
    m1k_20 affine2p banded-static med 1.22 0.74 8880 6550 3 5
    m1k_20 affine2p banded-static ultralow 0.88 0.99 9190 9084 7 11
    m1k_20 affine2p banded-adaptive cigar 0.80 0.89 1240 1108 23 30
    m1k_20 affine2p banded-adaptive score 0.74 0.78 1200 939 2 3
    m1k_20 affine2p banded-adaptive med 1.22 0.74 8840 6572 3 5
    m1k_20 affine2p banded-adaptive ultralow 0.88 1.00 9060 9065 7 11

    🤖 Generated with Claude Code

  3. cjfields commented on Oct 16, 2026

    @cjfields
    MemberAuthor

    Full EPYC 7713 sweep after #45 and #48 (hpcbioamd/compute-5-6, commit 8ccfb64, bench/sweep.py, 3 runs per case). Ratios are Rust / C, so below 1 means Rust is faster. "Before" is the quick run at aa091cf (1 run per case).

    Group (cases) before: time / C time / C (min–max) time / C-native (max) instructions / C peak RSS / C, max
    1k/10k cigar (20) 0.76 0.76 (0.64–0.88) 0.89 (1.04) 0.71 1.93
    1k/10k score (20) 0.74 0.74 (0.63–0.87) 0.87 (1.02) 0.69 1.00
    1k/10k med/low (40) 0.94 0.45 (0.30–0.79) 0.55 (0.87) 0.55 1.21
    1k/10k ultralow (20) 0.85 0.83 (0.73–1.15) 0.95 (2.35) 0.80 1.50
    span cigar/score (48) 0.61 0.61 (0.54–0.80) 0.57 (0.65) 0.62 1.89
    span med/low (24) 0.79 0.37 (0.27–0.70) 0.64 (0.77) 0.50 1.36
    heuristics (104) 0.87 0.83 (0.40–1.07) 0.92 (1.85) 0.83 1.55
    s150 med/low (10) 1.31 1.20 (1.12–1.31) 1.17 (1.29) 1.16 1.00
    s150 other (15) ~0.95 0.87–0.91 0.82–0.93 0.93–0.99 1.13
  4. cjfields commented on Oct 16, 2026

    @cjfields
    MemberAuthor

    Three-machine sweep at main after #45, #48 and #49. Generated by report.py (a draft generator, not yet in the repo) from three bench/sweep.py runs. The M2 run is c43fd1a, whose source is identical to fc51088; it's marked -dirty only because of an uncommitted build_c.sh change that gives C a cheap timer on macOS (see the note below). 3 runs per case.

    Rust time / C time: median (range) over the cases in each group. C is built without -march=native.

    Cases n Apple M2 Pro AMD EPYC 7713 Intel Xeon Gold 6150
    Score, high memory (1–10 kb) 20 0.97 (0.90–1.02) 0.71 (0.62–0.87) 0.85 (0.66–0.93)
    CIGAR, high memory (1–10 kb) 20 0.98 (0.93–1.02) 0.73 (0.64–0.87) 0.92 (0.78–1.04)
    CIGAR, med/low memory (1–10 kb) 40 0.52 (0.23–0.96) 0.44 (0.29–0.77) 0.51 (0.39–0.82)
    BiWFA, ultralow memory (1–10 kb) 20 0.95 (0.45–1.01) 0.75 (0.51–0.89) 0.79 (0.56–0.90)
    Ends-free / extension, high memory 48 0.84 (0.70–0.95) 0.60 (0.54–0.80) 0.71 (0.56–0.95)
    Ends-free / extension, med memory 24 0.45 (0.22–0.68) 0.38 (0.27–0.67) 0.46 (0.31–0.79)
    Heuristics (adaptive, X/Z-drop, banded) 104 0.90 (0.35–1.17) 0.78 (0.40–1.06) 0.89 (0.45–1.31)
    150 bp reads, high memory and BiWFA 15 0.90 (0.84–0.98) 0.89 (0.80–1.13) 0.97 (0.89–1.26)
    150 bp reads, med/low memory 10 1.24 (1.08–1.31) 1.20 (1.11–1.33) 1.34 (1.11–1.58)

    Overall

    Apple M2 Pro AMD EPYC 7713 Intel Xeon Gold 6150
    Cases faster than C 271 / 301 284 / 301 262 / 301
    Median Rust / C 0.87 0.69 0.80
    Median Rust / C -march=native 0.87 0.82 0.94
    Output identical to C 301 / 301 301 / 301 301 / 301

    Notes

    • Vectorize the BiWFA overlap scan (affine2p BiWFA ~2x faster) #49 on the Xeon: affine2p BiWFA went from 1.28–1.41× to 0.56–0.63× C, and from 2.5–3.0× to 1.07–1.29× C-native. C-native uses its AVX-512 scans there.
    • Skylake alignment effect: Vectorize the BiWFA overlap scan (affine2p BiWFA ~2x faster) #49 shifted the code layout. X-drop/Z-drop affine2p improved about 8%, while edit/indel/linear on m1k_5 didn't change. That's consistent with the effect depending on code layout (see README: note on the Skylake JCC erratum build flag #50).
    • Slower than C on all machines: 150 bp reads in med/low memory, at 1.1–1.6× (the Xeon is worst). This is the clearest remaining gap.
    • Also slower than C:
      • wfa-adaptive edit, at about 1.05–1.17×.
      • On the Xeon only, banded affine2p CIGAR on m1k_20, at about 1.3×.
    • macOS C timer: earlier M2 ratios were biased in Rust's favor whenever a run had many fast pairs. C's profiler_timer.c reads the clock through three Mach calls, about 1.3 µs per alignment inside the timed interval. That made s150 C look about 6× slower and m1k_5 fast modes about 20–25% slower. This M2 run uses clock_gettime instead, as on Linux.
  5. cjfields commented on Oct 16, 2026

    @cjfields
    MemberAuthor

    Xeon affine2p CIGAR gap: profiled (perf on compute-5-1 at fc51088).

    The extra cost of CIGAR over score-only is about 3× larger in Rust than in C:

    Rust CIGAR − score C CIGAR − score
    m1k_5 affine2p 6607 − 5059 = 1548 ms 6410 − 5940 = 470 ms
    m1k_20 banded affine2p 3530 − 1726 = 1804 ms 2690 − 2150 = 540 ms

    Counters, in millions:

    a2p CIGAR, Rust a2p CIGAR, C a2p score, Rust banded CIGAR, Rust banded CIGAR, C
    cycles 24,849 24,205 18,570 12,398 10,216
    instructions 53,936 81,744 51,342 20,536 27,866
    stall cycles 6,615 2,031 1,961 5,425 1,789
    store-buffer-full stall cycles (resource_stalls.sb) 4,751 415 566 4,457 504
    L2 misses 1,063 1,179 0 328 694
    • Cause: Rust's CIGAR mode is bound by the store buffer. Store-buffer-full cycles account for about two thirds of the extra CIGAR cycles. Score mode doesn't have the problem, because its wavefronts are reused in a ring and stay in cache.
    • Why: in high-memory mode, compute_affine2p_step allocates fresh slab wavefronts at every step for M, I1 and D1, and for I2/D2 unless both are inactive. All of them are written. C (wavefront_compute_allocate_output) gives each component whose inputs are null the shared victim wavefront, which stays hot in cache. So Rust writes many more cache lines of fresh memory per step. On Skylake each one is a read-for-ownership miss that holds a store-buffer entry. The EPYC absorbs this better.
    • This is probably also the CIGAR peak-memory gap, e.g. m1k_20 affine2p CIGAR at 64 vs 40 MB, and l10k_5 at 420 vs 222 MB.
    • Fix: implement the existing item "Per-component victim outputs for I1/D1 (and I2/D2)". Each output should use a victim per component, separate victims to avoid aliasing &mut outputs, whenever its inputs are null, as C does. That cuts the stores, the memory and the store-buffer stalls. It should also help gap-affine and other machines, to a lesser degree.
    • Second finding, Skylake front end: in all affine2p cases, including score mode, only 17–26% of Rust's micro-ops come from the uop cache (DSB), against 80–87% for C. The rest comes from the slower legacy decoder. Gap-affine is fine (81%). The inlined compute_affine2p_step seems too large or badly laid out for Skylake's uop cache. That likely explains part of why Rust slows more than C on the Xeon (1.47× vs 1.29× its EPYC time). Possible fixes: keep the kernel loops out of line, or rely on the JCC flag (README: note on the Skylake JCC erratum build flag #50).
  6. cjfields commented on Oct 16, 2026

    @cjfields
    MemberAuthor

    Correction to the Xeon CIGAR diagnosis above: victim outputs are not the cause. I implemented per-component victim outputs, where an I/D component with all-null inputs writes to its own victim and its slot is left null, as in C. Output was identical: the differential harness was clean and all tests passed. But it made no measurable difference on the M2 to time, instructions or peak memory, in any mode. I/D components with all-null inputs are rare after the first few steps. I've kept the change as a patch rather than merging it.

    What the data does show: the large CIGAR overhead (CIGAR time minus score-only time) is specific to the Xeon and to affine2p.

    CIGAR − score, Rust / C M2 EPYC Xeon
    m1k_5 affine2p +8% / +3% +10% / +5% +32% / +8%
    m1k_5 gap-affine +6% / +2% +7% / +5% 0% / −6%

    New hypothesis:

    • The AVX2 affine2p kernel issues 7 software prefetches per 8 diagonals while writing 5 output arrays. In CIGAR mode those outputs are fresh memory, so every write also misses the cache. That's about 12 concurrent miss streams.
    • Skylake has 10 L1 fill buffers, against about 2× that on Zen 3. The profile showed 2.5× C's fill-buffer-full events, plus the store-buffer stalls.
    • Gap-affine (4 prefetches, 3 outputs) stays under the limit, and score mode reuses outputs that are already in cache. The M2's NEON kernels don't prefetch.

    Test: an A/B on the Xeon comparing three builds: unchanged, no prefetches in the affine2p kernel, and no x86 software prefetches at all.

  7. cjfields commented on Oct 17, 2026

    @cjfields
    MemberAuthor

    Xeon affine2p CIGAR: narrowed to the kernel's output stores.

    Prefetch A/B, Xeon and EPYC:

    • Removing the x86 software prefetches does not change the store-buffer stalls (affine2p CIGAR: 4766 M → 4789 M cycles), so the prefetches aren't the cause.
    • Removing them is still 2–6% faster in most cases on both machines, with about 2.5% fewer instructions and identical output.

    Split stores: ruled out. Rust 12% of stores, C 10%.

    What does differ, m1k_5 affine2p CIGAR on the Xeon:

    Rust C
    stores 4.32 G 5.33 G
    RFOs (stores that need the line) 566 M 628 M
    RFOs that miss L2 160 M 114 M
    store-buffer stall cycles 4.67 G 0.38 G

    Where (perf record on those events, with source lines):

    • About 60% of Rust's store-buffer stalls are on the five output stores of the AVX2 affine2p kernel (affine2p.rs, the _mm256_storeu_si256 lines).
    • That's also where its L2-missing RFOs are: about 158 M, against about 88 M in C's kernel.
    • C's few stalls are mostly in init_ends and wavefront_init.

    So Rust's kernel writes to lines that are cold in L2 about 1.8× as often as C's. The allocation width isn't the reason: both use lo − (scope + 1) … hi + (scope + 1) merged with the historic range.

    Candidates, in a 4-way A/B on both x86 machines:

    • no software prefetches;
    • prefetchw on the 5 outputs, 4 lines ahead, to claim the lines before the stores arrive;
    • slab cache coloring off. Coloring offsets every allocation by 0–23 × 64 bytes, so consecutive wavefronts are never adjacent, and the L2 prefetcher may not be able to stream across them.
  8. cjfields commented on Oct 17, 2026

    @cjfields
    MemberAuthor

    Xeon affine2p CIGAR gap: fixed in #52. The affine2p kernel now issues PREFETCHW on its 5 outputs, 32 diagonals ahead, and the x86 input prefetches are removed. On the Xeon, store-buffer stalls drop to C's level (4.63 G → 0.43 G cycles). Against C, affine2p CIGAR goes from 1.01× to 0.89×, banded from 1.28× to 0.83×, and 10 kb from 0.97× to 0.87×. The EPYC is 2–6% faster on most affine2p cases, and output is identical. Of the distances tried (32, 64 and 128 diagonals), 32 was best on the Xeon; 64 and 128 gave smaller or no gains. Cache coloring isn't worth changing. The full A/B is in the PR.

  9. cjfields commented on Oct 18, 2026

    @cjfields
    MemberAuthor

    Begin-free ends-free costs about 3× per alignment (dada2-rs, HPCBio/dada2-rs#51)

    Measured on bef359a with Compat::Corrected. Honoring begin-free gaps (#17) makes each ends-free alignment about 3× slower, on long and short reads alike. It is now the main reason WFA is slower than NW in dada2-rs; the rest of the port is not slower than before.

    Pair benchmark

    7,980 pairs (top 20 × 400 most abundant uniques), match −5, mismatch 4, gap 8, ComputeAlignment, Apple Silicon. dada2-rs settings: BandedStatic over the NW band ([−lband, +rband], widened by the length difference), free ends bounded to band + |len1 − len2|, and "begin-free on" sets the begin-free values to that same bound. Times are the median of 2 runs.

    reads metric begin-free off begin-free on on / off
    PacBio HiFi 16S, ~1,450 bp, band 32 affine (o = 0) 2.66 s 9.89 s 3.7×
    linear 1.88 s 5.68 s 3.0×
    MiSeq V4 forward, 240 bp, band 16 affine (o = 0) 0.29 s 0.80 s 2.8×
    linear 0.19 s 0.53 s 2.8×

    The MiSeq reads are all 240 bp, so begin-free there is just 16. Even a begin-free span the size of the band triples the cost. Results are unchanged: dada2-rs's 10,000-pair parity sweep gives 0 disagreements with the NW optimum in both metrics.

    Bisect (long reads, on/off at each first-parent merge since 4cfeda8): the cost appears exactly at #17 (3.4 → 11.7 s), and #36 recovered part of it. With begin-free off, 5cdcea3 was 2.76 s against 4cfeda8's 3.13 s, so nothing else regressed. C 3282f83 shows the same multiplier (20 s → 62 s), so this is the algorithm's cost as C implements it, not a port quirk.

    End-to-end effect (dada2-rs learn-errors, wall time)

    data NW WFA affine, cap 50 WFA linear, cap 50
    PacBio fixture (2 samples), k=5 9.5 s 11.0 s 10.3 s
    PacBio fixture, k=7 2.6 s 2.2 s 2.0 s
    MiSeq SOP (20 samples), forward 15.7 s 44.2 s 32.4 s
    MiSeq SOP, reverse 8.6 s 24.1 s 16.1 s

    Before #17, WFA was about 6% faster than NW on MiSeq. Short reads now show the cost most, because NW is cheap there.

    A lead, not verified

    Seeding is already lazy (endsfree_begin_k: diagonal k enters at score k·|match|). My guess is that the cost comes afterwards. Once seeded, an off-diagonal begin-free cell stalls (it extends only until its first mismatch, which is short for non-repetitive sequence) but stays in the wavefront. Every later score step then computes the full [−band, +band] range, instead of the few diagonals a near-identical pair needs. Two checks would confirm or rule this out:

    • wavefront width per step, with begin-free on vs off, for a near-identical pair;
    • whether such a diagonal can be dropped once it can no longer beat the best end candidate under the match reward.

    Reproducing

    Pair lists: the PacBio pairs come from the dada2-rs concordance fixture (dev/concordance/data/pacbio, SRR8557463/4) after primer removal and filtering. The MiSeq pairs come from the MiSeq SOP forward reads after filterAndTrim (truncLen 240). Take the 400 most abundant uniques, one sequence per line. The benchmark loop:

    probe (Rust, depends on wfa2lib-rs at bef359a)
    fn band_adjust(l1: usize, l2: usize, b: i32) -> (i32, i32) {
        if l2 > l1 { (b, b + (l2 - l1) as i32) } else if l1 > l2 { (b + (l1 - l2) as i32, b) } else { (b, b) }
    }
    fn run<const N: usize>(mut a: WavefrontAligner<N>, seqs: &[Vec<u8>], band: i32, begin_free: bool) -> f64 {
        a.alignment_scope = AlignmentScope::ComputeAlignment;
        a.set_compat(Compat::Corrected);
        let t = std::time::Instant::now();
        for i in 0..20 { for j in 0..seqs.len() { if i == j { continue; }
            let (s1, s2) = (&seqs[i], &seqs[j]);
            let (lb, rb) = band_adjust(s1.len(), s2.len(), band);
            a.set_heuristic(HeuristicStrategy::BandedStatic { min_k: -lb, max_k: rb });
            let fe = band + (s1.len() as i32 - s2.len() as i32).abs();
            let (pf, tf) = (fe.min(s1.len() as i32), fe.min(s2.len() as i32));
            let (pb, tb) = if begin_free { (pf, tf) } else { (0, 0) };
            a.set_alignment_free_ends(pb, pf, tb, tf);
            a.align_endsfree(s1, s2);
        }}
        t.elapsed().as_secs_f64()
    }
    // affine: WavefrontAligner::<3>::new(WavefrontPenalties::new_affine(AffinePenalties { match_: -5, mismatch: 4, gap_opening: 0, gap_extension: 8 }))
    // linear: WavefrontAligner::<1>::new(WavefrontPenalties::new_linear(LinearPenalties { match_: -5, mismatch: 4, indel: 8 }))

    In dada2-rs, the same comparison end to end: branch perf/wfa-remeasure-51, cargo build --release --features wfa, then learn-errors ... --align-backend wfa2 --wfa-max-edits 50 against --align-backend nw, with DADA2RS_WFA_LINEAR=1 for the linear arm.

  10. cjfields commented on Oct 18, 2026

    @cjfields
    MemberAuthor

    Begin-free cost: the lead is confirmed, but the cost is per score step, not per cell

    Follow-up to the comment above (HPCBio/dada2-rs#51), on bef359a, Apple M2. Synthetic dada2-like pairs: 20 bases × 100 mutants each, ~20% of differences indels. Penalties: linear -5,4,8 and affine -5,4,0,8; corrected mode, banded-static ±band, free ends = band.

    1. Seeded diagonals do stay in the wavefront. Cells computed per alignment, counted with a temporary probe:

    pair (2% divergent) begin-free off begin-free on
    240 bp, band 16 35 cells, 10 steps, mean width 3.5 604 cells (17×), 43 steps, mean width 14.2
    1450 bp, band 32 2,039 cells, 116 steps, mean width 17.7 18,518 cells (9×), 368 steps, mean width 50 of 65

    2. Pruning them early doesn't help. I gave corrected mode's exact pruning (#42) an upper bound from the start: the cheaper of two gapless alignments (diagonal 0, and the end diagonal after a begin gap), computed in O(n) and used only when it's inside the band. Output was identical in every case, but time was 2–18% worse. With a match reward, a cell's lower bound grows by only |match| = 5 per diagonal of distance from the end diagonal, while each difference costs 18 (2·4 + 2·5). So the seeded diagonals can't be excluded until the search is nearly over. They're legitimate low-score cells that an exact search must visit, which is consistent with C showing the same 3×. Not merged.

    3. The time follows the number of steps. At 240 bp, 2% (linear), begin-free off is ~2.7 µs and 10 steps per alignment; on, it's ~8.3 µs and 43 steps. That's about 170 ns per extra step, or ~13 ns per extra cell, while the SIMD kernel's own work is under 1 ns per cell. Seeds make the wavefront non-null at many more scores, and each non-null step pays a fixed cost. Profile (sample), begin-free on, 240 bp:

    where share
    compute_linear_step itself (setup + kernel) ~33%
    NULL fills in it (memset_pattern16) ~8%
    slab allocation in it ~5%
    align_endsfree inlined: corrected extend, per-diagonal termination checks ~40%
    banded-static cutoff, every step ~6%
    slab reap per alignment ~3%

    That's also why banded NW wins at 240 bp despite ~13× more cells: it has no per-step setup.

    Next (new checklist item): reduce per-step overhead for narrow wavefronts, keeping output identical. Starting points:

    • skip the banded cutoff when the wavefront is already inside the band;
    • tighten the NULL-fill ranges;
    • streamline the corrected extend's termination checks.

    This helps every short-read alignment, with or without begin-free. I don't expect it to erase the 3×. With the ~10% divergence DADA2's screen allows, part of this is likely the inherent price of an exact WFA against banded NW.

  11. cjfields commented on Oct 18, 2026

    @cjfields
    MemberAuthor

    Global vs ends-free on primer-trimmed amplicons (HPCBio/dada2-rs#51)

    Context, from @cjfields: dada2-rs (like DADA2) denoises PCR-amplified data from a common region. The degenerate primers are trimmed before dereplication, screening, alignment and denoising. If trimming is done properly, the reads passing the screen are effectively anchored at both the 5' and 3' ends: every pair starts and ends at the same biological position. dada2-rs uses ends-free alignment to keep the full sequences of each pair, but the anchoring suggests a different alignment model.

    What anchoring implies:

    • 5' end: begin gaps are artefacts, from imperfect trimming or spurious "sliding" alignments on divergent pairs. The true alignment starts at (0, 0).
    • 3' end: with both primers trimmed (e.g. full-length PacBio 16S), length differences come from internal indels, so a global alignment is the right model. Reads truncated with truncLen (MiSeq) aren't anchored at 3'. There an internal indel leaves an overhang, so end-free only (begin 0) is the right model.
    • Scoring: match/mismatch/gap scores are unchanged; only which end gaps are free changes. In a global alignment every candidate consumes both sequences, so the match reward doesn't change the optimum. DADA2's 5/−4/−8 and match = 0 with WFA's internal 18/21 give the same alignments, and #102's ends-free suboptimality doesn't arise.

    Speed (b83c0ea, M2, ms per 2,000 synthetic dada2-like pairs, banded ±16 / ±32, median of 5; ends-free columns in corrected mode):

    data penalties global ends-free, begin 0 ends-free, begin = band (current dada2-rs)
    240 bp, 2% linear 4.4 4.9 14.7
    affine o = 0 6.6 7.1 22.6
    240 bp, 10% linear 12.0 14.4 41.8
    affine o = 0 18.4 20.6 61.5
    1450 bp, 2% linear 43.0 53.6 173.4
    affine o = 0 60.3 70.3 227.0
    1450 bp, 10% linear 113.6 152.7 445.8
    affine o = 0 161.7 198.3 586.6
    • Global WFA is 3.3–4× faster than the current configuration, and 10–25% faster than end-free only. Linear stays about 1.4× faster than affine o = 0.
    • The begin-free ~3× found above is avoided entirely: it came from seeding begin-free diagonals, and global alignment has none.

    Suggested for dada2-rs (no wfa2lib-rs change needed; global and end-free-only are already the fast paths, and #59 applies to both):

    • Make the alignment model a setting: global (WFA align_end2end; NW with end_gap_p = gap_p), end-free only, or the current fully free ends. Keep fully free ends as the default for parity with R DADA2, whose NW frees end gaps at both ends.
    • Before validating, histogram the leading and trailing gap lengths of the current NW optima on the MiSeq SOP and PacBio fixtures. That shows how many pairs actually use terminal gaps: those are the only ones whose alignment can change.
    • Then validate global or end-free-only against R DADA2 output (ASVs, abundances) on the concordance fixtures, and compare global NW with global WFA directly. Scores should agree for every pair; CIGARs may differ between co-optimal alignments.
  12. cjfields commented on Oct 18, 2026

    @cjfields
    MemberAuthor

    Correction to my earlier lead: the begin-free cost is a step, in both compat modes

    My earlier comment guessed that seeded begin-free diagonals stay in the wavefront and widen every later step to the full band. A sweep of the begin-free limit rules that out:

    Same 7,980 real pairs, bef359a, gap-affine/linear, match −5, BandedStatic over the NW band, end-free limits unchanged (band + |Δlen|), begin-free set to t on both sequences. Best of 2.

    begin-free t 0 1 2 3 8 band + |Δlen|
    PacBio linear, Compat::C 1.78 s 3.63 s 5.23 s 5.25 s 5.35 s 5.49 s
    PacBio linear, Corrected 1.93 s 3.89 s 5.60 s 5.66 s 5.60 s 5.88 s
    PacBio affine, Compat::C 2.63 s 5.78 s 9.17 s 9.35 s 9.47 s 10.00 s
    MiSeq linear, Compat::C 0.17 s 0.34 s 0.47 s — — 0.50 s
    • A step, not a slope. t = 1 doubles the cost, and t = 2 reaches the full cost. With t = 2 only four begin-free cells are added (diagonals ±1, ±2), yet it costs as much as a 32+ bound, so the width of the seeded range isn't what costs.
    • Not corrected mode. Compat::C shows the same step, so the cost lies in the begin-free path both modes share (endsfree_begin_k / endsfree_init seeding with match < 0), not in corrected mode's end handling. My second guess (corrected-mode accounting, in Opt-in alignment end model for primer-anchored amplicons (global / end-free-only) dada2-rs#289) is wrong as well.
    • Both read lengths. MiSeq (240 bp, band 16) has the same shape as PacBio (~1,450 bp, band 32).

    A per-step count of wavefront width (or cells computed) for one near-identical pair at t = 0, 1 and 2 should show where the extra work goes.

  13. cjfields commented on Oct 19, 2026

    @cjfields
    MemberAuthor

    Begin-free cost: cause found, mostly pruned exactly (perf/begin-free-dominance)

    Follow-up to the step-cost comment above. The per-step count it suggested shows where the extra work goes.

    Cause. Each begin-free seed k is placed at score k·|m| (5, 10, … with match −5). From there its cells advance only by mismatches and gaps (+18 / +21 here, the WFA costs of DADA2's 5/−4/−8). So each seed starts a lineage of wavefronts at scores the main lineage never reaches. With only two step costs, the main lineage is non-null at 33 of ~142 scores per pair (240 bp, 2%). The offsets of 5 and 10 fill most of the gaps: 94 non-null at a limit of 2, 120 at the full bound. That is why the cost is a step: once a couple of seeds exist, nearly every score already has a wavefront. It is the same in both compat modes because C does the same.

    Most of that work is dead. On the shifted diagonals almost nothing matches. The seed lineages trail cells the main path already reached more cheaply. For example, diagonal +1 sits at offset 2 at score 23, after the main path reached offset 10 there at score 21.

    Exact pruning. A cell is dominated if a cell of the same component reached at least as far on its diagonal at an earlier score. With non-negative WFA costs (match 0), any path on from a dominated cell can be redirected from the dominating one at no extra cost: its first gap run is one base shorter, or its first diagonal step is absorbed into it. The redirected path reaches the same end with the same skipped-base charge. So dropping dominated cells leaves the optimum unchanged. C and WFA in general keep them.

    Implementation. Corrected mode, with no heuristic or BandedStatic (dada2-rs's setting). After each compute step, dominated cells are trimmed from the ends of each M/I/D wavefront. The gate: always for gap-linear, and for affine/affine2p up to a begin-free limit of 4. Beyond 4, affine's three passes per score cost more than they save on short pairs. Fusing the check into the compute kernels would lift that limit.

    Correctness.

    • Output is byte-identical to main over 1,008 CLI configurations (linear/affine/affine2p, seven span shapes, high/med/low, with and without a band).
    • The DP-oracle tests pass, plus a new amplicon-like oracle test at begin-free limits 1–5.
    • dada2-rs at rev 05bc1a3 (Opt-in alignment end model for primer-anchored amplicons (global / end-free-only) dada2-rs#289) found no differences:
      • its NW parity sweeps pass;
      • the PacBio concordance tables and error models are identical;
      • scores agree on 7,980 real pairs per dataset at limits 0/1/2/3/8/full, and at the full bound they equal NW's optimum;
      • the max-steps fallback counts are unchanged.

    Speed (dada2-rs, real pairs, linear, s):

    begin-free 0 1 2 3 8 full
    PacBio b83c0ea 1.78 3.52 5.04 5.14 5.20 5.43
    PacBio 05bc1a3 1.79 2.14 2.26 2.23 2.37 2.74
    MiSeq b83c0ea 0.163 0.308 0.423 0.429 0.446 0.454
    MiSeq 05bc1a3 0.163 0.201 0.218 0.204 0.272 0.262

    Linear ends-free at the full bound (R DADA2's model) drops from 3.0× to 1.5–1.6× the begin-free-0 cost. A small tolerance (1–3 bases) now costs 1.2–1.3×. Synthetic pairs in .claude/PERFORMANCE_LOG.md show 0.47–0.95× for linear and up to 0.54× for affine at limits ≤ 4.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions