From b21c601052f85e392e48d903b8903a1a392fe786 Mon Sep 17 00:00:00 2001
From: Johannes Ranke <jranke@uni-bremen.de>
Date: Tue, 14 Apr 2015 19:50:57 +0200
Subject: Compile differential equation models with ccSolve package

If the ccSolve package is available, and time is not in the right hand
side of the equations (i.e. if only SFO and SFORB models are used),
the differential equation model is compiled from auto-generated C code.

Currently, one test (FOCUS 2006 D SFO_SFO) fails
---
 R/mkinmod.R     | 34 ++++++++++++++++++++++++++++++++++
 R/mkinpredict.R | 22 +++++++++++++---------
 2 files changed, 47 insertions(+), 9 deletions(-)

(limited to 'R')

diff --git a/R/mkinmod.R b/R/mkinmod.R
index d04f811f..37cd02b4 100644
--- a/R/mkinmod.R
+++ b/R/mkinmod.R
@@ -36,6 +36,9 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL)
   # differential equations (if supported), parameter names and a mapping from
   # model variables to observed variables. If possible, a matrix representation
   # of the differential equations is included
+  # Compiling the functions from the C code generated below using the ccSolve package
+  # only works if the implicit assumption about differential equations specified below
+  # is satisfied
   parms <- vector()
   # }}}
 
@@ -269,6 +272,37 @@ mkinmod <- function(..., use_of_ff = "min", speclist = NULL)
     model$coefmat <- m
   }#}}}
 
+  # Create a function compiled from C code if possible #{{{
+  if (require(ccSolve)) {
+    diffs.C <- paste(diffs, collapse = ";\n")
+    diffs.C <- paste0(diffs.C, ";")
+    for (i in seq_along(diffs)) {
+      obs_var <- names(diffs)[i]
+
+      # Replace d_... terms by f[i-1]
+      # First line
+      pattern <- paste0("^d_", obs_var)
+      replacement <- paste0("\nf[", i - 1, "]")
+      diffs.C <- gsub(pattern, replacement, diffs.C)
+      # Other lines
+      pattern <- paste0("\\nd_", obs_var)
+      replacement <- paste0("\nf[", i - 1, "]")
+      diffs.C <- gsub(pattern, replacement, diffs.C)
+
+      # Replace names of observed variables by y[i],
+      # making the implicit assumption that the observed variables only occur after "* "
+      pattern <- paste0("\\* ", obs_var)
+      replacement <- paste0("* y[", i - 1, "]")
+      diffs.C <- gsub(pattern, replacement, diffs.C)
+    }
+    # Unfortunately, the models with time occuring in the diffs do not compile
+    if (sum(sapply(spec, function(x) x$type %in% c("SFO", "SFORB"))) == length(spec)) {
+      message("Compiling differential equation model from auto-generated C code...")
+      model$compiled <- compile.ode(diffs.C, language = "C", parms = parms)
+    }
+  }
+  # }}}
+
   class(model) <- "mkinmod"
   return(model)
 }
diff --git a/R/mkinpredict.R b/R/mkinpredict.R
index da013d50..bed82441 100644
--- a/R/mkinpredict.R
+++ b/R/mkinpredict.R
@@ -87,17 +87,21 @@ mkinpredict <- function(mkinmod, odeparms, odeini,
     names(out) <- c("time", mod_vars)
   } 
   if (solution_type == "deSolve") {
-    mkindiff <- function(t, state, parms) {
+    if (!is.null(mkinmod$compiled)) {
+       mkindiff <- mkinmod$compiled
+     } else {
+       mkindiff <- function(t, state, parms) {
 
-      time <- t
-      diffs <- vector()
-      for (box in names(mkinmod$diffs))
-      {
-        diffname <- paste("d", box, sep="_")      
-        diffs[diffname] <- with(as.list(c(time, state, parms)),
-          eval(parse(text=mkinmod$diffs[[box]])))
+        time <- t
+        diffs <- vector()
+        for (box in names(mkinmod$diffs))
+        {
+          diffname <- paste("d", box, sep="_")      
+          diffs[diffname] <- with(as.list(c(time, state, parms)),
+            eval(parse(text=mkinmod$diffs[[box]])))
+        }
+        return(list(c(diffs)))
       }
-      return(list(c(diffs)))
     } 
     out <- ode(
       y = odeini,
-- 
cgit v1.2.1