Hi Dirk (and the rest of the list),
I think I have it working; attached is a patch that adds this
functionality to RcppDE. I also came along something which is probably a
bug: The passed environment is never used. In line 75 of DEoptim.R you
just create one. I changed this to
if(!hasArg(env)) env <- new.env()
Now you can use it like this:
#--------------------------
# RcppDE example
#--------------------------
library(inline)
library(RcppDE)
inc <- 'double rastrigin(SEXP xs, SEXP env) {
Rcpp::NumericVector x(xs);
Rcpp::Environment e(env);
std::vector<double> param = e["x"];
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]);
Rprintf("ap: %f\\n", param[0]);
}
return(sum);
}
'
## now via a class returning external pointer
src.xptr <- '
typedef double (*funcPtr)(SEXP, SEXP);
return(XPtr<funcPtr>(new funcPtr(&rastrigin)));
'
create_xptr <- cxxfunction(signature(), body=src.xptr, inc=inc,
plugin="Rcpp")
n <- 10
maxIt <- 100
env <- new.env()
env[["x"]] <- 1:4
res <- RcppDE::DEoptim(fn=create_xptr(), env=env, lower=rep(-25, n),
upper=rep(25, n),
control=list(NP=10*n, itermax=maxIt, trace=FALSE))
res$optim
#------------------------------------------------
Regards and thanks for the help,
Christoph
On 04/29/2013 02:41 PM, Dirk Eddelbuettel wrote:
On 29 April 2013 at 14:22, Christoph Bergmeir wrote:
| Hi,
|
| ok, that's exactly the file I was looking at right now. The problem is
| that my objective function is implemented in C++, and passed as an
| external pointer. So in evaluate.h, instead of in the EvalStandard case
| I'm in the EvalCompiled case, where the environment is not used. My
| question now is if I can include it there somehow..
Right. Just how EvalBase and EvalStandard extend the basic interface, we
"just" need to create something that lets an XPtr be evaluation along with an
environment.
Anyway, past bedtime here so no more from me...
Dirk
| Regards,
| Christoph
|
| On 04/29/2013 02:16 PM, Dirk Eddelbuettel wrote:
| >
| > On 29 April 2013 at 13:55, Christoph Bergmeir wrote:
| > | Hi,
| > |
| > | thanks Dirk for the prompt reply. Ok, this was also the solution I had
| > | in mind. I'll try to find the code in RcppDE or try to implement it
| > | myself..I'll keep you posted with any advances.
| >
| > I think this may do the trick: demo/environment.R in the RcppDE package:
| >
| >
| > suppressMessages(library(RcppDE))
| >
| > ## somewhat pathodological example with nuisance parameter mul
| > Rastrigin <- function(x) {
| > mul * (sum(x+2 - 10 * cos(2*pi*x)) + 20)
| > }
| >
| > ## create a new environment associated with the function
| > funenv <- environment(fun=Rastrigin)
| > assign("mul", 2, envir=funenv) ## set value
| >
| > out <- DEoptim(Rastrigin, -25, 25,
| > control=list(NP=10, trace=0),
| > env=funenv)
| > summary(out)
| >
| >
| >
| > Hth, Dirk
| >
| > |
| > | Regards,
| > | Christoph
| > |
| > | On 04/29/2013 01:44 PM, Dirk Eddelbuettel wrote:
| > | >
| > | > On 29 April 2013 at 12:37, Christoph Bergmeir wrote:
| > | > | 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
| > | >
| > | > I think it can. The DEoptim folks, particularly Josh, pointed this out
and
| > | > the best general way is to assign all you need in a new environment --
which
| > | > you can assign to from R and (thanks to Rcpp) from C++. Then pass that
down.
| > | >
| > | > I think I have an example of that in the package but I don't have time
right
| > | > now to chase this.
| > | >
| > | > But yes, this _is_ a very neat feature and something that needs broader
| > | > exposure.
| > | >
| > | > Maybe I can help in a few days.
| > | >
| > | > Dirk
| > | >
| > | > | 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.bergm...@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)
| > | > | _______________________________________________
| > | > | Rcpp-devel mailing list
| > | > | Rcpp-devel@lists.r-forge.r-project.org
| > | > |
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel
| > | >
| > |
| > | --
| > | Christoph Bergmeir
| > | e-mail: c.bergm...@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)
| >
|
| --
| Christoph Bergmeir
| e-mail: c.bergm...@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)
--
Christoph Bergmeir
e-mail: c.bergm...@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)
diff -rupN RcppDE_orig/R/DEoptim.R RcppDE/R/DEoptim.R
--- RcppDE_orig/R/DEoptim.R 2012-04-09 05:43:24.000000000 +1000
+++ RcppDE/R/DEoptim.R 2013-04-30 10:01:44.686949637 +1000
@@ -72,7 +72,7 @@ DEoptim <- function(fn, lower, upper, co
else
nam <- paste("par", 1:length(lower), sep = "")
- env <- new.env()
+ if(!hasArg(env)) env <- new.env()
ctrl <- do.call(DEoptim.control, as.list(control))
ctrl$npar <- length(lower)
Binary files RcppDE_orig/src/deoptim.o and RcppDE/src/deoptim.o differ
diff -rupN RcppDE_orig/src/devol.cpp RcppDE/src/devol.cpp
--- RcppDE_orig/src/devol.cpp 2012-04-09 06:24:02.000000000 +1000
+++ RcppDE/src/devol.cpp 2013-04-29 14:52:15.434839798 +1000
@@ -26,7 +26,7 @@ void devol(double VTR, double f_weight,
//ProfilerStart("/tmp/RcppDE.prof");
Rcpp::DE::EvalBase *ev = NULL; // pointer to abstract base class
if (TYPEOF(fcall) == EXTPTRSXP) { // non-standard mode: we are being passed an external pointer
- ev = new Rcpp::DE::EvalCompiled(fcall); // so assign a pointer using external pointer in fcall SEXP
+ ev = new Rcpp::DE::EvalCompiled(fcall, rho); // so assign a pointer using external pointer in fcall SEXP
} else { // standard mode: env_ is an env, fcall_ is a function
ev = new Rcpp::DE::EvalStandard(fcall, rho); // so assign R function and environment
}
Binary files RcppDE_orig/src/devol.o and RcppDE/src/devol.o differ
diff -rupN RcppDE_orig/src/evaluate.h RcppDE/src/evaluate.h
--- RcppDE_orig/src/evaluate.h 2012-04-09 06:24:02.000000000 +1000
+++ RcppDE/src/evaluate.h 2013-04-29 14:56:38.735580878 +1000
@@ -13,7 +13,7 @@
namespace Rcpp {
namespace DE {
- double genrose(SEXP xs) { // genrose function in C++
+ double genrose(SEXP xs, SEXP env) { // genrose function in C++
Rcpp::NumericVector x(xs);
int n = x.size();
double sum = 1.0;
@@ -23,7 +23,7 @@ namespace Rcpp {
return(sum);
}
- double wild(SEXP xs) { // wild function in C++
+ double wild(SEXP xs, SEXP env) { // wild function in C++
Rcpp::NumericVector x(xs);
int n = x.size();
double sum = 0.0;
@@ -35,7 +35,7 @@ namespace Rcpp {
return(sum);
}
- double rastrigin(SEXP xs) { // rastrigin function in C++
+ double rastrigin(SEXP xs, SEXP env) { // rastrigin function in C++
Rcpp::NumericVector x(xs);
int n = x.size();
double sum = 20.0;
@@ -73,22 +73,25 @@ namespace Rcpp {
}
};
- typedef double (*funcPtr)(SEXP);
+ typedef double (*funcPtr)(SEXP, SEXP);
class EvalCompiled : public EvalBase {
public:
- EvalCompiled( Rcpp::XPtr<funcPtr> xptr ) {
+ EvalCompiled( Rcpp::XPtr<funcPtr> xptr, SEXP env_ ) {
funptr = *(xptr);
+ env = env_;
};
- EvalCompiled( SEXP xps ) {
+ EvalCompiled( SEXP xps, SEXP env_ ) {
Rcpp::XPtr<funcPtr> xptr(xps);
funptr = *(xptr);
+ env = env_;
};
double eval(SEXP par) {
neval++;
- return funptr(par);
+ return funptr(par, env);
}
private:
funcPtr funptr;
+ SEXP env;
};
RcppExport SEXP putFunPtrInXPtr(SEXP funname) {
Binary files RcppDE_orig/src/permute.o and RcppDE/src/permute.o differ
Binary files RcppDE_orig/src/RcppDE.so and RcppDE/src/RcppDE.so differ
_______________________________________________
Rcpp-devel mailing list
Rcpp-devel@lists.r-forge.r-project.org
https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel