trap non-finite event-rate error

Author: kingaa

DESCRIPTION

@@ -1,14 +1,14 @@
 Package: phylopomp
 Type: Package
 Title: Phylodynamic Inference for POMP Models
-Version: 0.14.9.4
-Date: 2025-09-21
+Version: 0.14.10.0
+Date: 2025-12-16
 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"))
 	  )
 Description: Tools for phylodynamic analysis.
-biocViews:
 Depends: R(>= 4.5)
+Remotes: bioc::ggtree
 Imports: grid, ggtree(>= 3.16), yaml, ggplot2, scales, ape, dplyr, tibble, foreach, tidyr, pomp(>= 6.1), cowplot
 LinkingTo: pomp
 Suggests: tidyverse, doFuture, broom

src/popul_proc.h

@@ -138,16 +138,26 @@ class popul_proc_t  {
   void update_clocks (void) {
     double rate[nevent];
     double total_rate = event_rates(rate,nevent);
-    if (total_rate > 0) {
-      next = current+rexp(1/total_rate);
+    if (R_FINITE(total_rate)) {
+      if (total_rate > 0) {
+        next = current+rexp(1/total_rate);
+      } else {
+        next = R_PosInf;
+      }
     } else {
-      next = R_PosInf;
+      for (event = 0; event < nevent; event++) {
+        if (!R_FINITE(rate[event]))
+          Rprintf("in '%s' (%s line %d): invalid event rate[%zd]=%lg\n",
+                  __func__,__FILE__,__LINE__,event,rate[event]);
+      }
+      err("in '%s' (%s line %d): invalid total event rate=%lg",
+          __func__,__FILE__,__LINE__,total_rate);
     }
     double u = runif(0,total_rate);
     event = 0;
     while (u > rate[event] && event < nevent) {
       if (rate[event] < 0)
-        err("in '%s' (%s line %d): invalid negative rate[%zd]=%lg", // #nocov
+        err("in '%s' (%s line %d): invalid rate[%zd]=%lg", // #nocov
             __func__,__FILE__,__LINE__,event,rate[event]); // #nocov
       u -= rate[event];
       event++;

tests/twospecies2.R

@@ -0,0 +1,37 @@
+options(tidyverse.quiet=TRUE,digits=3)
+library(tidyverse)
+library(phylopomp)
+
+try(
+  simulate(
+    "TwoSpecies",
+    t0=0,time=3,
+    iota1=0,iota2=0,
+    Beta11=10,Beta12=0,
+    Beta21=0,Beta22=0,
+    gamma1=1,gamma2=1,
+    c1=1,c2=1,
+    psi1=1,psi2=1,
+    S1_0=1000,I1_0=2,R1_0=0,
+    S2_0=0,I2_0=0,R2_0=0,
+    b1=0,b2=0,
+    d1=0,d2=0,
+    omega1=0.5,omega2=0.5
+  )
+)
+
+simulate(
+  "TwoSpecies",
+  t0=0,time=3,
+  iota1=0,iota2=0,
+  Beta11=10,Beta12=0,
+  Beta21=0,Beta22=0,
+  gamma1=1,gamma2=1,
+  c1=1,c2=1,
+  psi1=1,psi2=1,
+  S1_0=1000,I1_0=2,R1_0=0,
+  S2_0=1,I2_0=0,R2_0=0,
+  b1=0,b2=0,
+  d1=0,d2=0,
+  omega1=0.5,omega2=0.5
+)

tests/twospecies2.Rout.save

@@ -0,0 +1,68 @@
+
+R version 4.5.2 (2025-10-31) -- "[Not] Part in a Rumble"
+Copyright (C) 2025 The R Foundation for Statistical Computing
+Platform: x86_64-pc-linux-gnu
+
+R is free software and comes with ABSOLUTELY NO WARRANTY.
+You are welcome to redistribute it under certain conditions.
+Type 'license()' or 'licence()' for distribution details.
+
+  Natural language support but running in an English locale
+
+R is a collaborative project with many contributors.
+Type 'contributors()' for more information and
+'citation()' on how to cite R or R packages in publications.
+
+Type 'demo()' for some demos, 'help()' for on-line help, or
+'help.start()' for an HTML browser interface to help.
+Type 'q()' to quit R.
+
+> options(tidyverse.quiet=TRUE,digits=3)
+> library(tidyverse)
+> library(phylopomp)
+
+Attaching package: 'phylopomp'
+
+The following object is masked from 'package:stats':
+
+    simulate
+
+> 
+> try(
++   simulate(
++     "TwoSpecies",
++     t0=0,time=3,
++     iota1=0,iota2=0,
++     Beta11=10,Beta12=0,
++     Beta21=0,Beta22=0,
++     gamma1=1,gamma2=1,
++     c1=1,c2=1,
++     psi1=1,psi2=1,
++     S1_0=1000,I1_0=2,R1_0=0,
++     S2_0=0,I2_0=0,R2_0=0,
++     b1=0,b2=0,
++     d1=0,d2=0,
++     omega1=0.5,omega2=0.5
++   )
++ )
+in 'update_clocks' (popul_proc.h line 151): invalid event rate[1]=-nan
+in 'update_clocks' (popul_proc.h line 151): invalid event rate[2]=-nan
+Error : in 'update_clocks' (popul_proc.h line 154): invalid total event rate=-nan
+> 
+> simulate(
++   "TwoSpecies",
++   t0=0,time=3,
++   iota1=0,iota2=0,
++   Beta11=10,Beta12=0,
++   Beta21=0,Beta22=0,
++   gamma1=1,gamma2=1,
++   c1=1,c2=1,
++   psi1=1,psi2=1,
++   S1_0=1000,I1_0=2,R1_0=0,
++   S2_0=1,I2_0=0,R2_0=0,
++   b1=0,b2=0,
++   d1=0,d2=0,
++   omega1=0.5,omega2=0.5
++ )
+<gpsim for TwoSpecies model>
+>