Corrected mode: Z-drop ends (#22), X-drop/Z-drop score with match == 0 (#55), BiWFA

.claude/FIXED_BUGS.md

@@ -483,8 +483,89 @@ known C bugs; the default `Compat::C` stays output-identical to C. It covers:
   source was pruned. The CIGAR merges them into one gap, so it rescores
   below the WFA score (never above): e.g. pair 298, affine `100,5`, 242 vs
   248. C's valid CIGARs show the same (affine2p `50,10`: 5/305 pairs).
-- **Tests:** `tests/corrected_xdrop.rs` (pair 201,
-  `tests/data/issue21_xdrop.seq`, all three gap metrics).
+- **Tests:** `tests/corrected_xdrop.rs` (all three gap metrics). First pair
+  201 with default penalties; since 14f changes X-drop's scores at
+  `match == 0`, pair 305 (`tests/data/issue21_xdrop.seq`) with match-reward
+  penalties and X-drop `20,2`.
+
+### 14e. Z-drop ends backtraced at the wrong score (HPCBio/wfa2lib-rs#22)
+
+- **C bug:** on a Z-drop, `wavefront_unialign_terminate` backtraces the best
+  cell's diagonal and offset (found at `max_wf_score`) starting at the
+  *current* score, where that cell needn't exist.
+  - High memory: the backtrace finds no source, stops, and pads with
+    unchecked `M`s plus leading indels. With default penalties that is
+    usually max-trimmed to empty; with `match < 0` the CIGARs are invalid
+    (e.g. `2I58M` on a 2 kb, 8%-error pair; up to 300/300 invalid on
+    `long.seq`, affine `-3,5,3,2`, Z-drop `100,5`).
+  - Med/low: the cell at the current score is often dead, and its
+    `bt_pcigar/bt_prev` were never initialized (leftovers of earlier
+    alignments in C; zeroes in Rust), so C's high and med results differ and
+    C's med depends on what it aligned before.
+- **Fix (corrected mode):** backtrace the best cell at `max_wf_score`
+  (`backtrace_zdrop_end`): high memory from the stored wavefronts; med/low
+  from the cell's `(pcigar, prev)` kept when it became the best
+  (`zdrop_best_bt`), pinned through buffer compaction (which now pins a
+  slice: the ends-free candidate and the Z-drop cell).
+- **Result:** Z-drop `20,1` / `50,1` / `50,10` / `100,5` × default and
+  `-3,5,3,2` × global / extension / ends-free × {`wfa.utest.seq`,
+  `long.seq`}: every partial CIGAR valid, and high = med except co-optimal
+  ties that score the same (an empty CIGAR vs `1M`/`2M`, score 0). C mode
+  is unchanged.
+- **Note:** with `match == 0`, C's Z-drop score makes the start the best
+  cell, so this fix alone leaves most such partial alignments empty; 14f
+  fixes that.
+- **BiWFA:** a base sub-alignment that Z-drops is backtraced the same way
+  in C (current score), and also starts in the sub-problem's end component,
+  which can be I/D; the best cell is an M cell. Corrected mode starts in M
+  at the best cell's score. BiWFA Z-drop `20,1` / `50,10` / `100,5` ×
+  default and `-3,5,3,2` × {`wfa.utest.seq`, `long.seq`}: no invalid
+  CIGARs (C: up to 64), and with default penalties far fewer empty results
+  (`long.seq` `50,10`: 0, C: 300).
+- **Tests:** `tests/corrected_zdrop.rs` (`long.seq` pair 1,
+  `tests/data/issue22_zdrop.seq`; a generated pair whose Z-drop comes after
+  two compaction probes, which fails without the pin; BiWFA on `long.seq`
+  pair 21, `tests/data/issue22_zdrop_biwfa.seq`).
+
+### 14f. X-drop/Z-drop with match == 0: sign-flipped SW score (HPCBio/wfa2lib-rs#55)
+
+- **C bug:** a regression from 51fb1aa (2023-03-19). X-drop and Z-drop score
+  cells with `(swg_match·(v+h) − wf_score)/2`. Before 51fb1aa, `match == 0`
+  used `(v+h) − wf_score` ("use match-score=1"); since then it is
+  `swg_match = -1`, so making progress lowers the score.
+  - X-drop prunes the diagonals that progress most: with default penalties
+    it abandons nearly every divergent pair (`long.seq`: 0–81 of 300
+    optimal, typically 0).
+  - Z-drop's best cell is the start, so its partial alignment is empty or
+    the leading matches (hidden in C by 14e).
+- **Fix (corrected mode):** `heuristic::swg_match` returns +1 for
+  `match == 0` (C: -1), passed to the heuristics in place of the match
+  penalty. The shared formula keeps the `/2`, so +1 has half the scale of
+  the pre-51fb1aa special case.
+- **Result:** X-drop with default penalties is optimal on 300/300 `long.seq`
+  pairs for affine/linear/affine2p at `50,10` and `100,5` (293–300 at
+  `50,1`); Z-drop gives real partial alignments (e.g. `50,10` on
+  `long.seq`: the full optimal alignment), high = med.
+- **Tests:** `test_compute_sw_scores_c_formula` (both stand-ins).
+
+### 14g. Corrected mode didn't reach BiWFA's sub-aligners (Rust bug)
+
+- **Symptom:** under BiWFA, corrected mode gave C's output: the #21 fix
+  (14d) and the #55 fix (14f) never applied. #54's "BiWFA identical to C"
+  in corrected mode was this, not parity. E.g. BiWFA X-drop `50,10`,
+  default penalties, `long.seq`: 107 invalid CIGARs and 2/300 optimal
+  scores in both modes.
+- **Root cause:** `take_biwfa_subs` hands out the forward, reverse and base
+  sub-aligners, which run the heuristics and the base backtrace, but they
+  were created with `Compat::C` and never given the parent's mode (only its
+  heuristic).
+- **Fix:** `take_biwfa_subs` sets their `compat` to the parent's on every
+  call (`set_compat` can change it between alignments).
+- **Result:** BiWFA X-drop `50,10` / `100,5` / `20,2`, default and
+  `-3,5,3,2`, `wfa.utest.seq` and `long.seq`: no invalid CIGARs in
+  corrected mode (C mode: up to 213/300); with default penalties 168–303
+  optimal (C mode: 0–240).
+- **Tests:** `tests/corrected_xdrop.rs` checks BiWFA too.
 
 ## 15. Empty pattern and text: CIGAR offset underflow (debug panic, HPCBio/wfa2lib-rs#41)
 

