fix bug in curtail()

Author: kingaa

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: phylopomp
 Type: Package
 Title: Phylodynamic Inference for POMP Models
-Version: 0.15.1.4
-Date: 2026-01-26
+Version: 0.15.2.0
+Date: 2026-01-29
 Authors@R: c(person(given=c("Aaron","A."),family="King",role=c("aut","cre"),email="kingaa@umich.edu",comment=c(ORCID="0000-0001-6159-3207")),
 	  person(given=c("Qianying"),family="Lin",role=c("aut"),comment=c(ORCID="0000-0001-8620-9910"))
 	  )

inst/NEWS

@@ -8,6 +8,10 @@ _C_h_a_n_g_e_s _i_n '_p_h_y_l_o_p_o_m_p' _v_e_r_s_i_o_n
         • New model which is identical to ‘TwoSpecies’ except in that
           lineages that terminate are replaced by _ghosts_.
 
+        • Bug fix in ‘curtail()’.
+
+        • Refactoring of ‘lbdp_exact()’.
+
 _C_h_a_n_g_e_s _i_n '_p_h_y_l_o_p_o_m_p' _v_e_r_s_i_o_n _0._1_4:
 
         • The ‘curtail()’ function now allows one to curtail a

inst/NEWS.Rd

@@ -5,6 +5,8 @@
     \item New model in which three distinct strains compete for susceptibles from a single pool.
     See \code{Strains}.
     \item New model which is identical to \code{TwoSpecies} except in that lineages that terminate are replaced by \emph{ghosts}.
