Skip to content
Back to formatted view

Raw Message

Message-ID: <6eca8baa-f016-5b75-9ea5-dfcf12583b99@tu-dresden.de>
Date: 2020-08-02T19:09:22Z
From: Thomas Petzoldt
Subject: [R-sig-dyn-mod] how to select the integrator in ode?
In-Reply-To: <d2bcd710-ca6c-fdaa-a82f-a1abde4cbfed@yeah.net>

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 at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-d