.claude/PLAN.md

@@ -43,6 +43,7 @@ tests/
   cli_output.rs                -- CLI output tests (expected outputs from C)
   corrected_endsfree.rs        -- Compat::Corrected ends-free vs a DP oracle
   corrected_xdrop.rs           -- Compat::Corrected X-drop CIGAR validity (#21)
+  corrected_zdrop.rs           -- Compat::Corrected Z-drop ends (#22)
   test_vs_reference.rs         -- library vs C reference outputs
   wfa_utest.sh                 -- stock C unit tests on the Rust CLI
   data/                        -- CLI test fixtures

.claude/PROGRESS.md

@@ -70,7 +70,13 @@ C changes after the port baseline that affect behaviour (HPCBio/wfa2lib-rs#10):
   modes. Root cause: in `med`/`low`, C backtraces from the end diagonal at
   the current score, often a dead cell whose chain comes from uninitialized
   memory, so C's output depends on earlier alignments. Rust reads such
-  cells as zeroes (deterministic; matches C on fresh memory).
+  cells as zeroes (deterministic; matches C on fresh memory). C's high mode
+  backtraces at the same wrong score (empty or invalid CIGARs). **Fixed in
+  corrected mode** (FIXED_BUGS 14e).
+- HPCBio/wfa2lib-rs#55: X-drop/Z-drop with `match == 0` use a sign-flipped
+  Smith-Waterman score (C regression from 51fb1aa): X-drop abandons most
+  divergent pairs, Z-drop's best cell is the start. **Fixed in corrected
+  mode** (FIXED_BUGS 14f).
 - HPCBio/wfa2lib-rs#23: score-only score for unreachable alignments (match < 0)
   is evaluated at C's unset end position (signed overflow); reproduced.
 - HPCBio/wfa2lib-rs#24: edit/indel + banded-static never terminates in C;
@@ -190,3 +196,13 @@ reports, inbound C ports) is tracked in HPCBio/wfa2lib-rs#19.
 - **2026-10-17:** HPCBio/wfa2lib-rs#21 fixed in corrected mode: a wavefront
   trimmed empty by a heuristic cutoff is null for the backtrace too, as for
   compute. X-drop with steps > 1 now gives valid CIGARs; scores unchanged.
+- **2026-10-17:** HPCBio/wfa2lib-rs#22 fixed in corrected mode: a Z-drop end
+  is backtraced from the best cell at the score it was found at (its
+  backtrace kept and pinned for med/low). All memory modes agree and every
+  partial CIGAR is valid; C's high mode had invalid ones with match < 0.
+- **2026-10-17:** HPCBio/wfa2lib-rs#55 (found while fixing #22): with
+  `match == 0`, C's X-drop/Z-drop score has its sign flipped since 51fb1aa.
+  Corrected mode uses +1: X-drop with default penalties goes from ~0/300 to
+  300/300 optimal on `long.seq`.
+- **2026-10-17:** corrected mode now reaches BiWFA's sub-aligners (a Rust bug
+  from #54, FIXED_BUGS 14g): the #21 and #55 fixes apply under BiWFA too.

README.md

@@ -189,8 +189,12 @@ By default output is identical to C, including results C gets wrong. With
 `--no-c-compat` (library: `set_compat(Compat::Corrected)`), known C bugs are
 fixed instead. Currently this covers ends-free alignment with a match reward
 (`match < 0`), where C can return a suboptimal alignment (corrected mode
-returns the optimal one), and X-drop with `steps_between_cutoffs > 1`, where
-C can return an invalid CIGAR. See `Compat` in the API docs.
+returns the optimal one), X-drop with `steps_between_cutoffs > 1`, where C
+can return an invalid CIGAR, and Z-drop, where C backtraces the partial
+alignment from the wrong score (empty or invalid CIGARs, and results that
+differ between memory modes). With `match == 0` (the default penalties), it
+also fixes the sign of the score X-drop and Z-drop use, which in C makes
+X-drop abandon most divergent pairs. See `Compat` in the API docs.
 
 ## Project structure
 

src/aligner.rs

@@ -81,6 +81,22 @@ pub enum MemoryMode {
 ///   are unchanged. This applies to any heuristic that trims a wavefront
 ///   to empty.
 ///
+/// - **Z-drop ends** (HPCBio/wfa2lib-rs#22): C backtraces the partial
+///   alignment from the best cell's diagonal and offset but at the current
+///   score, where that cell needn't exist. High memory then finds no
+///   source and returns an empty CIGAR or one with unchecked matches
+///   (invalid with `match < 0`); med/low read a dead cell whose backtrace
+///   C never initialized, so C's result depends on earlier alignments.
+///   Corrected mode backtraces the best cell at the score it was found at,
+///   keeping its backtrace for med/low, so all memory modes give a valid
+///   partial alignment with the same score.
+///
+/// - **X-drop and Z-drop with `match == 0`** (HPCBio/wfa2lib-rs#55): their
+///   Smith-Waterman score needs a stand-in match score, and C's is -1, so
+///   progress lowers the score: X-drop prunes the diagonals that progress
+///   most (abandoning most divergent pairs), and Z-drop's best cell is the
+///   start. Corrected mode uses +1, as C did before its commit 51fb1aa.
+///
 /// Otherwise, gap-linear/affine/affine2p with `match == 0` give identical
 /// results in both modes.
 #[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
@@ -172,6 +188,9 @@ pub struct WavefrontAligner<const N: usize> {
     pub compat: Compat,
     /// Corrected-mode ends-free end candidate (see [`Compat::Corrected`]).
     endsfree_best: Option<EndsfreeCandidate>,
+    /// Corrected-mode Z-drop best cell's piggyback backtrace `(pcigar, prev)`
+    /// (med/low), kept from the score it was found at (HPCBio/wfa2lib-rs#22).
+    zdrop_best_bt: (Pcigar, BtBlockIdx),
 
     // Heuristic
     pub heuristic: HeuristicStrategy,
@@ -423,6 +442,7 @@ impl<const N: usize> WavefrontAligner<N> {
             extension: false,
             compat: Compat::C,
             endsfree_best: None,
+            zdrop_best_bt: (0, BT_BLOCK_IDX_NULL),
             heuristic: HeuristicStrategy::None,
             heuristic_state: HeuristicState::default(),
             memory_mode: MemoryMode::High,
@@ -467,12 +487,16 @@ impl<const N: usize> WavefrontAligner<N> {
         if buf.get_size_used() <= self.max_memory_compact {
             return;
         }
-        let mut pinned = self.endsfree_best.map_or(BT_BLOCK_IDX_NULL, |b| b.bt.1);
+        let mut pinned = [
+            self.endsfree_best.map_or(BT_BLOCK_IDX_NULL, |b| b.bt.1),
+            self.zdrop_best_bt.1,
+        ];
         self.wf_components
             .compact_bt_buffer(score as usize, &mut pinned);
         if let Some(b) = self.endsfree_best.as_mut() {
-            b.bt.1 = pinned;
+            b.bt.1 = pinned[0];
         }
+        self.zdrop_best_bt.1 = pinned[1];
         let max_partial_compacts = if self.memory_mode == MemoryMode::Low {
             1
         } else {
@@ -629,14 +653,9 @@ impl<const N: usize> WavefrontAligner<N> {
             && self.alignment_end_pos.offset != OFFSET_NULL
         {
             // Z-dropped: partial alignment up to the best cell, max-trimmed
-            // (C's `wavefront_unialign_terminate`, which backtraces with the
-            // current score).
+            // (C's `wavefront_unialign_terminate`).
             if self.alignment_scope == AlignmentScope::ComputeAlignment {
-                self.backtrace_from(
-                    score,
-                    self.alignment_end_pos.k,
-                    self.alignment_end_pos.offset,
-                );
+                self.backtrace_zdrop_end(score);
                 self.maxtrim_cigar();
             }
             self.status = AlignStatus::Partial;
@@ -737,6 +756,8 @@ impl<const N: usize> WavefrontAligner<N> {
                         self.backtrace_from(score, k, offset);
                     }
                     self.cigar.score = best.total;
+                } else if unreachable {
+                    self.backtrace_zdrop_end(score);
                 } else {
                     self.backtrace_from(
                         score,
@@ -848,6 +869,7 @@ impl<const N: usize> WavefrontAligner<N> {
         };
         self.heuristic_state.reset(&self.heuristic);
         self.endsfree_best = None;
+        self.zdrop_best_bt = (0, BT_BLOCK_IDX_NULL);
 
         // Compute allocation limits based on sequence lengths.
         // k ranges from -plen to +tlen, so we need wavefronts that cover this
@@ -1164,6 +1186,24 @@ impl<const N: usize> WavefrontAligner<N> {
         }
     }
 
+    /// Backtrace a Z-drop end, the best cell in `alignment_end_pos`, at the
+    /// current `score`, as C does. The cell was found at an earlier score
+    /// and needn't exist at `score`: high mode then finds no source (an empty
+    /// CIGAR, or unchecked matches), and med/low read a dead cell
+    /// (HPCBio/wfa2lib-rs#22). Corrected mode backtraces it at the score it
+    /// was found at, from the wavefronts (high) or from its backtrace kept
+    /// then (med/low).
+    fn backtrace_zdrop_end(&mut self, score: i32) {
+        let WavefrontPos { score: best_score, k, offset } = self.alignment_end_pos;
+        if self.compat == Compat::C {
+            self.backtrace_from(score, k, offset);
+        } else if self.wf_components.bt_piggyback {
+            self.backtrace_piggyback_bt(self.zdrop_best_bt, best_score, k, offset);
+        } else {
+            self.backtrace_from(best_score, k, offset);
+        }
+    }
+
     /// Backtrace the CIGAR from M-wavefront cell (`k`, `offset`), starting at
     /// `score` (C's `wavefront_unialign_terminate`).
     fn backtrace_from(&mut self, score: i32, k: i32, offset: WfOffset) {
@@ -1306,6 +1346,7 @@ impl<const N: usize> WavefrontAligner<N> {
 
         let wf = unsafe { &mut *m_ptr };
         let (lo_base, hi_base) = (wf.lo, wf.hi);
+        let max_wf_score_base = self.heuristic_state.max_wf_score;
         let zdropped = heuristic::heuristic_cutoff(
             &self.heuristic,
             &mut self.heuristic_state,
@@ -1314,10 +1355,22 @@ impl<const N: usize> WavefrontAligner<N> {
             distances_base_k,
             plen,
             tlen,
-            self.penalties.match_,
+            heuristic::swg_match(self.penalties.match_, self.compat),
             score,
         );
 
+        // Corrected mode: keep the backtrace of a new Z-drop best cell, which
+        // the end backtrace needs after its wavefront is gone (med/low,
+        // HPCBio/wfa2lib-rs#22).
+        let st = &self.heuristic_state;
+        if self.compat == Compat::Corrected
+            && self.wf_components.bt_piggyback
+            && st.max_wf_score == score
+            && max_wf_score_base != score
+        {
+            self.zdrop_best_bt = wf.bt_or_zero(st.max_sw_score_k);
+        }
+
         if zdropped {
             // C's `wavefront_heuristic_zdrop`: the alignment ends at the best
             // cell seen; it is backtraced as a partial alignment.
@@ -2301,11 +2354,17 @@ impl<const N: usize> WavefrontAligner<N> {
             b.alignment_scope = AlignmentScope::ComputeAlignment;
             self.biwfa_base = Some(Box::new(b));
         }
-        (
+        let (mut fwd, mut rev, mut base) = (
             self.biwfa_fwd.take().unwrap(),
             self.biwfa_rev.take().unwrap(),
             self.biwfa_base.take().unwrap(),
-        )
+        );
+        // The sub-aligners run the heuristics and the base backtrace, so
+        // they follow this aligner's mode (set_compat can change it later).
+        fwd.compat = self.compat;
+        rev.compat = self.compat;
+        base.compat = self.compat;
+        (fwd, rev, base)
     }
 
     /// End-to-end alignment using a custom match function (Lambda mode).
@@ -2354,6 +2413,7 @@ impl<const N: usize> WavefrontAligner<N> {
             extension: false,
             compat: Compat::C,
             endsfree_best: None,
+            zdrop_best_bt: (0, BT_BLOCK_IDX_NULL),
             heuristic: HeuristicStrategy::None,
             heuristic_state: HeuristicState::default(),
             memory_mode: MemoryMode::High,
@@ -2744,7 +2804,14 @@ impl<const N: usize> WavefrontAligner<N> {
                 Ok(score)
             }
             AlignStatus::EndUnreachable if pos.offset != OFFSET_NULL => {
-                base.backtrace_components(score, pos.k, pos.offset, component_begin, component_end);
+                // C starts at the current score, in `component_end`;
+                // corrected mode at the best cell's score, in M (the cell is
+                // an M cell; HPCBio/wfa2lib-rs#22, as `backtrace_zdrop_end`).
+                if base.compat == Compat::C {
+                    base.backtrace_components(score, pos.k, pos.offset, component_begin, component_end);
+                } else {
+                    base.backtrace_components(pos.score, pos.k, pos.offset, component_begin, ComponentType::M);
+                }
                 base.maxtrim_cigar();
                 Err(BiStatus::Unattainable)
             }

src/components.rs

@@ -226,10 +226,11 @@ impl<const N: usize> WavefrontComponents<N> {
     /// C's `wavefront_components_compact_bt_buffer`: keep only the
     /// backtrace blocks reachable from live cells of the wavefronts still in
     /// scope (scores `score - max_score_scope + 1 ..= score`), compacting
-    /// the buffer and translating the cells' block indices. `pinned` is an
-    /// extra chain to keep (a corrected-mode ends-free end candidate, or
-    /// `BT_BLOCK_IDX_NULL`); it is translated too.
-    pub fn compact_bt_buffer(&mut self, score: usize, pinned: &mut BtBlockIdx) {
+    /// the buffer and translating the cells' block indices. `pinned` are
+    /// extra chains to keep (corrected mode: the ends-free end candidate and
+    /// the Z-drop best cell; `BT_BLOCK_IDX_NULL` if unused); they are
+    /// translated too.
+    pub fn compact_bt_buffer(&mut self, score: usize, pinned: &mut [BtBlockIdx]) {
         let Some(buf) = self.bt_buffer.as_mut() else {
             return;
         };
@@ -251,7 +252,9 @@ impl<const N: usize> WavefrontComponents<N> {
                 }
             }
         }
-        buf.mark_chain(*pinned, &mut bitmap);
+        for &p in pinned.iter() {
+            buf.mark_chain(p, &mut bitmap);
+        }
         bitmap.update_counters();
         let total = buf.compact_marked(&bitmap);
         // Translate the live cells' indices (C assumes they point past the
@@ -268,8 +271,10 @@ impl<const N: usize> WavefrontComponents<N> {
                 }
             }
         }
-        if *pinned != BT_BLOCK_IDX_NULL && *pinned >= compacted {
-            *pinned = compacted + bitmap.erank(*pinned) as BtBlockIdx;
+        for p in pinned.iter_mut() {
+            if *p != BT_BLOCK_IDX_NULL && *p >= compacted {
+                *p = compacted + bitmap.erank(*p) as BtBlockIdx;
+            }
         }
         buf.num_compacted_blocks = total;
     }