aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--DESCRIPTION2
-rw-r--r--R/mkinfit.R13
-rw-r--r--R/mkinpredict.R90
-rw-r--r--vignettes/web_only/benchmarks.html99
-rw-r--r--vignettes/web_only/mkin_benchmarks.rdabin1722 -> 1777 bytes
5 files changed, 104 insertions, 100 deletions
diff --git a/DESCRIPTION b/DESCRIPTION
index 2a3f231e..83a1a1c9 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -23,7 +23,7 @@ Description: Calculation routines based on the FOCUS Kinetics Report (2006,
note that no warranty is implied for correctness of results or fitness for a
particular purpose.
Depends: R (>= 2.15.1),
-Imports: stats, graphics, methods, parallel, deSolve, R6, inline (>= 0.3.19),
+Imports: stats, graphics, methods, parallel, deSolve (>= 1.35), R6, inline (>= 0.3.19),
numDeriv, lmtest, pkgbuild, nlme (>= 3.1-151), saemix (>= 3.1), rlang, vctrs
Suggests: knitr, rbenchmark, tikzDevice, testthat, rmarkdown, covr, vdiffr,
benchmarkme, tibble, stats4, readxl
diff --git a/R/mkinfit.R b/R/mkinfit.R
index 0d9246dd..6cca5616 100644
--- a/R/mkinfit.R
+++ b/R/mkinfit.R
@@ -501,12 +501,10 @@ mkinfit <- function(mkinmod, observed,
}
# Get native symbol before iterations info for speed
- call_lsoda <- getNativeSymbolInfo("call_lsoda", PACKAGE = "deSolve")
if (solution_type == "deSolve" & use_compiled[1] != FALSE) {
- mkinmod$diffs_address <- getNativeSymbolInfo("diffs",
- PACKAGE = mkinmod$dll_info[["name"]])$address
- mkinmod$initpar_address <- getNativeSymbolInfo("initpar",
- PACKAGE = mkinmod$dll_info[["name"]])$address
+ mkinmod[["symbols"]] <- deSolve::checkDLL(dllname = mkinmod$dll_info[["name"]],
+ func = "diffs", initfunc = "initpar",
+ jacfunc = NULL, nout = 0, outnames = NULL)
}
# Get the error model and the algorithm for fitting
@@ -903,9 +901,8 @@ mkinfit <- function(mkinmod, observed,
fit$time <- fit_time
# We also need the model and a model name for summary and plotting,
- # but without address info that will become invalid
- mkinmod$diffs_address <- NULL
- mkinmod$initpar_address <- NULL
+ # but without symbols because they could become invalid
+ fit$symbols <- NULL
fit$mkinmod <- mkinmod
fit$mkinmod$name <- mkinmod_name
fit$obs_vars <- obs_vars
diff --git a/R/mkinpredict.R b/R/mkinpredict.R
index 11a3b35b..f0d00319 100644
--- a/R/mkinpredict.R
+++ b/R/mkinpredict.R
@@ -24,6 +24,7 @@
#' parent compound.
#' @param method.ode The solution method passed via [mkinpredict] to [ode]] in
#' case the solution type is "deSolve" and we are not using compiled code.
+#' When using compiled code, only lsoda is supported.
#' @param use_compiled If set to \code{FALSE}, no compiled version of the
#' [mkinmod] model is used, even if is present.
#' @param atol Absolute error tolerance, passed to the ode solver.
@@ -34,7 +35,6 @@
#' FALSE). Setting this to FALSE has no effect for analytical solutions,
#' as these always return mapped output.
#' @param na_stop Should it be an error if [ode] returns NaN values
-#' @param call_lsoda The address of the compiled function "call_lsoda"
#' @param \dots Further arguments passed to the ode solver in case such a
#' solver is used.
#' @return A matrix with the numeric solution in wide format
@@ -56,11 +56,11 @@
#' mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
#' solution_type = "analytical")[21,]
#' mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
-#' method = "lsoda")[21,]
+#' method = "lsoda", use_compiled = FALSE)[21,]
#' mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
-#' method = "ode45")[21,]
+#' method = "ode45", use_compiled = FALSE)[21,]
#' mkinpredict(SFO, c(k_degradinol = 0.3), c(degradinol = 100), 0:20,
-#' method = "rk4")[21,]
+#' method = "rk4", use_compiled = FALSE)[21,]
#' # rk4 is not as precise here
#'
#' # The number of output times used to make a lot of difference until the
@@ -172,79 +172,19 @@ mkinpredict.mkinmod <- function(x,
if (solution_type == "deSolve") {
if (!is.null(x$cf) & use_compiled[1] != FALSE) {
-# out <- deSolve::ode(
-# y = odeini,
-# times = outtimes,
-# func = "diffs",
-# initfunc = "initpar",
-# dllname = x$dll_info[["name"]],
-# parms = odeparms[x$parms], # Order matters when using compiled models
-# method = method.ode,
-# atol = atol,
-# rtol = rtol,
-# maxsteps = maxsteps,
-# ...
-# )
-#
- # Prepare call to "call_lsoda"
- # Simplified code from deSolve::lsoda() adapted to our use case
- if (is.null(call_lsoda)) {
- call_lsoda <- getNativeSymbolInfo("call_lsoda", PACKAGE = "deSolve")
- }
- if (is.null(x$diffs_address)) {
- x$diffs_address <- getNativeSymbolInfo("diffs",
- PACKAGE = x$dll_info[["name"]])$address
- x$initpar_address <- getNativeSymbolInfo("initpar",
- PACKAGE = x$dll_info[["name"]])$address
- }
- rwork <- vector("double", 20)
- rwork[] <- 0.
- rwork[6] <- max(abs(diff(outtimes)))
-
- iwork <- vector("integer", 20)
- iwork[] <- 0
- iwork[6] <- maxsteps
-
- n <- length(odeini)
- lmat <- n^2 + 2 # from deSolve::lsoda(), for Jacobian type full, internal (2)
- # hard-coded default values of lsoda():
- maxordn <- 12L
- maxords <- 5L
- lrn <- 20 + n * (maxordn + 1) + 3 * n # length in case non-stiff method
- lrs <- 20 + n * (maxords + 1) + 3 * n + lmat # length in case stiff method
- lrw <- max(lrn, lrs) # actual length: max of both
- liw <- 20 + n
-
- on.exit(.C("unlock_solver", PACKAGE = "deSolve"))
-
- out_raw <- .Call(call_lsoda,
- as.double(odeini), # y
- as.double(outtimes), # times
- x$diffs_address, # derivfunc
- as.double(odeparms[x$parms]), # parms
- rtol, atol,
- NULL, NULL, # rho, tcrit
- NULL, # jacfunc
- x$initpar_address, # initfunc
- NULL, # eventfunc
- 0L, # verbose
- 1L, # iTask
- as.double(rwork), # rWork
- as.integer(iwork), # iWork
- 2L, # jT full Jacobian calculated internally
- 0L, # nOut
- as.integer(lrw), # lRw
- as.integer(liw), # lIw
- 1L, # Solver
- NULL, # rootfunc
- 0L, as.double(0), 0L, # nRoot, Rpar, Ipar
- 0L, # Type
- list(fmat = 0L, tmat = 0L, imat = 0L, ModelForc = NULL), # flist
- list(), # elist
- list(islag = 0L) # elag
+ out <- deSolve::lsoda(
+ y = odeini,
+ times = outtimes,
+ func = mkinmod[["symbols"]],
+ initfunc = "initpar",
+ dllname = x$dll_info[["name"]],
+ parms = odeparms[x$parms], # Order matters when using compiled models
+ atol = atol,
+ rtol = rtol,
+ maxsteps = maxsteps,
+ ...
)
- out <- t(out_raw)
colnames(out) <- c("time", mod_vars)
} else {
mkindiff <- function(t, state, parms) {
diff --git a/vignettes/web_only/benchmarks.html b/vignettes/web_only/benchmarks.html
index 451d8afe..fb799462 100644
--- a/vignettes/web_only/benchmarks.html
+++ b/vignettes/web_only/benchmarks.html
@@ -1592,7 +1592,7 @@ div.tocify {
<h1 class="title toc-ignore">Benchmark timings for mkin</h1>
<h4 class="author">Johannes Ranke</h4>
-<h4 class="date">Last change 14 July 2022 (rebuilt 2023-01-05)</h4>
+<h4 class="date">Last change 14 July 2022 (rebuilt 2023-02-13)</h4>
</div>
@@ -1636,11 +1636,31 @@ FOMC_SFO &lt;- mkinmod(
DFOP_SFO &lt;- mkinmod(
parent = mkinsub(&quot;FOMC&quot;, &quot;m1&quot;), # erroneously used FOMC twice, not fixed for consistency
m1 = mkinsub(&quot;SFO&quot;))
-t3 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D)))[[&quot;elapsed&quot;]]
-t4 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D),
- error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]
-t5 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D),
+t3 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D)))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t4 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D),
+ error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t5 &lt;- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D),
error_model = &quot;obs&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
<p>Two metabolites, synthetic data:</p>
<pre class="r"><code>m_synth_SFO_lin &lt;- mkinmod(parent = mkinsub(&quot;SFO&quot;, &quot;M1&quot;),
M1 = mkinsub(&quot;SFO&quot;, &quot;M2&quot;),
@@ -1656,18 +1676,36 @@ SFO_lin_a &lt;- synthetic_data_for_UBA_2014[[1]]$data
DFOP_par_c &lt;- synthetic_data_for_UBA_2014[[12]]$data
-t6 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a)))[[&quot;elapsed&quot;]]
-t7 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c)))[[&quot;elapsed&quot;]]
-
-t8 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a),
- error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]
-t9 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c),
- error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]
-
-t10 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a),
- error_model = &quot;obs&quot;))[[&quot;elapsed&quot;]]
-t11 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c),
+t6 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a)))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t7 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c)))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t8 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a),
+ error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t9 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c),
+ error_model = &quot;tc&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t10 &lt;- system.time(mmkin_bench(list(m_synth_SFO_lin), list(SFO_lin_a),
+ error_model = &quot;obs&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
+<pre class="r"><code>t11 &lt;- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c),
error_model = &quot;obs&quot;))[[&quot;elapsed&quot;]]</code></pre>
+<pre><code>## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable
+## Error in mkinmod[[&quot;symbols&quot;]] :
+## object of type &#39;closure&#39; is not subsettable</code></pre>
</div>
<div id="results" class="section level2">
<h2>Results</h2>
@@ -1858,6 +1896,14 @@ models fitted to two datasets, i.e. eight fits for each test.</p>
<td align="right">1.308</td>
<td align="right">1.793</td>
</tr>
+<tr class="odd">
+<td align="left">Linux</td>
+<td align="left">Ryzen 9 7950X 16-Core Processor</td>
+<td align="left">4.2.2</td>
+<td align="left">1.3.0</td>
+<td align="right">1.271</td>
+<td align="right">1.787</td>
+</tr>
</tbody>
</table>
</div>
@@ -2068,6 +2114,15 @@ for each test.</p>
<td align="right">2.364</td>
<td align="right">1.230</td>
</tr>
+<tr class="odd">
+<td align="left">Linux</td>
+<td align="left">Ryzen 9 7950X 16-Core Processor</td>
+<td align="left">4.2.2</td>
+<td align="left">1.3.0</td>
+<td align="right">0.526</td>
+<td align="right">0.623</td>
+<td align="right">0.631</td>
+</tr>
</tbody>
</table>
</div>
@@ -2344,6 +2399,18 @@ dataset, i.e. one fit for each test.</p>
<td align="right">0.801</td>
<td align="right">1.093</td>
</tr>
+<tr class="odd">
+<td align="left">Linux</td>
+<td align="left">Ryzen 9 7950X 16-Core Processor</td>
+<td align="left">4.2.2</td>
+<td align="left">1.3.0</td>
+<td align="right">0.242</td>
+<td align="right">0.237</td>
+<td align="right">0.238</td>
+<td align="right">0.239</td>
+<td align="right">0.237</td>
+<td align="right">0.237</td>
+</tr>
</tbody>
</table>
</div>
diff --git a/vignettes/web_only/mkin_benchmarks.rda b/vignettes/web_only/mkin_benchmarks.rda
index 64ac2680..93138541 100644
--- a/vignettes/web_only/mkin_benchmarks.rda
+++ b/vignettes/web_only/mkin_benchmarks.rda
Binary files differ

Contact - Imprint