How to solve this ode equations in R? -
i solve following ordinary differential equations in r runge kutta method.
dy1 <- (-51.33) * ((1-y[2]) / y[1]) dy2 <- 1.54 * y[1] * (1-y[2]) - 2.14 * y[2]
when y1
becomes zero, dy1
become infinity. avoid this, need write r code stating when y[1]
becomes less 0.001, stop y[1]
derivative , keeps zero. pasted r code below:
yini <- c(1,0) intabs <- function (t, y, parms) { ifelse (y[1] <= 0.01, dy1 <- 0, no) dy1 <- -51.33 * ((1-y[2]) / y[1]) dy2 <- 1.54 * y[1] * (1-y[2]) - 2.14 * y[2] list(c(dy1, dy2)) } times <- seq(from = 0, = 1, = .002) out <- ode (times = times, y=yini, func = intabs, parms = null, method = "rk4") head (out, n=50)
i used ifelse
statement denote y[1]
less or equal 0.001, keep dy1
zero. i'm getting results. seems made error, brings results of dy1
in negative values. i'm new in writing programmes. please if spot error..
such jump in ode function (if correctly encoded) poison step size controller. ode function has @ least smooth error order of integration method step size control work intended. still can fail if ode stiff.
a less invasive method de-singularize expression 1/y
use
y/(1e-6+y*y)
which y>1e-2
give approximation of 1/y
, smooth everywhere. change constant adapt region deviates 1/y
.
Comments
Post a Comment