I'm trying to call a C routine from the cubature package in a c++ function to perform multidimensional integration.
The basic R example I'm trying to reproduce is
library(cubature)
integrand <- function(x) sin(x)
adaptIntegrate(integrand, 0, pi)
I could just call this R function from Rcpp following this recipe from the gallery, but there would be some performance penalty in switching back and forth from c/c++ to R. It seems more sensible to directly call the C function from C++.
The C routine adapt_integrate
is exported from cubature
with
// R_RegisterCCallable("cubature", "adapt_integrate", (DL_FUNC) adapt_integrate);
I don't understand how to call it from c++, however. Here's my lame attempt,
sourceCpp(code = '
#include <Rcpp.h>
using namespace Rcpp;
// [[Rcpp::export]]
double integrand(double x){
return(sin(x));
}
// [[Rcpp::depends(cubature)]]
// [[Rcpp::export]]
Rcpp::List integratecpp(double llim, double ulim)
{
Rcpp::Function p_cubature = R_GetCCallable("cubature", "adapt_integrate");
Rcpp::List result = p_cubature(integrand, llim, ulim);
return(result);
}
'
)
integratecpp(0, pi)
This fails to compile; clearly I'm doing something very silly and missing some important steps to convert the output of R_GetCCallable
into an Rcpp::Function
(or call it directly?). I've read several related posts dealing with function pointers, but haven't seen an example using an external C function.
Unfortunately cubature
does not ship the headers in inst/include
, so you have to borrow that from them and do something like this in your code:
typedef void (*integrand) (unsigned ndim, const double *x, void *,
unsigned fdim, double *fval);
int adapt_integrate(
unsigned fdim, integrand f, void *fdata,
unsigned dim, const double *xmin, const double *xmax,
unsigned maxEval, double reqAbsError, double reqRelError,
double *val, double *err)
{
typedef int (*Fun)(unsigned,integrand,void*,unsigned,
const double*,const double*, unsigned, double, double, double*, double*) ;
Fun fun = (Fun) R_GetCCallable( "cubature", "adapt_integrate" ) ;
return fun(fdim,f,fdata,dim,xmin,xmax,maxEval,reqAbsError, reqRelError,val,err);
}
It might be a good idea to negociate with the maintainer of cubature
that he ships declarations in inst/include
so that you'd only have to use LinkingTo
.
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With