Prev

Next

Index-> contents reference index search external

Up-> CppAD AD ADValued atomic atomic_base atomic_get_started.cpp

ADValued-> Arithmetic unary_standard_math binary_math CondExp Discrete numeric_limits atomic

atomic-> checkpoint atomic_base

atomic_base-> atomic_ctor atomic_option atomic_afun atomic_forward atomic_reverse atomic_for_sparse_jac atomic_rev_sparse_jac atomic_for_sparse_hes atomic_rev_sparse_hes atomic_base_clear atomic_get_started.cpp atomic_norm_sq.cpp atomic_reciprocal.cpp atomic_set_sparsity.cpp atomic_tangent.cpp atomic_eigen_mat_mul.cpp atomic_eigen_mat_inv.cpp atomic_eigen_cholesky.cpp atomic_mat_mul.cpp

atomic_get_started.cpp

Headings-> Purpose Start Class Definition Constructor forward End Class Definition Use Atomic Function ---..Constructor ---..Recording ---..forward

# include <cppad/cppad.hpp>
namespace {          // isolate items below to this file
using CppAD::vector; // abbreviate as vector
class atomic_get_started : public CppAD::atomic_base<double> {

public:
     // constructor (could use const char* for name)
     atomic_get_started(const std::string& name) :
     // this example does not use any sparsity patterns
     CppAD::atomic_base<double>(name)
     { }
private:

     // forward mode routine called by CppAD
     virtual bool forward(
          size_t                    p ,
          size_t                    q ,
          const vector<bool>&      vx ,
                vector<bool>&      vy ,
          const vector<double>&    tx ,
                vector<double>&    ty
     )
     {
# ifndef NDEBUG
          size_t n = tx.size() / (q + 1);
          size_t m = ty.size() / (q + 1);
# endif
          assert( n == 1 );
          assert( m == 1 );

          // return flag
          bool ok = q == 0;
          if( ! ok )
               return ok;

          // check for defining variable information
          // This case must always be implemented
          if( vx.size() > 0 )
               vy[0] = vx[0];

          // Order zero forward mode.
          // This case must always be implemented
          // y^0 = f( x^0 ) = 1 / x^0
          double f = 1. / tx[0];
          if( p <= 0 )
               ty[0] = f;
          return ok;
     }

}; // End of atomic_get_started class
}  // End empty namespace

bool get_started(void)
{     bool ok = true;
     using CppAD::AD;
     using CppAD::NearEqual;
     double eps = 10. * CppAD::numeric_limits<double>::epsilon();

     // Create the atomic get_started object
     atomic_get_started afun("atomic_get_started");

     // Create the function f(x)
     //
     // domain space vector
     size_t n  = 1;
     double  x0 = 0.5;
     vector< AD<double> > ax(n);
     ax[0]     = x0;

     // declare independent variables and start tape recording
     CppAD::Independent(ax);

     // range space vector
     size_t m = 1;
     vector< AD<double> > ay(m);

     // call user function and store get_started(x) in au[0]
     vector< AD<double> > au(m);
     afun(ax, au);        // u = 1 / x

     // now use AD division to invert to invert the operation
     ay[0] = 1.0 / au[0]; // y = 1 / u = x

     // create f: x -> y and stop tape recording
     CppAD::ADFun<double> f;
     f.Dependent (ax, ay);  // f(x) = x

     // check function value
     double check = x0;
     ok &= NearEqual( Value(ay[0]) , check,  eps, eps);

     // check zero order forward mode
     size_t q;
     vector<double> x_q(n), y_q(m);
     q      = 0;
     x_q[0] = x0;
     y_q    = f.Forward(q, x_q);
     ok &= NearEqual(y_q[0] , check,  eps, eps);

     return ok;
}