From ed59a5129a79185b043e41ca12c2cb999ab9e004 Mon Sep 17 00:00:00 2001 From: R script Date: Sat, 4 Jul 2026 07:56:44 +0100 Subject: [PATCH] fix(sector): guard zero-Fitch-words empty-vector UB (EW analogue of PR#264) Any dataset that collapses to 0 Fitch blocks (total_words == 0) -- not just the all-hierarchy HSJ/xform case PR #264 fixed in wagner_tree/load_tip_states, but also plain equal-weights data where every character is constant or an autapomorphy -- leaves tip_states/prelim/edge_set genuinely empty std::vectors. Code that unconditionally took the address of element 0 of these empty vectors is undefined behaviour: it aborts under hardened libstdc++ assertions/ASan but is silent (out-of-bounds address, never dereferenced if the following block loop is also 0-trip) in a plain release build. Audited every zero-word-reachable site sharing this pattern and guarded each with the same have_words ternary wagner_tree already uses (never bare nullptr+0 arithmetic, which UBSan's pointer-overflow check would still flag): - ts_sector.cpp: compute_from_above_for_sector (called unconditionally on every build_reduced_dataset(), i.e. every sector reduction including the default rasStarts=1 -- the most severe of these, reachable with no non-default control), build_ras_sector (the RAS-restart start-tree builder, reachable via SearchControl(rasStarts>=2)), and build_reduced_dataset_collapsed (reachable via sectorCollapseTarget>0). - ts_prune_reinsert.cpp: build_reduced_dataset and expand_and_reinsert, reachable via pruneReinsertCycles>0 (auto-selected by the "large" preset at nTip>=120). - ts_temper.cpp: stochastic_tbr_phase (annealing), same "large"-preset reachability, needed its own top-level total_words==0 early return matching nni_search/tbr_search/ratchet_search/drift_search's existing style. - ts_tree.cpp: TreeState::save_node_state's prealloc fast path -- a shared primitive under TBR/ratchet/temper, guarded once at the source. - ts_wagner.cpp: wagner_goloboff_scores/wagner_entropy_scores (biased-Wagner tip scoring), also reachable via the "large" preset. Regression tests added to test-ts-sector.R (default control, rasStarts=3, sectorCollapseTarget=6, pruneReinsertCycles=1, annealCycles=1, all through MaximizeParsimony() on a 40-tip all-uninformative EW dataset) and test-ts-wagner.R (ts_wagner_bias_bench directly). Reachability of all 6 sites was confirmed empirically (temporary REprintf probes on each total_words==0 branch, fired under every new test, then removed). Verified the tests are non-tautological, not just reachability-proxy: a negative control reverting build_ras_sector's guard back to the unconditional &rd.data.tip_states[...] reproduces the exact predicted abort (Assertion '__n < this->size()' failed) under a local -D_GLIBCXX_ASSERTIONS -O0 -g build, with a gdb backtrace confirming build_ras_sector (ts_sector.cpp) -> search_sector -> xss_search -> run_single_replicate -> driven_search -> ts_driven_search, i.e. reachable from the public MaximizeParsimony() API. Restoring the guard reruns clean. Full targeted suite (sector, wagner, hsj, xform, anneal, driven, drift, ratchet, tbr, Concordance) passes under the same hardened build with no behavioural regression. Co-Authored-By: Claude Opus 4.8 --- src/ts_prune_reinsert.cpp | 46 +++++++++++++------ src/ts_sector.cpp | 56 ++++++++++++++++++----- src/ts_temper.cpp | 5 ++ src/ts_tree.cpp | 15 ++++-- src/ts_wagner.cpp | 8 ++++ tests/testthat/test-ts-sector.R | 81 +++++++++++++++++++++++++++++++++ tests/testthat/test-ts-wagner.R | 31 +++++++++++++ 7 files changed, 214 insertions(+), 28 deletions(-) diff --git a/src/ts_prune_reinsert.cpp b/src/ts_prune_reinsert.cpp index b2f99b8de..68a3ef092 100644 --- a/src/ts_prune_reinsert.cpp +++ b/src/ts_prune_reinsert.cpp @@ -323,11 +323,16 @@ DataSet build_reduced_dataset(const DataSet& ds, // Subset tip_states int tw = ds.total_words; red.tip_states.resize(static_cast(m) * tw); - for (int i = 0; i < m; ++i) { - int orig = tip_map[i]; - std::memcpy(&red.tip_states[static_cast(i) * tw], - &ds.tip_states[static_cast(orig) * tw], - tw * sizeof(uint64_t)); + // tw == 0 (zero Fitch blocks, e.g. every character constant/autapomorphic) + // leaves ds.tip_states empty; skip the memcpy rather than take the address + // of element 0 of an empty vector (undefined behaviour). + if (tw > 0) { + for (int i = 0; i < m; ++i) { + int orig = tip_map[i]; + std::memcpy(&red.tip_states[static_cast(i) * tw], + &ds.tip_states[static_cast(orig) * tw], + tw * sizeof(uint64_t)); + } } return red; @@ -428,18 +433,27 @@ void expand_and_reinsert( thread_local static int sa_delta_max = 0; #endif + // No Fitch words (e.g. every character constant/autapomorphic under equal + // weights): every insertion cost is identically zero. Skip the edge-set + // precompute and indirect length evaluation, which would otherwise take + // the address of element 0 of an empty vector (undefined behaviour; + // aborts under _GLIBCXX_ASSERTIONS) -- mirrors wagner_tree's guard. + const bool have_words = tw > 0; for (int tip : reinsert_order) { int new_internal = next_internal++; - const uint64_t* tip_prelim = - &ds.tip_states[static_cast(tip) * tw]; + const uint64_t* tip_prelim = have_words + ? &ds.tip_states[static_cast(tip) * tw] + : nullptr; // Exact insertion cost via directional edge sets: edge_set[D] = // combine(prelim[D], up[D]). Replaces the union-of-finals approximation // (final_[node] | final_[child]) that undercut insertion cost (~+30% Wagner // trees); mirrors the main Wagner builder. prelim is current here // (wagner_incremental_rescore maintains both prelim and final_). - compute_insertion_edge_sets(tree, ds, pr_edge_set, pr_up, pr_pre); + if (have_words) { + compute_insertion_edge_sets(tree, ds, pr_edge_set, pr_up, pr_pre); + } // Find best insertion edge via DFS from root int best_above = -1, best_below = -1; @@ -462,8 +476,11 @@ void expand_and_reinsert( if (lc < 0 || rc < 0) continue; // Evaluate edge (node, lc) - int extra = fitch_indirect_length_cached( - tip_prelim, &pr_edge_set[static_cast(lc) * tw], ds, best_extra); + int extra = have_words + ? fitch_indirect_length_cached( + tip_prelim, &pr_edge_set[static_cast(lc) * tw], ds, + best_extra) + : 0; if (extra < best_extra) { best_extra = extra; best_above = node; @@ -471,8 +488,11 @@ void expand_and_reinsert( } // Evaluate edge (node, rc) - extra = fitch_indirect_length_cached( - tip_prelim, &pr_edge_set[static_cast(rc) * tw], ds, best_extra); + extra = have_words + ? fitch_indirect_length_cached( + tip_prelim, &pr_edge_set[static_cast(rc) * tw], ds, + best_extra) + : 0; if (extra < best_extra) { best_extra = extra; best_above = node; @@ -496,7 +516,7 @@ void expand_and_reinsert( // is fully binary from root, so compute_insertion_edge_sets is exact and // safe. Tally Δ = exact_cost(E_bounded) − min_E exact_cost(E), per the // advisor: exact-suboptimality of the bounded choice, not raw edge flips. - if (best_below >= 0 && best_below != n_tip) { + if (have_words && best_below >= 0 && best_below != n_tip) { compute_insertion_edge_sets(tree, ds, sa_edge_set, sa_up, sa_pre); int exact_chosen = fitch_indirect_length_cached( tip_prelim, &sa_edge_set[static_cast(best_below) * tw], ds, INT_MAX); diff --git a/src/ts_sector.cpp b/src/ts_sector.cpp index c71b682b1..da2ac7466 100644 --- a/src/ts_sector.cpp +++ b/src/ts_sector.cpp @@ -72,8 +72,15 @@ static void compute_from_above_for_sector( // from_above[next] = fitch_join(from_above[node], prelim[sib]) // fitch_join: per-block, compute intersection; where empty, use union. - const uint64_t* sib_prelim = - &tree.prelim[static_cast(sib) * tw]; + // When tw == 0 (e.g. every character constant/autapomorphic under equal + // weights, leaving zero Fitch blocks), tree.prelim is empty and the + // block loop below never executes (ds.n_blocks == 0 too), so sib_prelim + // is never dereferenced -- but taking its address still needs guarding + // (address of element 0 of an empty vector is undefined behaviour; + // aborts under _GLIBCXX_ASSERTIONS). + const uint64_t* sib_prelim = (tw > 0) + ? &tree.prelim[static_cast(sib) * tw] + : nullptr; for (int b = 0; b < ds.n_blocks; ++b) { int off = ds.block_word_offset[b]; @@ -92,8 +99,10 @@ static void compute_from_above_for_sector( std::swap(from_above_cur, new_from_above); } - std::memcpy(from_above_out.data(), from_above_cur.data(), - tw * sizeof(uint64_t)); + if (tw > 0) { + std::memcpy(from_above_out.data(), from_above_cur.data(), + tw * sizeof(uint64_t)); + } } // ---- Conflict-guided sector selection ---- @@ -805,7 +814,11 @@ static ReducedDataset build_reduced_dataset_collapsed(const TreeState& tree, for (int i = 0; i < n_front; ++i) { // composite terminal states const int node = frontier[i]; const size_t dst = static_cast(i) * tw; - const uint64_t* src = (node < tree.n_tip) + // tw == 0 (zero Fitch blocks) leaves ds.tip_states/tree.prelim empty; + // guard the address-of (UB on an empty vector) -- the copy loop below + // is already a no-op (w < tw == 0), so src is never dereferenced. + const uint64_t* src = (tw == 0) ? nullptr + : (node < tree.n_tip) ? &ds.tip_states[static_cast(node) * tw] // real tip : &tree.prelim[static_cast(node) * tw]; // collapsed sub-clade for (int w = 0; w < tw; ++w) rd.data.tip_states[dst + w] = src[w]; @@ -921,12 +934,25 @@ static void build_ras_sector(ReducedDataset& rd, std::mt19937& rng) { // up-message buffer and preorder list each step. std::vector edge_set_up; std::vector edge_set_pre; + // When there are no Fitch words -- e.g. every character is constant or an + // autapomorphy, leaving the equal-weights dataset with zero blocks (the + // same all-uninformative case wagner_tree guards) -- rd.data.tip_states + // and edge_set are empty and every insertion cost is identically zero. + // Skip the edge-set precompute and the indirect length evaluation, which + // would otherwise take the address of element 0 of an empty vector + // (undefined behaviour; aborts under _GLIBCXX_ASSERTIONS/ASan). The DFS + // below still runs, so constraints are honoured and -- all costs being + // equal -- the first legal edge is chosen (mirrors wagner_tree's fix). + const bool have_words = rd.data.total_words > 0; for (int i = 2; i < n_real; ++i) { const int tip = order[i]; - const uint64_t* tip_prelim = - &rd.data.tip_states[static_cast(tip) * tw]; + const uint64_t* tip_prelim = have_words + ? &rd.data.tip_states[static_cast(tip) * tw] + : nullptr; - compute_insertion_edge_sets(t, rd.data, edge_set, edge_set_up, edge_set_pre); + if (have_words) { + compute_insertion_edge_sets(t, rd.data, edge_set, edge_set_up, edge_set_pre); + } int best_above = -1, best_below = -1, best_extra = INT_MAX; stack.clear(); @@ -939,14 +965,20 @@ static void build_ras_sector(ReducedDataset& rd, std::mt19937& rng) { int lc = t.left[ni]; int rc = t.right[ni]; if (lc >= 0) { - int extra = fitch_indirect_length_cached( - tip_prelim, &edge_set[static_cast(lc) * tw], rd.data, best_extra); + int extra = have_words + ? fitch_indirect_length_cached( + tip_prelim, &edge_set[static_cast(lc) * tw], + rd.data, best_extra) + : 0; if (extra < best_extra) { best_extra = extra; best_above = node; best_below = lc; } if (lc >= n_tip) stack.push_back(lc); } if (rc >= 0) { - int extra = fitch_indirect_length_cached( - tip_prelim, &edge_set[static_cast(rc) * tw], rd.data, best_extra); + int extra = have_words + ? fitch_indirect_length_cached( + tip_prelim, &edge_set[static_cast(rc) * tw], + rd.data, best_extra) + : 0; if (extra < best_extra) { best_extra = extra; best_above = node; best_below = rc; } if (rc >= n_tip) stack.push_back(rc); } diff --git a/src/ts_temper.cpp b/src/ts_temper.cpp index 41a903986..b995ae3d1 100644 --- a/src/ts_temper.cpp +++ b/src/ts_temper.cpp @@ -139,6 +139,11 @@ TemperResult stochastic_tbr_phase( double score = temper_full_rescore(tree, ds); double best_score = score; + // No informative characters: all trees have the same score. Skip the + // per-word clip/regraft scoring below, which would otherwise take the + // address of element 0 of the (empty, since total_words == 0) tree.prelim + // vector -- undefined behaviour (aborts under _GLIBCXX_ASSERTIONS). + if (ds.total_words == 0) return {best_score, score, 0, 0, 0}; const bool use_iw = std::isfinite(ds.concavity); const double eps = use_iw ? 1e-10 : 0.0; const double temperature = params.temperature; diff --git a/src/ts_tree.cpp b/src/ts_tree.cpp index 80357e218..a6a14aae4 100644 --- a/src/ts_tree.cpp +++ b/src/ts_tree.cpp @@ -52,9 +52,10 @@ void TreeState::load_tip_states(const DataSet& ds) { // T-262: bulk memcpy replaces per-element loop. Tip states occupy the // first n_tip * total_words entries of prelim/final_ (contiguous). size_t tip_bytes = static_cast(n_tip) * total_words * sizeof(uint64_t); - // All characters can be uninformative (e.g. every state a singleton), - // leaving zero blocks and a zero-length tip_states vector; its .data() - // is then permitted to be null, which memcpy's nonnull attribute forbids. + // All characters can be uninformative (e.g. every state constant or a + // singleton), leaving zero blocks and a zero-length tip_states vector; + // its .data() is then permitted to be null, which memcpy's nonnull + // attribute forbids. if (tip_bytes > 0) { std::memcpy(prelim.data(), ds.tip_states.data(), tip_bytes); std::memcpy(final_.data(), ds.tip_states.data(), tip_bytes); @@ -209,6 +210,14 @@ void TreeState::restore_prealloc_undo() { } void TreeState::save_node_state(int node) { + // No Fitch words (e.g. every character constant/autapomorphic under equal + // weights): there is no per-word state to save, so this is legitimately a + // no-op. Guard it explicitly -- the prealloc fast path below sizes its + // flat buffers to capacity * total_words, so total_words == 0 leaves them + // empty, and unconditionally memcpy'ing to/from element 0 of an empty + // vector is undefined behaviour (aborts under _GLIBCXX_ASSERTIONS). + if (total_words == 0) return; + // Fast path: use pre-allocated flat buffers (no heap allocation) if (prealloc_undo) { auto& u = *prealloc_undo; diff --git a/src/ts_wagner.cpp b/src/ts_wagner.cpp index 04c35cf6f..ada287f8b 100644 --- a/src/ts_wagner.cpp +++ b/src/ts_wagner.cpp @@ -611,6 +611,10 @@ std::vector wagner_goloboff_scores(const DataSet& ds) { int n_tip = ds.n_tips; int tw = ds.total_words; std::vector scores(n_tip, 0.0); + // No Fitch words (e.g. every character constant/autapomorphic): every tip + // is equally uninformative, so all scores are legitimately 0; skip the + // per-tip loop rather than take the address of an empty ds.tip_states. + if (tw == 0) return scores; for (int t = 0; t < n_tip; ++t) { double score = 0.0; @@ -644,6 +648,10 @@ std::vector wagner_entropy_scores(const DataSet& ds) { int n_tip = ds.n_tips; int tw = ds.total_words; std::vector scores(n_tip, 0.0); + // See wagner_goloboff_scores: no Fitch words means every score is + // legitimately 0; skip the loop rather than take the address of an empty + // ds.tip_states. + if (tw == 0) return scores; for (int t = 0; t < n_tip; ++t) { double score = 0.0; diff --git a/tests/testthat/test-ts-sector.R b/tests/testthat/test-ts-sector.R index f587bfdc9..892a00dc7 100644 --- a/tests/testthat/test-ts-sector.R +++ b/tests/testthat/test-ts-sector.R @@ -247,3 +247,84 @@ test_that("sector_diag with NA characters returns consistent scores", { expect_true(diag$clade_size >= 2) expect_true(diag$n_sector_tips == diag$clade_size + 1L) }) + +# ===== All-uninformative EW data: zero Fitch words (regression) ============= +# When every equal-weights character is constant or an autapomorphy (no state +# shared by >1 taxon), simplify_patterns removes every character, leaving +# DataSet::total_words == 0 and empty per-word state vectors (tip_states, +# tree.prelim, edge_set). Several sectorial code paths took the address of +# element 0 of these empty vectors -- undefined behaviour that aborts under +# hardened libstdc++ assertions / ASan -- including build_ras_sector() (the +# RAS-restart start-tree builder, only reached when rasStarts > 1) and +# compute_from_above_for_sector()/build_reduced_dataset_collapsed() (reached +# via ordinary rss_search/xss_search). 40 tips clears the sectorMinSize*2 +# gate (default sectorMinSize = 6) so a sector is always extracted. + +make_uninformative_ew <- function(n, n_char = 5L, seed) { + set.seed(seed) + mkchar <- function() { + states <- rep(1L, n) + # 6 tips get unique singleton states (2..7): no state is shared by more + # than one taxon, so every character is parsimony-uninformative. + states[sample(seq_len(n), 6L)] <- seq(2L, 7L) + states + } + mat <- vapply(seq_len(n_char), function(i) mkchar(), integer(n)) + rownames(mat) <- paste0("t", seq_len(n)) + MatrixToPhyDat(mat) +} + +test_that("Sectorial search handles all-uninformative EW data (zero Fitch words)", { + ds <- make_uninformative_ew(40L, seed = 1) + tree <- PectinateTree(names(ds)) + + res <- MaximizeParsimony(ds, tree = tree, maxReplicates = 2L, + targetHits = 1L, verbosity = 0L) + expect_s3_class(res[[1]], "phylo") + expect_equal(length(res[[1]]$tip.label), 40L) +}) + +test_that("build_ras_sector handles zero Fitch words (rasStarts > 1)", { + ds <- make_uninformative_ew(40L, seed = 1) + tree <- PectinateTree(names(ds)) + ctrl <- SearchControl(rasStarts = 3L) + + res <- MaximizeParsimony(ds, tree = tree, maxReplicates = 2L, + targetHits = 1L, verbosity = 0L, control = ctrl) + expect_s3_class(res[[1]], "phylo") + expect_equal(length(res[[1]]$tip.label), 40L) +}) + +test_that("build_reduced_dataset_collapsed handles zero Fitch words", { + ds <- make_uninformative_ew(40L, seed = 1) + tree <- PectinateTree(names(ds)) + ctrl <- SearchControl(sectorCollapseTarget = 6L) + + res <- MaximizeParsimony(ds, tree = tree, maxReplicates = 2L, + targetHits = 1L, verbosity = 0L, control = ctrl) + expect_s3_class(res[[1]], "phylo") + expect_equal(length(res[[1]]$tip.label), 40L) +}) + +test_that("expand_and_reinsert (prune-reinsert) handles zero Fitch words", { + ds <- make_uninformative_ew(40L, seed = 1) + tree <- PectinateTree(names(ds)) + ctrl <- SearchControl(pruneReinsertCycles = 1L) + + res <- MaximizeParsimony(ds, tree = tree, maxReplicates = 2L, + targetHits = 1L, verbosity = 0L, control = ctrl) + expect_s3_class(res[[1]], "phylo") + expect_equal(length(res[[1]]$tip.label), 40L) +}) + +test_that("stochastic_tbr_phase (annealing) handles zero Fitch words", { + ds <- make_uninformative_ew(40L, seed = 1) + tree <- PectinateTree(names(ds)) + ctrl <- SearchControl(annealCycles = 1L, annealPhases = 2L, + annealTStart = 5, annealTEnd = 0) + + res <- MaximizeParsimony(ds, tree = tree, maxReplicates = 2L, + targetHits = 1L, verbosity = 0L, control = ctrl) + expect_s3_class(res[[1]], "phylo") + expect_equal(length(res[[1]]$tip.label), 40L) +}) diff --git a/tests/testthat/test-ts-wagner.R b/tests/testthat/test-ts-wagner.R index d047b24cc..783e22c09 100644 --- a/tests/testthat/test-ts-wagner.R +++ b/tests/testthat/test-ts-wagner.R @@ -515,3 +515,34 @@ test_that("addition_order length/range/duplicate errors cleanly (T-323)", { "addition_order" ) }) + +# ===== All-uninformative data: zero Fitch words (regression) ================ +# When every character is constant or an autapomorphy (no state shared by +# more than one taxon), simplify_patterns removes every character, leaving +# DataSet::total_words == 0 and empty per-word state vectors. +# wagner_goloboff_scores()/wagner_entropy_scores() took the address of +# element 0 of the (empty) ds.tip_states vector -- undefined behaviour that +# aborts under hardened libstdc++ assertions / ASan. + +test_that("wagner_goloboff_scores/wagner_entropy_scores handle zero Fitch words", { + n <- 10L + set.seed(1) + # 3 characters, each a distinct permutation of 1..n: no state is shared by + # more than one taxon in any character, so every character is + # parsimony-uninformative (multi-character, avoiding the single-character + # vapply/t() degenerate case in prep_pd()). + mat <- vapply(seq_len(3L), function(i) as.character(sample(seq_len(n))), + character(n)) + rownames(mat) <- paste0("t", seq_len(n)) + pd <- MatrixToPhyDat(mat) + d <- prep_pd(pd) + + result <- TreeSearch:::ts_wagner_bias_bench( + d$contrast, d$tip_data, d$weight, d$levels, + min_steps = integer(0), concavity = -1, + bias = 1L, temperature = 1, n_reps = 1L, run_tbr = FALSE + ) + + expect_true(all(result$goloboff_scores == 0)) + expect_true(all(result$entropy_scores == 0)) +})