Prev Next cppad_poly.cpp

CppAD Speed: Second Derivative of a Polynomial

Operation Sequence
Note that the polynomial evaluation operation sequence does not depend on the argument to the polynomial. Hence we use the same ADFun object for all the argument values.

compute_poly
Routine that computes the second derivative of a polynomial using CppAD:
# include <cppad/cppad.hpp>
# include <cppad/speed/uniform_01.hpp>

void compute_poly(
     size_t                     size     , 
     size_t                     repeat   , 
     CppAD::vector<double>     &a        ,  // coefficients of polynomial
     CppAD::vector<double>     &z        ,  // polynomial argument value
     CppAD::vector<double>     &ddp      )  // second derivative w.r.t z  
{
     // -----------------------------------------------------
     // setup
     typedef CppAD::AD<double>     ADScalar; 
     typedef CppAD::vector<ADScalar> ADVector; 

     size_t i;      // temporary index
     size_t m = 1;  // number of dependent variables
     size_t n = 1;  // number of independent variables
     ADVector Z(n); // AD domain space vector
     ADVector P(m); // AD range space vector

     // choose the polynomial coefficients
     CppAD::uniform_01(size, a);

     // AD copy of the polynomial coefficients
     ADVector A(size);
     for(i = 0; i < size; i++)
          A[i] = a[i];

     // forward mode first and second differentials
     CppAD::vector<double> p(1), dp(1), dz(1), ddz(1);
     dz[0]  = 1.;
     ddz[0] = 0.;

     // choose an argument value
     CppAD::uniform_01(1, z);
     Z[0] = z[0];

     // declare independent variables
     Independent(Z);

     // AD computation of the function value 
     P[0] = CppAD::Poly(0, A, Z[0]);

     // create function object f : A -> detA
     CppAD::ADFun<double> f(Z, P);

     // ------------------------------------------------------
     while(repeat--)
     {    // get the next argument value
          CppAD::uniform_01(1, z);

          // evaluate the polynomial at the new argument value
          p = f.Forward(0, z);

          // evaluate first order Taylor coefficient
          dp = f.Forward(1, dz);

          // second derivative is twice second order Taylor coefficient
          ddp     = f.Forward(2, ddz);
          ddp[0] *= 2.;
     }
     return;
}

Input File: speed/cppad/poly.cpp