odeintr (version 1.7.0)

compile_sys_openmp: Compile ODE system with openmp multi-threading

Description

Generates an integrator using Rcpp and openmp

Usage

compile_sys_openmp(name, sys, pars = NULL, const = FALSE,
  method = "rk5_i", sys_dim = -1L, atol = 1e-06, rtol = 1e-06,
  globals = "", headers = "", footers = "", compile = TRUE,
  env = new.env(), ...)

Arguments

name

the name of the generated integration function

sys

a string containing C++ expressions

pars

a named vector of numbers or a vector of names or number of parameters

const

declare parameters const if true

method

a method string (see compile_sys)

sys_dim

length of the state vector

atol

absolute tolerance if using adaptive step size

rtol

relative tolerance if using adaptive step size

globals

a string with global C++ declarations

headers

code to appear before the odeintr namespace

footers

code to appear after the odeintr namespace

compile

if false, just return the code

env

install functions into this environment

...

passed to sourceCpp

Details

This functions behaves identically to compile_sys execpt that it does not allow one to override the default observer. In order to take advantage of openmp multi-threading, you must insert openmp pragmas into your system definition. See the examples. A special function laplace4 is defined and can be called from your system definition. It will compute a discrete 4-point Laplacian for use in solving PDE via the method of lines. The function takes x as its first argument, dxdt as its second argument and the diffusion coefficient D as its third parameter. This function uses the default openmp scheduling.

See Also

set_optimization, compile_sys

Examples

## Not run: ------------------------------------
# M = 200
# bistable = '
#  laplace4(x, dxdt, D);  // parallel 4-point discrete laplacian
#  #pragma omp parallel for
#  for (int i = 0; i < N; ++i)
#    dxdt[i] += a * x[i] * (1 - x[i]) * (x[i] - b);
# ' # bistable
# compile_sys_openmp("bistable", bistable, sys_dim = M * M,
#                    pars = c(D = 0.1, a = 1.0, b = 1/2),
#                    const = TRUE)
# at = 10 ^ (0:3)
# inic = rbinom(M * M, 1, 1/2)
# system.time({x = bistable_at(inic, at)})
# par(mfrow = rep(2, 2), mar = rep(1, 4), oma = rep(1, 4))
# for (i in 1:4){
#   image(matrix(unlist(x[i, -1]), M, M),
#         asp = 1, col = c("black", "lightgray"),
#         axes = FALSE)
#   title(main=paste("Time =", x[i, 1]))}
## ---------------------------------------------