In preparation for a course I came across a model that required that the points where the ODE is evaluated be dependent on parameters, necessitating their gradients. Currently we require that the ODE time points be data, but it’s straightforward to promote them to parameters.

The corresponding Jacobian terms are simply

dz_{i} (t) / dt |*{t = T} = dz*{i}/dt(t = T, theta, z(t = T)).

In other words, the derivative of the states with respect to a given time point is just ODE RHS evaluated at that time point.

Operationally this would require

- promoting the time input from
`vector<double>`

to`vector<T>`

- using meta programs to grab the value of the times where they are currently used in the code
- passing the vector of ODE RHS values at each time point (n_times * n_states values) into
`decouple_states`

so that they can be included in the`precomputed_gradients`

call

Does this all sounds good? Charles and I would like to get this in within the next week so that we can use it as a demonstration in the aforementioned course so prompt comments are appreciated.