Smitom Borah 2/21/2021
This example for the ODIN package is taken from the README markdown file of the ODIN package. The example is based on the Lorenz system which provides a simple mathematical model of atmospheric advection.
The equations are given below:
here, (\sigma), R and (\beta) are three parameter and y1, y2 and y3 are the three coordinates respectively.
The ODIN package essentially has three parts:
lorenz <- odin::odin({
## Derivatives
deriv(y1) <- sigma * (y2 - y1)
deriv(y2) <- R * y1 - y2 - y1 * y3
deriv(y3) <- -b * y3 + y1 * y2
## Initial conditions
initial(y1) <- 10.0
initial(y2) <- 1.0
initial(y3) <- 1.0
## parameters
sigma <- 10.0
R <- 28.0
b <- 8.0 / 3.0
})
mod <- lorenz()
t <- seq(0, 100, length.out = 50000)
y <- mod$run(t)
# This generates an object that can be used to integrate the set of differential equations, by default starting at the initial conditions specified above (though custom initial conditions can be given). The equations are translated into C, compiled, loaded, and bundled into an object. lorenz here is a function that generates an instance of the model.
# to view the last 10 rows of the y deSolve
tail(y, n=10L)## t y1 y2 y3
## [49991,] 99.982 -1.2005969 1.172437 23.85095
## [49992,] 99.984 -1.1537287 1.160185 23.72135
## [49993,] 99.986 -1.1080299 1.148055 23.59256
## [49994,] 99.988 -1.0634749 1.136063 23.46460
## [49995,] 99.990 -1.0200382 1.124222 23.33743
## [49996,] 99.992 -0.9776944 1.112545 23.21107
## [49997,] 99.994 -0.9364186 1.101047 23.08549
## [49998,] 99.996 -0.8961860 1.089736 22.96069
## [49999,] 99.998 -0.8569720 1.078626 22.83666
## [50000,] 100.000 -0.8187525 1.067725 22.71340
So the values of y1, y2 and y3 will be -0.82, 1.07 and 22.71, respectively, at t=50000.
Sometimes, instead of giving the values for the initial condition explicitly, an user() is used. When used as, say sigma<-user(10.0), it means we allow a user parameter called sigma to be specified, but if omitted it has a default value of 1. When used as user(), it means that there is no default value and sigma must be provided.
For more details, go through the following link: