---
title: "odin functions"
author: "Rich FitzJohn"
date: "`r Sys.Date()`"
output:
  rmarkdown::html_vignette:
    toc: true
vignette: >
  %\VignetteIndexEntry{odin functions}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

``` {r include = FALSE}
path <- system.file("functions.json", package = "odin", mustWork = TRUE)
functions <- jsonlite::fromJSON(path, FALSE)
vcapply <- odin:::vcapply

render <- function(x) {
  render1 <- function(el) {
    if (!is.null(el$description)) {
      ret <- sprintf("**%s**: %s", el$title, el$description)
    } else {
      ret <- sprintf("**%s**", el$title)
    }
    if (!is.null(el$examples)) {
      examples <- vcapply(el$examples, function(ex) {
        sprintf("`%s` &rarr; `%s`", ex[[1]], ex[[2]])
      })
      ret <- sprintf("%s (e.g., %s)", ret, paste(examples, collapse = "; "))
    }
    ret
  }
  value <- vapply(x, render1, "", USE.NAMES = FALSE)
  paste(sprintf("* `%s` -- %s\n", names(x), value), collapse = "\n")
}
```

`odin` supports many functions that you'd expect to see for
constructing differential equation models; primarily mathematical
functions available through R's "Rmath" library.  These include all
mathematical operations, and many more obscure mathematical
functions.  Special support is provided for working with arrays.  A
further set of functions is available for working discrete time
stochastic models.

## Basic operators

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$basic))
```

Because general programming is not supported in `odin` and because
every line must contain an assignment, instead of writing

```r
if (mycondition) {
  a <- true_value
} else {
  a <- false_value
}
```

instead write

```r
a <- if (mycondition) true_value else false_value
```

(this works in normal R too!)

## Array support

There are a group of functions for interacting with arrays

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$arrays))
```

When working with arrays, use generally implies a "for loop" in the
generated C code.  For example, in the example in [the main package
vignette](odin.html) the derivatives are computed as

```r
deriv(y[]) <- r[i] * y[i] * (1 - sum(ay[i, ]))
```

The indexes on the right hand side can be one of `i`, `j`, `k`, `l`
`i5`, `i6`, `i7` or `i8` corresponding to the index on the *left hand
side* being iterated over (`odin` supports arrays up to 8 dimensions).
The left-hand-side here contains no explicit entry (`y[]`) which is
equivalent to `y[1:length(y)]`, which expands (approximately) to the
"for loop"

```r
for (i in 1:length(y)) {
  deriv(y[i]) <- r[i] * y[i] * (1 - sum(ay[i, ]))
}
```

(except slightly different, and in C).

Similarly, the expression

```r
ay[, ] <- a[i, j] * y[j]
```

involves loops over two dimensions (`ay[, ]` becomes `ay[1:dim(ay,
1), 1:dim(ay, 2)]` and so the loop becomes

```r
for (i in 1:dim(ay, 1)) {
  for (j in 1:dim(ay, 2)) {
    ay[i, j] <- a[i, j] * y[j]
  }
}
```

Due to constraints with using C, few things can be used as an index; in particular the following will not work:

```r
idx <- if (t > 5) 2 else 1
x <- vec[idx]
```

(or where `idx` is some general odin variable as the result of a different assignment). You must use `as.integer` to cast this to integer immediately before indexing:

```r
idx <- if (t > 5) 2 else 1
x <- vec[as.integer(idx)]
```

This will *truncate* the value (same behaviour as `truncate`) so be warned if passing in things that may be approximately integer - you may want to use `as.integer(round(x))` in that case.

The interpolation functions are described in more detail in the
[main package vignette](odin.html)

## Operators

A number of logical-returning operators exist, primarily to support
the `if` statement; all the usual comparison operators exist
(though not vectorised `|` or `&`).

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$operators))
```

Be wary of strict equality with `==` or `!=` as numbers may be
floating point numbers, which have some surprising properties for
the uninitiated, for example
``` {r }
sqrt(3)^2 == 3
```

## Mathematical operators

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$maths))
```

The exact for `%%` and `%/%` for floating point numbers and signed
numbers are complicated - please see `?Arithmetic`.  The rules for
operators in `odin` are exactly those in R as the same underlying
functions are used.

Similarly, for the differences between `round`, `floor`, `ceiling`
and `truncate`, see the help page `?round`.  Note that R's
behaviour for rounding away from 0.5 is exactly followed and that
this slightly changed behaviour at version 4.0.0

All the usual trig functions are also available:

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$trig))
```

## Stochastic models

For discrete time stochastic models, all of R's normal stochastic
distribution functions are available:

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$stochastic))
```

With random number functions we can write:

```r
x <- runif(10, 20)
```

which will generate a random number from the uniform distribution.  If you write:

```r
x[] <- runif(10, 20)
```

then each element of `x` will be filled with a *different* random number drawn from this distribution (which is generally what you want). Random numbers are considered to be *time varying* which means they will automatically generate each time step, so if you write

```r
x <- rnorm(0, 10)
update(y[]) <- y[i] + x
```

then at each time step, each element of `y` will be updated by the same random number from a normal distribution with a mean of zero and a standard deviation of 10 - the number will change each time step but be the same for each element of `y` in the example above.

In addition, two functions that are vector returning and require
some care to use:

``` {r echo = FALSE, results = "asis"}
writeLines(render(functions$inplace))
```

Both these functions require a vector input (of probabilities for `rmultinom` and of counts for `rmhyper`) and return a vector the same length.  So the expression

```r
y[] <- rmultinom(10, p)
```

will produce a vector `y` of samples from the multinomial distribution with parameters `size = 10` (so after wards `sum(y)` is 10) and probabilities `p`. It is very important that `y` and `p` have the same size.

At the moment it is not possible to use expressions like

```r
y[1, ] <- rmultinom(10, p[i, ])
```

but this is planned for implementation in the future.  A full example of using `rmultinom` is given in the [discrete models](discrete.html) vignette.