3 // Copyright (C) 2007, 2008, 2009, 2010 Free Software Foundation, Inc.
5 // This file is part of the GNU ISO C++ Library. This library is free
6 // software; you can redistribute it and/or modify it under the terms
7 // of the GNU General Public License as published by the Free Software
8 // Foundation; either version 3, or (at your option) any later
11 // This library is distributed in the hope that it will be useful, but
12 // WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // General Public License for more details.
16 // Under Section 7 of GPL version 3, you are granted additional
17 // permissions described in the GCC Runtime Library Exception, version
18 // 3.1, as published by the Free Software Foundation.
20 // You should have received a copy of the GNU General Public License and
21 // a copy of the GCC Runtime Library Exception along with this program;
22 // see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
23 // <http://www.gnu.org/licenses/>.
25 /** @file parallel/multiseq_selection.h
26 * @brief Functions to find elements of a certain global __rank in
27 * multiple sorted sequences. Also serves for splitting such
30 * The algorithm description can be found in
32 * P. J. Varman, S. D. Scheufler, B. R. Iyer, and G. R. Ricard.
33 * Merging Multiple Lists on Hierarchical-Memory Multiprocessors.
34 * Journal of Parallel and Distributed Computing, 12(2):171–177, 1991.
36 * This file is a GNU parallel extension to the Standard C++ Library.
39 // Written by Johannes Singler.
41 #ifndef _GLIBCXX_PARALLEL_MULTISEQ_SELECTION_H
42 #define _GLIBCXX_PARALLEL_MULTISEQ_SELECTION_H 1
47 #include <bits/stl_algo.h>
49 #include <parallel/sort.h>
51 namespace __gnu_parallel
53 /** @brief Compare __a pair of types lexicographically, ascending. */
54 template<typename _T1, typename _T2, typename _Compare>
56 : public std::binary_function<std::pair<_T1, _T2>,
57 std::pair<_T1, _T2>, bool>
63 _Lexicographic(_Compare& __comp) : _M_comp(__comp) { }
66 operator()(const std::pair<_T1, _T2>& __p1,
67 const std::pair<_T1, _T2>& __p2) const
69 if (_M_comp(__p1.first, __p2.first))
72 if (_M_comp(__p2.first, __p1.first))
76 return __p1.second < __p2.second;
80 /** @brief Compare __a pair of types lexicographically, descending. */
81 template<typename _T1, typename _T2, typename _Compare>
82 class _LexicographicReverse : public std::binary_function<_T1, _T2, bool>
88 _LexicographicReverse(_Compare& __comp) : _M_comp(__comp) { }
91 operator()(const std::pair<_T1, _T2>& __p1,
92 const std::pair<_T1, _T2>& __p2) const
94 if (_M_comp(__p2.first, __p1.first))
97 if (_M_comp(__p1.first, __p2.first))
101 return __p2.second < __p1.second;
106 * @brief Splits several sorted sequences at a certain global __rank,
107 * resulting in a splitting point for each sequence.
108 * The sequences are passed via a sequence of random-access
109 * iterator pairs, none of the sequences may be empty. If there
110 * are several equal elements across the split, the ones on the
111 * __left side will be chosen from sequences with smaller number.
112 * @param __begin_seqs Begin of the sequence of iterator pairs.
113 * @param __end_seqs End of the sequence of iterator pairs.
114 * @param __rank The global rank to partition at.
115 * @param __begin_offsets A random-access __sequence __begin where the
116 * __result will be stored in. Each element of the sequence is an
117 * iterator that points to the first element on the greater part of
118 * the respective __sequence.
119 * @param __comp The ordering functor, defaults to std::less<_Tp>.
121 template<typename _RanSeqs, typename _RankType, typename _RankIterator,
124 multiseq_partition(_RanSeqs __begin_seqs, _RanSeqs __end_seqs,
126 _RankIterator __begin_offsets,
127 _Compare __comp = std::less<
128 typename std::iterator_traits<typename
129 std::iterator_traits<_RanSeqs>::value_type::
130 first_type>::value_type>()) // std::less<_Tp>
132 _GLIBCXX_CALL(__end_seqs - __begin_seqs)
134 typedef typename std::iterator_traits<_RanSeqs>::value_type::first_type
136 typedef typename std::iterator_traits<_RanSeqs>::difference_type
138 typedef typename std::iterator_traits<_It>::difference_type
140 typedef typename std::iterator_traits<_It>::value_type _ValueType;
142 _Lexicographic<_ValueType, _SeqNumber, _Compare> __lcomp(__comp);
143 _LexicographicReverse<_ValueType, _SeqNumber, _Compare> __lrcomp(__comp);
145 // Number of sequences, number of elements in total (possibly
146 // including padding).
147 _DifferenceType __m = std::distance(__begin_seqs, __end_seqs), __nn = 0,
150 for (_SeqNumber __i = 0; __i < __m; __i++)
152 __nn += std::distance(__begin_seqs[__i].first,
153 __begin_seqs[__i].second);
154 _GLIBCXX_PARALLEL_ASSERT(
155 std::distance(__begin_seqs[__i].first,
156 __begin_seqs[__i].second) > 0);
161 for (_SeqNumber __i = 0; __i < __m; __i++)
162 __begin_offsets[__i] = __begin_seqs[__i].second; // Very end.
167 _GLIBCXX_PARALLEL_ASSERT(__m != 0);
168 _GLIBCXX_PARALLEL_ASSERT(__nn != 0);
169 _GLIBCXX_PARALLEL_ASSERT(__rank >= 0);
170 _GLIBCXX_PARALLEL_ASSERT(__rank < __nn);
172 _DifferenceType* __ns = new _DifferenceType[__m];
173 _DifferenceType* __a = new _DifferenceType[__m];
174 _DifferenceType* __b = new _DifferenceType[__m];
177 __ns[0] = std::distance(__begin_seqs[0].first, __begin_seqs[0].second);
179 for (_SeqNumber __i = 0; __i < __m; __i++)
181 __ns[__i] = std::distance(__begin_seqs[__i].first,
182 __begin_seqs[__i].second);
183 __nmax = std::max(__nmax, __ns[__i]);
186 __r = __rd_log2(__nmax) + 1;
188 // Pad all lists to this length, at least as long as any ns[__i],
189 // equality iff __nmax = 2^__k - 1.
190 __l = (1ULL << __r) - 1;
192 for (_SeqNumber __i = 0; __i < __m; __i++)
200 // 0 <= __a[__i] <= __ns[__i], 0 <= __b[__i] <= __l
202 #define __S(__i) (__begin_seqs[__i].first)
204 // Initial partition.
205 std::vector<std::pair<_ValueType, _SeqNumber> > __sample;
207 for (_SeqNumber __i = 0; __i < __m; __i++)
208 if (__n < __ns[__i]) //__sequence long enough
209 __sample.push_back(std::make_pair(__S(__i)[__n], __i));
210 __gnu_sequential::sort(__sample.begin(), __sample.end(), __lcomp);
212 for (_SeqNumber __i = 0; __i < __m; __i++) //conceptual infinity
213 if (__n >= __ns[__i]) //__sequence too short, conceptual infinity
215 std::make_pair(__S(__i)[0] /*__dummy element*/, __i));
217 _DifferenceType __localrank = __rank / __l;
221 __j < __localrank && ((__n + 1) <= __ns[__sample[__j].second]);
223 __a[__sample[__j].second] += __n + 1;
224 for (; __j < __m; __j++)
225 __b[__sample[__j].second] -= __n + 1;
227 // Further refinement.
232 _SeqNumber __lmax_seq = -1; // to avoid warning
233 const _ValueType* __lmax = 0; // impossible to avoid the warning?
234 for (_SeqNumber __i = 0; __i < __m; __i++)
240 __lmax = &(__S(__i)[__a[__i] - 1]);
245 // Max, favor rear sequences.
246 if (!__comp(__S(__i)[__a[__i] - 1], *__lmax))
248 __lmax = &(__S(__i)[__a[__i] - 1]);
256 for (__i = 0; __i < __m; __i++)
258 _DifferenceType __middle = (__b[__i] + __a[__i]) / 2;
259 if (__lmax && __middle < __ns[__i] &&
260 __lcomp(std::make_pair(__S(__i)[__middle], __i),
261 std::make_pair(*__lmax, __lmax_seq)))
262 __a[__i] = std::min(__a[__i] + __n + 1, __ns[__i]);
267 _DifferenceType __leftsize = 0;
268 for (_SeqNumber __i = 0; __i < __m; __i++)
269 __leftsize += __a[__i] / (__n + 1);
271 _DifferenceType __skew = __rank / (__n + 1) - __leftsize;
275 // Move to the left, find smallest.
276 std::priority_queue<std::pair<_ValueType, _SeqNumber>,
277 std::vector<std::pair<_ValueType, _SeqNumber> >,
278 _LexicographicReverse<_ValueType, _SeqNumber, _Compare> >
281 for (_SeqNumber __i = 0; __i < __m; __i++)
282 if (__b[__i] < __ns[__i])
283 __pq.push(std::make_pair(__S(__i)[__b[__i]], __i));
285 for (; __skew != 0 && !__pq.empty(); --__skew)
287 _SeqNumber __source = __pq.top().second;
291 = std::min(__a[__source] + __n + 1, __ns[__source]);
292 __b[__source] += __n + 1;
294 if (__b[__source] < __ns[__source])
296 std::make_pair(__S(__source)[__b[__source]], __source));
301 // Move to the right, find greatest.
302 std::priority_queue<std::pair<_ValueType, _SeqNumber>,
303 std::vector<std::pair<_ValueType, _SeqNumber> >,
304 _Lexicographic<_ValueType, _SeqNumber, _Compare> >
307 for (_SeqNumber __i = 0; __i < __m; __i++)
309 __pq.push(std::make_pair(__S(__i)[__a[__i] - 1], __i));
311 for (; __skew != 0; ++__skew)
313 _SeqNumber __source = __pq.top().second;
316 __a[__source] -= __n + 1;
317 __b[__source] -= __n + 1;
319 if (__a[__source] > 0)
320 __pq.push(std::make_pair(
321 __S(__source)[__a[__source] - 1], __source));
327 // __a[__i] == __b[__i] in most cases, except when __a[__i] has been
328 // clamped because of having reached the boundary
330 // Now return the result, calculate the offset.
332 // Compare the keys on both edges of the border.
334 // Maximum of left edge, minimum of right edge.
335 _ValueType* __maxleft = 0;
336 _ValueType* __minright = 0;
337 for (_SeqNumber __i = 0; __i < __m; __i++)
342 __maxleft = &(__S(__i)[__a[__i] - 1]);
345 // Max, favor rear sequences.
346 if (!__comp(__S(__i)[__a[__i] - 1], *__maxleft))
347 __maxleft = &(__S(__i)[__a[__i] - 1]);
350 if (__b[__i] < __ns[__i])
353 __minright = &(__S(__i)[__b[__i]]);
356 // Min, favor fore sequences.
357 if (__comp(__S(__i)[__b[__i]], *__minright))
358 __minright = &(__S(__i)[__b[__i]]);
363 _SeqNumber __seq = 0;
364 for (_SeqNumber __i = 0; __i < __m; __i++)
365 __begin_offsets[__i] = __S(__i) + __a[__i];
374 * @brief Selects the element at a certain global __rank from several
377 * The sequences are passed via a sequence of random-access
378 * iterator pairs, none of the sequences may be empty.
379 * @param __begin_seqs Begin of the sequence of iterator pairs.
380 * @param __end_seqs End of the sequence of iterator pairs.
381 * @param __rank The global rank to partition at.
382 * @param __offset The rank of the selected element in the global
383 * subsequence of elements equal to the selected element. If the
384 * selected element is unique, this number is 0.
385 * @param __comp The ordering functor, defaults to std::less.
387 template<typename _Tp, typename _RanSeqs, typename _RankType,
390 multiseq_selection(_RanSeqs __begin_seqs, _RanSeqs __end_seqs,
392 _RankType& __offset, _Compare __comp = std::less<_Tp>())
394 _GLIBCXX_CALL(__end_seqs - __begin_seqs)
396 typedef typename std::iterator_traits<_RanSeqs>::value_type::first_type
398 typedef typename std::iterator_traits<_RanSeqs>::difference_type
400 typedef typename std::iterator_traits<_It>::difference_type
403 _Lexicographic<_Tp, _SeqNumber, _Compare> __lcomp(__comp);
404 _LexicographicReverse<_Tp, _SeqNumber, _Compare> __lrcomp(__comp);
406 // Number of sequences, number of elements in total (possibly
407 // including padding).
408 _DifferenceType __m = std::distance(__begin_seqs, __end_seqs);
409 _DifferenceType __nn = 0;
410 _DifferenceType __nmax, __n, __r;
412 for (_SeqNumber __i = 0; __i < __m; __i++)
413 __nn += std::distance(__begin_seqs[__i].first,
414 __begin_seqs[__i].second);
416 if (__m == 0 || __nn == 0 || __rank < 0 || __rank >= __nn)
418 // result undefined if there is no data or __rank is outside bounds
419 throw std::exception();
423 _DifferenceType* __ns = new _DifferenceType[__m];
424 _DifferenceType* __a = new _DifferenceType[__m];
425 _DifferenceType* __b = new _DifferenceType[__m];
428 __ns[0] = std::distance(__begin_seqs[0].first, __begin_seqs[0].second);
430 for (_SeqNumber __i = 0; __i < __m; ++__i)
432 __ns[__i] = std::distance(__begin_seqs[__i].first,
433 __begin_seqs[__i].second);
434 __nmax = std::max(__nmax, __ns[__i]);
437 __r = __rd_log2(__nmax) + 1;
439 // Pad all lists to this length, at least as long as any ns[__i],
440 // equality iff __nmax = 2^__k - 1
441 __l = __round_up_to_pow2(__r) - 1;
443 for (_SeqNumber __i = 0; __i < __m; ++__i)
451 // 0 <= __a[__i] <= __ns[__i], 0 <= __b[__i] <= __l
453 #define __S(__i) (__begin_seqs[__i].first)
455 // Initial partition.
456 std::vector<std::pair<_Tp, _SeqNumber> > __sample;
458 for (_SeqNumber __i = 0; __i < __m; __i++)
460 __sample.push_back(std::make_pair(__S(__i)[__n], __i));
461 __gnu_sequential::sort(__sample.begin(), __sample.end(),
462 __lcomp, sequential_tag());
464 // Conceptual infinity.
465 for (_SeqNumber __i = 0; __i < __m; __i++)
466 if (__n >= __ns[__i])
468 std::make_pair(__S(__i)[0] /*__dummy element*/, __i));
470 _DifferenceType __localrank = __rank / __l;
474 __j < __localrank && ((__n + 1) <= __ns[__sample[__j].second]);
476 __a[__sample[__j].second] += __n + 1;
477 for (; __j < __m; ++__j)
478 __b[__sample[__j].second] -= __n + 1;
480 // Further refinement.
485 const _Tp* __lmax = 0;
486 for (_SeqNumber __i = 0; __i < __m; ++__i)
491 __lmax = &(__S(__i)[__a[__i] - 1]);
494 if (__comp(*__lmax, __S(__i)[__a[__i] - 1])) //max
495 __lmax = &(__S(__i)[__a[__i] - 1]);
501 for (__i = 0; __i < __m; __i++)
503 _DifferenceType __middle = (__b[__i] + __a[__i]) / 2;
504 if (__lmax && __middle < __ns[__i]
505 && __comp(__S(__i)[__middle], *__lmax))
506 __a[__i] = std::min(__a[__i] + __n + 1, __ns[__i]);
511 _DifferenceType __leftsize = 0;
512 for (_SeqNumber __i = 0; __i < __m; ++__i)
513 __leftsize += __a[__i] / (__n + 1);
515 _DifferenceType __skew = __rank / (__n + 1) - __leftsize;
519 // Move to the left, find smallest.
520 std::priority_queue<std::pair<_Tp, _SeqNumber>,
521 std::vector<std::pair<_Tp, _SeqNumber> >,
522 _LexicographicReverse<_Tp, _SeqNumber, _Compare> >
525 for (_SeqNumber __i = 0; __i < __m; ++__i)
526 if (__b[__i] < __ns[__i])
527 __pq.push(std::make_pair(__S(__i)[__b[__i]], __i));
529 for (; __skew != 0 && !__pq.empty(); --__skew)
531 _SeqNumber __source = __pq.top().second;
535 = std::min(__a[__source] + __n + 1, __ns[__source]);
536 __b[__source] += __n + 1;
538 if (__b[__source] < __ns[__source])
540 std::make_pair(__S(__source)[__b[__source]], __source));
545 // Move to the right, find greatest.
546 std::priority_queue<std::pair<_Tp, _SeqNumber>,
547 std::vector<std::pair<_Tp, _SeqNumber> >,
548 _Lexicographic<_Tp, _SeqNumber, _Compare> > __pq(__lcomp);
550 for (_SeqNumber __i = 0; __i < __m; ++__i)
552 __pq.push(std::make_pair(__S(__i)[__a[__i] - 1], __i));
554 for (; __skew != 0; ++__skew)
556 _SeqNumber __source = __pq.top().second;
559 __a[__source] -= __n + 1;
560 __b[__source] -= __n + 1;
562 if (__a[__source] > 0)
563 __pq.push(std::make_pair(
564 __S(__source)[__a[__source] - 1], __source));
570 // __a[__i] == __b[__i] in most cases, except when __a[__i] has been
571 // clamped because of having reached the boundary
573 // Now return the result, calculate the offset.
575 // Compare the keys on both edges of the border.
577 // Maximum of left edge, minimum of right edge.
578 bool __maxleftset = false, __minrightset = false;
580 // Impossible to avoid the warning?
581 _Tp __maxleft, __minright;
582 for (_SeqNumber __i = 0; __i < __m; ++__i)
588 __maxleft = __S(__i)[__a[__i] - 1];
594 if (__comp(__maxleft, __S(__i)[__a[__i] - 1]))
595 __maxleft = __S(__i)[__a[__i] - 1];
598 if (__b[__i] < __ns[__i])
602 __minright = __S(__i)[__b[__i]];
603 __minrightset = true;
608 if (__comp(__S(__i)[__b[__i]], __minright))
609 __minright = __S(__i)[__b[__i]];
614 // Minright is the __splitter, in any case.
616 if (!__maxleftset || __comp(__minright, __maxleft))
618 // Good luck, everything is split unambiguously.
623 // We have to calculate an offset.
626 for (_SeqNumber __i = 0; __i < __m; ++__i)
629 = std::lower_bound(__S(__i), __S(__i) + __ns[__i],
632 __offset += __a[__i] - lb;
646 #endif /* _GLIBCXX_PARALLEL_MULTISEQ_SELECTION_H */