diff options
Diffstat (limited to 'custom_lsoda_call_edited.patch')
| -rw-r--r-- | custom_lsoda_call_edited.patch | 116 |
1 files changed, 116 insertions, 0 deletions
diff --git a/custom_lsoda_call_edited.patch b/custom_lsoda_call_edited.patch new file mode 100644 index 00000000..a79cbbcd --- /dev/null +++ b/custom_lsoda_call_edited.patch @@ -0,0 +1,116 @@ +--- a/R/mkinfit.R ++++ b/R/mkinfit.R + ++ # 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 ++ } ++ + +@@ -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, ...) + +--- a/R/mkinpredict.R ++++ b/R/mkinpredict.R + +@@ -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, +- ... ++ # 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) { + |