Merge pull request #182 from mrc-ide/release-0.1.12

Release 0.1.12

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: individual
 Title: Framework for Specifying and Simulating Individual Based Models
-Version: 0.1.11
+Version: 0.1.12
 Authors@R: c(
   person(
     given = "Giovanni",
@@ -27,6 +27,7 @@ Authors@R: c(
     given = "Paul",
     family = "Liétar",
     role = c('aut'),
+    comment = c(ORCID = "0009-0000-3813-6227"),
     email = '[email protected]'
   ),
   person(
@@ -64,7 +65,7 @@ Suggests:
     testthat (>= 2.1.0),
     xml2,
     bench
-RoxygenNote: 7.2.1.9000
+RoxygenNote: 7.2.3
 VignetteBuilder: knitr
 LinkingTo:
   Rcpp,

NEWS.md

@@ -1,3 +1,7 @@
+# individual 0.1.12
+
+  * Simulation state can be saved and restored, allowing the simulation to be resumed.
+
 # individual 0.1.11
 
   * Optimised rendering memory usage and speed

R/RcppExports.R

@@ -173,8 +173,16 @@ create_targeted_event <- function(size) {
     .Call(`_individual_create_targeted_event`, size)
 }
 
-event_tick <- function(event) {
-    invisible(.Call(`_individual_event_tick`, event))
+event_base_tick <- function(event) {
+    invisible(.Call(`_individual_event_base_tick`, event))
+}
+
+event_base_get_timestep <- function(event) {
+    .Call(`_individual_event_base_get_timestep`, event)
+}
+
+event_base_should_trigger <- function(event) {
+    .Call(`_individual_event_base_should_trigger`, event)
 }
 
 event_schedule <- function(event, delays) {
@@ -185,6 +193,14 @@ event_clear_schedule <- function(event) {
     invisible(.Call(`_individual_event_clear_schedule`, event))
 }
 
+event_checkpoint <- function(event) {
+    .Call(`_individual_event_checkpoint`, event)
+}
+
+event_restore <- function(event, time, schedule) {
+    invisible(.Call(`_individual_event_restore`, event, time, schedule))
+}
+
 targeted_event_clear_schedule_vector <- function(event, target) {
     invisible(.Call(`_individual_targeted_event_clear_schedule_vector`, event, target))
 }
@@ -229,14 +245,6 @@ targeted_event_schedule_multi_delay_vector <- function(event, target, delay) {
     invisible(.Call(`_individual_targeted_event_schedule_multi_delay_vector`, event, target, delay))
 }
 
-event_get_timestep <- function(event) {
-    .Call(`_individual_event_get_timestep`, event)
-}
-
-event_should_trigger <- function(event) {
-    .Call(`_individual_event_should_trigger`, event)
-}
-
 targeted_event_get_target <- function(event) {
     .Call(`_individual_targeted_event_get_target`, event)
 }
@@ -245,6 +253,14 @@ targeted_event_resize <- function(event) {
     invisible(.Call(`_individual_targeted_event_resize`, event))
 }
 
+targeted_event_checkpoint <- function(event) {
+    .Call(`_individual_targeted_event_checkpoint`, event)
+}
+
+targeted_event_restore <- function(event, time, state) {
+    invisible(.Call(`_individual_targeted_event_restore`, event, time, state))
+}
+
 process_listener <- function(event, listener) {
     invisible(.Call(`_individual_process_listener`, event, listener))
 }

R/categorical_variable.R

@@ -38,7 +38,7 @@ CategoricalVariable <- R6Class(
     #' @description return a character vector of possible values.
     #' Note that the order of the returned vector may not be the same order
     #' that was given when the variable was intitialized, due to the underlying
-    #' unordered storage type. 
+    #' unordered storage type.
     get_categories = function() {
       categorical_variable_get_categories(self$.variable)
     },
@@ -94,6 +94,23 @@ CategoricalVariable <- R6Class(
     size = function() variable_get_size(self$.variable),
 
     .update = function() variable_update(self$.variable),
-    .resize = function() variable_resize(self$.variable)
+    .resize = function() variable_resize(self$.variable),
+
+    .checkpoint = function() {
+      categories <- self$get_categories()
+      values <- lapply(categories, function(c) self$get_index_of(c)$to_vector())
+      names(values) <- categories
+      values
+    },
+
+    .restore = function(values) {
+      stopifnot(names(values) == self$get_categories())
+      stopifnot(sum(sapply(values, length)) == categorical_variable_get_size(self$.variable))
+
+      for (c in names(values)) {
+        self$queue_update(c, values[[c]])
+      }
+      self$.update()
+    }
   )
 )

R/double_variable.R

@@ -145,6 +145,13 @@ DoubleVariable <- R6Class(
     size = function() variable_get_size(self$.variable),
 
     .update = function() variable_update(self$.variable),
-    .resize = function() variable_resize(self$.variable)
+    .resize = function() variable_resize(self$.variable),
+
+    .checkpoint = function() self$get_values(),
+    .restore = function(values) {
+      stopifnot(length(values) == variable_get_size(self$.variable))
+      self$queue_update(values)
+      self$.update()
+    }
   )
 )

R/event.R

@@ -1,61 +1,78 @@
-#' @title Event Class
-#' @description Describes a general event in the simulation.
+#' @title EventBase Class
+#' @description Common functionality shared between simple and targeted events.
 #' @importFrom R6 R6Class
-#' @export
-Event <- R6Class(
-  'Event',
+EventBase <- R6Class(
+  'EventBase',
   public = list(
-
     .event = NULL,
     .listeners = list(),
 
-    #' @description Initialise an Event.
-    initialize = function() {
-      self$.event <- create_event()
-    },
-
     #' @description Add an event listener.
     #' @param listener the function to be executed on the event, which takes a single
     #' argument giving the time step when this event is triggered.
     add_listener = function(listener) {
       self$.listeners <- c(self$.listeners, listener)
     },
 
-    #' @description Schedule this event to occur in the future.
-    #' @param delay the number of time steps to wait before triggering the event,
-    #' can be a scalar or a vector of values for events that should be triggered
-    #' multiple times.
-    schedule = function(delay) event_schedule(self$.event, delay),
+    .timestep = function() event_base_get_timestep(self$.event),
 
-    #' @description Stop a future event from triggering.
-    clear_schedule = function() event_clear_schedule(self$.event),
-
-    .tick = function() event_tick(self$.event),
+    .tick = function() event_base_tick(self$.event),
 
     .process = function() {
       for (listener in self$.listeners) {
-        if (event_should_trigger(self$.event)) {
+        if (event_base_should_trigger(self$.event)) {
           if (inherits(listener, "externalptr")) {
             self$.process_listener_cpp(listener)
           } else {
             self$.process_listener(listener)
           }
         }
       }
+    }
+  )
+)
+
+#' @title Event Class
+#' @description Describes a general event in the simulation.
+#' @importFrom R6 R6Class
+#' @export
+Event <- R6Class(
+  'Event',
+  inherit=EventBase,
+  public = list(
+    #' @description Initialise an Event.
+    initialize = function() {
+      self$.event <- create_event()
     },
 
+    #' @description Schedule this event to occur in the future.
+    #' @param delay the number of time steps to wait before triggering the event,
+    #' can be a scalar or a vector of values for events that should be triggered
+    #' multiple times.
+    schedule = function(delay) event_schedule(self$.event, delay),
+
+    #' @description Stop a future event from triggering.
+    clear_schedule = function() event_clear_schedule(self$.event),
+
     .process_listener = function(listener) {
-      listener(event_get_timestep(self$.event))
+      listener(self$.timestep())
     },
 
-    .process_listener_cpp = function(listener){
+    .process_listener_cpp = function(listener) {
       process_listener(
         event = self$.event,
         listener = listener
       )
     },
 
     # NOTE: intentionally empty
-    .resize = function() {}
+    .resize = function() {},
+
+    .checkpoint = function() {
+      event_checkpoint(self$.event)
+    },
+    .restore = function(time, schedule) {
+      event_restore(self$.event, time, schedule)
+    }
   )
 )

R/integer_variable.R

@@ -188,6 +188,13 @@ IntegerVariable <- R6Class(
     size = function() variable_get_size(self$.variable),
 
     .update = function() variable_update(self$.variable),
-    .resize = function() variable_resize(self$.variable)
+    .resize = function() variable_resize(self$.variable),
+
+    .checkpoint = function() self$get_values(),
+    .restore = function(values) {
+      stopifnot(length(values) == variable_get_size(self$.variable))
+      self$queue_update(values)
+      self$.update()
+    }
   )
 )

R/ragged_double.R

@@ -150,6 +150,13 @@ RaggedDouble <- R6Class(
     size = function() variable_get_size(self$.variable),
 
     .update = function() variable_update(self$.variable),
-    .resize = function() variable_resize(self$.variable)
+    .resize = function() variable_resize(self$.variable),
+
+    .checkpoint = function() self$get_values(),
+    .restore = function(values) {
+      stopifnot(length(values) == variable_get_size(self$.variable))
+      self$queue_update(values)
+      self$.update()
+    }
   )
 )

R/ragged_integer.R

@@ -150,6 +150,13 @@ RaggedInteger <- R6Class(
     size = function() variable_get_size(self$.variable),
 
     .update = function() variable_update(self$.variable),
-    .resize = function() variable_resize(self$.variable)
+    .resize = function() variable_resize(self$.variable),
+
+    .checkpoint = function() self$get_values(),
+    .restore = function(values) {
+      stopifnot(length(values) == variable_get_size(self$.variable))
+      self$queue_update(values)
+      self$.update()
+    }
   )
 )

R/simulation.R

@@ -4,7 +4,9 @@
 #' @param variables a list of Variables
 #' @param events a list of Events
 #' @param processes a list of processes to execute on each timestep
-#' @param timesteps the number of timesteps to simulate
+#' @param timesteps the end timestep of the simulation. If `state` is not NULL, timesteps must be greater than `state$timestep`
+#' @param state a checkpoint from which to resume the simulation
+#' @param restore_random_state if TRUE, restore R's global random number generator's state from the checkpoint.
 #' @examples
 #' population <- 4
 #' timesteps <- 5
@@ -35,12 +37,23 @@ simulation_loop <- function(
   variables = list(),
   events = list(),
   processes = list(),
-  timesteps
+  timesteps,
+  state = NULL,
+  restore_random_state = FALSE
   ) {
   if (timesteps <= 0) {
     stop('End timestep must be > 0')
   }
-  for (t in seq_len(timesteps)) {
+
+  start <- 1
+  if (!is.null(state)) {
+    start <- restore_state(state, variables, events, restore_random_state)
+    if (start > timesteps) {
+      stop("Restored state is already longer than timesteps")
+    }
+  }
+
+  for (t in seq(start, timesteps)) {
     for (process in processes) {
       execute_any_process(process, t)
     }
@@ -60,6 +73,57 @@ simulation_loop <- function(
       event$.tick()
     }
   }
+
+  invisible(checkpoint_state(timesteps, variables, events))
+}
+
+#' @title Save the simulation state
+#' @description Save the simulation state in an R object, allowing it to be
+#' resumed later using \code{\link[individual]{restore_state}}.
+#' @param timesteps <- the number of time steps that have already been simulated
+#' @param variables the list of Variables
+#' @param events the list of Events
+checkpoint_state <- function(timesteps, variables, events) {
+  random_state <- .GlobalEnv$.Random.seed
+  list(
+    variables=lapply(variables, function(v) v$.checkpoint()),
+    events=lapply(events, function(e) e$.checkpoint()),
+    timesteps=timesteps,
+    random_state=random_state
+  )
+}
+
+#' @title Restore the simulation state
+#' @description Restore the simulation state from a previous checkpoint.
+#' The state of passed events and variables is overwritten to match the state they
+#' had when the simulation was checkpointed. Returns the time step at which the
+#' simulation should resume.
+#' @param state the simulation state to restore, as returned by \code{\link[individual]{restore_state}}.
+#' @param variables the list of Variables
+#' @param events the list of Events
+#' @param restore_random_state if TRUE, restore R's global random number generator's state from the checkpoint.
+restore_state <- function(state, variables, events, restore_random_state) {
+  timesteps <- state$timesteps + 1
+
+  if (length(variables) != length(state$variables)) {
+    stop("Checkpoint's variables do not match simulation's")
+  }
+  for (i in seq_along(variables)) {
+    variables[[i]]$.restore(state$variables[[i]])
+  }
+
+  if (length(events) != length(state$events)) {
+    stop("Checkpoint's events do not match simulation's")
+  }
+  for (i in seq_along(events)) {
+    events[[i]]$.restore(timesteps, state$events[[i]])
+  }
+
+  if (restore_random_state) {
+    .GlobalEnv$.Random.seed <- state$random_state
+  }
+
+  timesteps
 }
 
 #' @title Execute a C++ or R process in the simulation