adcomp
adcomp copied to clipboard
Namespace declarations to enable Eigen/Tensor in TMB models
This PR makes a few small header edits to enable including the Eigen/Tensor module in TMB previously discussed on the TMB mailing list here: https://groups.google.com/g/tmb-users/c/lIsIJ0njMOA
It make the following changes:
-
Explicitly reference
tmbutils::arrayintmbutils/density.hppandtmbutils/R_inla.hpp.This avoids a namespace clash between
Eigen::andtmbutils::when Eigen/Tensor headers are included:> TMB/include/tmbutils/density.hpp:106:3: error: > reference to 'array' is ambiguous -
Add
#define R_NO_REMAPbefore#include <Rinternals.h>in TMB.hpp. This avoids a compiler error when including<TMB.hpp>is included before<Tensor>(compiler Apple clang version 13.0.0):
> unsupported/Eigen/CXX11/src/Tensor/TensorTraits.h:122:8: error:
> explicit specialization of non-template struct 'Rf_eval'
-
Specify
Rf_prefix in a few locations intmb_core.hppanddynamic_data.hppwhere it was omitted (Rf_findVar,Rf_install,Rf_setAttrib,Rf_mkString,Rf_ScalarLogical). This was necessitated by the addition of#define R_NO_REMAP. I think these were just oversights because theRf_prefix is used in most calls to those functions in those files. -
Update
inst/include/unsupported/headers to Eigen 3.4.0 for consistency with updated headers ininst/include/Eigen.
I have aimed to be relatively light touch with edits because my understanding of the internals is not very deep. But let me know if there is anything more that I can do to help or clarify.
Thanks for your consideration, and, as ever, many thanks for the brilliant software.
Thanks, Jeff
I see in the test failures that adding #define R_NO_REMAP has knock on effects for glmmTMB:
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/glmmTMB.cpp#L170
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/glmmTMB.cpp#L171
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/glmmTMB.cpp#L176
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/glmmTMB.cpp#L177
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/distrib.h#L182
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/distrib.h#L276
- https://github.com/glmmTMB/glmmTMB/blob/master/glmmTMB/src/distrib.h#L298
and maybe elsewhere. Curious your thoughts about making PRs to add Rf_ prefix there also, or another solution than adding #define R_NO_REMAP here.
@jeffeaton I agree that getting Eigen tensors to work in TMB seems like a relevant request. Obviously, I can't merge this PR because it breaks too much stuff, but this is a good starting point - thanks! I'll experiment a bit with this myself for now.
Curious your thoughts about making PRs to add Rf_ prefix there also
It's my impression that error and warning are so widely used that it would be a pain to change. Another solution might be to set R_NO_REMAP and then add a few remaps manually:
#define error Rf_error
#define warning Rf_warning
Many thanks! Let me know if anything further I can do to help.
@jeffeaton
I updated Eigen/unsupported to version 3.4 and can now do (without further changes):
tensor.cpp
#include <TMB.hpp>
#undef eval // Workaround
#include <unsupported/Eigen/CXX11/Tensor>
template<class Type>
Type objective_function<Type>::operator() ()
{
PARAMETER(p);
Eigen::Tensor<Type, 2> a(2, 3);
a.setValues({{p,p+1,p+2},{p+3,p+4,p+5}});
Rcout << a << "\n\n";
Eigen::Tensor<Type, 0> b = a.sum();
return b(0);
}
tensor.R
library(TMB)
## Compile and load the model
compile("tensor.cpp")
dyn.load(dynlib("tensor"))
## Data and parameters
data <- list()
parameters <- list(p=0)
## Make a function object
obj <- MakeADFun(data, parameters, DLL="tensor")
Thanks so much for the easy solution on this.
Just leaving a note that to get the example above working on Windows (with Rtools 4.2), I had to add two more #undef statements:
tensor.cpp
#include <TMB.hpp>
#undef eval // Workaround
#undef Realloc
#undef Free
#include <unsupported/Eigen/CXX11/Tensor>
template<class Type>
Type objective_function<Type>::operator() ()
{
PARAMETER(p);
Eigen::Tensor<Type, 2> a(2, 3);
a.setValues({{p,p+1,p+2},{p+3,p+4,p+5}});
Rcout << a << "\n\n";
Eigen::Tensor<Type, 0> b = a.sum();
return b(0);
}
The resolution was derived from this SO post: https://stackoverflow.com/questions/11588765/using-rcpp-with-windows-specific-includes