[Rcpp-commits] r2429 - in pkg/RcppDE: scripts src
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Nov 12 04:00:56 CET 2010
Author: edd
Date: 2010-11-12 04:00:53 +0100 (Fri, 12 Nov 2010)
New Revision: 2429
Modified:
pkg/RcppDE/scripts/openmp.r
pkg/RcppDE/src/Makevars
pkg/RcppDE/src/devol.cpp
pkg/RcppDE/src/devolMP.cpp
pkg/RcppDE/src/evaluateMP.cpp
Log:
checking some pending OpenMP changes before splitting this into RcppDE and RcppParDE
Modified: pkg/RcppDE/scripts/openmp.r
===================================================================
--- pkg/RcppDE/scripts/openmp.r 2010-11-11 16:31:10 UTC (rev 2428)
+++ pkg/RcppDE/scripts/openmp.r 2010-11-12 03:00:53 UTC (rev 2429)
@@ -4,25 +4,76 @@
suppressMessages(library(DEoptim)) # the original, currently 2.0.7
suppressMessages(library(RcppDE)) # the contender
+suppressMessages(library(inline)) # helps with
-Genrose <- function(x) { ## One generalization of the Rosenbrock banana valley function (n parameters)
+Genrose <- function(x) { ## One generalization of the Rosenbrock banana valley function (n parameters)
n <- length(x)
1.0 + sum (100 * (x[-n]^2 - x[-1])^2 + (x[-1] - 1)^2)
}
-maxIt <- 500 # not excessive but so that we get some run-time on simple problems
-n <- 40
+#maxIt <- 500 # not excessive but so that we get some run-time on simple problems
+#n <- 40
+maxIt <- 500
+n <- 20
+
+inc <- 'double genrose(SEXP xs) {
+ Rcpp::NumericVector x(xs);
+ int n = x.size();
+ double sum = 1.0;
+ for (int i=1; i<n; i++) {
+ sum += 100*( pow(x[i-1]*x[i-1] - x[i], 2)) + (x[i] - 1)*(x[i] - 1);
+ }
+ return(sum);
+ }
+
+ double wild(SEXP xs) {
+ Rcpp::NumericVector x(xs);
+ int n = x.size();
+ double sum = 0.0;
+ for (int i=0; i<n; i++) {
+ sum += 10 * sin(0.3 * x[i]) * sin(1.3 * x[i]*x[i]) + 0.00001 * x[i]*x[i]*x[i]*x[i] + 0.2 * x[i] + 80;
+ }
+ sum /= n;
+ return(sum);
+ }
+
+ double rastrigin(SEXP xs) {
+ 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);
+ }
+
+ '
+
+## now via a class returning external pointer
+src.xptr <- 'std::string fstr = Rcpp::as<std::string>(funname);
+ typedef double (*funcPtr)(SEXP);
+ if (fstr == "genrose")
+ return(XPtr<funcPtr>(new funcPtr(&genrose)));
+ else if (fstr == "wild")
+ return(XPtr<funcPtr>(new funcPtr(&wild)));
+ else
+ return(XPtr<funcPtr>(new funcPtr(&rastrigin)));
+ '
+
+create_xptr <- cxxfunction(signature(funname="character"), body=src.xptr, inc=inc, plugin="Rcpp")
+xptr <- create_xptr("genrose")
+
set.seed(42)
print(system.time( {
- res <- RcppDE::DEoptim(fn=Genrose, lower=rep(-25, n), upper=rep(25, n), control=list(NP=10*n, itermax=maxIt, trace=FALSE))
- print(res[[1]])
+ res <- RcppDE::DEoptim(fn=xptr, lower=rep(-25, n), upper=rep(25, n), control=list(NP=10*n, itermax=maxIt, trace=FALSE))
+ print(res[[1]][["bestval"]])
}))
set.seed(42)
print(system.time( {
res <- DEoptim::DEoptim(fn=Genrose, lower=rep(-25, n), upper=rep(25, n), control=list(NP=10*n, itermax=maxIt, trace=FALSE))
- print(res[[1]])
+ print(res[[1]][["bestval"]])
}))
Modified: pkg/RcppDE/src/Makevars
===================================================================
--- pkg/RcppDE/src/Makevars 2010-11-11 16:31:10 UTC (rev 2428)
+++ pkg/RcppDE/src/Makevars 2010-11-12 03:00:53 UTC (rev 2429)
@@ -1,11 +1,11 @@
## Hey Emacs make this a -*- mode: makefile; -*- file
##
## -- for OpenMP (with -D macro to switch to OpenMP enabled source file)
-#PKG_CXXFLAGS=-fopenmp -DUSE_OPENMP
-#PKG_LIBS= -fopenmp -lgomp $(shell $(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()") $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)
+PKG_CXXFLAGS=-fopenmp -DUSE_OPENMP
+PKG_LIBS= -fopenmp -lgomp $(shell $(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()") $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)
##
## -- for Google Perftools profiling
## PKG_LIBS= $(shell $(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()") $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS) -lprofiler
##
## -- default
-PKG_LIBS= $(shell $(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()") $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)
+#PKG_LIBS= $(shell $(R_HOME)/bin/Rscript -e "Rcpp:::LdFlags()") $(LAPACK_LIBS) $(BLAS_LIBS) $(FLIBS)
Modified: pkg/RcppDE/src/devol.cpp
===================================================================
--- pkg/RcppDE/src/devol.cpp 2010-11-11 16:31:10 UTC (rev 2428)
+++ pkg/RcppDE/src/devol.cpp 2010-11-12 03:00:53 UTC (rev 2429)
@@ -8,8 +8,9 @@
// (http://www.icsi.berkeley.edu/~storn/DeWin.zip)
#ifndef USE_OPENMP
-#include <RcppArmadillo.h>
-#include "evaluate.h"
+#include <RcppArmadillo.h> // declarations for both Rcpp and RcppArmadillo offering Armadillo classes
+#include <omp.h> // OpenMP for compiler-generated multithreading
+#include "evaluate.h" // simple function evaluation framework
void permute(int ia_urn2[], int i_urn2_depth, int i_NP, int i_avoid, int ia_urntmp[]);
@@ -23,14 +24,13 @@
arma::mat &d_pop, Rcpp::List &d_storepop, arma::mat & d_bestmemit, arma::colvec & d_bestvalit,
int & i_iterations, double i_pPct, long & l_nfeval) {
+ //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
} else { // standard mode: env_ is an env, fcall_ is a function
ev = new Rcpp::DE::EvalStandard(fcall, rho); // so assign R function and environment
}
-
- //ProfilerStart("/tmp/RcppDE.prof");
const int urn_depth = 5; // 4 + one index to avoid
Rcpp::NumericVector par(i_D); // initialize parameter vector to pass to evaluate function
arma::icolvec::fixed<urn_depth> ia_urn2; // fixed-size vector for urn draws
@@ -191,7 +191,7 @@
// ------Trial mutation now in t_tmpP-----------------
memcpy(REAL(par), t_tmpP.memptr(), Rf_nrows(par) * sizeof(double));
- double t_tmpC = ev->eval(par);
+ double t_tmpC = ev->eval(par); // Evaluate mutant in t_tmpP
if (t_tmpC <= ta_oldC[i] || i_bs_flag) { // i_bs_flag means will choose best NP later
ta_newP.col(i) = t_tmpP; // replace target with mutant
ta_newC[i] = t_tmpC;
@@ -247,7 +247,7 @@
if (same && i_iter > 1) {
i_xav++;
memcpy(REAL(par), t_bestP.memptr(), Rf_nrows(par) * sizeof(double));
- double tmp_best = ev->eval(par);
+ double tmp_best = ev->eval(par);// if re-evaluation of winner
if (i_av_winner) // possibly letting the winner be the average of all past generations
t_bestC = ((1/(double)i_xav) * t_bestC) + ((1/(double)i_xav) * tmp_best) +
(d_bestvalit[i_iter-1] * ((double)(i_xav - 2))/(double)i_xav);
Modified: pkg/RcppDE/src/devolMP.cpp
===================================================================
--- pkg/RcppDE/src/devolMP.cpp 2010-11-11 16:31:10 UTC (rev 2428)
+++ pkg/RcppDE/src/devolMP.cpp 2010-11-12 03:00:53 UTC (rev 2429)
@@ -10,10 +10,26 @@
#ifdef USE_OPENMP
#include <RcppArmadillo.h> // declarations for both Rcpp and RcppArmadillo offering Armadillo classes
#include <omp.h> // OpenMP for compiler-generated multithreading
+#include "evaluate.h" // simple function evaluation framework
void permute(int ia_urn2[], int i_urn2_depth, int i_NP, int i_avoid, int ia_urntmp[]);
-double evaluate(long &l_nfeval, const double *param, SEXP parS, SEXP fcall, SEXP env);
+inline double drndu(void) {
+#ifdef USE_OPENMP
+ return arma::randu();
+#else
+ return ::unif_rand();
+#endif
+}
+
+inline int irndu(const double val) {
+#ifdef USE_OPENMP
+ return static_cast<int>(arma::randu() * val);
+#else
+ return static_cast<int>(::unif_rand() * val);
+#endif
+}
+
void devol(double VTR, double f_weight, double f_cross, int i_bs_flag,
arma::colvec & fa_minbound, arma::colvec & fa_maxbound, SEXP fcall, SEXP rho, int i_trace,
int i_strategy, int i_D, int i_NP, int i_itermax, arma::mat & initialpopm,
@@ -25,6 +41,12 @@
int & i_iterations, double i_pPct, long & l_nfeval) {
//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
+ } else { // standard mode: env_ is an env, fcall_ is a function
+ ev = new Rcpp::DE::EvalStandard(fcall, rho); // so assign R function and environment
+ }
const int urn_depth = 5; // 4 + one index to avoid
Rcpp::NumericVector par(i_D); // initialize parameter vector to pass to evaluate function
arma::icolvec::fixed<urn_depth> ia_urn2; // fixed-size vector for urn draws
@@ -53,13 +75,13 @@
for (int i = 0; i < i_NP; i++) { // ------Initialization-----------------------------
if (i_specinitialpop <= 0) { // random initial member
for (int j = 0; j < i_D; j++) {
- //ta_popP.at(j,i) = fa_minbound[j] + ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
- ta_popP.at(j,i) = fa_minbound[j] + arma::randu() * (fa_maxbound[j] - fa_minbound[j]);
+ ta_popP.at(j,i) = fa_minbound[j] + drndu() * (fa_maxbound[j] - fa_minbound[j]);
}
} else { // or user-specified initial member
ta_popP.col(i) = initialpop.col(i);
}
- ta_popC[i] = evaluate(l_nfeval, ta_popP.colptr(i), par, fcall, rho);
+ memcpy(REAL(par), ta_popP.colptr(i), i_D * sizeof(double));
+ ta_popC[i] = ev->eval(par);
if (i == 0 || ta_popC[i] <= t_bestC) {
t_bestC = ta_popC[i];
t_bestP = ta_popP.unsafe_col(i);
@@ -73,7 +95,7 @@
int popcnt = 0;
int i_xav = 1;
- for (i_iter=0; (i_iter < i_itermax) && (t_bestC > VTR); i_iter++) { // main loop ====================================
+ for (i_iter=0; (i_iter < i_itermax) && (t_bestC > VTR); i_iter++) { // main loop =============================
if (i_iter % i_storepopfreq == 0 && i_iter >= i_storepopfrom) { // store intermediate populations
d_storepop[popcnt++] = Rcpp::wrap( trans(ta_oldP) );
} // end store pop
@@ -83,15 +105,14 @@
t_bestitP = t_bestP;
//i_iter++; // increase iteration counter
- //double f_dither = f_weight + ::unif_rand() * (1.0 - f_weight); // ----computer dithering factor -----------------
- double f_dither = f_weight + arma::randu() * (1.0 - f_weight); // ----computer dithering factor -----------------
+ double f_dither = f_weight + drndu() * (1.0 - f_weight); // ----computer dithering factor ----------
- if (i_strategy == 6) { // ---DE/current-to-p-best/1 -----------------------------------------------------
- arma::colvec temp_oldC = ta_oldC; // create a copy of ta_oldC to avoid changing it
- rsort_with_index( temp_oldC.memptr(), sortIndex.begin(), i_NP ); // sort temp_oldC to use sortIndex later
+ if (i_strategy == 6) { // ---DE/current-to-p-best/1 --------------------------------------
+ arma::colvec temp_oldC = ta_oldC; // copy of ta_oldC to not change it
+ rsort_with_index( temp_oldC.memptr(), sortIndex.begin(), i_NP ); // to use sortIndex later
}
-#pragma omp parallel for shared(ia_urn2,ta_oldP,ta_newP,ta_newC,t_tmpP,ia_urntmp) schedule(static)
+#pragma omp parallel for shared(ia_urn2,ta_oldP,ta_newP,ta_newC,t_tmpP,ia_urntmp) schedule(static,1)
for (int i = 0; i < i_NP; i++) { // ----start of loop through ensemble------------------------
t_tmpP = ta_oldP.col(i); // t_tmpP is the vector to mutate and eventually select
@@ -102,107 +123,94 @@
// ===Choice of strategy=======================================================
switch (i_strategy) {
- case 1: { // ---classical strategy DE/rand/1/bin-----------------------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 1: { // ---classical strategy DE/rand/1/bin-----------------------------
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target
- t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + f_weight * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
+ t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + f_weight *
+ (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- case 2: { // ---DE/local-to-best/1/bin---------------------------------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 2: { // ---DE/local-to-best/1/bin---------------------------------------
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target
- t_tmpP[j] = t_tmpP[j] + f_weight * (t_bestitP[j] - t_tmpP[j]) + f_weight * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
+ t_tmpP[j] = t_tmpP[j] + f_weight * (t_bestitP[j] - t_tmpP[j]) + f_weight *
+ (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- case 3: { // ---DE/best/1/bin with jitter------------------------------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 3: { // ---DE/best/1/bin with jitter------------------------------------
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target
- //double f_jitter = 0.0001 * ::unif_rand() + f_weight;
- double f_jitter = 0.0001 * arma::randu() + f_weight;
+ double f_jitter = 0.0001 * drndu() + f_weight;
t_tmpP[j] = t_bestitP[j] + f_jitter * (ta_oldP.at(j,ia_urn2[1]) - ta_oldP.at(j,ia_urn2[2]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- case 4: { // ---DE/rand/1/bin with per-vector-dither-------------------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 4: { // ---DE/rand/1/bin with per-vector-dither-------------------------
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target *
- //t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + (f_weight + ::unif_rand()*(1.0 - f_weight))* (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
- t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + (f_weight + arma::randu()*(1.0 - f_weight))* (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
+ t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + (f_weight + drndu()*(1.0 - f_weight))
+ * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- case 5: { // ---DE/rand/1/bin with per-generation-dither---------------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 5: { // ---DE/rand/1/bin with per-generation-dither---------------------
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target
- t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + f_dither * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
+ t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + f_dither
+ * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- case 6: { // ---DE/current-to-p-best/1 (JADE)--------------------------------------------
- //int i_pbest = sortIndex[static_cast<int>(::unif_rand() * p_NP)]; // select from [0, 1, 2, ..., (pNP-1)]
- int i_pbest = sortIndex[static_cast<int>(arma::randu() * p_NP)]; // select from [0, 1, 2, ..., (pNP-1)]
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
+ case 6: { // ---DE/current-to-p-best/1 (JADE)--------------------------------
+ int i_pbest = sortIndex[irndu(p_NP)]; // select from [0, 1, 2, ..., (pNP-1)]
+ int j = irndu(i_D); // random parameter
do { // add fluctuation to random target
- t_tmpP[j] = ta_oldP.at(j,i) + f_weight * (ta_oldP.at(j,i_pbest) - ta_oldP.at(j,i)) + f_weight * (ta_oldP.at(j,ia_urn2[1]) - ta_oldP.at(j,ia_urn2[2]));
+ t_tmpP[j] = ta_oldP.at(j,i) + f_weight * (ta_oldP.at(j,i_pbest) - ta_oldP.at(j,i)) +
+ f_weight * (ta_oldP.at(j,ia_urn2[1]) - ta_oldP.at(j,ia_urn2[2]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
break;
}
- default: { // ---variation to DE/rand/1/bin: either-or-algorithm--------------------------
- //int j = static_cast<int>(::unif_rand() * i_D); // random parameter
- int j = static_cast<int>(arma::randu() * i_D); // random parameter
- //if (::unif_rand() < 0.5) { // differential mutation, Pmu = 0.5
- if (arma::randu() < 0.5) { // differential mutation, Pmu = 0.5
+ default: { // ---variation to DE/rand/1/bin: either-or-algorithm--------------
+ int j = irndu(i_D); // random parameter
+ if (drndu() < 0.5) { // differential mutation, Pmu = 0.5
do { // add fluctuation to random target */
- t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + f_weight * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
+ t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) +
+ f_weight * (ta_oldP.at(j,ia_urn2[2]) - ta_oldP.at(j,ia_urn2[3]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
} else { // recombination with K = 0.5*(F+1) -. F-K-Rule
do { // add fluctuation to random target */
- t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + 0.5 * (f_weight + 1.0) * (ta_oldP.at(j,ia_urn2[2]) + ta_oldP.at(j,ia_urn2[3]) - 2 * ta_oldP.at(j,ia_urn2[1]));
+ t_tmpP[j] = ta_oldP.at(j,ia_urn2[1]) + 0.5 * (f_weight + 1.0) *
+ (ta_oldP.at(j,ia_urn2[2]) + ta_oldP.at(j,ia_urn2[3]) - 2 * ta_oldP.at(j,ia_urn2[1]));
j = (j + 1) % i_D;
- //} while ((::unif_rand() < f_cross) && (++k < i_D));
- } while ((arma::randu() < f_cross) && (++k < i_D));
+ } while ((drndu() < f_cross) && (++k < i_D));
}
break;
}
} // end switch (i_strategy) ...
- for (int j = 0; j < i_D; j++) { // ----boundary constraints, bounce-back method was not enforcing bounds correctly
+ for (int j = 0; j < i_D; j++) { // boundary constr, bounce-back method was not enforcing bounds
if (t_tmpP[j] < fa_minbound[j]) {
- //t_tmpP[j] = fa_minbound[j] + ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
- t_tmpP[j] = fa_minbound[j] + arma::randu() * (fa_maxbound[j] - fa_minbound[j]);
+ t_tmpP[j] = fa_minbound[j] + drndu() * (fa_maxbound[j] - fa_minbound[j]);
}
if (t_tmpP[j] > fa_maxbound[j]) {
- //t_tmpP[j] = fa_maxbound[j] - ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
- t_tmpP[j] = fa_maxbound[j] - arma::randu() * (fa_maxbound[j] - fa_minbound[j]);
+ t_tmpP[j] = fa_maxbound[j] - drndu() * (fa_maxbound[j] - fa_minbound[j]);
}
}
// ------Trial mutation now in t_tmpP-----------------
- double t_tmpC = evaluate(l_nfeval, t_tmpP.memptr(), par, fcall, rho); // Evaluate mutant in t_tmpP[]
- if (t_tmpC <= ta_oldC[i] || i_bs_flag) { // i_bs_flag means that we will choose best NP vectors from old and new population later
+ memcpy(REAL(par), t_tmpP.memptr(), i_D * sizeof(double));
+ double t_tmpC = ev->eval(par); // Evaluate mutant in t_tmpP
+ if (t_tmpC <= ta_oldC[i] || i_bs_flag) { // i_bs_flag means will choose best NP later
ta_newP.col(i) = t_tmpP; // replace target with mutant
ta_newC[i] = t_tmpC;
if (t_tmpC <= t_bestC) {
@@ -213,7 +221,11 @@
ta_newP.col(i) = ta_oldP.col(i);
ta_newC[i] = ta_oldC[i];
}
+ //#ifdef USE_OPENMP
+ //std::cerr << "Thread " << omp_get_thread_num() << "\n";
+ //#endif
} // End mutation loop through pop., ie the "for (i = 0; i < i_NP; i++)"
+#pragma omp barrier
if (i_bs_flag) { // examine old and new pop. and take the best NP members into next generation
@@ -256,9 +268,9 @@
}
if (same && i_iter > 1) {
i_xav++;
- double tmp_best = evaluate(l_nfeval, t_bestP.memptr(), par, fcall, rho); // if re-evaluation of winner
-
- if (i_av_winner) // possibly letting the winner be the average of all past generations
+ memcpy(REAL(par), t_bestP.memptr(), i_D * sizeof(double));
+ double tmp_best = ev->eval(par);// if re-evaluation of winner
+ if (i_av_winner) // possibly letting winner be the average of all past generations
t_bestC = ((1/(double)i_xav) * t_bestC) + ((1/(double)i_xav) * tmp_best) +
(d_bestvalit[i_iter-1] * ((double)(i_xav - 2))/(double)i_xav);
else
@@ -278,8 +290,8 @@
} // end loop through generations
d_pop = ta_oldP;
- i_iterations = i_iter;
-
+ i_iterations = i_iter;
+ l_nfeval = ev->getNbEvals();
//PutRNGstate();
// ProfilerStop();
}
Modified: pkg/RcppDE/src/evaluateMP.cpp
===================================================================
--- pkg/RcppDE/src/evaluateMP.cpp 2010-11-11 16:31:10 UTC (rev 2428)
+++ pkg/RcppDE/src/evaluateMP.cpp 2010-11-12 03:00:53 UTC (rev 2429)
@@ -4,7 +4,7 @@
// Copyright (C) 2010 Dirk Eddelbuettel <edd at debian.org>
//
// DEoptim is Copyright (C) 2009 David Ardia and Katharine Mullen
-
+#if 0
#ifdef USE_OPENMP
#include <RcppArmadillo.h>
@@ -22,4 +22,4 @@
return(sum);
}
#endif
-
+#endif
More information about the Rcpp-commits
mailing list