[CppAD] filling Eigen sparse matrix in presence of CppAD
Braun, Michael
braunm at mail.smu.edu
Thu Dec 18 18:42:40 EST 2014
The problem with deleting typedef Scalar value_type is that when I try to use the CppAD functionality, I get a
/opt/cppad/include/cppad/local/independent.hpp:159:30: error: no type named
'value_type' in 'Eigen::Matrix<CppAD::AD<double>, -1, 1, 0, -1, 1>'
{ typedef typename VectorAD::value_type ADBase;
Here is another MWE that uses a tape. It does not compile. The plugin files are below the main trips.cpp file. In short, removing the typedef means that there is no value_type. Including it gives me the “reserve is ambiguous” error.
// trips.cpp
#define EIGEN_MATRIX_PLUGIN <eigen_plugin.h>
#define EIGEN_SPARSEMATRIX_PLUGIN <eigen_sparse_plugin.h>
#include <cppad/cppad.hpp>
#include <Eigen/Eigen>
#include <Eigen/Sparse>
#include <cppad/example/cppad_eigen.hpp>
void triptestAD () {
using Eigen::Matrix;
using Eigen::Dynamic;
using Eigen::Triplet;
using Eigen::VectorXd;
typedef CppAD::AD<double> AScalar;
typedef Matrix<AScalar, Dynamic, Dynamic> MatrixXA;
typedef Matrix<AScalar, Dynamic, 1> VectorXA;
typedef Triplet<double> TT;
CppAD::ADFun<double> tape;
int n = 5;
VectorXA f(1); // to hold result
VectorXd Pd(5);
Pd << 1, 3, 5, 7, 9;
VectorXA P = Pd.cast<AScalar>();
std::vector<TT> trips;
Eigen::SparseMatrix<AScalar> M(n,n);
for (int i=0; i<n; i++) {
trips.push_back (TT(i,i,i*i));
}
M.setFromTriplets(trips.begin(), trips.end());
CppAD::Independent(P);
f(0) = (M * P).sum();
tape.Dependent(P, f);
tape.optimize();
VectorXd val(1);
std::cout << "f = " << val(0) << "\n\n";
val = tape.Forward(0, Pd);
VectorXd grad = tape.Jacobian(Pd);
std::cout << "grad = " << grad << "\n";
}
int main (void)
{
triptestAD();
return 0;
}
// eigen_plugin.h
typedef Scalar value_type;
// eigen_sparse_plugin.h (following suggestion in Eigen support forum).
// typedef Scalar value_type;
template<typename S, int R, int C, int Opt, int MR, int MC>
void reserve(const Matrix<S,R,Opt,MR,MC> &reserveSizes)
{
reserveInnerVectors(reserveSizes);
}
On Dec 18, 2014, at 6:25 AM, Brad Bell <bradbell at seanet.com<mailto:bradbell at seanet.com>> wrote:
Using g++ 4.8.3, if you delete
typedef Scalar value_type;
at line 24 of
https://projects.coin-or.org/CppAD/browser/trunk/cppad/example/eigen_plugin.hpp
the error in compiling the example below goes away.
CppAD expects
Vector::value_type
to be the type of the elements of a vector; see
http://www.coin-or.org/CppAD/Doc/simplevector.xml#Value%20Type
This is the same as for the C++ standard vector; as can be seen by searching for value_type on
http://www.cplusplus.com/reference/vector/vector/
I wanted to ask about this on the eigen user's list
http://listengine.tuxfamily.org/lists.tuxfamily.org/eigen/2014/12/maillist.html
but I did not see your posting there.
Brad
On 12/17/2014 11:33 AM, Braun, Michael wrote:
I posted this same question on the Eigen forum because I do not know if it is really an Eigen issue or a CppAD issue.
As part of a larger application that uses CppAD, I am trying to fill an Eigen sparse matrix. The following MWE does not use CppAD, and when the first line (#define EIGEN_MATRIX_PLUGIN ... ) is commented out, this code runs fine.
#define EIGEN_MATRIX_PLUGIN <cppad/example/eigen_plugin.hpp>
#include <Eigen/Eigen>
#include <Eigen/Sparse>
#include <cppad/cppad.hpp>
#include <cppad/example/cppad_eigen.hpp>
void triptest () {
typedef Eigen::Triplet<double> TT;
int n = 5;
std::vector<TT> trips;
Eigen::SparseMatrix<double> M(n,n);
for (int i=0; i<n; i++) {
trips.push_back (TT(i,i,i*i));
}
M.setFromTriplets(trips.begin(), trips.end());
std::cout << M << "\n";
}
int main (void)
{
triptest();
return 0;
}
But when left in, I get the following compile error:
$ clang++ -I/opt/cppad/include -I$R_HOME/library/RcppEigen/include trips.cpp -o trips.so
In file included from trips.cpp:4:
In file included from .../Eigen/Sparse:19:
In file included from .../Eigen/SparseCore:40:
/Library/Frameworks/R.framework/Resources/library/RcppEigen/include/Eigen/src/SparseCore/SparseMatrix.h:957:11: error:
call to member function 'reserve' is ambiguous
trMat.reserve(wi);
~~~~~~^~~~~~~
.../Eigen/src/SparseCore/SparseMatrix.h:1013:13: note:
in instantiation of function template specialization
'Eigen::internal::set_from_triplets<std::__1::__wrap_iter<Eigen::Triplet<double, int> *>,
Eigen::SparseMatrix<double, 0, int> >' requested here
internal::set_from_triplets(begin, end, *this);
^
trips.cpp:18:5: note: in instantiation of function template specialization 'Eigen::SparseMatrix<double, 0,
int>::setFromTriplets<std::__1::__wrap_iter<Eigen::Triplet<double, int> *> >' requested here
M.setFromTriplets(trips.begin(), trips.end());
^
.../Eigen/src/SparseCore/SparseMatrix.h:270:17: note:
candidate function [with SizesType = Eigen::Matrix<int, -1, 1, 0, -1, 1>]
inline void reserve(const SizesType& reserveSizes, const typename SizesType::value_type& enableif = typenam...
^
.../Eigen/src/SparseCore/SparseMatrix.h:276:17: note:
candidate function [with SizesType = Eigen::Matrix<int, -1, 1, 0, -1, 1>]
inline void reserve(const SizesType& reserveSizes, const typename SizesType::Scalar& enableif =
^
.../Eigen/src/SparseCore/SparseMatrix.h:256:17: note:
candidate function not viable: no known conversion from 'Matrix<Index, Dynamic, 1>' to 'Index' (aka 'int') for 1st
argument
inline void reserve(Index reserveSize)
^
Any idea what's going on? I have been able to narrow the problem directly to that first #define. Also, this appears to happen only with EIGEN_MATRIX_PLUGIN, and not ARRAY or SPARSEMATRIX.
Thanks,
Michael
--------------------------
Michael Braun
Associate Professor of Marketing
Cox School of Business
Southern Methodist University
Dallas, TX 75275
braunm at smu.edu<mailto:braunm at smu.edu>
_______________________________________________
CppAD mailing list
CppAD at list.coin-or.org<mailto:CppAD at list.coin-or.org>
http://list.coin-or.org/mailman/listinfo/cppad
--------------------------------------------
Michael Braun, Ph.D.
Associate Professor of Marketing
Cox School of Business
Southern Methodist University
braunm at smu.edu<mailto:braunm at smu.edu>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://list.coin-or.org/pipermail/cppad/attachments/20141218/a083e2b9/attachment.html>
More information about the CppAD
mailing list