12.2 Coding an algebraic system

A system of algebraic equations is coded directly in Stan as a function with a strictly specified signature. For example, the nonlinear system given above can be coded using the following function in Stan (see the user-defined functions section for more information on coding user-defined functions).

vector system(vector y,        // unknowns
              vector theta,    // parameters
              data real[] x_r,      // data (real)
              int[] x_i) {     // data (integer)
  vector[2] z;
  z[1] = y[1] - theta[1];
  z[2] = y[1] * y[2] - theta[2];
  return z;

The function takes the unknowns we wish to solve for in y (a vector), the system parameters in theta (a vector), the real data in x_r (a real array) and the integer data in x_i (an integer array). The system function returns the value of the function (a vector), for which we want to compute the roots. Our example does not use real or integer data. Nevertheless, these unused arguments must be included in the system function with exactly the signature above.

The body of the system function here could also be coded using a row vector constructor and transposition,

return [ y[1] - theta[1],
         y[1] * y[2] - theta[2] ]';

As systems get more complicated, naming the intermediate expressions goes a long way toward readability.

Strict signature

The function defining the system must have exactly these argument types and return type. This may require passing in zero-length arrays for data or a zero-length vector for parameters if the system does not involve data or parameters.