[R-sig-dyn-mod] how to select the integrator in ode?

[Jinsong Zhao]

Thanks a lot for the suggestions and the deSolve package.

The help package of "forcings" solve my question. And it also solve the 
question I encountered in your book "A practical guide to ecological 
modelling":

outtimes <- seq(1,365,by=1)
Carbon <- meanDepo*optpar$par[2]/optpar$par[1]
names(Carbon) <-"Carbon"
out <- as.data.frame(ode(Carbon,outtimes,minmod, optpar$par))
names(out) <- c("time","Carbon","minrate")

when I print out, I found the last row of the out is NA.

 > tail(out)
    time    Carbon   minrate
360  360  99.51982 0.9179177
361  361  99.80277 0.9205274
362  362 100.14148 0.9236516
363  363 100.53546 0.9272854
364  364 100.98419 0.9314243
365  365        NA        NA

I learned that "euler" is fixed time step, and I am familiar with the 
method because it always used as an demo in most introduction materials 
of numerical analysis. The reason that I asked this question is 
sometimes the "euler" method could be replaced by "lsoda", and produce 
same results, and sometimes they can not be used interchangeably.

Thanks again for the explanation and suggestions.

Best,
Jinsong


On 2020/8/3 3:09, Thomas Petzoldt wrote:
> Hi,
> 
> thanks for considering deSolve. Before answering your question, I want 
> to emphasize that an ordinary differential equation system (ODE) is 
> deterministic by definition. Never use a stochastic component like 
> runif() inside of the model function of an ODE, especially not when 
> using an automatic step-size solver like "lsoda". If you want to 
> simulate stochastic forcing, use the "forcing" mechanism instead, see 
> the ?forcings help page.
> 
> So first of all: remove the runif() component from the model!
> 
> The question which solver to use is best understood by reading a book 
> about differential equations or even better by attending a course in 
> numerical mathematics. Methods "euler" and the default solver "lsoda" 
> are fundamentally different, so (besides your random element in the 
> model) differences are completely natural.
> 
> In short: with a relatively large time step, "euler" solves a difference 
> equation using a given fixed time step while "lsoda" approximates a 
> differential equation with variable internal time step.
> 
> ThPe
> 
> 
> 
> Am 02.08.2020 um 16:59 schrieb Jinsong Zhao:
>> Hi there,
>>
>> I try to solve a simple ODE model:
>>
>> logistic <- function(time, state, params) {
>>    with(as.list(c(state, params)), {
>>            K <- K.mean + runif(1, -10, 10)
>>            dPOP <- rate * POP * (1 - POP / K)
>>            list(c(dPOP),
>>                 c(K = K))
>>    })
>> }
>>
>> params <- c(rate = 0.1,
>>             K.mean = 100)
>> state <- c(POP = 10)
>> time <- seq(0,120)
>> out <- ode(state, time, logistic, params, method = "euler")
>> plot(out)
>>
>> However, if I remove `method = "euler"`, i.e., using the default 
>> integrator, the output is very different.
>>
>> I don't know when I should use "euler" or any other integrator. Any 
>> suggestion? Thanks a lot.
>>
>> Best,
>> Jinsong
>>
>> _______________________________________________
>> R-sig-dynamic-models mailing list
>> R-sig-dynamic-models using r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-d
>