Benaim: Model 'Benaim'

In CharlotteJana/pdmppoly: Moment Approximation for Polynomial PDMPs

Description

This Model was first introduced by Benaim et al in 2012 (example 1.3). It is an interesting example for the importance of the switching rates.

Usage

Benaim

polyBenaim

Format

Benaim is an object of class pdmpModel,
polyBenaim is an object of class polyPdmpModel.

Simulation

The simulations in PROM were done with slot times set to

• from = 0, to = 1000, by = 0.1.

The following parameter sets were simulated:

• b = 2, beta = 1.6

• b = 2, beta = 1.4

• b = 2, beta = 1.3

• b = 2, beta = 0.8

• b = 2, beta = 0.5

• b = 2, beta = 0.3

Source

The model is introduced in [BenaimCo2012a] as example 1.3.

Examples

library(spray)
#------ code to generate the pdmpModel version -----

Benaim <- new("pdmpModel",
    descr = "Modell B",
    parms = list(b = 2, β = 1.4),
    init = c(f1 = 10, f2 = 10, d = 1), 
    discStates = list(d = 0:1),
    dynfunc = function(t, x, parms) {
      df <- with(as.list(c(x, parms)), 
                 c(-f1, -f2) + switch(d+1, c(2*b*f2, 0), 
                                       c(0, 2*b*f1)))
      return(c(df, 0))
    }, 
    ratefunc = function(t, x, parms) {
      return(with(as.list(c(x, parms)), β/2))
    }, 
    jumpfunc = function(t, x, parms, jtype) {
      c(x[1:2], 1 - x[3])
    }, 
    times = c(from = 0, to = 100, by = 0.1),
    solver = "lsodar")

#------ code to generate the polyPdmpModel version -----

polyBenaim <- new("polyPdmpModel",
    descr = "Modell B",
    parms = list(b = 2, β = 1.4),
    init = c(f1 = 10, f2 = 10, d = 1), 
    discStates = list(d = 0:1),
    dynpolys = quote(list(
      list(overall = -lone(1,3), specific = list(2*b*lone(2,3), 0)),
      list(overall = -lone(2,3), specific = list(0, 2*b*lone(1,3)))
    )), 
    ratepolys = quote(list(
      list(β/2, β/2)
    )), 
    jumpfunc = function(t, x, parms, jtype) {
      c(x[1:2], 1 - x[3])
    }, 
    times = c(from = 0, to = 100, by = 0.1),
    solver = "lsodar")

#------- comparison of the models --------------

identical(sim(Benaim, outSlot = FALSE, seed = 10),
          sim(polyBenaim, outSlot = FALSE, seed = 10))

CharlotteJana/pdmppoly documentation built on Sept. 4, 2019, 4:40 p.m.