// -*- C++ -*-
// Copyright (C) 2007-2015 Free Software Foundation, Inc.
//
// This file is part of the GNU ISO C++ Library. This library is free
// software; you can redistribute it and/or modify it under the terms
// of the GNU General Public License as published by the Free Software
// Foundation; either version 3, or (at your option) any later
// version.
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// General Public License for more details.
// Under Section 7 of GPL version 3, you are granted additional
// permissions described in the GCC Runtime Library Exception, version
// 3.1, as published by the Free Software Foundation.
// You should have received a copy of the GNU General Public License and
// a copy of the GCC Runtime Library Exception along with this program;
// see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
// <http://www.gnu.org/licenses/>.
/** @file parallel/random_shuffle.h
* @brief Parallel implementation of std::random_shuffle().
* This file is a GNU parallel extension to the Standard C++ Library.
*/
// Written by Johannes Singler.
#ifndef _GLIBCXX_PARALLEL_RANDOM_SHUFFLE_H
#define _GLIBCXX_PARALLEL_RANDOM_SHUFFLE_H 1
#include <limits>
#include <bits/stl_numeric.h>
#include <parallel/parallel.h>
#include <parallel/random_number.h>
namespace __gnu_parallel
{
/** @brief Type to hold the index of a bin.
*
* Since many variables of this type are allocated, it should be
* chosen as small as possible.
*/
typedef unsigned short _BinIndex;
/** @brief Data known to every thread participating in
__gnu_parallel::__parallel_random_shuffle(). */
template<typename _RAIter>
struct _DRandomShufflingGlobalData
{
typedef std::iterator_traits<_RAIter> _TraitsType;
typedef typename _TraitsType::value_type _ValueType;
typedef typename _TraitsType::difference_type _DifferenceType;
/** @brief Begin iterator of the __source. */
_RAIter& _M_source;
/** @brief Temporary arrays for each thread. */
_ValueType** _M_temporaries;
/** @brief Two-dimensional array to hold the thread-bin distribution.
*
* Dimensions (_M_num_threads + 1) __x (_M_num_bins + 1). */
_DifferenceType** _M_dist;
/** @brief Start indexes of the threads' __chunks. */
_DifferenceType* _M_starts;
/** @brief Number of the thread that will further process the
corresponding bin. */
_ThreadIndex* _M_bin_proc;
/** @brief Number of bins to distribute to. */
int _M_num_bins;
/** @brief Number of bits needed to address the bins. */
int _M_num_bits;
/** @brief Constructor. */
_DRandomShufflingGlobalData(_RAIter& __source)
: _M_source(__source) { }
};
/** @brief Local data for a thread participating in
__gnu_parallel::__parallel_random_shuffle().
*/
template<typename _RAIter, typename _RandomNumberGenerator>
struct _DRSSorterPU
{
/** @brief Number of threads participating in total. */
int _M_num_threads;
/** @brief Begin index for bins taken care of by this thread. */
_BinIndex _M_bins_begin;
/** @brief End index for bins taken care of by this thread. */
_BinIndex __bins_end;
/** @brief Random _M_seed for this thread. */
uint32_t _M_seed;
/** @brief Pointer to global data. */
_DRandomShufflingGlobalData<_RAIter>* _M_sd;
};
/** @brief Generate a random number in @c [0,2^__logp).
* @param __logp Logarithm (basis 2) of the upper range __bound.
* @param __rng Random number generator to use.
*/
template<typename _RandomNumberGenerator>
inline int
__random_number_pow2(int __logp, _RandomNumberGenerator& __rng)
{ return __rng.__genrand_bits(__logp); }
/** @brief Random shuffle code executed by each thread.
* @param __pus Array of thread-local data records. */
template<typename _RAIter, typename _RandomNumberGenerator>
void
__parallel_random_shuffle_drs_pu(_DRSSorterPU<_RAIter,
_RandomNumberGenerator>* __pus)
{
typedef std::iterator_traits<_RAIter> _TraitsType;
typedef typename _TraitsType::value_type _ValueType;
typedef typename _TraitsType::difference_type _DifferenceType;
_ThreadIndex __iam = omp_get_thread_num();
_DRSSorterPU<_RAIter, _RandomNumberGenerator>* __d = &__pus[__iam];
_DRandomShufflingGlobalData<_RAIter>* __sd = __d->_M_sd;
// Indexing: _M_dist[bin][processor]
_DifferenceType __length = (__sd->_M_starts[__iam + 1]
- __sd->_M_starts[__iam]);
_BinIndex* __oracles = new _BinIndex[__length];
_DifferenceType* __dist = new _DifferenceType[__sd->_M_num_bins + 1];
_BinIndex* __bin_proc = new _BinIndex[__sd->_M_num_bins];
_ValueType** __temporaries = new _ValueType*[__d->_M_num_threads];
// Compute oracles and count appearances.
for (_BinIndex __b = 0; __b < __sd->_M_num_bins + 1; ++__b)
__dist[__b] = 0;
int __num_bits = __sd->_M_num_bits;
_RandomNumber __rng(__d->_M_seed);
// First main loop.
for (_DifferenceType __i = 0; __i < __length; ++__i)
{
_BinIndex __oracle = __random_number_pow2(__num_bits, __rng);
__oracles[__i] = __oracle;
// To allow prefix (partial) sum.
++(__dist[__oracle + 1]);
}
for (_BinIndex __b = 0; __b < __sd->_M_num_bins + 1; ++__b)
__sd->_M_dist[__b][__iam + 1] = __dist[__b];
# pragma omp barrier
# pragma omp single
{
// Sum up bins, __sd->_M_dist[__s + 1][__d->_M_num_threads] now
// contains the total number of items in bin __s
for (_BinIndex __s = 0; __s < __sd->_M_num_bins; ++__s)
__gnu_sequential::partial_sum(__sd->_M_dist[__s + 1],
__sd->_M_dist[__s + 1]
+ __d->_M_num_threads + 1,
__sd->_M_dist[__s + 1]);
}
# pragma omp barrier
_SequenceIndex __offset = 0, __global_offset = 0;
for (_BinIndex __s = 0; __s < __d->_M_bins_begin; ++__s)
__global_offset += __sd->_M_dist[__s + 1][__d->_M_num_threads];
# pragma omp barrier
for (_BinIndex __s = __d->_M_bins_begin; __s < __d->__bins_end; ++__s)
{
for (int __t = 0; __t < __d->_M_num_threads + 1; ++__t)
__sd->_M_dist[__s + 1][__t] += __offset;
__offset = __sd->_M_dist[__s + 1][__d->_M_num_threads];
}
__sd->_M_temporaries[__iam] = static_cast<_ValueType*>
(::operator new(sizeof(_ValueType) * __offset));
# pragma omp barrier
// Draw local copies to avoid false sharing.
for (_BinIndex __b = 0; __b < __sd->_M_num_bins + 1; ++__b)
__dist[__b] = __sd->_M_dist[__b][__iam];
for (_BinIndex __b = 0; __b < __sd->_M_num_bins; ++__b)
__bin_proc[__b] = __sd->_M_bin_proc[__b];
for (_ThreadIndex __t = 0; __t < __d->_M_num_threads; ++__t)
__temporaries[__t] = __sd->_M_temporaries[__t];
_RAIter __source = __sd->_M_source;
_DifferenceType __start = __sd->_M_starts[__iam];
// Distribute according to oracles, second main loop.
for (_DifferenceType __i = 0; __i < __length; ++__i)
{
_BinIndex __target_bin = __oracles[__i];
_ThreadIndex __target_p = __bin_proc[__target_bin];
// Last column [__d->_M_num_threads] stays unchanged.
::new(&(__temporaries[__target_p][__dist[__target_bin + 1]++]))
_ValueType(*(__source + __i + __start));
}
delete[] __oracles;
delete[] __dist;
delete[] __bin_proc;
delete[] __temporaries;
# pragma omp barrier
// Shuffle bins internally.
for (_BinIndex __b = __d->_M_bins_begin; __b < __d->__bins_end; ++__b)
{
_ValueType* __begin =
(__sd->_M_temporaries[__iam]
+ (__b == __d->_M_bins_begin
? 0 : __sd->_M_dist[__b][__d->_M_num_threads])),
*__end = (__sd->_M_temporaries[__iam]
+ __sd->_M_dist[__b + 1][__d->_M_num_threads]);
__sequential_random_shuffle(__begin, __end, __rng);
std::copy(__begin, __end, __sd->_M_source + __global_offset
+ (__b == __d->_M_bins_begin
? 0 : __sd->_M_dist[__b][__d->_M_num_threads]));
}
for (_SequenceIndex __i = 0; __i < __offset; ++__i)
__sd->_M_temporaries[__iam][__i].~_ValueType();
::operator delete(__sd->_M_temporaries[__iam]);
}
/** @brief Round up to the next greater power of 2.
* @param __x _Integer to round up */
template<typename _Tp>
_Tp
__round_up_to_pow2(_Tp __x)
{
if (__x <= 1)
return 1;
else
return (_Tp)1 << (__rd_log2(__x - 1) + 1);
}
/** @brief Main parallel random shuffle step.
* @param __begin Begin iterator of sequence.
* @param __end End iterator of sequence.
* @param __n Length of sequence.
* @param __num_threads Number of threads to use.
* @param __rng Random number generator to use.
*/
template<typename _RAIter, typename _RandomNumberGenerator>
void
__parallel_random_shuffle_drs(_RAIter __begin, _RAIter __end,
typename std::iterator_traits
<_RAIter>::difference_type __n,
_ThreadIndex __num_threads,
_RandomNumberGenerator& __rng)
{
typedef std::iterator_traits<_RAIter> _TraitsType;
typedef typename _TraitsType::value_type _ValueType;
typedef typename _TraitsType::difference_type _DifferenceType;
_GLIBCXX_CALL(__n)
const _Settings& __s = _Settings::get();
if (__num_threads > __n)
__num_threads = static_cast<_ThreadIndex>(__n);
_BinIndex __num_bins, __num_bins_cache;
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_L1
// Try the L1 cache first.
// Must fit into L1.
__num_bins_cache =
std::max<_DifferenceType>(1, __n / (__s.L1_cache_size_lb
/ sizeof(_ValueType)));
__num_bins_cache = __round_up_to_pow2(__num_bins_cache);
// No more buckets than TLB entries, power of 2
// Power of 2 and at least one element per bin, at most the TLB size.
__num_bins = std::min<_DifferenceType>(__n, __num_bins_cache);
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_TLB
// 2 TLB entries needed per bin.
__num_bins = std::min<_DifferenceType>(__s.TLB_size / 2, __num_bins);
#endif
__num_bins = __round_up_to_pow2(__num_bins);
if (__num_bins < __num_bins_cache)
{
#endif
// Now try the L2 cache
// Must fit into L2
__num_bins_cache = static_cast<_BinIndex>
(std::max<_DifferenceType>(1, __n / (__s.L2_cache_size
/ sizeof(_ValueType))));
__num_bins_cache = __round_up_to_pow2(__num_bins_cache);
// No more buckets than TLB entries, power of 2.
__num_bins = static_cast<_BinIndex>
(std::min(__n, static_cast<_DifferenceType>(__num_bins_cache)));
// Power of 2 and at least one element per bin, at most the TLB size.
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_TLB
// 2 TLB entries needed per bin.
__num_bins = std::min(static_cast<_DifferenceType>(__s.TLB_size / 2),
__num_bins);
#endif
__num_bins = __round_up_to_pow2(__num_bins);
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_L1
}
#endif
__num_bins = __round_up_to_pow2(
std::max<_BinIndex>(__num_threads, __num_bins));
if (__num_threads <= 1)
{
_RandomNumber __derived_rng(
__rng(std::numeric_limits<uint32_t>::max()));
__sequential_random_shuffle(__begin, __end, __derived_rng);
return;
}
_DRandomShufflingGlobalData<_RAIter> __sd(__begin);
_DRSSorterPU<_RAIter, _RandomNumber >* __pus;
_DifferenceType* __starts;
# pragma omp parallel num_threads(__num_threads)
{
_ThreadIndex __num_threads = omp_get_num_threads();
# pragma omp single
{
__pus = new _DRSSorterPU<_RAIter, _RandomNumber>[__num_threads];
__sd._M_temporaries = new _ValueType*[__num_threads];
__sd._M_dist = new _DifferenceType*[__num_bins + 1];
__sd._M_bin_proc = new _ThreadIndex[__num_bins];
for (_BinIndex __b = 0; __b < __num_bins + 1; ++__b)
__sd._M_dist[__b] = new _DifferenceType[__num_threads + 1];
for (_BinIndex __b = 0; __b < (__num_bins + 1); ++__b)
{
__sd._M_dist[0][0] = 0;
__sd._M_dist[__b][0] = 0;
}
__starts = __sd._M_starts = new _DifferenceType[__num_threads + 1];
int __bin_cursor = 0;
__sd._M_num_bins = __num_bins;
__sd._M_num_bits = __rd_log2(__num_bins);
_DifferenceType __chunk_length = __n / __num_threads,
__split = __n % __num_threads,
__start = 0;
_DifferenceType __bin_chunk_length = __num_bins / __num_threads,
__bin_split = __num_bins % __num_threads;
for (_ThreadIndex __i = 0; __i < __num_threads; ++__i)
{
__starts[__i] = __start;
__start += (__i < __split
? (__chunk_length + 1) : __chunk_length);
int __j = __pus[__i]._M_bins_begin = __bin_cursor;
// Range of bins for this processor.
__bin_cursor += (__i < __bin_split
? (__bin_chunk_length + 1)
: __bin_chunk_length);
__pus[__i].__bins_end = __bin_cursor;
for (; __j < __bin_cursor; ++__j)
__sd._M_bin_proc[__j] = __i;
__pus[__i]._M_num_threads = __num_threads;
__pus[__i]._M_seed = __rng(std::numeric_limits<uint32_t>::max());
__pus[__i]._M_sd = &__sd;
}
__starts[__num_threads] = __start;
} //single
// Now shuffle in parallel.
__parallel_random_shuffle_drs_pu(__pus);
} // parallel
delete[] __starts;
delete[] __sd._M_bin_proc;
for (int __s = 0; __s < (__num_bins + 1); ++__s)
delete[] __sd._M_dist[__s];
delete[] __sd._M_dist;
delete[] __sd._M_temporaries;
delete[] __pus;
}
/** @brief Sequential cache-efficient random shuffle.
* @param __begin Begin iterator of sequence.
* @param __end End iterator of sequence.
* @param __rng Random number generator to use.
*/
template<typename _RAIter, typename _RandomNumberGenerator>
void
__sequential_random_shuffle(_RAIter __begin, _RAIter __end,
_RandomNumberGenerator& __rng)
{
typedef std::iterator_traits<_RAIter> _TraitsType;
typedef typename _TraitsType::value_type _ValueType;
typedef typename _TraitsType::difference_type _DifferenceType;
_DifferenceType __n = __end - __begin;
const _Settings& __s = _Settings::get();
_BinIndex __num_bins, __num_bins_cache;
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_L1
// Try the L1 cache first, must fit into L1.
__num_bins_cache = std::max<_DifferenceType>
(1, __n / (__s.L1_cache_size_lb / sizeof(_ValueType)));
__num_bins_cache = __round_up_to_pow2(__num_bins_cache);
// No more buckets than TLB entries, power of 2
// Power of 2 and at least one element per bin, at most the TLB size
__num_bins = std::min(__n, (_DifferenceType)__num_bins_cache);
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_TLB
// 2 TLB entries needed per bin
__num_bins = std::min((_DifferenceType)__s.TLB_size / 2, __num_bins);
#endif
__num_bins = __round_up_to_pow2(__num_bins);
if (__num_bins < __num_bins_cache)
{
#endif
// Now try the L2 cache, must fit into L2.
__num_bins_cache = static_cast<_BinIndex>
(std::max<_DifferenceType>(1, __n / (__s.L2_cache_size
/ sizeof(_ValueType))));
__num_bins_cache = __round_up_to_pow2(__num_bins_cache);
// No more buckets than TLB entries, power of 2
// Power of 2 and at least one element per bin, at most the TLB size.
__num_bins = static_cast<_BinIndex>
(std::min(__n, static_cast<_DifferenceType>(__num_bins_cache)));
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_TLB
// 2 TLB entries needed per bin
__num_bins = std::min<_DifferenceType>(__s.TLB_size / 2, __num_bins);
#endif
__num_bins = __round_up_to_pow2(__num_bins);
#if _GLIBCXX_RANDOM_SHUFFLE_CONSIDER_L1
}
#endif
int __num_bits = __rd_log2(__num_bins);
if (__num_bins > 1)
{
_ValueType* __target =
static_cast<_ValueType*>(::operator new(sizeof(_ValueType) * __n));
_BinIndex* __oracles = new _BinIndex[__n];
_DifferenceType* __dist0 = new _DifferenceType[__num_bins + 1],
* __dist1 = new _DifferenceType[__num_bins + 1];
for (int __b = 0; __b < __num_bins + 1; ++__b)
__dist0[__b] = 0;
_RandomNumber __bitrng(__rng(0xFFFFFFFF));
for (_DifferenceType __i = 0; __i < __n; ++__i)
{
_BinIndex __oracle = __random_number_pow2(__num_bits, __bitrng);
__oracles[__i] = __oracle;
// To allow prefix (partial) sum.
++(__dist0[__oracle + 1]);
}
// Sum up bins.
__gnu_sequential::partial_sum(__dist0, __dist0 + __num_bins + 1,
__dist0);
for (int __b = 0; __b < __num_bins + 1; ++__b)
__dist1[__b] = __dist0[__b];
// Distribute according to oracles.
for (_DifferenceType __i = 0; __i < __n; ++__i)
::new(&(__target[(__dist0[__oracles[__i]])++]))
_ValueType(*(__begin + __i));
for (int __b = 0; __b < __num_bins; ++__b)
__sequential_random_shuffle(__target + __dist1[__b],
__target + __dist1[__b + 1], __rng);
// Copy elements back.
std::copy(__target, __target + __n, __begin);
delete[] __dist0;
delete[] __dist1;
delete[] __oracles;
for (_DifferenceType __i = 0; __i < __n; ++__i)
__target[__i].~_ValueType();
::operator delete(__target);
}
else
__gnu_sequential::random_shuffle(__begin, __end, __rng);
}
/** @brief Parallel random public call.
* @param __begin Begin iterator of sequence.
* @param __end End iterator of sequence.
* @param __rng Random number generator to use.
*/
template<typename _RAIter, typename _RandomNumberGenerator>
inline void
__parallel_random_shuffle(_RAIter __begin, _RAIter __end,
_RandomNumberGenerator __rng = _RandomNumber())
{
typedef std::iterator_traits<_RAIter> _TraitsType;
typedef typename _TraitsType::difference_type _DifferenceType;
_DifferenceType __n = __end - __begin;
__parallel_random_shuffle_drs(__begin, __end, __n,
__get_max_threads(), __rng);
}
}
#endif /* _GLIBCXX_PARALLEL_RANDOM_SHUFFLE_H */