1 #ifndef STAN_MATH_REV_MAT_FUNCTOR_ALGEBRA_SYSTEM_HPP 2 #define STAN_MATH_REV_MAT_FUNCTOR_ALGEBRA_SYSTEM_HPP 24 template <
typename F,
typename T0,
typename T1,
bool x_is_iv>
29 Eigen::Matrix<T0, Eigen::Dynamic, 1>
x_;
31 Eigen::Matrix<T1, Eigen::Dynamic, 1>
y_;
42 const Eigen::Matrix<T1, Eigen::Dynamic, 1>& y,
43 const std::vector<double>& dat,
44 const std::vector<int>& dat_int, std::ostream* msgs)
45 : f_(f), x_(x), y_(y), dat_(dat), dat_int_(dat_int), msgs_(msgs) {}
57 const Eigen::Matrix<T, Eigen::Dynamic, 1>& iv)
const {
59 return f_(iv, y_, dat_, dat_int_, msgs_);
61 return f_(x_, iv, dat_, dat_int_, msgs_);
72 template <
typename T,
int NX = Eigen::Dynamic,
int NY = Eigen::Dynamic>
78 nlo_functor(
int inputs,
int values) : m_inputs(inputs), m_values(values) {}
80 int inputs()
const {
return m_inputs; }
81 int values()
const {
return m_values; }
93 template <
typename S,
typename F,
typename T0,
typename T1>
103 const Eigen::Matrix<T0, Eigen::Dynamic, 1>& x,
104 const Eigen::Matrix<T1, Eigen::Dynamic, 1>& y,
105 const std::vector<double>& dat,
106 const std::vector<int>& dat_int, std::ostream* msgs)
107 : fs_(f, x, y, dat, dat_int, msgs), x_size_(x.
size()) {}
116 int operator()(
const Eigen::VectorXd& iv, Eigen::VectorXd& fvec) {
127 int df(
const Eigen::VectorXd& iv, Eigen::MatrixXd& fjac)
const {
139 Eigen::VectorXd fvec;
150 Eigen::VectorXd
get_value(
const Eigen::VectorXd& iv)
const {
return fs_(iv); }
std::vector< double > dat_
real data
Eigen::Matrix< T0, Eigen::Dynamic, 1 > x_
unknowns
int df(const Eigen::VectorXd &iv, Eigen::MatrixXd &fjac) const
Assign the Jacobian to fjac (signature required by Eigen).
int operator()(const Eigen::VectorXd &iv, Eigen::VectorXd &fvec)
Computes the value the algebraic function, f, when pluging in the independent variables, and the Jacobian w.r.t unknowns.
A functor that allows us to treat either x or y as the independent variable.
Eigen::VectorXd get_value(const Eigen::VectorXd &iv) const
Performs the same task as df(), but returns the value of algebraic function, instead of saving it ins...
Eigen::MatrixXd J_
Jacobian of algebraic function wrt unknowns.
void jacobian(const F &f, const Eigen::Matrix< T, Eigen::Dynamic, 1 > &x, Eigen::Matrix< T, Eigen::Dynamic, 1 > &fx, Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &J)
nlo_functor(int inputs, int values)
A functor with the required operators to call Eigen's algebraic solver.
F f_
algebraic system functor
hybrj_functor_solver(const S &fs, const F &f, const Eigen::Matrix< T0, Eigen::Dynamic, 1 > &x, const Eigen::Matrix< T1, Eigen::Dynamic, 1 > &y, const std::vector< double > &dat, const std::vector< int > &dat_int, std::ostream *msgs)
system_functor(const F &f, const Eigen::Matrix< T0, Eigen::Dynamic, 1 > &x, const Eigen::Matrix< T1, Eigen::Dynamic, 1 > &y, const std::vector< double > &dat, const std::vector< int > &dat_int, std::ostream *msgs)
int x_size_
number of unknowns
int size(const std::vector< T > &x)
Return the size of the specified standard vector.
Eigen::MatrixXd get_jacobian(const Eigen::VectorXd &iv)
Performs the same task as the operator(), but returns the Jacobian, instead of saving it inside an ar...
S fs_
Wrapper around algebraic system.
Eigen::Matrix< T, Eigen::Dynamic, 1 > operator()(const Eigen::Matrix< T, Eigen::Dynamic, 1 > &iv) const
An operator that takes in an independent variable.
Eigen::Matrix< T1, Eigen::Dynamic, 1 > y_
auxiliary parameters
std::ostream * msgs_
stream message
std::vector< int > dat_int_
integer data
A structure which gets passed to Eigen's dogleg algebraic solver.