diff --git a/man/figures/Stemwardness.png b/man/figures/Stemwardness.png index 3b49ff717..d2954f608 100644 Binary files a/man/figures/Stemwardness.png and b/man/figures/Stemwardness.png differ diff --git a/src/splits_to_tree.cpp b/src/splits_to_tree.cpp index 9984170e3..cbd3e5aba 100644 --- a/src/splits_to_tree.cpp +++ b/src/splits_to_tree.cpp @@ -7,8 +7,8 @@ using namespace Rcpp; using namespace TreeTools; inline void insert_ancestor(const int16 tip, const int16 *next_node, - std::array& parent, - std::array& patriarch) { + int16* parent, + int16* patriarch) { if (patriarch[tip]) { parent[patriarch[tip]] = *next_node; } else { @@ -19,9 +19,6 @@ inline void insert_ancestor(const int16 tip, const int16 *next_node, // [[Rcpp::export]] IntegerMatrix splits_to_edge(const RawMatrix splits, const IntegerVector nTip) { - if (double(nTip[0]) > 2048) { - Rcpp::stop("This many leaves are not (yet) supported."); - } const int16 n_tip = int16(nTip[0]); if (splits.nrow() == 0) { IntegerMatrix ret(n_tip, 2); @@ -32,12 +29,38 @@ IntegerMatrix splits_to_edge(const RawMatrix splits, const IntegerVector nTip) { return ret; } const SplitList x(splits); - alignas(64) std::array parent{}; - alignas(64) std::array patriarch{}; - - std::array split_order; - std::iota(split_order.begin(), split_order.begin() + x.n_splits, 0); - std::sort(split_order.begin(), split_order.begin() + x.n_splits, + + // Decide whether to use stack or heap allocation based on tree size + const bool use_heap = (n_tip > SL_MAX_TIPS) || (x.n_splits > SL_MAX_SPLITS); + + // Stack allocation for small trees (fast path) + alignas(64) std::array stack_parent{}; + alignas(64) std::array stack_patriarch{}; + + // Heap allocation for large trees + std::vector heap_parent; + std::vector heap_patriarch; + + // Pointers to active storage + int16* parent; + int16* patriarch; + + if (use_heap) { + const size_t parent_size = static_cast(n_tip) + + static_cast(x.n_splits); + heap_parent.resize(parent_size, 0); + heap_patriarch.resize(n_tip, 0); + parent = heap_parent.data(); + patriarch = heap_patriarch.data(); + } else { + parent = stack_parent.data(); + patriarch = stack_patriarch.data(); + } + + // Allocate split_order appropriately + std::vector split_order(x.n_splits); + std::iota(split_order.begin(), split_order.end(), 0); + std::sort(split_order.begin(), split_order.end(), [&in_split = x.in_split](int16 a, int16 b) { return in_split[a] > in_split[b]; }); diff --git a/tests/testthat/test-Splits.R b/tests/testthat/test-Splits.R index ad130c6cb..974a02e3e 100644 --- a/tests/testthat/test-Splits.R +++ b/tests/testthat/test-Splits.R @@ -512,13 +512,15 @@ test_that("Split combination", { #TODO: Fully test splits with large (> 8 tip) trees }) -test_that("as.phylo.Splits() fails gracefully", { - expect_no_error(as.Splits(BalancedTree(3000))) +test_that("as.phylo.Splits() supports large trees", { + tree3000 <- BalancedTree(3000) + expect_no_error(splits3000 <- as.Splits(tree3000)) + result <- as.phylo(splits3000) - expect_error( - as.phylo(as.Splits(BalancedTree(3000))), - "many leaves are not .*supported" - ) + # Verify it's a valid phylo object + expect_s3_class(result, "phylo") + expect_equal(length(result$tip.label), 3000) + expect_equal(NTip(result), 3000) }) test_that("as.phylo.Splits()", {