Automatic Differentiation
 
Loading...
Searching...
No Matches

◆ ode_adjoint_impl() [1/2]

template<typename F , typename T_y0 , typename T_t0 , typename T_ts , typename T_abs_tol_fwd , typename T_abs_tol_bwd , typename... T_Args, require_all_eigen_col_vector_t< T_y0, T_abs_tol_fwd, T_abs_tol_bwd > * = nullptr, require_any_not_st_arithmetic< T_y0, T_t0, T_ts, T_Args... > * = nullptr>
auto stan::math::ode_adjoint_impl ( const char *  function_name,
F &&  f,
const T_y0 &  y0,
const T_t0 &  t0,
const std::vector< T_ts > &  ts,
double  relative_tolerance_forward,
const T_abs_tol_fwd &  absolute_tolerance_forward,
double  relative_tolerance_backward,
const T_abs_tol_bwd &  absolute_tolerance_backward,
double  relative_tolerance_quadrature,
double  absolute_tolerance_quadrature,
long int  max_num_steps,
long int  num_steps_between_checkpoints,
int  interpolation_polynomial,
int  solver_forward,
int  solver_backward,
std::ostream *  msgs,
const T_Args &...  args 
)

Solve the ODE initial value problem y' = f(t, y), y(t0) = y0 at a set of times, { t1, t2, t3, ... } using the stiff backward differentiation formula BDF solver or the non-stiff Adams solver from CVODES.

The ODE system is integrated using the adjoint sensitivity approach of CVODES.

f must define an operator() with the signature as: template<typename T_t, typename T_y, typename... T_Args> Eigen::Matrix<stan::return_type_t<T_t, T_y, T_Args...>, Eigen::Dynamic, 1> operator()(const T_t& t, const Eigen::Matrix<T_y, Eigen::Dynamic, 1>& y, std::ostream* msgs, const T_Args&... args);

t is the time, y is the state, msgs is a stream for error messages, and args are optional arguments passed to the ODE solve function (which are passed through to f without modification).

Template Parameters
FType of ODE right hand side
T_y0Type of initial state
T_t0Type of initial time
T_tsType of output times
T_ArgsTypes of pass-through parameters
Parameters
function_nameCalling function name (for printing debugging messages)
fRight hand side of the ODE
y0Initial state
t0Initial time
tsTimes at which to solve the ODE at. All values must be sorted and not less than t0.
relative_tolerance_forwardRelative tolerance for forward problem passed to CVODES
absolute_tolerance_forwardAbsolute tolerance per ODE state for forward problem passed to CVODES
relative_tolerance_backwardRelative tolerance for backward problem passed to CVODES
absolute_tolerance_backwardAbsolute tolerance per ODE state for backward problem passed to CVODES
relative_tolerance_quadratureRelative tolerance for quadrature problem passed to CVODES
absolute_tolerance_quadratureAbsolute tolerance for quadrature problem passed to CVODES
max_num_stepsUpper limit on the number of integration steps to take between each output (error if exceeded)
num_steps_between_checkpointsNumber of integrator steps after which a checkpoint is stored for the backward pass
interpolation_polynomialtype of polynomial used for interpolation
solver_forwardsolver used for forward pass
solver_backwardsolver used for backward pass
[in,out]msgsthe print stream for warning messages
argsExtra arguments passed unmodified through to ODE right hand side
Returns
An std::vector of Eigen column vectors with scalars equal to the least upper bound of T_y0, T_t0, T_ts, and the lambda's arguments. This represents the solution to ODE at times ts

Definition at line 71 of file ode_adjoint.hpp.