library(rODE)
library(ggplot2)
setClass("EulerError", slots = c(
stack = "environment" # environment object inside the class
),
contains = c("ODE")
)
setMethod("initialize", "EulerError", function(.Object, ...) {
.Object@stack$n <- 0 # "n" belongs to the class environment
.Object@state <- vector("numeric", 1)
return(.Object)
})
setMethod("getExactSolution", "EulerError", function(object, t, ...) {
# analytical solution
return(exp(t) - t - 1)
})
setMethod("getState", "EulerError", function(object, ...) {
object@state
})
setMethod("getRate", "EulerError", function(object, state, ...) {
object@rate[1] <- state[1] + state[2] # x + y
object@rate[2] <- 1 # dx/dx
object@stack$n <- object@stack$n + 1 # add 1 to the rate count
object@rate
})
# constructor
EulerError <- function(y) {
.EulerError <- new("EulerError")
.EulerError@state[1] = y # y
.EulerError@state[2] = 0 # x = t
return(.EulerError)
}
# class implementation
EulerErrorApp <- function(stepSize = 0.1) {
initial_y <- 0
xmax <- 1
stepSize <- stepSize
n_steps <- as.integer((xmax + stepSize / 2) / stepSize)
ode <- EulerError(initial_y)
ode_solver <- Euler(ode)
ode_solver <- setStepSize(ode_solver, stepSize)
steps <- 0
rowVector <- vector("list")
i <- 1
while (steps < n_steps) {
ode_solver <- step(ode_solver)
state <- getState(ode_solver@ode)
steps <- ode_solver@ode@stack$n
rowVector[[i]] <- list(
x = state[2], # x = t
y = state[1], # y
TrueY = getExactSolution(ode_solver@ode, state[2]),
steps = steps)
i <- i + 1
}
data.table::rbindlist(rowVector)
}
