diff options
| author | Johannes Ranke <jranke@uni-bremen.de> | 2022-11-16 09:15:36 +0100 |
|---|---|---|
| committer | Johannes Ranke <jranke@uni-bremen.de> | 2022-11-16 09:15:36 +0100 |
| commit | 51d63256a7b3020ee11931d61b4db97b9ded02c0 (patch) | |
| tree | cb6d628211c99cb6dd1938428a18ef4dd5a997dc | |
| parent | 679cf716192cdfd91dfd232578cbd4e30d7eac12 (diff) | |
We get about 25% performance gain
with the custom lsoda call, avoiding repeated getNativeSymbolInfo calls.
It's just that we should not be calling foreign functions from different
packages, because the may change without notice. Using
getNativeSymbolInfo for "call_lsoda" avoids the CRAN note, and a similar
call could probably be used for "unlock_solver", avoiding the NOTE in
checks for cran, but we should not do this in a CRAN package.
| -rw-r--r-- | DESCRIPTION | 8 | ||||
| -rw-r--r-- | NAMESPACE | 1 | ||||
| -rw-r--r-- | R/mkinfit.R | 17 | ||||
| -rw-r--r-- | R/mkinmod.R | 4 | ||||
| -rw-r--r-- | R/mkinpredict.R | 107 | ||||
| -rw-r--r-- | R/nlme.mmkin.R | 14 | ||||
| -rw-r--r-- | R/saem.R | 22 | ||||
| -rw-r--r-- | log/build.log | 2 | ||||
| -rw-r--r-- | log/check.log | 23 | ||||
| -rw-r--r-- | log/test.log | 34 | ||||
| -rw-r--r-- | man/mkinmod.Rd | 4 | ||||
| -rw-r--r-- | man/mkinpredict.Rd | 21 | ||||
| -rw-r--r-- | man/saem.Rd | 10 | ||||
| -rw-r--r-- | tests/testthat/test_dmta.R | 3 | ||||
| -rw-r--r-- | vignettes/mkin.html | 2 | ||||
| -rw-r--r-- | vignettes/web_only/benchmarks.html | 91 | ||||
| -rw-r--r-- | vignettes/web_only/mkin_benchmarks.rda | bin | 1586 -> 1641 bytes | |||
| -rw-r--r-- | vignettes/web_only/saem_benchmarks.html | 59 | ||||
| -rw-r--r-- | vignettes/web_only/saem_benchmarks.rda | bin | 359 -> 471 bytes |
19 files changed, 328 insertions, 94 deletions
diff --git a/DESCRIPTION b/DESCRIPTION index bb9ae3a5..a006a446 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,8 +1,8 @@ Package: mkin Type: Package Title: Kinetic Evaluation of Chemical Degradation Data -Version: 1.2.0 -Date: 2022-11-01 +Version: 1.3.0 +Date: 2022-11-15 Authors@R: c( person("Johannes", "Ranke", role = c("aut", "cre", "cph"), email = "johannes.ranke@jrwb.de", @@ -22,8 +22,8 @@ Description: Calculation routines based on the FOCUS Kinetics Report (2006, Ranke et al. (2021) <doi:10.3390/environments8080071>. Please 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), +Depends: R (>= 2.15.1), deSolve +Imports: stats, graphics, methods, parallel, 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 @@ -147,7 +147,6 @@ export(status) export(tex_listing) export(transform_odeparms) export(which.best) -import(deSolve) import(graphics) import(nlme) importFrom(R6,R6Class) diff --git a/R/mkinfit.R b/R/mkinfit.R index 693778fd..0d9246dd 100644 --- a/R/mkinfit.R +++ b/R/mkinfit.R @@ -500,6 +500,15 @@ 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 + } + # Get the error model and the algorithm for fitting err_mod <- match.arg(error_model) error_model_algorithm = match.arg(error_model_algorithm) @@ -610,7 +619,8 @@ mkinfit <- function(mkinmod, observed, solution_type = solution_type, use_compiled = use_compiled, method.ode = method.ode, - atol = atol, rtol = rtol, ...) + atol = atol, rtol = rtol, + call_lsoda = call_lsoda, ...) observed_index <- cbind(as.character(observed$time), as.character(observed$name)) observed$predicted <- out[observed_index] @@ -892,7 +902,10 @@ mkinfit <- function(mkinmod, observed, fit$calls <- calls fit$time <- fit_time - # We also need the model and a model name for summary and plotting + # 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 fit$mkinmod <- mkinmod fit$mkinmod$name <- mkinmod_name fit$obs_vars <- obs_vars diff --git a/R/mkinmod.R b/R/mkinmod.R index d8740aed..b1fb57cb 100644 --- a/R/mkinmod.R +++ b/R/mkinmod.R @@ -45,7 +45,9 @@ #' @param dll_dir Directory where an DLL object, if generated internally by #' [inline::cfunction()], should be saved. The DLL will only be stored in a #' permanent location for use in future sessions, if 'dll_dir' and 'name' -#' are specified. +#' are specified. This is helpful if fit objects are cached e.g. by knitr, +#' as the cache remains functional across sessions if the DLL is stored in +#' a user defined location. #' @param unload If a DLL from the target location in 'dll_dir' is already #' loaded, should that be unloaded first? #' @param overwrite If a file exists at the target DLL location in 'dll_dir', diff --git a/R/mkinpredict.R b/R/mkinpredict.R index 0dc9cf55..11a3b35b 100644 --- a/R/mkinpredict.R +++ b/R/mkinpredict.R @@ -19,26 +19,24 @@ #' @param solution_type The method that should be used for producing the #' predictions. This should generally be "analytical" if there is only one #' observed variable, and usually "deSolve" in the case of several observed -#' variables. The third possibility "eigen" is faster but not applicable to -#' some models e.g. using FOMC for the parent compound. +#' variables. The third possibility "eigen" is fast in comparison to uncompiled +#' ODE models, but not applicable to some models, e.g. using FOMC for the +#' parent compound. #' @param method.ode The solution method passed via [mkinpredict] to [ode]] in -#' case the solution type is "deSolve". The default "lsoda" is performant, but -#' sometimes fails to converge. +#' case the solution type is "deSolve" and we are not using compiled code. #' @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 [ode]. Default is 1e-8, -#' lower than in [lsoda]. -#' @param rtol Absolute error tolerance, passed to [ode]. Default is 1e-10, -#' much lower than in [lsoda]. -#' @param maxsteps Maximum number of steps, passed to [ode]. +#' @param atol Absolute error tolerance, passed to the ode solver. +#' @param rtol Absolute error tolerance, passed to the ode solver. +#' @param maxsteps Maximum number of steps, passed to the ode solver. #' @param map_output Boolean to specify if the output should list values for #' the observed variables (default) or for all state variables (if set to #' 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. -#' @import deSolve #' @return A matrix with the numeric solution in wide format #' @author Johannes Ranke #' @examples @@ -116,9 +114,10 @@ mkinpredict.mkinmod <- function(x, outtimes = seq(0, 120, by = 0.1), solution_type = "deSolve", use_compiled = "auto", - method.ode = "lsoda", atol = 1e-8, rtol = 1e-10, maxsteps = 20000, + method.ode = "lsoda", atol = 1e-8, rtol = 1e-10, maxsteps = 20000L, map_output = TRUE, na_stop = TRUE, + call_lsoda = NULL, ...) { @@ -173,20 +172,80 @@ 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 = if (is.null(x$dll_info)) inline::getDynLib(x$cf)[["name"]] - else x$dll_info[["name"]], - parms = odeparms[x$parms], # Order matters when using compiled models - method = method.ode, - atol = atol, - rtol = rtol, - maxsteps = maxsteps, - ... +# 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 <- t(out_raw) + colnames(out) <- c("time", mod_vars) } else { mkindiff <- function(t, state, parms) { diff --git a/R/nlme.mmkin.R b/R/nlme.mmkin.R index 09cb84b8..e193e5e3 100644 --- a/R/nlme.mmkin.R +++ b/R/nlme.mmkin.R @@ -186,10 +186,24 @@ nlme.mmkin <- function(model, data = "auto", thisCall[["control"]] <- control } + # Provide the address of call_lsoda to the fitting function + call_lsoda <- getNativeSymbolInfo("call_lsoda", PACKAGE = "deSolve") + if (model[[1]]$solution_type == "deSolve" & !is.null(model[[1]]$mkinmod$cf)) { + # The mkinmod stored in the first fit will be used by nlme + model[[1]]$mkinmod$diffs_address <- getNativeSymbolInfo("diffs", + PACKAGE = model[[1]]$mkinmod$dll_info[["name"]])$address + model[[1]]$mkinmod$initpar_address <- getNativeSymbolInfo("initpar", + PACKAGE = model[[1]]$mkinmod$dll_info[["name"]])$address + } + fit_time <- system.time(val <- do.call("nlme.formula", thisCall)) val$time <- fit_time val$mkinmod <- model[[1]]$mkinmod + # Don't return addresses that will become invalid + val$mkinmod$diffs_address <- NULL + val$mkinmod$initpar_address <- NULL + val$data <- thisCall[["data"]] val$mmkin <- model if (is.list(start)) val$mean_dp_start <- start$fixed @@ -120,12 +120,12 @@ utils::globalVariables(c("predicted", "std")) #' summary(f_saem_dfop_sfo, data = TRUE) #' #' # The following takes about 6 minutes -#' #f_saem_dfop_sfo_deSolve <- saem(f_mmkin["DFOP-SFO", ], solution_type = "deSolve", -#' # control = list(nbiter.saemix = c(200, 80), nbdisplay = 10)) +#' f_saem_dfop_sfo_deSolve <- saem(f_mmkin["DFOP-SFO", ], solution_type = "deSolve", +#' nbiter.saemix = c(200, 80)) #' -#' #saemix::compare.saemix(list( -#' # f_saem_dfop_sfo$so, -#' # f_saem_dfop_sfo_deSolve$so)) +#' #anova( +#' # f_saem_dfop_sfo, +#' # f_saem_dfop_sfo_deSolve)) #' #' # If the model supports it, we can also use eigenvalue based solutions, which #' # take a similar amount of time @@ -580,6 +580,15 @@ saemix_model <- function(object, solution_type = "auto", transform_rates <- object[[1]]$transform_rates transform_fractions <- object[[1]]$transform_fractions + # Get native symbol info for speed + call_lsoda <- getNativeSymbolInfo("call_lsoda", PACKAGE = "deSolve") + if (solution_type == "deSolve" & !is.null(mkin_model$cf)) { + mkin_model$diffs_address <- getNativeSymbolInfo("diffs", + PACKAGE = mkin_model$dll_info[["name"]])$address + mkin_model$initpar_address <- getNativeSymbolInfo("initpar", + PACKAGE = mkin_model$dll_info[["name"]])$address + } + # Define the model function model_function <- function(psi, id, xidep) { @@ -613,7 +622,8 @@ saemix_model <- function(object, solution_type = "auto", odeparms = odeparms, odeini = odeini, solution_type = solution_type, outtimes = sort(unique(i_time)), - na_stop = FALSE + na_stop = FALSE, + call_lsoda = call_lsoda ) out_index <- cbind(as.character(i_time), as.character(i_name)) diff --git a/log/build.log b/log/build.log index c4f9b8a2..245bd205 100644 --- a/log/build.log +++ b/log/build.log @@ -5,5 +5,5 @@ * creating vignettes ... OK * checking for LF line-endings in source and make files and shell scripts * checking for empty or unneeded directories -* building ‘mkin_1.2.0.tar.gz’ +* building ‘mkin_1.3.0.tar.gz’ diff --git a/log/check.log b/log/check.log index 3fea2ec6..e7d4d327 100644 --- a/log/check.log +++ b/log/check.log @@ -5,7 +5,7 @@ * using options ‘--no-tests --as-cran’ * checking for file ‘mkin/DESCRIPTION’ ... OK * checking extension type ... Package -* this is package ‘mkin’ version ‘1.2.0’ +* this is package ‘mkin’ version ‘1.3.0’ * package encoding: UTF-8 * checking CRAN incoming feasibility ... Note_to_CRAN_maintainers Maintainer: ‘Johannes Ranke <johannes.ranke@jrwb.de>’ @@ -37,11 +37,18 @@ Maintainer: ‘Johannes Ranke <johannes.ranke@jrwb.de>’ * checking whether the namespace can be unloaded cleanly ... OK * checking loading without being on the library search path ... OK * checking use of S3 registration ... OK -* checking dependencies in R code ... OK +* checking dependencies in R code ... NOTE +Package in Depends field not imported from: ‘deSolve’ + These packages need to be imported from (in the NAMESPACE file) + for when this namespace is loaded but not attached. * checking S3 generic/method consistency ... OK * checking replacement functions ... OK -* checking foreign function calls ... OK -* checking R code for possible problems ... [14s/14s] OK +* checking foreign function calls ... NOTE +Foreign function call to a different package: + .C("unlock_solver", ..., PACKAGE = "deSolve") +See chapter ‘System and foreign language interfaces’ in the ‘Writing R +Extensions’ manual. +* checking R code for possible problems ... [17s/17s] OK * checking Rd files ... OK * checking Rd metadata ... OK * checking Rd line widths ... OK @@ -57,7 +64,7 @@ Maintainer: ‘Johannes Ranke <johannes.ranke@jrwb.de>’ * checking data for ASCII and uncompressed saves ... OK * checking installed files from ‘inst/doc’ ... OK * checking files in ‘vignettes’ ... OK -* checking examples ... [15s/15s] OK +* checking examples ... [20s/20s] OK * checking for unstated dependencies in ‘tests’ ... OK * checking tests ... SKIPPED * checking for unstated dependencies in vignettes ... OK @@ -69,5 +76,9 @@ Maintainer: ‘Johannes Ranke <johannes.ranke@jrwb.de>’ * checking for detritus in the temp directory ... OK * DONE -Status: OK +Status: 2 NOTEs +See + ‘/home/jranke/git/mkin/mkin.Rcheck/00check.log’ +for details. + diff --git a/log/test.log b/log/test.log index 7005ac37..10c0aa76 100644 --- a/log/test.log +++ b/log/test.log @@ -1,22 +1,22 @@ ℹ Testing mkin ✔ | F W S OK | Context ✔ | 5 | AIC calculation -✔ | 5 | Analytical solutions for coupled models [3.2s] +✔ | 5 | Analytical solutions for coupled models [2.9s] ✔ | 5 | Calculation of Akaike weights ✔ | 3 | Export dataset for reading into CAKE ✔ | 12 | Confidence intervals and p-values [1.0s] -✔ | 1 12 | Dimethenamid data from 2018 [32.4s] +✔ | 1 12 | Dimethenamid data from 2018 [29.5s] ──────────────────────────────────────────────────────────────────────────────── -Skip ('test_dmta.R:98'): Different backends get consistent results for SFO-SFO3+, dimethenamid data +Skip ('test_dmta.R:99'): Different backends get consistent results for SFO-SFO3+, dimethenamid data Reason: Fitting this ODE model with saemix takes about 15 minutes on my system ──────────────────────────────────────────────────────────────────────────────── ✔ | 14 | Error model fitting [4.9s] ✔ | 5 | Time step normalisation ✔ | 4 | Calculation of FOCUS chi2 error levels [0.6s] -✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.8s] -✔ | 4 | Test fitting the decline of metabolites from their maximum [0.4s] +✔ | 14 | Results for FOCUS D established in expertise for UBA (Ranke 2014) [0.7s] +✔ | 4 | Test fitting the decline of metabolites from their maximum [0.3s] ✔ | 1 | Fitting the logistic model [0.2s] -✔ | 10 | Batch fitting and diagnosing hierarchical kinetic models [24.2s] +✔ | 10 | Batch fitting and diagnosing hierarchical kinetic models [24.7s] ✔ | 1 12 | Nonlinear mixed-effects models [0.3s] ──────────────────────────────────────────────────────────────────────────────── Skip ('test_mixed.R:74'): saemix results are reproducible for biphasic fits @@ -26,28 +26,28 @@ Reason: Fitting with saemix takes around 10 minutes when using deSolve ✔ | 10 | Special cases of mkinfit calls [0.5s] ✔ | 3 | mkinfit features [0.7s] ✔ | 8 | mkinmod model generation and printing [0.2s] -✔ | 3 | Model predictions with mkinpredict [0.3s] -✔ | 7 | Multistart method for saem.mmkin models [36.9s] -✔ | 16 | Evaluations according to 2015 NAFTA guidance [2.3s] -✔ | 9 | Nonlinear mixed-effects models with nlme [8.6s] -✔ | 16 | Plotting [9.8s] +✔ | 3 | Model predictions with mkinpredict [0.4s] +✔ | 7 | Multistart method for saem.mmkin models [36.8s] +✔ | 16 | Evaluations according to 2015 NAFTA guidance [2.5s] +✔ | 9 | Nonlinear mixed-effects models with nlme [8.8s] +✔ | 16 | Plotting [10.0s] ✔ | 4 | Residuals extracted from mkinfit models -✔ | 1 36 | saemix parent models [66.3s] +✔ | 1 36 | saemix parent models [66.1s] ──────────────────────────────────────────────────────────────────────────────── Skip ('test_saemix_parent.R:143'): We can also use mkin solution methods for saem Reason: This still takes almost 2.5 minutes although we do not solve ODEs ──────────────────────────────────────────────────────────────────────────────── -✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [1.4s] +✔ | 2 | Complex test case from Schaefer et al. (2007) Piacenza paper [1.1s] ✔ | 11 | Processing of residue series -✔ | 10 | Fitting the SFORB model [3.7s] +✔ | 10 | Fitting the SFORB model [3.4s] ✔ | 1 | Summaries of old mkinfit objects ✔ | 5 | Summary [0.2s] -✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [2.2s] -✔ | 9 | Hypothesis tests [8.0s] +✔ | 4 | Results for synthetic data established in expertise for UBA (Ranke 2014) [1.6s] +✔ | 9 | Hypothesis tests [6.1s] ✔ | 4 | Calculation of maximum time weighted average concentrations (TWAs) [2.2s] ══ Results ═════════════════════════════════════════════════════════════════════ -Duration: 211.7 s +Duration: 206.2 s ── Skipped tests ────────────────────────────────────────────────────────────── • Fitting this ODE model with saemix takes about 15 minutes on my system (1) diff --git a/man/mkinmod.Rd b/man/mkinmod.Rd index 87ce9016..612c3c2b 100644 --- a/man/mkinmod.Rd +++ b/man/mkinmod.Rd @@ -58,7 +58,9 @@ applicable to give detailed information about the C function being built.} \item{dll_dir}{Directory where an DLL object, if generated internally by \code{\link[inline:cfunction]{inline::cfunction()}}, should be saved. The DLL will only be stored in a permanent location for use in future sessions, if 'dll_dir' and 'name' -are specified.} +are specified. This is helpful if fit objects are cached e.g. by knitr, +as the cache remains functional across sessions if the DLL is stored in +a user defined location.} \item{unload}{If a DLL from the target location in 'dll_dir' is already loaded, should that be unloaded first?} diff --git a/man/mkinpredict.Rd b/man/mkinpredict.Rd index 0797f259..d93c0753 100644 --- a/man/mkinpredict.Rd +++ b/man/mkinpredict.Rd @@ -18,9 +18,10 @@ mkinpredict(x, odeparms, odeini, outtimes, ...) method.ode = "lsoda", atol = 1e-08, rtol = 1e-10, - maxsteps = 20000, + maxsteps = 20000L, map_output = TRUE, na_stop = TRUE, + call_lsoda = NULL, ... ) @@ -60,23 +61,21 @@ solver is used.} \item{solution_type}{The method that should be used for producing the predictions. This should generally be "analytical" if there is only one observed variable, and usually "deSolve" in the case of several observed -variables. The third possibility "eigen" is faster but not applicable to -some models e.g. using FOMC for the parent compound.} +variables. The third possibility "eigen" is fast in comparison to uncompiled +ODE models, but not applicable to some models, e.g. using FOMC for the +parent compound.} \item{use_compiled}{If set to \code{FALSE}, no compiled version of the \link{mkinmod} model is used, even if is present.} \item{method.ode}{The solution method passed via \link{mkinpredict} to \link{ode}] in -case the solution type is "deSolve". The default "lsoda" is performant, but -sometimes fails to converge.} +case the solution type is "deSolve" and we are not using compiled code.} -\item{atol}{Absolute error tolerance, passed to \link{ode}. Default is 1e-8, -lower than in \link{lsoda}.} +\item{atol}{Absolute error tolerance, passed to the ode solver.} -\item{rtol}{Absolute error tolerance, passed to \link{ode}. Default is 1e-10, -much lower than in \link{lsoda}.} +\item{rtol}{Absolute error tolerance, passed to the ode solver.} -\item{maxsteps}{Maximum number of steps, passed to \link{ode}.} +\item{maxsteps}{Maximum number of steps, passed to the ode solver.} \item{map_output}{Boolean to specify if the output should list values for the observed variables (default) or for all state variables (if set to @@ -84,6 +83,8 @@ FALSE). Setting this to FALSE has no effect for analytical solutions, as these always return mapped output.} \item{na_stop}{Should it be an error if \link{ode} returns NaN values} + +\item{call_lsoda}{The address of the compiled function "call_lsoda"} } \value{ A matrix with the numeric solution in wide format diff --git a/man/saem.Rd b/man/saem.Rd index 11463351..3a5abada 100644 --- a/man/saem.Rd +++ b/man/saem.Rd @@ -206,12 +206,12 @@ plot(f_saem_dfop_sfo) summary(f_saem_dfop_sfo, data = TRUE) # The following takes about 6 minutes -#f_saem_dfop_sfo_deSolve <- saem(f_mmkin["DFOP-SFO", ], solution_type = "deSolve", -# control = list(nbiter.saemix = c(200, 80), nbdisplay = 10)) +f_saem_dfop_sfo_deSolve <- saem(f_mmkin["DFOP-SFO", ], solution_type = "deSolve", + nbiter.saemix = c(200, 80)) -#saemix::compare.saemix(list( -# f_saem_dfop_sfo$so, -# f_saem_dfop_sfo_deSolve$so)) +#anova( +# f_saem_dfop_sfo, +# f_saem_dfop_sfo_deSolve)) # If the model supports it, we can also use eigenvalue based solutions, which # take a similar amount of time diff --git a/tests/testthat/test_dmta.R b/tests/testthat/test_dmta.R index 30c5d7c4..5cfc61d2 100644 --- a/tests/testthat/test_dmta.R +++ b/tests/testthat/test_dmta.R @@ -85,7 +85,8 @@ sfo_sfo3p <- mkinmod( ) dmta_sfo_sfo3p_tc <- mmkin(list("SFO-SFO3+" = sfo_sfo3p), - dmta_ds, error_model = "tc", quiet = TRUE, cores = n_cores) + dmta_ds, error_model = "tc", quiet = TRUE, + cores = n_cores) test_that("Different backends get consistent results for SFO-SFO3+, dimethenamid data", { diff --git a/vignettes/mkin.html b/vignettes/mkin.html index 0d5ed6fc..38c44a0f 100644 --- a/vignettes/mkin.html +++ b/vignettes/mkin.html @@ -1599,7 +1599,7 @@ div.tocify { <h1 class="title toc-ignore">Introduction to mkin</h1> <h4 class="author">Johannes Ranke</h4> -<h4 class="date">Last change 15 February 2021 (rebuilt 2022-07-12)</h4> +<h4 class="date">Last change 15 February 2021 (rebuilt 2022-11-15)</h4> </div> diff --git a/vignettes/web_only/benchmarks.html b/vignettes/web_only/benchmarks.html index 5f81f39f..5376c1f5 100644 --- a/vignettes/web_only/benchmarks.html +++ b/vignettes/web_only/benchmarks.html @@ -1599,7 +1599,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 2022-07-14)</h4> +<h4 class="date">Last change 14 July 2022 (rebuilt 2022-11-15)</h4> </div> @@ -1635,7 +1635,7 @@ FOMC_SFO <- mkinmod( parent = mkinsub("FOMC", "m1"), m1 = mkinsub("SFO")) DFOP_SFO <- mkinmod( - parent = mkinsub("FOMC", "m1"), + parent = mkinsub("FOMC", "m1"), # erroneously used FOMC twice, not fixed for consistency m1 = mkinsub("SFO")) t3 <- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D)))[["elapsed"]] t4 <- system.time(mmkin_bench(list(SFO_SFO, FOMC_SFO, DFOP_SFO), list(FOCUS_D), @@ -1816,6 +1816,30 @@ t11 <- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c), <td align="right">1.770</td> <td align="right">3.377</td> </tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.1</td> +<td align="left">1.1.2</td> +<td align="right">1.957</td> +<td align="right">3.633</td> +</tr> +<tr class="even"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.2.0</td> +<td align="right">2.129</td> +<td align="right">3.784</td> +</tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.3.0</td> +<td align="right">2.046</td> +<td align="right">3.693</td> +</tr> </tbody> </table> </div> @@ -1979,6 +2003,33 @@ t11 <- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c), <td align="right">5.758</td> <td align="right">2.558</td> </tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.1</td> +<td align="left">1.1.2</td> +<td align="right">1.503</td> +<td align="right">6.147</td> +<td align="right">2.803</td> +</tr> +<tr class="even"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.2.0</td> +<td align="right">1.559</td> +<td align="right">6.097</td> +<td align="right">2.841</td> +</tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.3.0</td> +<td align="right">1.230</td> +<td align="right">4.333</td> +<td align="right">2.187</td> +</tr> </tbody> </table> </div> @@ -2193,6 +2244,42 @@ t11 <- system.time(mmkin_bench(list(m_synth_DFOP_par), list(DFOP_par_c), <td align="right">1.744</td> <td align="right">2.566</td> </tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.1</td> +<td align="left">1.1.2</td> +<td align="right">0.861</td> +<td align="right">1.295</td> +<td align="right">1.507</td> +<td align="right">3.102</td> +<td align="right">1.961</td> +<td align="right">2.852</td> +</tr> +<tr class="even"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.2.0</td> +<td align="right">0.911</td> +<td align="right">1.328</td> +<td align="right">1.519</td> +<td align="right">2.986</td> +<td align="right">1.957</td> +<td align="right">2.769</td> +</tr> +<tr class="odd"> +<td align="left">Linux</td> +<td align="left">Ryzen 7 1700</td> +<td align="left">4.2.2</td> +<td align="left">1.3.0</td> +<td align="right">0.693</td> +<td align="right">0.996</td> +<td align="right">1.121</td> +<td align="right">2.174</td> +<td align="right">1.427</td> +<td align="right">2.026</td> +</tr> </tbody> </table> </div> diff --git a/vignettes/web_only/mkin_benchmarks.rda b/vignettes/web_only/mkin_benchmarks.rda Binary files differindex 2d3deb26..c07def65 100644 --- a/vignettes/web_only/mkin_benchmarks.rda +++ b/vignettes/web_only/mkin_benchmarks.rda diff --git a/vignettes/web_only/saem_benchmarks.html b/vignettes/web_only/saem_benchmarks.html index 4875bb1b..714dc1ff 100644 --- a/vignettes/web_only/saem_benchmarks.html +++ b/vignettes/web_only/saem_benchmarks.html @@ -1599,7 +1599,7 @@ div.tocify { <h1 class="title toc-ignore">Benchmark timings for saem.mmkin</h1> <h4 class="author">Johannes Ranke</h4> -<h4 class="date">Last change 14 November 2022 (rebuilt 2022-11-14)</h4> +<h4 class="date">Last change 14 November 2022 (rebuilt 2022-11-15)</h4> </div> @@ -1781,10 +1781,20 @@ t11 <- system.time(sforb_sfo3_plus_const <- saem(three_met_sep_tc["SF <td align="left">Linux</td> <td align="left">1.2.0</td> <td align="left">3.2</td> -<td align="right">2.996</td> -<td align="right">5.207</td> -<td align="right">5.317</td> -<td align="right">5.171</td> +<td align="right">2.110</td> +<td align="right">4.632</td> +<td align="right">4.264</td> +<td align="right">4.930</td> +</tr> +<tr class="even"> +<td align="left">Ryzen 7 1700</td> +<td align="left">Linux</td> +<td align="left">1.3.0</td> +<td align="left">3.2</td> +<td align="right">2.394</td> +<td align="right">4.748</td> +<td align="right">4.883</td> +<td align="right">4.937</td> </tr> </tbody> </table> @@ -1808,10 +1818,20 @@ t11 <- system.time(sforb_sfo3_plus_const <- saem(three_met_sep_tc["SF <td align="left">Linux</td> <td align="left">1.2.0</td> <td align="left">3.2</td> -<td align="right">5.671</td> -<td align="right">7.696</td> -<td align="right">8.166</td> -<td align="right">8.168</td> +<td align="right">5.602</td> +<td align="right">7.373</td> +<td align="right">7.815</td> +<td align="right">7.831</td> +</tr> +<tr class="even"> +<td align="left">Ryzen 7 1700</td> +<td align="left">Linux</td> +<td align="left">1.3.0</td> +<td align="left">3.2</td> +<td align="right">5.622</td> +<td align="right">7.445</td> +<td align="right">8.297</td> +<td align="right">7.740</td> </tr> </tbody> </table> @@ -1836,8 +1856,16 @@ t11 <- system.time(sforb_sfo3_plus_const <- saem(three_met_sep_tc["SF <td align="left">Linux</td> <td align="left">1.2.0</td> <td align="left">3.2</td> -<td align="right">24.883</td> -<td align="right">818.157</td> +<td align="right">24.014</td> +<td align="right">749.699</td> +</tr> +<tr class="even"> +<td align="left">Ryzen 7 1700</td> +<td align="left">Linux</td> +<td align="left">1.3.0</td> +<td align="left">3.2</td> +<td align="right">24.480</td> +<td align="right">519.087</td> </tr> </tbody> </table> @@ -1861,7 +1889,14 @@ t11 <- system.time(sforb_sfo3_plus_const <- saem(three_met_sep_tc["SF <td align="left">Linux</td> <td align="left">1.2.0</td> <td align="left">3.2</td> -<td align="right">1355.036</td> +<td align="right">1249.834</td> +</tr> +<tr class="even"> +<td align="left">Ryzen 7 1700</td> +<td align="left">Linux</td> +<td align="left">1.3.0</td> +<td align="left">3.2</td> +<td align="right">944.471</td> </tr> </tbody> </table> diff --git a/vignettes/web_only/saem_benchmarks.rda b/vignettes/web_only/saem_benchmarks.rda Binary files differindex e8b139da..992b58cf 100644 --- a/vignettes/web_only/saem_benchmarks.rda +++ b/vignettes/web_only/saem_benchmarks.rda |