Prev Next Index-> contents reference index search external Up-> CppAD Appendix deprecated cppad_ipopt_nlp ipopt_nlp_ode ipopt_nlp_ode_problem ipopt_nlp_ode_problem.hpp cppad_ipopt_nlp-> ipopt_nlp_get_started.cpp ipopt_nlp_ode ipopt_ode_speed.cpp ipopt_nlp_ode-> ipopt_nlp_ode_problem ipopt_nlp_ode_simple ipopt_nlp_ode_fast ipopt_nlp_ode_run.hpp ipopt_nlp_ode_check.cpp ipopt_nlp_ode_problem-> ipopt_nlp_ode_problem.hpp ipopt_nlp_ode_problem.hpp Headings

$\newcommand{\W}[1]{ \; #1 \; } \newcommand{\R}[1]{ {\rm #1} } \newcommand{\B}[1]{ {\bf #1} } \newcommand{\D}[2]{ \frac{\partial #1}{\partial #2} } \newcommand{\DD}[3]{ \frac{\partial^2 #1}{\partial #2 \partial #3} } \newcommand{\Dpow}[2]{ \frac{\partial^{#1}}{\partial {#2}^{#1}} } \newcommand{\dpow}[2]{ \frac{ {\rm d}^{#1}}{{\rm d}\, {#2}^{#1}} }$
ODE Inverse Problem Definitions: Source Code
# include "../src/cppad_ipopt_nlp.hpp" namespace { //------------------------------------------------------------------ typedef Ipopt::Number Number; Number a0 = 1.; // simulation value for a[0] Number a1 = 2.; // simulation value for a[1] Number a2 = 1.; // simulatioln value for a[2] // function used to simulate data Number y_one(Number t) { Number y_1 = a0*a1 * (exp(-a2*t) - exp(-a1*t)) / (a1 - a2); return y_1; } // time points were we have data (no data at first point) double s[] = { 0.0, 0.5, 1.0, 1.5, 2.0 }; // Simulated data for case with no noise (first point is not used) double z[] = { 0.0, y_one(0.5), y_one(1.0), y_one(1.5), y_one(2.0) }; // Number of measurement values size_t Nz = sizeof(z) / sizeof(z[0]) - 1; // Number of components in the function y(t, a) size_t Ny = 2; // Number of components in the vectro a size_t Na = 3; // Initial Condition function, F(a) = y(t, a) at t = 0 // (for this particular example) template <class Vector> Vector eval_F(const Vector &a) { Vector F(Ny); // y_0 (t) = a[0]*exp(-a[1] * t) F[0] = a[0]; // y_1 (t) = // a[0]*a[1]*(exp(-a[2] * t) - exp(-a[1] * t))/(a[1] - a[2]) F[1] = 0.; return F; } // G(y, a) = \partial_t y(t, a); i.e. the differential equation // (for this particular example) template <class Vector> Vector eval_G(const Vector &y , const Vector &a) { Vector G(Ny); // y_0 (t) = a[0]*exp(-a[1] * t) G[0] = -a[1] * y[0]; // y_1 (t) = // a[0]*a[1]*(exp(-a[2] * t) - exp(-a[1] * t))/(a[1] - a[2]) G[1] = +a[1] * y[0] - a[2] * y[1]; return G; } // H(i, y, a) = contribution to objective at i-th data point // (for this particular example) template <class Scalar, class Vector> Scalar eval_H(size_t i, const Vector &y, const Vector &a) { // This particular H is for a case where y_1 (t) is measured Scalar diff = z[i] - y[1]; return diff * diff; } // function used to count the number of calls to eval_r size_t count_eval_r(void) { static size_t count = 0; ++count; return count; } }