[Rcpp-devel] Additional parameters for an objective function, e.g. in RcppDE

Christoph Bergmeir c.bergmeir at decsai.ugr.es
Mon Apr 29 04:37:07 CEST 2013


Dear list,

I'm looking for some advice on a specific problem. Using RcppDE there is 
the possibility to give the optimizer directly an external pointer to 
the C++ function it will use as the objective function. I found this 
mechanism pretty useful as it may speed up things quite a lot (I have a 
problem where the speedup is from 17 minutes to some seconds), so that I 
use the same mechanism as RcppDE in our package Rmalschains and in the 
Rdonlp2 package, which is available from Rmetrics on Rforge.

The problem that this mechanism has is that it cannot handle additional 
parameters to the objective function. Having additional parameters is 
often essential, because if you fit a model to data you need the data 
available in the target function. I illustrate the problem with an 
example I took from the RcppDE tests:

#-----------------------------------------

library(inline)
library(RcppDE)

inc <- 'double rastrigin(SEXP xs) { //here I want to give it an 
additional parameter: SEXP additional_parameter

   //Do something with the parameter, e.g. use it for result 
calculation. Here we just want to print it
   //double my_additional_parameter = 
Rcpp::as<double>(additional_parameter);
   //Rprintf("ap: %f\\n", my_additional_parameter);

   Rcpp::NumericVector x(xs);
   int n = x.size();
   double sum = 20.0;
   for (int i=0; i<n; i++) {
   sum += x[i]+2 - 10*cos(2*M_PI*x[i]);

}
return(sum);
}
'

src.xptr <- '
     typedef double (*funcPtr)(SEXP);
     return(XPtr<funcPtr>(new funcPtr(&rastrigin)));
     '
create_xptr <- cxxfunction(signature(), body=src.xptr, inc=inc, 
plugin="Rcpp")

n <- 10
maxIt <- 100

res <- RcppDE::DEoptim(fn=create_xptr(), lower=rep(-25, n), 
upper=rep(25, n),
       control=list(NP=10*n, itermax=maxIt, trace=FALSE)) #, 
additional_paramater=25)

res$optim

#-----------------------------------------

I currently get around this by having a global singleton object which 
holds these parameters. This works but of course is not very nice when 
it comes to parallelization. The code is more or less like this:

//----------------------------------------------
class TargetFunction {

   private:

   static TargetFunction *TargetFunctionSingleton;
   std::vector<double> param;
   double objval;

   public:

   void eval(const double* x, int n) {
     double sum = 20.0;
     for (int i=0; i<n; i++) {
       sum += x[i]+2 - 10*cos(2*M_PI*x[i]);
     };

//here I can use the parameter now!!
     Rprintf("ap: %f\\n", param[0]);

     this->objval = sum;
   };

   void init(std::vector<double> & p_param) {
	  this->param = p_param;
   };

   static TargetFunction* getTargetFunctionSingleton() {
	  if( TargetFunctionSingleton == 0 )
		  TargetFunctionSingleton = new TargetFunction();
	  return TargetFunctionSingleton;
   };

   static void deleteTargetFunctionSingleton(){
	  if( TargetFunctionSingleton == 0 ) return;
	  else {
		  delete TargetFunctionSingleton;
		  TargetFunctionSingleton = 0;
	  }
	  return;
   };

   double getObjVal() {
     return(objval);
   };


};

TargetFunction* TargetFunction::TargetFunctionSingleton = 0;

RcppExport SEXP targetFunction(SEXP p_par)
{
	Rcpp::NumericVector par(p_par);

	TargetFunction* sp = TargetFunction::getTargetFunctionSingleton();

	sp->eval(par.begin(), par.size());

	return Rcpp::wrap(sp->getObjVal());

}

RcppExport SEXP targetFunctionInit(SEXP p_param) {

	TargetFunction::deleteTargetFunctionSingleton();

	TargetFunction* sp = TargetFunction::getTargetFunctionSingleton();

   std::vector<double> param = Rcpp::as< std::vector<double> >(p_param);

	sp->init(param);

	return R_NilValue;

}

RcppExport SEXP GetTargetFunctionPtr() {

	typedef SEXP (*funcPtr)(SEXP);

	return (Rcpp::XPtr<funcPtr>(new funcPtr(&targetFunction)));
}
//-----------------------------------------------------

Now, before doing the optimization, I call targetFunctionInit and set 
the additional parameters. Afterwards, everything is as in the example 
above, and I have the additional parameters available in the target 
function. Now the question is how I could solve this more elegantly, or 
more R like. The first thing that comes to mind is to use an R 
environment instead of the singleton.  However, how can I do this? I 
could have a singleton list of objects and then use the address of the R 
environment as a hash to find the right object in the list. But this is 
probably not really the way R environments should be used, and I wonder 
if this will cause any trouble.

Any advise is highly appreciated.

Regards,
Christoph

-- 
Christoph Bergmeir
e-mail: c.bergmeir at decsai.ugr.es
Grupo SCI2S, DiCITS Lab          (http://sci2s.ugr.es/DiCITS)
Dpto. de Ciencias de la Computacion e Inteligencia Artificial
E.T.S. Ingenierias de Informatica y Telecomunicacion
Universidad de Granada
18071 - GRANADA (Spain)


More information about the Rcpp-devel mailing list