## 19.1 Basic functions

Here’s an example of a skeletal Stan program with a user-defined relative difference function employed in the generated quantities block to compute a relative differences between two parameters.

```
functions {
real relative_diff(real x, real y) {
real abs_diff;
real avg_scale;
abs_diff = fabs(x - y);
avg_scale = (fabs(x) + fabs(y)) / 2;
return abs_diff / avg_scale;
}
}
// ...
generated quantities {
real rdiff;
rdiff = relative_diff(alpha, beta);
}
```

The function is named `relative_diff`

, and is declared to have
two real-valued arguments and return a real-valued result. It is
used the same way a built-in function would be used in the generated
quantities block.

### User-defined functions block

All functions are defined in their own block, which is labeled
`functions`

and must appear before all other program blocks. The
user-defined functions block is optional.

### Function bodies

The body (the part between the curly braces) contains ordinary Stan code, including local variables. The new function is used in the generated quantities block just as any of Stan’s built-in functions would be used.

### Return statements

Return statements, such as the one on the last line of the definition
of `relative_diff`

above, are only allowed in the bodies of
function definitions. Return statements may appear anywhere in a
function, but functions with non-void return types must end in a
return statement.

### Reject statements

The Stan reject statement provides a mechanism to report errors or problematic values encountered during program execution. It accepts any number of quoted string literals or Stan expressions as arguments. This statement is typically embedded in a conditional statement in order to detect bad or illegal outcomes of some processing step.

#### Catching errors

Rejection is used to flag errors that arise in inputs or in program state. It is far better to fail early with a localized informative error message than to run into problems much further downstream (as in rejecting a state or failing to compute a derivative).

The most common errors that are coded is to test that all of the arguments to a function are legal. The following function takes a square root of its input, so requires non-negative inputs; it is coded to guard against illegal inputs.

```
real dbl_sqrt(real x) {
if (!(x >= 0)) {
reject("dblsqrt(x): x must be positive; found x = ", x);
}
return 2 * sqrt(x);
}
```

The negation of the positive test is important, because it also
catches the case where `x`

is a not-a-number value. If the
condition had been coded as `(x < 0)`

it would not catch the
not-a-number case, though it could be written as
`(x < 0 || is_nan(x))`

. The positive infinite case is allowed
through, but could also be checked with the `is_inf(x)`

function.
The square root function does not itself reject, but some downstream
consumer of `dbl_sqrt(-2)`

would be likely to raise an error, at
which point the origin of the illegal input requires detective work.
Or even worse, as Matt Simpson pointed out in the GitHub comments, the
function could go into an infinite loop if it starts with an infinite
value and tries to reduce it by arithmetic, likely consuming all
available memory and crashing an interface. Much better to catch
errors early and report on their origin.

The effect of rejection depends on the program block in which the
rejection is executed. In transformed data, rejections cause the
program to fail to load. In transformed parameters or in the model
block, rejections cause the current state to be rejected in the
Metropolis sense.^{28}

In generated quantities, rejections cause execution to halt because there is no way to recover and generate the remaining parameters, so extra care should be taken in calling functions in the generated quantities block.

### Type declarations for functions

Function argument and return types for vector and matrix types are not declared with their sizes, unlike type declarations for variables. Function argument type declarations may not be declared with constraints, either lower or upper bounds or structured constraints like forming a simplex or correlation matrix, (as is also the case for local variables); see the table of types in the reference manual for full details.

For example, here’s a function to compute the entropy of a categorical
distribution with simplex parameter `theta`

.

```
real entropy(vector theta) {
return sum(theta .* log(theta));
}
```

Although `theta`

must be a simplex, only the type `vector`

is used.^{29}

Upper or lower bounds on values or constrained types are not allowed as return types or argument types in function declarations.

### Array types for function declarations

Array arguments have their own syntax, which follows that used in this manual for function signatures. For example, a function that operates on a two-dimensional array to produce a one-dimensional array might be declared as follows.

`array[] real baz(array[,] real x);`

The notation `[ ]`

is used for one-dimensional arrays (as in the
return above), `[ , ]`

for two-dimensional arrays,
`[ , , ]`

for three-dimensional arrays, and so on.

Functions support arrays of any type, including matrix and vector types. As with other types, no constraints are allowed.

### Data-only function arguments

A function argument which is a real-valued type or
a container of a real-valued type,
i.e., not an integer type or integer array type,
can be qualified using the prefix qualifier `data`

.
The following is an example of a data-only function argument.

```
real foo(real y, data real mu) {
return -0.5 * (y - mu)^2;
}
```

This qualifier restricts this argument to being invoked with expressions which consist only of data variables, transformed data variables, literals, and function calls. A data-only function argument cannot involve real variables declared in the parameters, transformed parameters, or model block. Attempts to invoke a function using an expression which contains parameter, transformed parameters, or model block variables as a data-only argument will result in an error message from the parser.

Use of the `data`

qualifier must be consistent between the
forward declaration and the definition of a functions.

This qualifier should be used when writing functions that call the
built-in ordinary differential equation (ODE) solvers, algebraic
solvers, or map functions. These higher-order functions have strictly
specified signatures where some arguments of are data only
expressions. (See the ODE solver chapter for
more usage details and the functions reference manual for full
definitions.) When writing a function which calls the ODE or
algebraic solver, arguments to that function which are passed into the
call to the solver, either directly or indirectly, should have the
`data`

prefix qualifier. This allows for compile-time type checking
and increases overall program understandability.

Just because this makes it possible to code a rejection sampler does not make it a good idea. Rejections break differentiability and the smooth exploration of the posterior. In Hamiltonian Monte Carlo, it can cause the sampler to be reduced to a diffusive random walk.↩︎

A range of built-in validation routines is coming to Stan soon! Alternatively, the

`reject`

statement can be used to check constraints on the simplex.↩︎