diff --git a/R/RcppExports.R b/R/RcppExports.R index f5d3db8..18b51bd 100644 --- a/R/RcppExports.R +++ b/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)) } diff --git a/R/event.R b/R/event.R index 9667691..c7d3443 100644 --- a/R/event.R +++ b/R/event.R @@ -1,19 +1,12 @@ -#' @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. @@ -21,20 +14,13 @@ Event <- R6Class( 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 { @@ -42,13 +28,37 @@ Event <- R6Class( } } } + } + ) +) + +#' @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 @@ -56,6 +66,13 @@ Event <- R6Class( }, # NOTE: intentionally empty - .resize = function() {} + .resize = function() {}, + + .checkpoint = function() { + event_checkpoint(self$.event) + }, + .restore = function(time, schedule) { + event_restore(self$.event, time, schedule) + } ) ) diff --git a/R/simulation.R b/R/simulation.R index f2bf320..0229c00 100644 --- a/R/simulation.R +++ b/R/simulation.R @@ -85,6 +85,7 @@ 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 ) @@ -99,19 +100,25 @@ checkpoint_state <- function(timesteps, variables, events) { #' @param variables the list of Variables #' @param events the list of Events restore_state <- function(state, variables, events) { + 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) > 0) { - stop("Events cannot be restored yet") + + 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]]) } .GlobalEnv$.Random.seed <- state$random_state - state$timesteps + 1 + timesteps } #' @title Execute a C++ or R process in the simulation diff --git a/R/targeted_event.R b/R/targeted_event.R index 446557d..e265b27 100644 --- a/R/targeted_event.R +++ b/R/targeted_event.R @@ -5,7 +5,7 @@ #' @export TargetedEvent <- R6Class( 'TargetedEvent', - inherit = Event, + inherit = EventBase, public = list( #' @description Initialise a TargetedEvent. @@ -103,7 +103,7 @@ TargetedEvent <- R6Class( .process_listener = function(listener) { listener( - event_get_timestep(self$.event), + self$.timestep(), Bitset$new(from=targeted_event_get_target(self$.event)) ) }, @@ -116,6 +116,13 @@ TargetedEvent <- R6Class( ) }, - .resize = function() targeted_event_resize(self$.event) + .resize = function() targeted_event_resize(self$.event), + + .checkpoint = function() { + targeted_event_checkpoint(self$.event) + }, + .restore = function(time, schedule) { + targeted_event_restore(self$.event, time, schedule) + } ) ) diff --git a/inst/include/Event.h b/inst/include/Event.h index 2305a24..4f08aa2 100644 --- a/inst/include/Event.h +++ b/inst/include/Event.h @@ -12,6 +12,7 @@ #include #include #include +#include #include #include #include @@ -40,6 +41,7 @@ inline std::vector round_delay(const std::vector& delay) { //' @title abstract base class for events class EventBase { +protected: size_t t = 1; public: virtual void tick(); @@ -78,7 +80,9 @@ class Event : public EventBase { virtual void schedule(std::vector delays); virtual void clear_schedule(); - + + virtual std::vector checkpoint(); + virtual void restore(size_t time, std::vector schedule); }; //' @title process an event by calling a listener @@ -109,6 +113,17 @@ inline void Event::clear_schedule() { simple_schedule.clear(); } +//' @title save this event's state +inline std::vector Event::checkpoint() { + return {simple_schedule.begin(), simple_schedule.end()}; +} + +//' @title restore this event's state from a previous checkpoint +inline void Event::restore(size_t time, std::vector schedule) { + t = time; + simple_schedule.clear(); + simple_schedule.insert(schedule.begin(), schedule.end()); +} //' @title a targeted event in the simulation //' @description This class provides functionality for targeted events which are @@ -155,6 +170,8 @@ class TargetedEvent : public EventBase { virtual void clear_schedule(const individual_index_t&); virtual individual_index_t get_scheduled() const; + virtual std::vector> checkpoint() const; + virtual void restore(size_t time, std::vector> schedule); }; inline TargetedEvent::TargetedEvent(size_t size) @@ -364,4 +381,20 @@ inline void TargetedEvent::resize() { } } +//' @title save this event's state +inline std::vector> +TargetedEvent::checkpoint() const { + return {targeted_schedule.begin(), targeted_schedule.end()}; +} + +//' @title restore this event's state from a previous checkpoint +inline void TargetedEvent::restore( + size_t time, + std::vector> schedule +) { + t = time; + targeted_schedule.clear(); + targeted_schedule.insert(schedule.begin(), schedule.end()); +} + #endif /* INST_INCLUDE_EVENT_H_ */ diff --git a/man/Event.Rd b/man/Event.Rd index 6577508..6d23ec4 100644 --- a/man/Event.Rd +++ b/man/Event.Rd @@ -6,21 +6,33 @@ \description{ Describes a general event in the simulation. } +\section{Super class}{ +\code{individual::EventBase} -> \code{Event} +} \section{Methods}{ \subsection{Public methods}{ \itemize{ \item \href{#method-Event-new}{\code{Event$new()}} -\item \href{#method-Event-add_listener}{\code{Event$add_listener()}} \item \href{#method-Event-schedule}{\code{Event$schedule()}} \item \href{#method-Event-clear_schedule}{\code{Event$clear_schedule()}} -\item \href{#method-Event-.tick}{\code{Event$.tick()}} -\item \href{#method-Event-.process}{\code{Event$.process()}} \item \href{#method-Event-.process_listener}{\code{Event$.process_listener()}} \item \href{#method-Event-.process_listener_cpp}{\code{Event$.process_listener_cpp()}} \item \href{#method-Event-.resize}{\code{Event$.resize()}} +\item \href{#method-Event-.checkpoint}{\code{Event$.checkpoint()}} +\item \href{#method-Event-.restore}{\code{Event$.restore()}} \item \href{#method-Event-clone}{\code{Event$clone()}} } } +\if{html}{\out{ +
Inherited methods + +
+}} \if{html}{\out{
}} \if{html}{\out{}} \if{latex}{\out{\hypertarget{method-Event-new}{}}} @@ -30,24 +42,6 @@ Initialise an Event. \if{html}{\out{
}}\preformatted{Event$new()}\if{html}{\out{
}} } -} -\if{html}{\out{
}} -\if{html}{\out{}} -\if{latex}{\out{\hypertarget{method-Event-add_listener}{}}} -\subsection{Method \code{add_listener()}}{ -Add an event listener. -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{Event$add_listener(listener)}\if{html}{\out{
}} -} - -\subsection{Arguments}{ -\if{html}{\out{
}} -\describe{ -\item{\code{listener}}{the function to be executed on the event, which takes a single -argument giving the time step when this event is triggered.} -} -\if{html}{\out{
}} -} } \if{html}{\out{
}} \if{html}{\out{}} @@ -77,24 +71,6 @@ Stop a future event from triggering. \if{html}{\out{
}}\preformatted{Event$clear_schedule()}\if{html}{\out{
}} } -} -\if{html}{\out{
}} -\if{html}{\out{}} -\if{latex}{\out{\hypertarget{method-Event-.tick}{}}} -\subsection{Method \code{.tick()}}{ -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{Event$.tick()}\if{html}{\out{
}} -} - -} -\if{html}{\out{
}} -\if{html}{\out{}} -\if{latex}{\out{\hypertarget{method-Event-.process}{}}} -\subsection{Method \code{.process()}}{ -\subsection{Usage}{ -\if{html}{\out{
}}\preformatted{Event$.process()}\if{html}{\out{
}} -} - } \if{html}{\out{
}} \if{html}{\out{}} @@ -122,6 +98,24 @@ Stop a future event from triggering. \if{html}{\out{
}}\preformatted{Event$.resize()}\if{html}{\out{
}} } +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-Event-.checkpoint}{}}} +\subsection{Method \code{.checkpoint()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{Event$.checkpoint()}\if{html}{\out{
}} +} + +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-Event-.restore}{}}} +\subsection{Method \code{.restore()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{Event$.restore(time, schedule)}\if{html}{\out{
}} +} + } \if{html}{\out{
}} \if{html}{\out{}} diff --git a/man/EventBase.Rd b/man/EventBase.Rd new file mode 100644 index 0000000..267c3d2 --- /dev/null +++ b/man/EventBase.Rd @@ -0,0 +1,81 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/event.R +\name{EventBase} +\alias{EventBase} +\title{EventBase Class} +\description{ +Common functionality shared between simple and targeted events. +} +\section{Methods}{ +\subsection{Public methods}{ +\itemize{ +\item \href{#method-EventBase-add_listener}{\code{EventBase$add_listener()}} +\item \href{#method-EventBase-.timestep}{\code{EventBase$.timestep()}} +\item \href{#method-EventBase-.tick}{\code{EventBase$.tick()}} +\item \href{#method-EventBase-.process}{\code{EventBase$.process()}} +\item \href{#method-EventBase-clone}{\code{EventBase$clone()}} +} +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-EventBase-add_listener}{}}} +\subsection{Method \code{add_listener()}}{ +Add an event listener. +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{EventBase$add_listener(listener)}\if{html}{\out{
}} +} + +\subsection{Arguments}{ +\if{html}{\out{
}} +\describe{ +\item{\code{listener}}{the function to be executed on the event, which takes a single +argument giving the time step when this event is triggered.} +} +\if{html}{\out{
}} +} +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-EventBase-.timestep}{}}} +\subsection{Method \code{.timestep()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{EventBase$.timestep()}\if{html}{\out{
}} +} + +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-EventBase-.tick}{}}} +\subsection{Method \code{.tick()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{EventBase$.tick()}\if{html}{\out{
}} +} + +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-EventBase-.process}{}}} +\subsection{Method \code{.process()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{EventBase$.process()}\if{html}{\out{
}} +} + +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-EventBase-clone}{}}} +\subsection{Method \code{clone()}}{ +The objects of this class are cloneable with this method. +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{EventBase$clone(deep = FALSE)}\if{html}{\out{
}} +} + +\subsection{Arguments}{ +\if{html}{\out{
}} +\describe{ +\item{\code{deep}}{Whether to make a deep clone.} +} +\if{html}{\out{
}} +} +} +} diff --git a/man/TargetedEvent.Rd b/man/TargetedEvent.Rd index bb65d91..b691293 100644 --- a/man/TargetedEvent.Rd +++ b/man/TargetedEvent.Rd @@ -8,7 +8,7 @@ Describes a targeted event in the simulation. This is useful for events which are triggered for a sub-population. } \section{Super class}{ -\code{\link[individual:Event]{individual::Event}} -> \code{TargetedEvent} +\code{individual::EventBase} -> \code{TargetedEvent} } \section{Methods}{ \subsection{Public methods}{ @@ -23,15 +23,18 @@ This is useful for events which are triggered for a sub-population. \item \href{#method-TargetedEvent-.process_listener}{\code{TargetedEvent$.process_listener()}} \item \href{#method-TargetedEvent-.process_listener_cpp}{\code{TargetedEvent$.process_listener_cpp()}} \item \href{#method-TargetedEvent-.resize}{\code{TargetedEvent$.resize()}} +\item \href{#method-TargetedEvent-.checkpoint}{\code{TargetedEvent$.checkpoint()}} +\item \href{#method-TargetedEvent-.restore}{\code{TargetedEvent$.restore()}} \item \href{#method-TargetedEvent-clone}{\code{TargetedEvent$clone()}} } } \if{html}{\out{
Inherited methods
}} @@ -180,6 +183,24 @@ Shrink the TargetedEvent. \if{html}{\out{
}}\preformatted{TargetedEvent$.resize()}\if{html}{\out{
}} } +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-TargetedEvent-.checkpoint}{}}} +\subsection{Method \code{.checkpoint()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{TargetedEvent$.checkpoint()}\if{html}{\out{
}} +} + +} +\if{html}{\out{
}} +\if{html}{\out{}} +\if{latex}{\out{\hypertarget{method-TargetedEvent-.restore}{}}} +\subsection{Method \code{.restore()}}{ +\subsection{Usage}{ +\if{html}{\out{
}}\preformatted{TargetedEvent$.restore(time, schedule)}\if{html}{\out{
}} +} + } \if{html}{\out{
}} \if{html}{\out{}} diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index 6c40980..7cb6374 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -512,7 +512,7 @@ BEGIN_RCPP END_RCPP } // create_event -Rcpp::XPtr create_event(); +Rcpp::XPtr create_event(); RcppExport SEXP _individual_create_event() { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; @@ -522,7 +522,7 @@ BEGIN_RCPP END_RCPP } // create_targeted_event -Rcpp::XPtr create_targeted_event(size_t size); +Rcpp::XPtr create_targeted_event(size_t size); RcppExport SEXP _individual_create_targeted_event(SEXP sizeSEXP) { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; @@ -532,16 +532,38 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } -// event_tick -void event_tick(const Rcpp::XPtr event); -RcppExport SEXP _individual_event_tick(SEXP eventSEXP) { +// event_base_tick +void event_base_tick(const Rcpp::XPtr event); +RcppExport SEXP _individual_event_base_tick(SEXP eventSEXP) { BEGIN_RCPP Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); - event_tick(event); + event_base_tick(event); return R_NilValue; END_RCPP } +// event_base_get_timestep +size_t event_base_get_timestep(const Rcpp::XPtr event); +RcppExport SEXP _individual_event_base_get_timestep(SEXP eventSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + rcpp_result_gen = Rcpp::wrap(event_base_get_timestep(event)); + return rcpp_result_gen; +END_RCPP +} +// event_base_should_trigger +bool event_base_should_trigger(const Rcpp::XPtr event); +RcppExport SEXP _individual_event_base_should_trigger(SEXP eventSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + rcpp_result_gen = Rcpp::wrap(event_base_should_trigger(event)); + return rcpp_result_gen; +END_RCPP +} // event_schedule void event_schedule(const Rcpp::XPtr event, std::vector delays); RcppExport SEXP _individual_event_schedule(SEXP eventSEXP, SEXP delaysSEXP) { @@ -563,6 +585,29 @@ BEGIN_RCPP return R_NilValue; END_RCPP } +// event_checkpoint +std::vector event_checkpoint(const Rcpp::XPtr event); +RcppExport SEXP _individual_event_checkpoint(SEXP eventSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + rcpp_result_gen = Rcpp::wrap(event_checkpoint(event)); + return rcpp_result_gen; +END_RCPP +} +// event_restore +void event_restore(const Rcpp::XPtr event, size_t time, std::vector schedule); +RcppExport SEXP _individual_event_restore(SEXP eventSEXP, SEXP timeSEXP, SEXP scheduleSEXP) { +BEGIN_RCPP + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + Rcpp::traits::input_parameter< size_t >::type time(timeSEXP); + Rcpp::traits::input_parameter< std::vector >::type schedule(scheduleSEXP); + event_restore(event, time, schedule); + return R_NilValue; +END_RCPP +} // targeted_event_clear_schedule_vector void targeted_event_clear_schedule_vector(const Rcpp::XPtr event, std::vector target); RcppExport SEXP _individual_targeted_event_clear_schedule_vector(SEXP eventSEXP, SEXP targetSEXP) { @@ -688,46 +733,47 @@ BEGIN_RCPP return R_NilValue; END_RCPP } -// event_get_timestep -size_t event_get_timestep(const Rcpp::XPtr event); -RcppExport SEXP _individual_event_get_timestep(SEXP eventSEXP) { +// targeted_event_get_target +Rcpp::XPtr targeted_event_get_target(const Rcpp::XPtr event); +RcppExport SEXP _individual_targeted_event_get_target(SEXP eventSEXP) { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; - Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); - rcpp_result_gen = Rcpp::wrap(event_get_timestep(event)); + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + rcpp_result_gen = Rcpp::wrap(targeted_event_get_target(event)); return rcpp_result_gen; END_RCPP } -// event_should_trigger -bool event_should_trigger(const Rcpp::XPtr event); -RcppExport SEXP _individual_event_should_trigger(SEXP eventSEXP) { +// targeted_event_resize +void targeted_event_resize(const Rcpp::XPtr event); +RcppExport SEXP _individual_targeted_event_resize(SEXP eventSEXP) { BEGIN_RCPP - Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; - Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); - rcpp_result_gen = Rcpp::wrap(event_should_trigger(event)); - return rcpp_result_gen; + Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); + targeted_event_resize(event); + return R_NilValue; END_RCPP } -// targeted_event_get_target -Rcpp::XPtr targeted_event_get_target(const Rcpp::XPtr event); -RcppExport SEXP _individual_targeted_event_get_target(SEXP eventSEXP) { +// targeted_event_checkpoint +Rcpp::List targeted_event_checkpoint(const Rcpp::XPtr event); +RcppExport SEXP _individual_targeted_event_checkpoint(SEXP eventSEXP) { BEGIN_RCPP Rcpp::RObject rcpp_result_gen; Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); - rcpp_result_gen = Rcpp::wrap(targeted_event_get_target(event)); + rcpp_result_gen = Rcpp::wrap(targeted_event_checkpoint(event)); return rcpp_result_gen; END_RCPP } -// targeted_event_resize -void targeted_event_resize(const Rcpp::XPtr event); -RcppExport SEXP _individual_targeted_event_resize(SEXP eventSEXP) { +// targeted_event_restore +void targeted_event_restore(const Rcpp::XPtr event, size_t time, Rcpp::List state); +RcppExport SEXP _individual_targeted_event_restore(SEXP eventSEXP, SEXP timeSEXP, SEXP stateSEXP) { BEGIN_RCPP Rcpp::RNGScope rcpp_rngScope_gen; Rcpp::traits::input_parameter< const Rcpp::XPtr >::type event(eventSEXP); - targeted_event_resize(event); + Rcpp::traits::input_parameter< size_t >::type time(timeSEXP); + Rcpp::traits::input_parameter< Rcpp::List >::type state(stateSEXP); + targeted_event_restore(event, time, state); return R_NilValue; END_RCPP } @@ -1443,9 +1489,13 @@ static const R_CallMethodDef CallEntries[] = { {"_individual_double_variable_queue_shrink_bitset", (DL_FUNC) &_individual_double_variable_queue_shrink_bitset, 2}, {"_individual_create_event", (DL_FUNC) &_individual_create_event, 0}, {"_individual_create_targeted_event", (DL_FUNC) &_individual_create_targeted_event, 1}, - {"_individual_event_tick", (DL_FUNC) &_individual_event_tick, 1}, + {"_individual_event_base_tick", (DL_FUNC) &_individual_event_base_tick, 1}, + {"_individual_event_base_get_timestep", (DL_FUNC) &_individual_event_base_get_timestep, 1}, + {"_individual_event_base_should_trigger", (DL_FUNC) &_individual_event_base_should_trigger, 1}, {"_individual_event_schedule", (DL_FUNC) &_individual_event_schedule, 2}, {"_individual_event_clear_schedule", (DL_FUNC) &_individual_event_clear_schedule, 1}, + {"_individual_event_checkpoint", (DL_FUNC) &_individual_event_checkpoint, 1}, + {"_individual_event_restore", (DL_FUNC) &_individual_event_restore, 3}, {"_individual_targeted_event_clear_schedule_vector", (DL_FUNC) &_individual_targeted_event_clear_schedule_vector, 2}, {"_individual_targeted_event_clear_schedule", (DL_FUNC) &_individual_targeted_event_clear_schedule, 2}, {"_individual_targeted_event_get_scheduled", (DL_FUNC) &_individual_targeted_event_get_scheduled, 1}, @@ -1457,10 +1507,10 @@ static const R_CallMethodDef CallEntries[] = { {"_individual_targeted_event_schedule_vector", (DL_FUNC) &_individual_targeted_event_schedule_vector, 3}, {"_individual_targeted_event_schedule_multi_delay", (DL_FUNC) &_individual_targeted_event_schedule_multi_delay, 3}, {"_individual_targeted_event_schedule_multi_delay_vector", (DL_FUNC) &_individual_targeted_event_schedule_multi_delay_vector, 3}, - {"_individual_event_get_timestep", (DL_FUNC) &_individual_event_get_timestep, 1}, - {"_individual_event_should_trigger", (DL_FUNC) &_individual_event_should_trigger, 1}, {"_individual_targeted_event_get_target", (DL_FUNC) &_individual_targeted_event_get_target, 1}, {"_individual_targeted_event_resize", (DL_FUNC) &_individual_targeted_event_resize, 1}, + {"_individual_targeted_event_checkpoint", (DL_FUNC) &_individual_targeted_event_checkpoint, 1}, + {"_individual_targeted_event_restore", (DL_FUNC) &_individual_targeted_event_restore, 3}, {"_individual_process_listener", (DL_FUNC) &_individual_process_listener, 2}, {"_individual_process_targeted_listener", (DL_FUNC) &_individual_process_targeted_listener, 3}, {"_individual_create_integer_variable", (DL_FUNC) &_individual_create_integer_variable, 1}, diff --git a/src/event.cpp b/src/event.cpp index 25b9352..522291f 100644 --- a/src/event.cpp +++ b/src/event.cpp @@ -9,20 +9,30 @@ #include "utils.h" //[[Rcpp::export]] -Rcpp::XPtr create_event() { - return Rcpp::XPtr(new Event(), true); +Rcpp::XPtr create_event() { + return Rcpp::XPtr(new Event(), true); } //[[Rcpp::export]] -Rcpp::XPtr create_targeted_event(size_t size) { - return Rcpp::XPtr(new TargetedEvent(size), true); +Rcpp::XPtr create_targeted_event(size_t size) { + return Rcpp::XPtr(new TargetedEvent(size), true); } //[[Rcpp::export]] -void event_tick(const Rcpp::XPtr event) { +void event_base_tick(const Rcpp::XPtr event) { event->tick(); } +//[[Rcpp::export]] +size_t event_base_get_timestep(const Rcpp::XPtr event) { + return event->get_time(); +} + +//[[Rcpp::export]] +bool event_base_should_trigger(const Rcpp::XPtr event) { + return event->should_trigger(); +} + //[[Rcpp::export]] void event_schedule(const Rcpp::XPtr event, std::vector delays) { event->schedule(delays); @@ -33,6 +43,21 @@ void event_clear_schedule(const Rcpp::XPtr event) { event->clear_schedule(); } +//[[Rcpp::export]] +std::vector event_checkpoint(const Rcpp::XPtr event) { + return event->checkpoint(); +} + +//[[Rcpp::export]] +void event_restore(const Rcpp::XPtr event, size_t time, std::vector schedule) { + for (size_t event_ts: schedule) { + if (event_ts < time) { + Rcpp::stop("schedule is in the past"); + } + } + event->restore(time, schedule); +} + //[[Rcpp::export]] void targeted_event_clear_schedule_vector( const Rcpp::XPtr event, @@ -147,16 +172,6 @@ void targeted_event_schedule_multi_delay_vector( event->schedule(target, delay); } -//[[Rcpp::export]] -size_t event_get_timestep(const Rcpp::XPtr event) { - return event->get_time(); -} - -//[[Rcpp::export]] -bool event_should_trigger(const Rcpp::XPtr event) { - return event->should_trigger(); -} - //[[Rcpp::export]] Rcpp::XPtr targeted_event_get_target(const Rcpp::XPtr event) { return Rcpp::XPtr( @@ -170,6 +185,45 @@ void targeted_event_resize(const Rcpp::XPtr event) { event->resize(); } +//[[Rcpp::export]] +Rcpp::List targeted_event_checkpoint(const Rcpp::XPtr event) { + std::vector timesteps; + std::vector> targets; + + const auto schedule = event->checkpoint(); + for (const auto& it: schedule) { + timesteps.push_back(it.first); + targets.push_back(bitset_to_vector_internal(it.second)); + } + + return Rcpp::List::create( + Rcpp::Named("timesteps") = timesteps, + Rcpp::Named("targets") = targets); +} + +//[[Rcpp::export]] +void targeted_event_restore(const Rcpp::XPtr event, size_t time, Rcpp::List state) { + std::vector timesteps = state["timesteps"]; + std::vector> targets = state["targets"];; + + if (timesteps.size() != targets.size()) { + Rcpp::stop("length mismatch between timesteps and targets"); + } + + std::vector> schedule; + for (size_t i = 0; i < timesteps.size(); i++) { + if (timesteps[i] < time) { + Rcpp::stop("schedule is in the past"); + } + decrement(targets[i]); + auto bitmap = individual_index_t(event->size()); + bitmap.insert_safe(targets[i].cbegin(), targets[i].cend()); + schedule.push_back({timesteps[i], bitmap}); + } + + event->restore(time, schedule); +} + // [[Rcpp::export]] void process_listener( const Rcpp::XPtr event, diff --git a/tests/testthat/test-checkpoint.R b/tests/testthat/test-checkpoint.R index 6126311..2755eac 100644 --- a/tests/testthat/test-checkpoint.R +++ b/tests/testthat/test-checkpoint.R @@ -45,6 +45,66 @@ test_that("deterministic state model is resumable", { expect_mapequal(second_phase$data[6:10,], expected[6:10,]) }) +test_that("deterministic model with events is resumable", { + simulation <- function(timesteps, state=NULL) { + population <- 10 + health <- CategoricalVariable$new(c('S', 'I', 'R'), rep('S', population)) + render <- Render$new(timesteps) + + infection <- TargetedEvent$new(population) + recovery <- TargetedEvent$new(population) + + infection$add_listener(function(simulation, target) { + health$queue_update('I', target) + }) + + recovery$add_listener(function(simulation, target) { + health$queue_update('R', target) + }) + + delayed_shift_generator <- function(from, event, delay, rate) { + function(t) { + from_state <- health$get_index_of(from) + # remove the already scheduled individuals + from_state$and(event$get_scheduled()$not(TRUE)) + target <- from_state$to_vector()[seq_len(min(rate,from_state$size()))] + event$schedule(target, delay); + } + } + + processes <- list( + delayed_shift_generator('S', infection, 1, 2), + delayed_shift_generator('I', recovery, 2, 1), + categorical_count_renderer_process(render, health, c('S', 'I', 'R')) + ) + + new_state <- simulation_loop( + variables = list(health), + events = list(infection, recovery), + processes = processes, + timesteps = timesteps, + state = state + ) + list(state=new_state, data=render$to_dataframe()) + } + + first_phase <- simulation(5, state=NULL) + second_phase <- simulation(15, state=first_phase$state) + + expected <- data.frame( + timestep = c(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15), + S_count = c(10, 10, 8, 6, 4, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0), + I_count = c(0, 0, 2, 4, 6, 7, 8, 7, 6, 5, 4, 3, 2, 1, 0), + R_count = c(0, 0, 0, 0, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10) + ) + expected_na <- data.frame(timestep=1:5, S_count=NA_real_, I_count=NA_real_, R_count=NA_real_) + + expect_mapequal(first_phase$data, expected[1:5,]) + expect_mapequal(second_phase$data[1:5,], expected_na) + expect_mapequal(second_phase$data[6:15,], expected[6:15,]) +}) + + test_that("stochastic simulation is repeatable", { simulation <- function(timesteps, state = NULL) { render <- Render$new(timesteps) @@ -66,20 +126,6 @@ test_that("stochastic simulation is repeatable", { expect_equal(first, second) }) -test_that("events are not resumable", { - simulation <- function(timesteps, state = NULL) { - event <- Event$new() - simulation_loop( - events = list(event), - timesteps = timesteps, - state = state - ) - } - - initial <- simulation(5) - expect_error(simulation(10, state=initial), "Events cannot be restored yet") -}) - test_that("cannot add nor remove variables when resuming", { make_variables <- function(count) { lapply(seq_len(count), function(i) DoubleVariable$new(1:10)) diff --git a/tests/testthat/test-events.R b/tests/testthat/test-events.R index b4beb2a..c674ec1 100644 --- a/tests/testthat/test-events.R +++ b/tests/testthat/test-events.R @@ -15,23 +15,6 @@ test_that("first event is triggered at t=1", { event$.tick() }) -test_that("first event is triggered at t=1", { - event <- Event$new() - listener <- mockery::mock() - event$add_listener(listener) - event$schedule(c(0, 1)) - - #time = 1 - event$.process() - mockery::expect_args(listener, 1, t = 1) - event$.tick() - - #time = 2 - event$.process() - mockery::expect_args(listener, 2, t = 2) - event$.tick() -}) - test_that("events can be scheduled for the future", { event <- Event$new() listener <- mockery::mock() @@ -117,3 +100,78 @@ test_that("events can be scheduled and canceled", { mockery::expect_called(listener, 1) mockery::expect_args(listener, 1, t = 2) }) + +test_that("events can be saved and restored", { + listener <- mockery::mock() + + old_event <- Event$new() + old_event$add_listener(listener) + + # This schedules at t=2 and t=4 + old_event$schedule(c(1, 3)) + + #time = 1 + old_event$.process() + mockery::expect_called(listener, 0) + old_event$.tick() + + #time = 2 + old_event$.process() + mockery::expect_called(listener, 1) + mockery::expect_args(listener, 1, t = 2) + old_event$.tick() + + new_event <- Event$new() + new_event$add_listener(listener) + new_event$.restore( + old_event$.timestep(), + old_event$.checkpoint()) + + #time = 3 + new_event$.process() + mockery::expect_called(listener, 1) + new_event$.tick() + + #time = 4 + new_event$.process() + mockery::expect_called(listener, 2) + mockery::expect_args(listener, 2, t = 4) +}) + +test_that("events are cleared when restored", { + listener <- mockery::mock() + + old_event <- Event$new() + old_event$schedule(3) # t=4 + + new_event <- Event$new() + new_event$add_listener(listener) + + # Schedule at t=2. This will be cleared and overridden when restoring, + # replaced by the earlier t=4 schedule. + new_event$schedule(1) + new_event$.restore( + old_event$.timestep(), + old_event$.checkpoint()) + + #time=1 + new_event$.process() + mockery::expect_called(listener, 0) + new_event$.tick() + + #time=2 + new_event$.process() + mockery::expect_called(listener, 0) + new_event$.tick() + + #time=2 + new_event$.process() + mockery::expect_called(listener, 0) + new_event$.tick() + + #time=4 + new_event$.process() + mockery::expect_called(listener, 1) + mockery::expect_args(listener, 1, t = 4) + new_event$.tick() +}) diff --git a/tests/testthat/test-targetedevent.R b/tests/testthat/test-targetedevent.R index 7cdeb44..a85b4e1 100644 --- a/tests/testthat/test-targetedevent.R +++ b/tests/testthat/test-targetedevent.R @@ -526,3 +526,79 @@ test_that("targeted events work for vector delay, bitset target", { expect_error(event$schedule(target = target,delay = delay)) }) + +test_that("targeted events can be saved and restored", { + listener <- mockery::mock() + old_event <- TargetedEvent$new(10) + old_event$add_listener(listener) + old_event$schedule(c(2, 4), 2) + old_event$schedule(c(5, 7), 4) + + #time = 1 + old_event$.process() + mockery::expect_called(listener, 0) + old_event$.tick() + + #time = 2 + old_event$.process() + mockery::expect_called(listener, 0) + old_event$.tick() + + #time = 3 + old_event$.process() + mockery::expect_called(listener, 1) + expect_targeted_listener(listener, 1, t = 3, target = c(2, 4)) + old_event$.tick() + + new_event <- TargetedEvent$new(10) + new_event$add_listener(listener) + new_event$.restore( + old_event$.timestep(), + old_event$.checkpoint()) + + #time = 4 + new_event$.process() + mockery::expect_called(listener, 1) + new_event$.tick() + + #time = 5 + new_event$.process() + mockery::expect_called(listener, 2) + expect_targeted_listener(listener, 2, t = 5, target = c(5, 7)) +}) + +test_that("targeted events are cleared when restored", { + listener <- mockery::mock() + old_event <- TargetedEvent$new(5) + old_event$schedule(c(2, 4), 2) + + new_event <- TargetedEvent$new(10) + new_event$add_listener(listener) + new_event$schedule(c(1, 3), 2) + new_event$schedule(c(5), 3) + + new_event$.restore( + old_event$.timestep(), + old_event$.checkpoint()) + + #time=1 + new_event$.process() + mockery::expect_called(listener, 0) + new_event$.tick() + + #time=2 + new_event$.process() + mockery::expect_called(listener, 0) + new_event$.tick() + + #time=3 + new_event$.process() + mockery::expect_called(listener, 1) + # Only the individuals from the first round of scheduling are targeted. + expect_targeted_listener(listener, 1, t = 3, target = c(2, 4)) + new_event$.tick() + + #time=4 + new_event$.process() + mockery::expect_called(listener, 1) +}) diff --git a/vignettes/Checkpoint.Rmd b/vignettes/Checkpoint.Rmd index 08a657d..70d1ff6 100644 --- a/vignettes/Checkpoint.Rmd +++ b/vignettes/Checkpoint.Rmd @@ -35,7 +35,7 @@ When modeling the impact of an intervention on a disease, it is common to have a Individual allows the user to run a simulation for a number of time steps, save the state of the simulation and resume it multiple times, with different parameters each time. This way, the initial phase before the intervention only needs to be simulated once. -## Usage {usage} +## Usage {#usage} The typical way to use this feature is to define a simulation function which creates all the relevant simulation data and then calls `simulation_loop`. The function we define takes in an optional `state` parameter that is passed through to `simulation_loop`. @@ -57,7 +57,7 @@ state <- run_simulation(timesteps = 50) ``` Finally, the simulation is resumed with a larger number of time steps, passing in the state object as an argument. The `timesteps` argument refers to the total number of time steps, including both the original simulation run and the new one. In this case, `run_simulation` will only simulate 50 extra steps. Before running the actual simulation, `simulation_loop` will reload the simulation state from its argument, overwriting any values we had set when initializing the variables. -```{r} +```{r, results='hide'} run_simulation(timesteps = 100, state=state) ``` @@ -99,16 +99,20 @@ make_infection_process <- function(health, vaccinated, N, beta, vaccine_efficacy } ``` -At the start of the simulation no vaccination takes place. Only after a number of time steps, determined by the `vaccination_ts` parameter, does the intervention begin. At each time step after that, an individual has a fixed probability of becoming vaccinated. If `vaccination_ts` is `NULL`, the intervention never begins. +For a while at the start of the simulation no vaccination takes place. Only after a number of time steps, determined by the `vaccination_start` parameter, does the intervention begin. Periodically, an event will fire and every individual has a fixed probability of becoming vaccinated. If `vaccination_start` is `NULL`, the intervention never begins. ```{r} -make_vaccination_process <- function(vaccinated, vaccination_ts, vaccination_rate) { - function(t) { - if (!is.null(vaccination_ts) && t >= vaccination_ts) { - vaccinated$queue_update(value = "Y", - vaccinated$get_index_of("N")$sample(vaccination_rate)) - } +make_vaccination_event <- function(vaccinated, vaccination_start, vaccination_interval, vaccination_rate) { + e <- Event$new() + e$add_listener(function(t) { + vaccinated$queue_update(value = "Y", + vaccinated$get_index_of("N")$sample(vaccination_rate)) + e$schedule(vaccination_interval) + }) + if (!is.null(vaccination_start)) { + e$schedule(vaccination_start - 1) } + e } ``` @@ -122,9 +126,10 @@ run_simulation <- function( beta = 0.1, # S -> I gamma = 0.05, # I -> R xi = 0.03, # R -> S - vaccination_ts = NULL, + vaccination_start = NULL, + vaccination_interval = 10, + vaccination_rate = 0.05, # N -> Y vaccine_efficacy = 1, - vaccination_rate = 0.005, # N -> Y state = NULL) { variables <- make_variables(N, I0) @@ -134,8 +139,8 @@ run_simulation <- function( N, beta, vaccine_efficacy) recovery_process <- bernoulli_process(variables$health, "I", "R", gamma) return_process <- bernoulli_process(variables$health, "R", "S", xi) - vaccination_process <- make_vaccination_process( - variables$vaccinated, vaccination_ts, vaccination_rate) + vaccination_event <- make_vaccination_event( + variables$vaccinated, vaccination_start, vaccination_interval, vaccination_rate) renderer <- Render$new(timesteps = steps) health_render_process <- categorical_count_renderer_process( @@ -148,11 +153,11 @@ run_simulation <- function( infection_process, recovery_process, return_process, - vaccination_process, health_render_process) final_state <- simulation_loop( variables=variables, + events=list(vaccination_event), processes=processes, timesteps=steps, state=state @@ -185,7 +190,7 @@ We see that the simulation takes some time to settle from its initial parameters We will now enable the vaccine intervention, but only starting at a point after the simulation has settled, for example at `t=500`. ```{r} -data <- run_simulation(steps=1500, vaccination_ts = 500, vaccine_efficacy = 1)$result +data <- run_simulation(steps=1500, vaccination_start = 500, vaccine_efficacy = 1)$result colours <- c("royalblue3","firebrick3","darkorchid3") matplot( x=data["timestep"], @@ -202,19 +207,19 @@ legend( ) ``` -The simulation above clearly shows the effect of the vaccination campaign, starting at `ts=500`. However, it made the optimistic assumption of a 100% vaccine efficacy. We wish to run the simulation again but with varying levels of efficacy, in order the compare its impact. +The simulation above clearly shows the effect of the vaccination campaign, starting at `t=500`. However, it made the optimistic assumption of a 100% vaccine efficacy. We wish to run the simulation again but with varying levels of efficacy, in order the compare its impact. -While we could run the code above many times over, each simulation would repeat the first 499 timesteps, despite the result being identical each time. Instead we start by running only these timesteps, and saving the result. We omit any intervention parameters, since these are irrelevant anyway. +While we could run the code above many times over, each simulation would repeat the first 499 timesteps, despite the result being identical each time. Instead we start by running only these timesteps, and saving the result. We do need to specify the start of the intervention, as it is necessary to schedule the first vaccination event. However the details of the intervention (ie. `vaccine_efficacy`) are irrelevant and can be omitted. ```{r} -initial <- run_simulation(steps=499) +initial <- run_simulation(steps=1, vaccination_start = 500) ``` From this initial result, we can resume the simulation, but using different values of vaccine efficacy each time. We also include a control simulation, in which no vaccination takes place. Each of these simulation will skip the first 499 steps and only run the next 1001 time steps. ```{r} -control <- run_simulation(steps=1500, state=initial$state)$result -vaccine30 <- run_simulation(steps=1500, vaccination_ts = 500, vaccine_efficacy=0.3, state=initial$state)$result -vaccine50 <- run_simulation(steps=1500, vaccination_ts = 500, vaccine_efficacy=0.5, state=initial$state)$result -vaccine100 <- run_simulation(steps=1500, vaccination_ts = 500, vaccine_efficacy=1.0, state=initial$state)$result +control <- run_simulation(steps=1500, vaccination_start = 500, state=initial$state)$result +vaccine30 <- run_simulation(steps=1500, vaccination_start = 500, vaccine_efficacy=0.3, state=initial$state)$result +vaccine50 <- run_simulation(steps=1500, vaccination_start = 500, vaccine_efficacy=0.5, state=initial$state)$result +vaccine100 <- run_simulation(steps=1500, vaccination_start = 500, vaccine_efficacy=1.0, state=initial$state)$result ``` Finally we aggregate and plot the results from all these simulations. We also need to include the data from our initial run, which we will plot the same colour as our control simulation. @@ -249,7 +254,7 @@ Saving and restoring the simulation state comes with a number of caveats. - The state object's structure is not stable and is expected to change. One should not expect to serialize the state to disk and have it work with future versions of the individual package. - The simulation must be re-created in an identical way. Variables and events may not be added or removed, variable sizes must remain constant, the list of categories in a `CategoricalVariable` cannot be modified, etc. The order of variables and events passed to the `run_simulation` function must remain stable. - Restoring a simulation state also restores R's global random number generator state. This can have side effects on other parts of a program. -- Events are not yet supported, although this is planned. +- If an event is scheduled before the checkpoint, the time at which it will execute cannot be changed when resuming, even if that time is in the future. For example in the SIRS model above, we would not be able to resume the simulation with different values for `vaccination_start`; changing that parameter would have no effect. While parameters of the simulation can be changed between the initial run and the subsequent runs (as demonstrated with the `vaccine_efficacy` parameter above), in general you should not modify parameters that would have been already had an impact on the first part of the simulation. Doing so would produce results that can only be produced through checkpoint and resume, and not as a single simulation.