1 #ifndef STAN_MATH_PRIM_ARR_FUNCTOR_INTEGRATE_ODE_HPP
2 #define STAN_MATH_PRIM_ARR_FUNCTOR_INTEGRATE_ODE_HPP
12 #include <boost/numeric/odeint.hpp>
58 template <
typename F,
typename T1,
typename T2>
59 std::vector<std::vector<typename stan::return_type<T1, T2>::type> >
61 const std::vector<T1> y0,
63 const std::vector<double>& ts,
64 const std::vector<T2>& theta,
65 const std::vector<double>& x,
66 const std::vector<int>& x_int,
68 using boost::numeric::odeint::integrate_times;
69 using boost::numeric::odeint::make_dense_output;
70 using boost::numeric::odeint::runge_kutta_dopri5;
83 const double absolute_tolerance = 1
e-6;
84 const double relative_tolerance = 1
e-6;
85 const double step_size = 0.1;
89 coupled_system(f, y0, theta, x, x_int, msgs);
92 std::vector<double> ts_vec(ts.size() + 1);
94 for (
size_t n = 0; n < ts.size(); n++)
97 std::vector<std::vector<double> > y_coupled(ts_vec.size());
101 std::vector<double> initial_coupled_state
102 = coupled_system.initial_state();
104 integrate_times(make_dense_output(absolute_tolerance,
106 runge_kutta_dopri5<std::vector<double>,
111 initial_coupled_state,
112 boost::begin(ts_vec), boost::end(ts_vec),
117 y_coupled.erase(y_coupled.begin());
120 return coupled_system.decouple_states(y_coupled);
bool check_less(const char *function, const char *name, const T_y &y, const T_high &high)
Return true if y is strictly less than high.
bool check_nonzero_size(const char *function, const char *name, const T_y &y)
Return true if the specified matrix/vector is of non-zero size.
Observer for the coupled states.
bool check_ordered(const char *function, const char *name, const Eigen::Matrix< T_y, Eigen::Dynamic, 1 > &y)
Return true if the specified vector is sorted into strictly increasing order.
double e()
Return the base of the natural logarithm.
bool check_finite(const char *function, const char *name, const T_y &y)
Return true if y is finite.
Base template class for a coupled ordinary differential equation system, which adds sensitivities to ...
std::vector< std::vector< typename stan::return_type< T1, T2 >::type > > integrate_ode(const F &f, const std::vector< T1 > y0, const double t0, const std::vector< double > &ts, const std::vector< T2 > &theta, const std::vector< double > &x, const std::vector< int > &x_int, std::ostream *msgs)
Return the solutions for the specified system of ordinary differential equations given the specified ...