Clarify ode methods & continue when min step size reached: Document c…

…ode and vignettes (#413)

This PR better documents the ODE stepper used in plant, both in code and
via a new article; and enable system to continue when  minimum step size is reached. 
In addition,

- move relevant files for ode_solver into a subdirectory
`inst/include/plant/ode_solver`
- rename some items for clarity
   - `advance` -> `advance_adaptive`
   - `node_schedule_ode_times` -> `ode_times`
   - rename file `plant_tools` -> `individual_tools`
- clarify how to control time steps (closes #411, needed for TF24 Soil
water milestone)

Gives rise to #415

R/RcppExports.R

@@ -37,8 +37,8 @@ OdeRunner___Lorenz__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___Lorenz__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___Lorenz__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___Lorenz__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___Lorenz__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___Lorenz__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___Lorenz__advance_fixed <- function(obj_, time) {
@@ -81,8 +81,8 @@ OdeRunner___OdeR__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___OdeR__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___OdeR__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___OdeR__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___OdeR__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___OdeR__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___OdeR__advance_fixed <- function(obj_, time) {
@@ -125,8 +125,8 @@ OdeRunner___FF16__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___FF16__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___FF16__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___FF16__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___FF16__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___FF16__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___FF16__advance_fixed <- function(obj_, time) {
@@ -169,8 +169,8 @@ OdeRunner___TF24__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___TF24__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___TF24__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___TF24__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___TF24__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___TF24__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___TF24__advance_fixed <- function(obj_, time) {
@@ -213,8 +213,8 @@ OdeRunner___FF16w__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___FF16w__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___FF16w__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___FF16w__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___FF16w__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___FF16w__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___FF16w__advance_fixed <- function(obj_, time) {
@@ -257,8 +257,8 @@ OdeRunner___FF16r__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___FF16r__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___FF16r__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___FF16r__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___FF16r__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___FF16r__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___FF16r__advance_fixed <- function(obj_, time) {
@@ -301,8 +301,8 @@ OdeRunner___K93__ctor <- function(obj, control) {
     .Call('_plant_OdeRunner___K93__ctor', PACKAGE = 'plant', obj, control)
 }
 
-OdeRunner___K93__advance <- function(obj_, time) {
-    invisible(.Call('_plant_OdeRunner___K93__advance', PACKAGE = 'plant', obj_, time))
+OdeRunner___K93__advance_adaptive <- function(obj_, time) {
+    invisible(.Call('_plant_OdeRunner___K93__advance_adaptive', PACKAGE = 'plant', obj_, time))
 }
 
 OdeRunner___K93__advance_fixed <- function(obj_, time) {

R/RcppR6.R

@@ -1,6 +1,6 @@
 ## Generated by RcppR6: do not edit by hand
 ## Version: 0.2.4
-## Hash: 0238ec3fc83ef4c3da5199ad96b782b0
+## Hash: 86d2dab10f9081823f16a63c8035a2ae
 
 ##' @importFrom Rcpp evalCpp
 ##' @importFrom R6 R6Class
@@ -107,8 +107,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___Lorenz__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___Lorenz__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___Lorenz__advance_fixed(self, time)
@@ -169,8 +169,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___OdeR__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___OdeR__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___OdeR__advance_fixed(self, time)
@@ -231,8 +231,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___FF16__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___FF16__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___FF16__advance_fixed(self, time)
@@ -293,8 +293,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___TF24__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___TF24__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___TF24__advance_fixed(self, time)
@@ -355,8 +355,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___FF16w__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___FF16w__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___FF16w__advance_fixed(self, time)
@@ -417,8 +417,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___FF16r__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___FF16r__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___FF16r__advance_fixed(self, time)
@@ -479,8 +479,8 @@ OdeRunner <- function(T) {
       initialize = function(ptr) {
         self$.ptr <- ptr
       },
-      advance = function(time) {
-        OdeRunner___K93__advance(self, time)
+      advance_adaptive = function(time) {
+        OdeRunner___K93__advance_adaptive(self, time)
       },
       advance_fixed = function(time) {
         OdeRunner___K93__advance_fixed(self, time)

R/build_schedule.R

@@ -1,4 +1,4 @@
-##' Build an appropriately refined schedule.
+##' Build an appropriately refined schedule for node introduction.
 ##'
 ##' There are control options (within the \code{Parameters} object)
 ##' that affect how this function runs, in particular
@@ -51,7 +51,7 @@ build_schedule <- function(p, env = make_environment(parameters = p),
     plant_log_debug(msg, routine="schedule", event="split", round=i)
   }
 
-  p$node_schedule_ode_times <- res$ode_times
+  p$ode_times <- res$ode_times
   ## Useful to record the last offspring produced:
 
   attr(p, "offspring_production") <- offspring_production

inst/RcppR6_classes.yml

@@ -36,7 +36,7 @@ OdeRunner:
     times: {type: "std::vector<double>", access: member}
     object: {type: T, access: member}
   methods:
-    advance: {return_type: void, args: [time: double]}
+    advance_adaptive: {return_type: void, args: [time: double]}
     advance_fixed: {return_type: void, args: [time: "std::vector<double>"]}
     step: {return_type: void}
     step_to: {return_type: void, args: [time: double]}
@@ -354,7 +354,7 @@ Parameters:
     - strategy_default: "T"
     - node_schedule_times_default: "std::vector<double>"
     - node_schedule_times: "std::vector<std::vector<double> >"
-    - node_schedule_ode_times: "std::vector<double>"
+    - ode_times: "std::vector<double>"
 
 Node:
   name_cpp: "plant::Node<T,E>"

inst/include/plant.h

@@ -9,10 +9,10 @@
 #include <plant/interpolator.h>
 #include <plant/adaptive_interpolator.h>
 
-#include <plant/ode_control.h>
-#include <plant/ode_step.h>
-#include <plant/ode_solver.h>
-#include <plant/ode_runner.h>
+#include <plant/ode_solver/ode_control.h>
+#include <plant/ode_solver/ode_step.h>
+#include <plant/ode_solver/ode_solver.h>
+#include <plant/ode_solver/ode_runner.h>
 
 #include <plant/environment.h>
 #include <plant/resource_spline.h>
@@ -48,10 +48,10 @@
 #include <plant/stochastic_patch.h>
 #include <plant/stochastic_patch_runner.h>
 
-#include <plant/plant_runner.h>
+#include <plant/individual_runner.h>
 
 // Purely for testing
-#include <plant/lorenz.h>
+#include <plant/ode_solver/lorenz.h>
 
 // Include this early on.  It can be either after classes have been
 // declared (but before Rcpp has been loaded) or first.  This file will
@@ -65,14 +65,13 @@
 // them earlier up.
 
 #include <Rcpp.h>
-#include <plant/ode_r.h>
+#include <plant/ode_solver/ode_r.h>
 
 // This line can safely be the last line in the file, but may go any
 // point after RcppR6_pre.hpp is included.
 #include <plant/RcppR6_post.hpp>
 #include <plant/util_post_rcpp.h>
 #include <plant/get_state.h>
 #include <plant/get_aux.h>
-#include <plant/individual_tools.h>
 
 #endif

inst/include/plant/RcppR6_post.hpp

@@ -596,7 +596,7 @@ template <> inline SEXP wrap(const plant::Parameters<plant::FF16_Strategy,plant:
   ret["strategy_default"] = Rcpp::wrap(x.strategy_default);
   ret["node_schedule_times_default"] = Rcpp::wrap(x.node_schedule_times_default);
   ret["node_schedule_times"] = Rcpp::wrap(x.node_schedule_times);
-  ret["node_schedule_ode_times"] = Rcpp::wrap(x.node_schedule_ode_times);
+  ret["ode_times"] = Rcpp::wrap(x.ode_times);
   ret.attr("class") = Rcpp::CharacterVector::create("Parameters<FF16,FF16_Env>", "Parameters");
   return ret;
 }
@@ -625,8 +625,8 @@ template <> inline plant::Parameters<plant::FF16_Strategy,plant::FF16_Environmen
   ret.node_schedule_times_default = Rcpp::as<std::vector<double> >(xl["node_schedule_times_default"]);
   // ret.node_schedule_times = Rcpp::as<decltype(retnode_schedule_times) >(xl["node_schedule_times"]);
   ret.node_schedule_times = Rcpp::as<std::vector<std::vector<double> > >(xl["node_schedule_times"]);
-  // ret.node_schedule_ode_times = Rcpp::as<decltype(retnode_schedule_ode_times) >(xl["node_schedule_ode_times"]);
-  ret.node_schedule_ode_times = Rcpp::as<std::vector<double> >(xl["node_schedule_ode_times"]);
+  // ret.ode_times = Rcpp::as<decltype(retode_times) >(xl["ode_times"]);
+  ret.ode_times = Rcpp::as<std::vector<double> >(xl["ode_times"]);
   ret.validate();
   return ret;
 }
@@ -641,7 +641,7 @@ template <> inline SEXP wrap(const plant::Parameters<plant::TF24_Strategy,plant:
   ret["strategy_default"] = Rcpp::wrap(x.strategy_default);
   ret["node_schedule_times_default"] = Rcpp::wrap(x.node_schedule_times_default);
   ret["node_schedule_times"] = Rcpp::wrap(x.node_schedule_times);
-  ret["node_schedule_ode_times"] = Rcpp::wrap(x.node_schedule_ode_times);
+  ret["ode_times"] = Rcpp::wrap(x.ode_times);
   ret.attr("class") = Rcpp::CharacterVector::create("Parameters<TF24,TF24_Env>", "Parameters");
   return ret;
 }
@@ -670,8 +670,8 @@ template <> inline plant::Parameters<plant::TF24_Strategy,plant::TF24_Environmen
   ret.node_schedule_times_default = Rcpp::as<std::vector<double> >(xl["node_schedule_times_default"]);
   // ret.node_schedule_times = Rcpp::as<decltype(retnode_schedule_times) >(xl["node_schedule_times"]);
   ret.node_schedule_times = Rcpp::as<std::vector<std::vector<double> > >(xl["node_schedule_times"]);
-  // ret.node_schedule_ode_times = Rcpp::as<decltype(retnode_schedule_ode_times) >(xl["node_schedule_ode_times"]);
-  ret.node_schedule_ode_times = Rcpp::as<std::vector<double> >(xl["node_schedule_ode_times"]);
+  // ret.ode_times = Rcpp::as<decltype(retode_times) >(xl["ode_times"]);
+  ret.ode_times = Rcpp::as<std::vector<double> >(xl["ode_times"]);
   ret.validate();
   return ret;
 }
@@ -686,7 +686,7 @@ template <> inline SEXP wrap(const plant::Parameters<plant::FF16w_Strategy,plant
   ret["strategy_default"] = Rcpp::wrap(x.strategy_default);
   ret["node_schedule_times_default"] = Rcpp::wrap(x.node_schedule_times_default);
   ret["node_schedule_times"] = Rcpp::wrap(x.node_schedule_times);
-  ret["node_schedule_ode_times"] = Rcpp::wrap(x.node_schedule_ode_times);
+  ret["ode_times"] = Rcpp::wrap(x.ode_times);
   ret.attr("class") = Rcpp::CharacterVector::create("Parameters<FF16w,FF16_Env>", "Parameters");
   return ret;
 }
@@ -715,8 +715,8 @@ template <> inline plant::Parameters<plant::FF16w_Strategy,plant::FF16_Environme
   ret.node_schedule_times_default = Rcpp::as<std::vector<double> >(xl["node_schedule_times_default"]);
   // ret.node_schedule_times = Rcpp::as<decltype(retnode_schedule_times) >(xl["node_schedule_times"]);
   ret.node_schedule_times = Rcpp::as<std::vector<std::vector<double> > >(xl["node_schedule_times"]);
-  // ret.node_schedule_ode_times = Rcpp::as<decltype(retnode_schedule_ode_times) >(xl["node_schedule_ode_times"]);
-  ret.node_schedule_ode_times = Rcpp::as<std::vector<double> >(xl["node_schedule_ode_times"]);
+  // ret.ode_times = Rcpp::as<decltype(retode_times) >(xl["ode_times"]);
+  ret.ode_times = Rcpp::as<std::vector<double> >(xl["ode_times"]);
   ret.validate();
   return ret;
 }
@@ -731,7 +731,7 @@ template <> inline SEXP wrap(const plant::Parameters<plant::FF16r_Strategy,plant
   ret["strategy_default"] = Rcpp::wrap(x.strategy_default);
   ret["node_schedule_times_default"] = Rcpp::wrap(x.node_schedule_times_default);
   ret["node_schedule_times"] = Rcpp::wrap(x.node_schedule_times);
-  ret["node_schedule_ode_times"] = Rcpp::wrap(x.node_schedule_ode_times);
+  ret["ode_times"] = Rcpp::wrap(x.ode_times);
   ret.attr("class") = Rcpp::CharacterVector::create("Parameters<FF16r,FF16_Env>", "Parameters");
   return ret;
 }
@@ -760,8 +760,8 @@ template <> inline plant::Parameters<plant::FF16r_Strategy,plant::FF16_Environme
   ret.node_schedule_times_default = Rcpp::as<std::vector<double> >(xl["node_schedule_times_default"]);
   // ret.node_schedule_times = Rcpp::as<decltype(retnode_schedule_times) >(xl["node_schedule_times"]);
   ret.node_schedule_times = Rcpp::as<std::vector<std::vector<double> > >(xl["node_schedule_times"]);
-  // ret.node_schedule_ode_times = Rcpp::as<decltype(retnode_schedule_ode_times) >(xl["node_schedule_ode_times"]);
-  ret.node_schedule_ode_times = Rcpp::as<std::vector<double> >(xl["node_schedule_ode_times"]);
+  // ret.ode_times = Rcpp::as<decltype(retode_times) >(xl["ode_times"]);
+  ret.ode_times = Rcpp::as<std::vector<double> >(xl["ode_times"]);
   ret.validate();
   return ret;
 }
@@ -776,7 +776,7 @@ template <> inline SEXP wrap(const plant::Parameters<plant::K93_Strategy,plant::
   ret["strategy_default"] = Rcpp::wrap(x.strategy_default);
   ret["node_schedule_times_default"] = Rcpp::wrap(x.node_schedule_times_default);
   ret["node_schedule_times"] = Rcpp::wrap(x.node_schedule_times);
-  ret["node_schedule_ode_times"] = Rcpp::wrap(x.node_schedule_ode_times);
+  ret["ode_times"] = Rcpp::wrap(x.ode_times);
   ret.attr("class") = Rcpp::CharacterVector::create("Parameters<K93,K93_Env>", "Parameters");
   return ret;
 }
@@ -805,8 +805,8 @@ template <> inline plant::Parameters<plant::K93_Strategy,plant::K93_Environment>
   ret.node_schedule_times_default = Rcpp::as<std::vector<double> >(xl["node_schedule_times_default"]);
   // ret.node_schedule_times = Rcpp::as<decltype(retnode_schedule_times) >(xl["node_schedule_times"]);
   ret.node_schedule_times = Rcpp::as<std::vector<std::vector<double> > >(xl["node_schedule_times"]);
-  // ret.node_schedule_ode_times = Rcpp::as<decltype(retnode_schedule_ode_times) >(xl["node_schedule_ode_times"]);
-  ret.node_schedule_ode_times = Rcpp::as<std::vector<double> >(xl["node_schedule_ode_times"]);
+  // ret.ode_times = Rcpp::as<decltype(retode_times) >(xl["ode_times"]);
+  ret.ode_times = Rcpp::as<std::vector<double> >(xl["ode_times"]);
   ret.validate();
   return ret;
 }

inst/include/plant/control.h

@@ -3,7 +3,7 @@
 #define PLANT_PLANT_CONTROL_H_
 
 #include <plant/qag.h>
-#include <plant/ode_control.h>
+#include <plant/ode_solver/ode_control.h>
 #include <string>
 
 // The `Control` object holds all the non-biological control

inst/include/plant/environment.h

@@ -5,7 +5,7 @@
 #include <plant/control.h>
 #include <plant/interpolator.h>
 #include <plant/adaptive_interpolator.h>
-#include <plant/ode_interface.h>
+#include <plant/ode_solver/ode_interface.h>
 #include <plant/internals.h>
 #include <plant/util.h>
 #include <unordered_map>

inst/include/plant/individual.h

@@ -3,7 +3,7 @@
 #define PLANT_PLANT_PLANT_MINIMAL_H_
 
 #include <memory> // std::shared_ptr
-#include <plant/ode_interface.h>
+#include <plant/ode_solver/ode_interface.h>
 #include <vector>
 #include <plant/internals.h>
 #include <plant/uniroot.h>

inst/include/plant/internals.h

@@ -7,7 +7,7 @@
 #define FECUNDITY_INDEX 2
 
 #include <memory> // std::shared_ptr
-#include <plant/ode_interface.h>
+#include <plant/ode_solver/ode_interface.h>
 #include <vector>
 // #include <plant/plant_internals.h>
 

inst/include/plant/node.h

@@ -4,7 +4,7 @@
 
 #include <plant/environment.h>
 #include <plant/gradient.h>
-#include <plant/ode_interface.h>
+#include <plant/ode_solver/ode_interface.h>
 
 namespace plant {
 

inst/include/plant/lorenz.h → inst/include/plant/ode_solver/lorenz.h

@@ -2,7 +2,7 @@
 #ifndef PLANT_PLANT_LORENZ_H_
 #define PLANT_PLANT_LORENZ_H_
 
-#include <plant/ode_interface.h>
+#include <plant/ode_solver/ode_interface.h>
 
 namespace plant {
 namespace ode {