## Defining the ODE System

The routines solve the general n-dimensional first-order system,

    dy_i(t)/dt = f_i(t, y_1(t), ..., y_n(t))

for i = 1, \dots, n. The stepping functions rely on the vector of derivatives f_i and the Jacobian matrix, J_{ij} = df_i(t,y(t)) / dy_j. A system of equations is defined using the `gsl_odeiv_system` datatype.

Data Type: gsl_odeiv_system
  This data type defines a general ODE system with arbitrary parameters.
  `int (* function) (double t, const double y[], double dydt[], void * params)`
    This function should store the vector elements f_i(t,y,params) in the array dydt, for arguments (t,y) and parameters params
  `int (* jacobian) (double t, const double y[], double * dfdy, double dfdt[], void * params);`
    This function should store the vector elements df_i(t,y,params)/dt in the array dfdt and the Jacobian matrix J_{ij} in the array dfdy regarded as a row-ordered matrix `J(i,j) = dfdy[i * dim + j]` where `dim` is the dimension of the system. Some of the simpler solver algorithms do not make use of the Jacobian matrix, so it is not always strictly necessary to provide it (this element of the struct can be replace by a null pointer). However, it is useful to provide the Jacobian to allow the solver algorithms to be interchanged -- the best algorithms make use of the Jacobian.
  `size_t dimension;`
    This is the dimension of the system of equations
  `void * params`
    This is a pointer to the arbitrary parameters of the system.