+    \item Bug fix in \code{curtail()}.
+    \item Refactoring of \code{lbdp_exact()}.
   }
 }
 \section{Changes in \pkg{phylopomp} version 0.14}{

src/genealogy.h

@@ -540,7 +540,7 @@ class genealogy_t : public nodeseq_t {
     if (!empty() && troot > timezero()) {
       node_t *p = front();
       node_t *q;
-      while (!empty() && p->slate <= troot) {
+      while (!empty() && p->slate < troot) {
         ball_t *b;
         assert(p->holds_own());
         while (p->size() > 1) {
@@ -552,12 +552,16 @@ class genealogy_t : public nodeseq_t {
           case black:
             q = make_node(b->deme());
             q->slate = troot;
-            b->holder() = q;
             q->insert(b); p->erase(b);
+            b->holder() = q;
             push_back(q);
             break;
           case green:
             q = b->child();
+            if (q == p) {
+              b = p->first_ball();
+              q = b->child();
+            }
             if (q->slate < troot) {
               q->insert(b); p->erase(b);
               b->holder() = q;

src/nodeseq.h

@@ -92,10 +92,12 @@ class nodeseq_t : public std::list<node_t*> {
 private:
 
   //! Order relation among nodes.
-  //! Nodes should be ordered by time, then by unique name.
+  //! An ancestor node should always come before its descendants.
+  //! Nodes should be ordered by time, then arbitrarily.
   static bool compare (node_t* p, node_t* q) {
     return (p->slate < q->slate) ||
-      ((p->slate == q->slate) && (p->uniq < q->uniq));
+      ((p->slate == q->slate) &&
+       ((p==q->green_ball()->holder()) || (p->uniq < q->uniq)));
   };
 
 public:

src/pocket.h

@@ -120,6 +120,10 @@ class pocket_t : public std::set<ball_t*,ball_order> {
     }
     return result;
   };
+  //! retrieve the first ball
+  ball_t* first_ball (void) const {
+    return *cbegin();
+  };
   //! retrieve the last ball
   ball_t* last_ball (void) const {
     return *crbegin();

tests/parse.Rout.save

@@ -136,7 +136,7 @@ Error : in 'scan_tree': premature end of Newick string.
 +   r"{((o_9_1:1.000000,b_3_2:0.500000,o_1_3:1.000000)m_0_0:0.000000);}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 620): invalid Newick label: expected one of 'b','g','m', or 'o', got ';'.
+Error : in 'scan_label' (genealogy.h line 624): invalid Newick label: expected one of 'b','g','m', or 'o', got ';'.
 > 
 > try(
 +   r"{(o_9_1:1.000000,b_3_2:0.500000,o_1_3:1.000000)m_0_0:0.000000);}" |>
@@ -167,13 +167,13 @@ Error : in 'parse': invalid Newick string.
 +   r"{(o_9_1:1.000000)_g_1_3:1.000000)m_0_0:0.000000;}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 620): invalid Newick label: expected one of 'b','g','m', or 'o', got '_'.
+Error : in 'scan_label' (genealogy.h line 624): invalid Newick label: expected one of 'b','g','m', or 'o', got '_'.
 > 
 > try(
 +   r"{(o_9_1:1.000000)z_1_3:1.000000)m_0_0:0.000000;}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 620): invalid Newick label: expected one of 'b','g','m', or 'o', got 'z'.
+Error : in 'scan_label' (genealogy.h line 624): invalid Newick label: expected one of 'b','g','m', or 'o', got 'z'.
 > 
 > r"{((,o_9_1:1.000000,)g_1_3:1.000000)m_0_0:0.000000;}" |>
 +    parse_newick() |>
@@ -212,7 +212,7 @@ Error : in 'parse': invalid Newick string: unbalanced parentheses.
 +   r"{((((o_9_1:1.000000)g_1_3:1.000000)m_0_0:0.000000)g_3:22)()))));}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 620): invalid Newick label: expected one of 'b','g','m', or 'o', got '('.
+Error : in 'scan_label' (genealogy.h line 624): invalid Newick label: expected one of 'b','g','m', or 'o', got '('.
 > 
 > try(
 +   r"{(o_9_1:1.000000,b_1_3)m_0_0:0.000000;}" |>
@@ -258,7 +258,7 @@ $newick
 +   r"{(o_a9_1:1.000000,b_0_0:3)m_0_0:0.000000,b_____2:17;}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 637): invalid Newick format: deme should be indicated with an integer.
+Error : in 'scan_label' (genealogy.h line 641): invalid Newick format: deme should be indicated with an integer.
 > 
 > try(
 +   r"{(o_9_1:1.000000,b_0_0:3)m_0_0(:0.500000,b_____2:17;}" |>
@@ -270,7 +270,7 @@ Error : in 'scan_label': invalid Newick format: missing or invalid branch length
 +   r"{(o_9_1:1.000000,b_0_0:3)m_)0_0:0.500000,b_____2:17;}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 627): invalid Newick format.
+Error : in 'scan_label' (genealogy.h line 631): invalid Newick format.
 > 
 > try(
 +   r"{(o_9_1:1.000000,b_0_0:3)m_0_0:(0.000000,b_____2:17;}" |>
@@ -321,7 +321,7 @@ $newick
 +   r"{(o_9_1:1.000000,b_0_0:3))m_0_0:0.000000,b_____2:17;}" |>
 +      parse_newick()
 + )
-Error : in 'scan_label' (genealogy.h line 620): invalid Newick label: expected one of 'b','g','m', or 'o', got ')'.
+Error : in 'scan_label' (genealogy.h line 624): invalid Newick label: expected one of 'b','g','m', or 'o', got ')'.
 > 
 > try(
 +   r"{(o_9_1:1.000000,b_0_0:3)m_0_0:0.000000,b_____2:17))));}" |>
@@ -359,7 +359,7 @@ Error : in 'scan_label': invalid Newick format: missing or invalid branch length
 +      parse_newick() |>
 +      newick(prune=FALSE,obscure=FALSE)
 + )
-Error : in 'scan_label' (genealogy.h line 606): invalid Newick format: empty label.
+Error : in 'scan_label' (genealogy.h line 610): invalid Newick format: empty label.
 > 
 > try(
 +   r"{(o_9_1:  1.000000,b_0:3,)m_0_0:-0.00001,b_2_45:17}" |>