diff --git a/similarity_search/apps/bench_distfunc.cc b/similarity_search/apps/bench_distfunc.cc index 18de9b6..8011d88 100644 --- a/similarity_search/apps/bench_distfunc.cc +++ b/similarity_search/apps/bench_distfunc.cc @@ -1964,6 +1964,82 @@ void TestSparseJaccardSimilarity(const string& dataFile, size_t N, size_t Rep) { } +template +void TestRenyiDivSlow(size_t N, size_t dim, size_t Rep, T alpha) { + T* pArr = new T[N * dim]; + + T *p = pArr; + for (size_t i = 0; i < N; ++i, p+= dim) { + GenRandVect(p, dim, T(RANGE_SMALL), T(1.0), true /* norm. for regular KL */); + } + + WallClockTimer t; + + t.reset(); + + T DiffSum = 0; + + T fract = T(1)/N; + + for (size_t i = 0; i < Rep; ++i) { + for (size_t j = 1; j < N; ++j) { + DiffSum += 0.01f * renyiDivergenceSlow(pArr + j*dim, pArr + (j-1)*dim, dim, alpha) / N; + } + /* + * Multiplying by 0.01 and dividing the sum by N is to prevent Intel from "cheating": + * + * http://searchivarius.org/blog/problem-previous-version-intels-library-benchmark + */ + DiffSum *= fract; + } + + uint64_t tDiff = t.split(); + + LOG(LIB_INFO) << "Ignore: " << DiffSum; + LOG(LIB_INFO) << typeid(T).name() << " " << "Elapsed: " << tDiff / 1e3 << " ms " << " # of slow Renyi-div. (alpha=" << alpha << ") per second: " << (1e6/tDiff) * N * Rep ; + + delete [] pArr; + +} + +template +void TestRenyiDivFast(size_t N, size_t dim, size_t Rep, T alpha) { + T* pArr = new T[N * dim]; + + T *p = pArr; + for (size_t i = 0; i < N; ++i, p+= dim) { + GenRandVect(p, dim, T(RANGE_SMALL), T(1.0), true /* norm. for regular KL */); + } + + WallClockTimer t; + + t.reset(); + + T DiffSum = 0; + + T fract = T(1)/N; + + for (size_t i = 0; i < Rep; ++i) { + for (size_t j = 1; j < N; ++j) { + DiffSum += 0.01f * renyiDivergenceFast(pArr + j*dim, pArr + (j-1)*dim, dim, alpha) / N; + } + /* + * Multiplying by 0.01 and dividing the sum by N is to prevent Intel from "cheating": + * + * http://searchivarius.org/blog/problem-previous-version-intels-library-benchmark + */ + DiffSum *= fract; + } + + uint64_t tDiff = t.split(); + + LOG(LIB_INFO) << "Ignore: " << DiffSum; + LOG(LIB_INFO) << typeid(T).name() << " " << "Elapsed: " << tDiff / 1e3 << " ms " << " # of fast Renyi-div. (alpha=" << alpha << ") per second: " << (1e6/tDiff) * N * Rep ; + + delete [] pArr; + +} + } // namespace similarity using namespace similarity; @@ -1977,6 +2053,34 @@ int main(int argc, char* argv[]) { int dim = 128; + float alphaStepSlow = 1.0f/4; + for (float alpha = alphaStepSlow; alpha <= 2; alpha+= alphaStepSlow) { + if (fabs(alpha - 1) < 1e-6) continue; + nTest++; + TestRenyiDivSlow(1024, dim, 100, alpha); + } +#if TEST_SPEED_DOUBLE + for (double alpha = alphaStepSlow; alpha <= 2; alpha+= alphaStepSlow) { + if (fabs(alpha - 1) < 1e-6) continue; + nTest++; + TestRenyiDivSlow(1024, dim, 100, alpha); + } +#endif + + float alphaStepFast = 1.0f/32; + for (float alpha = alphaStepFast; alpha <= 2; alpha+= alphaStepFast) { + if (fabs(alpha - 1) < 1e-6) continue; + nTest++; + TestRenyiDivFast(1024, dim, 100, alpha); + } +#if TEST_SPEED_DOUBLE + for (double alpha = alphaStepFast; alpha <= 2; alpha+= alphaStepFast) { + if (fabs(alpha - 1) < 1e-6) continue; + nTest++; + TestRenyiDivFast(1024, dim, 100, alpha); + } +#endif + #if defined(WITH_EXTRAS) nTest++; TestSQFDMinus(2000, 50); diff --git a/similarity_search/apps/tune_vptree.cc b/similarity_search/apps/tune_vptree.cc index 4dd4755..b5d2da9 100644 --- a/similarity_search/apps/tune_vptree.cc +++ b/similarity_search/apps/tune_vptree.cc @@ -113,8 +113,7 @@ void RunExper(unsigned AddRestartQty, LOG(LIB_INFO) << "We are going to tune parameters for " << MethodName; - static std::random_device rd; - static std::mt19937 engine(rd()); + static thread_local auto& engine(GET_RANDOM_GENERATOR); static std::normal_distribution<> normGen(0.0f, log(FullFactor)); AnyParamManager pmgr(IndexParams); diff --git a/similarity_search/include/utils.h b/similarity_search/include/utils.h index 977950b..7cfd6b6 100644 --- a/similarity_search/include/utils.h +++ b/similarity_search/include/utils.h @@ -65,13 +65,13 @@ typedef SSIZE_T ssize_t; #define FIELD_DELIMITER ':' +#define GET_RANDOM_GENERATOR getRandomGenerator() + namespace similarity { using std::string; using std::vector; using std::stringstream; -using std::random_device; -using std::mt19937; using namespace std; @@ -85,6 +85,19 @@ bool DoesFileExist(const char *filename); inline bool DoesFileExist(const string &filename) { return DoesFileExist(filename.c_str()); } extern int randomSeed; + +/* + * Random number generation is thread safe when respective + * objects are not shared among threads. So, we will keep one + * random number generator per thread. + */ +template +inline RandGenType & getRandomGenerator() { + static thread_local RandGenType gen(randomSeed); + + return gen; +} + // random 32-bit integer number inline int32_t RandomInt() { /* @@ -93,10 +106,9 @@ inline int32_t RandomInt() { * random number generator per thread. */ // thread_local is static by default, but let's keep it static for clarity - static thread_local mt19937 gen(randomSeed); static thread_local std::uniform_int_distribution distr(0, std::numeric_limits::max()); - return distr(gen); + return distr(GET_RANDOM_GENERATOR); } template @@ -108,10 +120,9 @@ inline T RandomReal() { * random number generator per thread. */ // thread_local is static by default, but let's keep it static for clarity - static thread_local mt19937 gen(randomSeed); static thread_local std::uniform_real_distribution distr(0, 1); - return distr(gen); + return distr(GET_RANDOM_GENERATOR); } void RStrip(char* str); diff --git a/similarity_search/src/randproj_util.cc b/similarity_search/src/randproj_util.cc index 1f99e2b..a01dcf7 100644 --- a/similarity_search/src/randproj_util.cc +++ b/similarity_search/src/randproj_util.cc @@ -23,6 +23,7 @@ #include "randproj_util.h" #include "distcomp.h" #include "logging.h" +#include "utils.h" namespace similarity { @@ -32,8 +33,7 @@ template void initRandProj(size_t nSrcDim, size_t nDstDim, bool bDoOrth, vector>& projMatr) { // Static is thread-safe in C++-11 - static std::random_device rd; - static std::mt19937 engine(rd()); + static thread_local auto& randGen(GET_RANDOM_GENERATOR); static std::normal_distribution<> normGen(0.0f, 1.0f); // 1. Create normally distributed vectors @@ -41,7 +41,7 @@ template void initRandProj(size_t nSrcDim, size_t nDstDim, for (size_t i = 0; i < nDstDim; ++i) { projMatr[i].resize(nSrcDim); for (size_t j = 0; j < nSrcDim; ++j) - projMatr[i][j] = normGen(engine); + projMatr[i][j] = normGen(randGen); } /* * 2. If bDoOrth == true, normalize the basis using the numerically stable diff --git a/similarity_search/src/searchoracle.cc b/similarity_search/src/searchoracle.cc index 12e7c87..4dd4435 100644 --- a/similarity_search/src/searchoracle.cc +++ b/similarity_search/src/searchoracle.cc @@ -205,8 +205,8 @@ void PolynomialPruner::SetIndexTimeParams(AnyParamManager& pmgr) { float recall = 0, time_best = 0, impr_best = -1, alpha_left = 0, alpha_right = 0; unsigned exp_left = 0, exp_right = 0; - static std::random_device rd; - static std::mt19937 engine(rd()); + + static thread_local auto& randGen(GET_RANDOM_GENERATOR); static std::normal_distribution<> normGen(0.0f, log(fullFactor)); @@ -218,8 +218,8 @@ void PolynomialPruner::SetIndexTimeParams(AnyParamManager& pmgr) { if (k > 0) { // Let's do some random normal fun - alpha_left_loc = exp(normGen(engine)); - alpha_right_loc = exp(normGen(engine)); + alpha_left_loc = exp(normGen(randGen)); + alpha_right_loc = exp(normGen(randGen)); LOG(LIB_INFO) << " RANDOM STARTING POINTS: " << alpha_left_loc << " " << alpha_right_loc; }