[Rcpp-commits] r2381 - in pkg/RcppDE: . src
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Oct 31 23:54:09 CET 2010
Author: edd
Date: 2010-10-31 23:54:09 +0100 (Sun, 31 Oct 2010)
New Revision: 2381
Modified:
pkg/RcppDE/openmp.r
pkg/RcppDE/src/Makevars
pkg/RcppDE/src/devolMP.cpp
pkg/RcppDE/src/permuteMP.cpp
Log:
OpenMP ready by removing all calls to libRmath / unif_rand / *RNGstate
Modified: pkg/RcppDE/openmp.r
===================================================================
--- pkg/RcppDE/openmp.r 2010-10-31 22:52:58 UTC (rev 2380)
+++ pkg/RcppDE/openmp.r 2010-10-31 22:54:09 UTC (rev 2381)
@@ -12,7 +12,7 @@
maxIt <- 500 # not excessive but so that we get some run-time on simple problems
-n <- 20
+n <- 40
set.seed(42)
print(system.time( {
Modified: pkg/RcppDE/src/Makevars
===================================================================
--- pkg/RcppDE/src/Makevars 2010-10-31 22:52:58 UTC (rev 2380)
+++ pkg/RcppDE/src/Makevars 2010-10-31 22:54:09 UTC (rev 2381)
@@ -1,9 +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/devolMP.cpp
===================================================================
--- pkg/RcppDE/src/devolMP.cpp 2010-10-31 22:52:58 UTC (rev 2380)
+++ pkg/RcppDE/src/devolMP.cpp 2010-10-31 22:54:09 UTC (rev 2381)
@@ -38,7 +38,7 @@
for (int i = 0; i < i_NP; i++)
sortIndex[i] = i;
}
- GetRNGstate();
+ //GetRNGstate();
initialpop.zeros(); // initialize initial popuplation
d_bestmemit.zeros(); // initialize best members
@@ -53,7 +53,8 @@
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] + ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
+ ta_popP.at(j,i) = fa_minbound[j] + arma::randu() * (fa_maxbound[j] - fa_minbound[j]);
}
} else { // or user-specified initial member
ta_popP.col(i) = initialpop.col(i);
@@ -82,14 +83,15 @@
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 + ::unif_rand() * (1.0 - f_weight); // ----computer dithering factor -----------------
+ double f_dither = f_weight + arma::randu() * (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
}
-#pragma omp parallel for shared(ta_oldP,ta_newP,ta_newC,t_tmpP) private(i) schedule(dynamic)
+#pragma omp parallel for shared(ia_urn2,ta_oldP,ta_newP,ta_newC,t_tmpP,ia_urntmp) schedule(static)
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
@@ -101,68 +103,87 @@
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>(::unif_rand() * i_D); // random parameter
+ int j = static_cast<int>(arma::randu() * 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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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>(::unif_rand() * i_D); // random parameter
+ int j = static_cast<int>(arma::randu() * 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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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>(::unif_rand() * i_D); // random parameter
+ int j = static_cast<int>(arma::randu() * i_D); // random parameter
do { // add fluctuation to random target
- double f_jitter = 0.0001 * ::unif_rand() + f_weight;
+ //double f_jitter = 0.0001 * ::unif_rand() + f_weight;
+ double f_jitter = 0.0001 * arma::randu() + 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 ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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>(::unif_rand() * i_D); // random parameter
+ int j = static_cast<int>(arma::randu() * 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 + ::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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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>(::unif_rand() * i_D); // random parameter
+ int j = static_cast<int>(arma::randu() * 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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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 j = static_cast<int>(::unif_rand() * i_D); // random parameter
+ //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
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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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
- if (::unif_rand() < 0.5) { // differential mutation, Pmu = 0.5
+ //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
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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < 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]));
j = (j + 1) % i_D;
- } while ((::unif_rand() < f_cross) && (++k < i_D));
+ //} while ((::unif_rand() < f_cross) && (++k < i_D));
+ } while ((arma::randu() < f_cross) && (++k < i_D));
}
break;
}
@@ -170,10 +191,12 @@
for (int j = 0; j < i_D; j++) { // ----boundary constraints, bounce-back method was not enforcing bounds correctly
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] + ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
+ t_tmpP[j] = fa_minbound[j] + arma::randu() * (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] - ::unif_rand() * (fa_maxbound[j] - fa_minbound[j]);
+ t_tmpP[j] = fa_maxbound[j] - arma::randu() * (fa_maxbound[j] - fa_minbound[j]);
}
}
@@ -257,7 +280,7 @@
d_pop = ta_oldP;
i_iterations = i_iter;
- PutRNGstate();
+ //PutRNGstate();
// ProfilerStop();
}
#endif
Modified: pkg/RcppDE/src/permuteMP.cpp
===================================================================
--- pkg/RcppDE/src/permuteMP.cpp 2010-10-31 22:52:58 UTC (rev 2380)
+++ pkg/RcppDE/src/permuteMP.cpp 2010-10-31 22:54:09 UTC (rev 2381)
@@ -30,7 +30,7 @@
// indices from ia_urn2[i], i = 1, 2, 3, ..., i_urn2_depth.
// Return Value : -
void permute(int ia_urn2[], int i_urn2_depth, int i_NP, int i_avoid, int ia_urn1[]) {
- GetRNGstate();
+ //GetRNGstate();
int i_urn1 = 0, i_urn2 = 0, k = i_NP;
for (int i = 0; i < i_NP; i++)
ia_urn1[i] = i; /* initialize urn1 */
@@ -43,8 +43,8 @@
//k = k - 1; /* reduce number of accessible indices */
i_urn2 = i_urn2 + 1; /* next position in urn2 */
//i_urn1 = static_cast<int>(::unif_rand() * k); /* choose a random index */
- i_urn1 = static_cast<int>(::unif_rand() * (k-1)); /* choose a random index */
+ i_urn1 = static_cast<int>(arma::randu() * (k-1)); /* choose a random index */
}
- PutRNGstate();
+ //PutRNGstate();
}
#endif
More information about the Rcpp-commits
mailing list