Error in result

[drjdn]

Hi,

The following code produces the wrong result:

#include <cppad/example/cppad_eigen.hpp>
#include <Eigen/Dense>
#include <iostream>

using CppAD::AD;
using CppAD::ADFun;
using CppAD::Independent;
using Eigen::Matrix;
using Eigen::Dynamic;
using Eigen::VectorXd;
using Eigen::MatrixXd;


int main() {

  typedef Matrix< AD<double> , Dynamic, Dynamic > a_MatrixXd;
  typedef Matrix< AD<double> , Dynamic , 1>       a_VectorXd;

  int size = 3, i , j;

  a_VectorXd a_x(size);
  VectorXd x(size);

  x << 1, 2, 3;
  for(i = 0; i < size; i++) {
    a_x[i] = x[i];
  }

  a_MatrixXd a_X(size, size);
  a_X << 1, 2, 3, 4, 5, 6, 7, 8, 9;

  Independent(a_x);
  AD<double> dd = 2.0;
  a_VectorXd a_y = (a_x.transpose() * a_X * a_x)/dd;


  // create f: x -> y and stop tape recording
  ADFun<double> f(a_x, a_y);

  VectorXd jac = f.Jacobian(x);
  VectorXd H = f.Hessian(x,0);

  MatrixXd Hss(size,size);
  Hss = H;
  Hss.resize(size,size);

  std::cout << "True Jacobian/gradient =\n" << a_X * a_x << "\n\n";
  std::cout << "CppAD Jacobian/gradient =\n" << jac << "  [FAIL]\n\n";
  std::cout << "True Hessian =\n" << a_X << "\n\n";
  std::cout << "CppAD Hessian =\n" << Hss << "  [FAIL]\n\n";

  return 0;
}

I'm running Ubuntu 15.10 on an x86-64 machine using CppAD master and the most recent stable release of Eigen (3.2.7). If I set a_X and a_b using for loops I get the correct result.

Thanks,
Jason

[bradbell]

Would you please simplify your example. Just output the value you get and the value you expect for the simplest case that you consider wrong. For example, if it is the Jacobian, you could output
std::cout << "CppAD Jacobian = " << jac << "\n\n"
std::cout << "Expected Jacobian = " << ?? << "\n\n";

[drjdn]

I've simplified as much as I could above and provided better info in the output as suggested. It just computes a quadratic form and calculates its Jacobian/gradient and Hessian via CppAD. I could also remove the Hessian code if you prefer (under the assumption that the same issues is causing the problem in both cases).

Thanks,
Jason

[bradbell]

I have reproduced and checked in a copy of your bug report; see
27158c6
When I run the script I get
True f'(x) = 7 10
CppAD f'(x) = 6 10.5

Thanks.

[bradbell]

The matrix Q is not symmetric, hence the true derivative is more complicated than your formula; see
7823e24
Now when I run the script I get
True f'(x) = 6 10.5
CppAD f'(x) = 6 10.5

Please close this bug if you are satisfied with this result

[drjdn]

Arrrgh! Very sorry for the noise. I work almost exclusively with symmetric positive definite Q matrices and hence my oversight.... Apologies.

Jason