CPPAD_DISCRETE_FUNCTION does not work through base2ad

[perrydv]

Hi @bradbell, it looks like CPPAD_DISCRETE_FUNCTION does not work correctly through use of base2ad. Reproducible example below. Is it supposed to work? Am I trying to use it correctly? Thanks!

In the code below, there is a flag to use a simple alternative calculation in the AD tape, as a check that the code works except for the CPPAD_DISCRETE_FUNCTION issue.

#include <cppad/cppad.hpp>
using CppAD::AD;
using CppAD::ADFun;
using std::cout;
using std::endl;

bool do_case_of_interest_or_alternative_for_checking;

double identity(const double &x) {
  return x;
}
CPPAD_DISCRETE_FUNCTION(double, identity);

void record_tape(double x, double y, ADFun<double> &adtape) {
  // f(x, y) = identity(x) * y, where identity should set derivatives for x to 0.
  std::vector<AD<double> > indVars(2);
  std::vector<AD<double> > depVars(1);

  indVars[0] = x;
  indVars[1] = y;
  CppAD::Independent(indVars);
  if(do_case_of_interest_or_alternative_for_checking) {
    depVars[0] = identity(indVars[0]) * indVars[1]; // this line fails in the double tape. CPPAD_DISCRETE_FUNCTION seems like it can't be double-taped.
  } else {
    depVars[0] = indVars[0] * indVars[1]; // this line works in double taping.
  }
  adtape.Dependent(indVars, depVars);
}

void record_double_tape(double x, double y, ADFun<AD<double>, double> &innertape, ADFun<double> &adtape) {
  // y = x * x * x
  std::vector<AD<double> > indVars(2);
  std::vector<AD<double> > depVars(1);

  indVars[0] = x;
  indVars[1] = y;
  CppAD::Independent(indVars);
  depVars = innertape.Forward(0, indVars);
  adtape.Dependent(indVars, depVars);
}

void discrete_test() {
  ADFun<double> adtape;
  double xRec = 1.0;
  double yRec = 2.0;
  record_tape(xRec, yRec, adtape);

  std::vector<double> newInd(2);
  std::vector<double> res(2);
  std::vector<double> w(1);

  newInd[0] = 10.;
  newInd[1] = 20.;

  adtape.Forward(0, newInd);
  w[0] = 1.;
  res = adtape.Reverse(1, w);
  std::cout<<"single tape results"<<std::endl;
  if(do_case_of_interest_or_alternative_for_checking) {
    std::cout<<"df/dx = "<<res[0]<<" (should be 0) "<<" df/dy = "<<res[1]<<" (should be 10)"<<std::endl;
  } else {
    std::cout<<"df/dx = "<<res[0]<<" (should be 20) "<<" df/dy = "<<res[1]<<" (should be 10)"<<std::endl;
  }
  ADFun<double> addoubletape;
  ADFun<AD<double>, double > adtape2;
  adtape2 = adtape.base2ad();
  record_double_tape(xRec, yRec, adtape2, addoubletape);

  addoubletape.Forward(0, newInd);
  res = addoubletape.Reverse(1, w);
  std::cout<<"double tape results"<<std::endl;
  if(do_case_of_interest_or_alternative_for_checking) {
    std::cout<<"df/dx = "<<res[0]<<" (should be 0) "<<" df/dy = "<<res[1]<<" (should be 10)"<<std::endl;
  } else {
    std::cout<<"df/dx = "<<res[0]<<" (should be 20) "<<" df/dy = "<<res[1]<<" (should be 10)"<<std::endl;
  }
}

int main(void) {
  do_case_of_interest_or_alternative_for_checking = false;
  cout<<"\nchecking alternative case to show code makes sense (doing this first because next case errors out)"<<endl;
  discrete_test();

  do_case_of_interest_or_alternative_for_checking = true;
  cout<<"\nchecking case of interest (crash expected before double-taping results)"<<endl;
  discrete_test();
}

[bradbell]

@perrydv Thanks for the bug report. I think the following commit fixes this bug: https://github.com/coin-or/CppAD/commit/a7fe6c0fe50014b34b75448199af1212e09d5d6f

Please try the current master and see if it works for you. If so, then please close the bug.

Thanks again.

[perrydv]

@bradbell Sorry for the delay and yes I confirm this fixes the test case I was using. Thanks. I will close now.