Chombo + EB  3.0
sort_utils.H
Go to the documentation of this file.
1 /***********************************************************
2  * This file is part of DEMSort, a high-performance *
3  * distributed-memory external-memory sorting program *
4  * *
5  * Copyright (C) 2010 Mirko Rahn and Johannes Singler *
6  * Karlsruhe Institute of Technology (KIT), Germany *
7  * *
8  * All rights reserved. Free for academic use. *
9  ***********************************************************/
10 
11 #ifndef _SORT_UTILS_H_
12 #define _SORT_UTILS_H_
13 
14 // to mark unused parameters explicitly
15 #ifdef __GNUC__
16 #define PARAM_UNUSED(x) x __attribute__ ((unused))
17 #else
18 #define PARAM_UNUSED(x) x
19 #endif
20 
21 typedef unsigned long long key_type;
22 typedef long long index_type;
23 
24 template<class S, class Comparator>
26 {
27  S* data;
28  Comparator& less;
29 
30  IndirectLess(S* data, Comparator& less) : data(data), less(less)
31  {
32  }
33 
34  inline bool operator () (const index_type& a, const index_type& b)
35  {
36  return (less(data[a], data[b])
37  ? true
38  : (less(data[b], data[a])
39  ? false
40  : (a < b) //break ties on position in sequence
41  )
42  );
43  }
44 };
45 
46 template<class S, class Comparator>
48 {
49  S* data;
50  Comparator& less;
51 
52  IndirectInvertedLess(S* data, Comparator& less) : data(data), less(less)
53  {
54  }
55 
56  inline bool operator () (const index_type& a, const index_type& b)
57  {
58  return (less(data[b], data[a])
59  ? true
60  : (less(data[a], data[b])
61  ? false
62  : (b < a) //break ties on position in sequence
63  )
64  );
65  }
66 };
67 
68 template<class InIt, class Comparator>
70 {
71  InIt* data;
72  Comparator& less;
73 
74  IndirectIndirectLess(InIt* data, Comparator& less) : data(data), less(less)
75  {
76  }
77 
78  inline bool operator () (const unsigned& a, const unsigned& b)
79  {
80  return (less(*(data[a]), *(data[b]))
81  ? true
82  : (less(*(data[b]), *(data[a]))
83  ? false
84  : (a < b) //break ties on position in sequence
85  )
86  );
87  }
88 };
89 
90 template<class InIt, class Comparator>
92 {
93  InIt* data;
94  Comparator& less;
95 
96  IndirectIndirectInvertedLess(InIt* data, Comparator& less) : data(data), less(less)
97  {
98  }
99 
100  inline bool operator () (const unsigned& a, const unsigned& b)
101  {
102  return (less(*(data[b]), *(data[a]))
103  ? true
104  : (less(*(data[a]), *(data[b]))
105  ? false
106  : (b < a) //break ties on position in sequence
107  )
108  );
109  }
110 };
111 
112 #ifdef CH_MPI
113 
114 #define MPI_INDEX_TYPE MPI_LONG_LONG_INT
115 #define MPI_KEY_TYPE MPI_UNSIGNED_LONG_LONG
116 
117 // given by MPI_Standard
118 typedef int PE_NUM;
119 
120 //template <class T>
121 class SortComm
122 {
123  public:
124  SortComm( MPI_Comm a_comm ) : m_comm(a_comm)
125  {
126  MPI_Comm_size(a_comm,&m_numpes);
127  MPI_Comm_rank(a_comm,&m_mype);
128  }
129  MPI_Comm m_comm;
130  int m_mype;
131  int m_numpes;
132 };
133 
134 // ** print utils in sort_env.h ******************************************** //
135 
136 std::string sort_format_pe_number(int i)
137 {
138  std::ostringstream ss;
139 
140  ss << std::setw(4) << i;
141 
142  return ss.str();
143 }
144 
145 // ************************************************************************* //
146 
147 static std::string sort_format_iam(const SortComm* env)
148 {
149  return sort_format_pe_number(env->m_mype) + ": ";
150 }
151 
152 // ************************************************************************* //
153 
154 std::ostream& sort_pe_out(const SortComm* env)
155 {
156  return std::cout << sort_format_iam(env);
157 }
158 
159 std::ostream& sort_pe_err(const SortComm* env)
160 {
161  return std::cerr << sort_format_iam(env) << "ERROR: ";
162 }
163 
164 std::ostream& sort_nullstream() // eats everything, prints nothing
165 {
166  static std::ostream rc(NULL);
167 
168  return rc;
169 }
170 
171 std::ostream& sort_root_out(const SortComm* env)
172 {
173  return (env->m_mype == 0) ? sort_pe_out(env) : sort_nullstream();
174 }
175 
176 std::ostream& sort_root_err(const SortComm* env)
177 {
178  return (env->m_mype == 0) ? sort_pe_err(env) : sort_nullstream();
179 }
180 
181 // ************************************************************************* //
182 
183 void sort_pe_msg(const SortComm* env, const std::string msg)
184 {
185  sort_pe_out(env) << msg << std::endl;
186 }
187 
188 void sort_root_msg(const SortComm* env, const std::string msg)
189 {
190  sort_root_out(env) << msg << std::endl;
191 }
192 
193 PE_NUM pe_right(PE_NUM id, const SortComm* env)
194 {
195  return (id + 1) % env->m_numpes;
196 }
197 
198 PE_NUM pe_left(PE_NUM id, const SortComm* env)
199 {
200  return (id - 1 + env->m_numpes) % env->m_numpes;
201 }
202 
203 enum Tag
204 {
205 // in makeGlobalRuns, in case of block with filler
206  TAG_ALLTOALLV_MANUAL=10001
207 // in alltoall
208  , TAG_EXTERNAL_ALLTOALL=10002
209 
210 // for msr
211  , TAG_NEED=10003,
212  TAG_QUERY=10004,
213  TAG_GIVE=20005,
214  TAG_ANSWER=10006,
215  TAG_DONE=10007,
216  TAG_STOP=10008,
217  TAG_SHIFT_RIGHT_DONE=10009,
218  TAG_START_SHIFT_LEFT=10010
219 
220  , TAG_CHECK=10011
221 
222  , TAG_LEFT_MAX=10012,
223  TAG_RIGHT_MIN=10013
224 
225  , TAG_GRD=10014
226 
227  , TAG_PART_DIST=10015
228 };
229 
230 /* ************************************************************************* */
231 // functions to post-process the result of multiway_select
232 #include <limits.h>
233 struct MultiwaySelectResult
234 {
235  int* m_sdispl;
236  int* m_scount;
237  int* m_rcount;
238  int* m_rdispl;
239  index_type *m_lborder;
240 
241  MultiwaySelectResult(index_type* a_lborder,
242  int a_len,
243  SortComm* a_env
244  )
245  : m_lborder(a_lborder)
246  {
247  m_sdispl = new int[a_env->m_numpes + 1];
248  m_scount = new int[a_env->m_numpes];
249  m_rcount = new int[a_env->m_numpes];
250  m_rdispl = new int[a_env->m_numpes];
251  const int P = a_env->m_numpes;
252  {
253  int *lborder= new int[P];
254  for (int ii=0;ii<P; ii++)
255  {
256  CH_assert(a_lborder[ii]<INT_MAX);
257  lborder[ii] = a_lborder[ii];
258  }
259 
260  MPI_Alltoall( lborder, 1, MPI_INT,
261  m_sdispl, 1, MPI_INT, a_env->m_comm );
262 
263  delete [] lborder;
264  }
265  m_sdispl[P] = a_len;
266  //if (a_env->m_mype==0)std::cout << "\t[" << a_env->m_mype <<"] m_sdispl = " << m_sdispl[0] << " " << m_sdispl[1] << " " << m_sdispl[2] << " " << m_sdispl[3] << std::endl;
267 
268  for (int ii = 0; ii < P; ++ii)
269  {
270  m_scount[ii] = m_sdispl[ii + 1] - m_sdispl[ii];
271  }
272  //std::cout << "\t[" << a_env->m_mype <<"] m_scount = " << m_scount[0] << " " << m_scount[1] << " " << m_scount[2] << " " << m_scount[3] << std::endl;
273 
274  MPI_Alltoall(m_scount, 1, MPI_INT,
275  m_rcount, 1, MPI_INT, a_env->m_comm);
276  //std::cout << "\t[" << a_env->m_mype <<"] m_rcount = " << m_rcount[0] << " " << m_rcount[1] << " " << m_rcount[2] << " " << m_rcount[3] << std::endl;
277 
278  m_rdispl[0] = 0;
279  for (int ii = 1; ii < P; ++ii)
280  {
281  m_rdispl[ii] = m_rdispl[ii - 1] + m_rcount[ii - 1];
282  }
283  //std::cout << "\t[" << a_env->m_mype <<"] m_rdispl = " << m_rdispl[0] << " " << m_rdispl[1] << " " << m_rdispl[2] << " " << m_rdispl[3] << std::endl;
284  }
285 
286  // constructor with m_sdispl
287  MultiwaySelectResult(int* a_sdisp,
288  SortComm* env
289  )
290  {
291  m_scount = new int[env->m_numpes];
292  m_rcount = new int[env->m_numpes];
293  m_rdispl = new int[env->m_numpes];
294  m_sdispl = a_sdisp;
295 
296  // MPI_Alltoall(lborder, 1, MPI_INT,
297  // m_sdispl, 1, MPI_INT, env->m_comm);
298 
299  // m_sdispl[env->m_numpes] = a_len;
300  std::cout << "\t[" << env->m_mype <<"]MultiwaySelectResult m_sdispl = " << m_sdispl[0] << " " << m_sdispl[1] << " " << m_sdispl[2] << " " << m_sdispl[3] << " " << m_sdispl[4] << std::endl;
301 
302  for (int i = 0; i < env->m_numpes; ++i)
303  {
304  m_scount[i] = m_sdispl[i + 1] - m_sdispl[i];
305  }
306  //MPI_Barrier(env->m_comm);
307  //std::cout << "\t[" << env->m_mype <<"]MultiwaySelectResult m_scount = " << m_scount[0] << " " << m_scount[1] << " " << m_scount[2] << " " << m_scount[3] << std::endl;
308  MPI_Alltoall(m_scount, 1, MPI_INT,
309  m_rcount, 1, MPI_INT, env->m_comm);
310 
311  std::cout << "\t[" << env->m_mype <<"]MultiwaySelectResult m_rcount = " << m_rcount[0] << " " << m_rcount[1] << " " << m_rcount[2] << " " << m_rcount[3] << std::endl;
312  m_rdispl[0] = 0;
313 
314  for (int i = 1; i < env->m_numpes; ++i)
315  {
316  m_rdispl[i] = m_rdispl[i - 1] + m_rcount[i - 1];
317  }
318  MPI_Barrier(env->m_comm);
319  std::cout << "\t[" << env->m_mype <<"]MultiwaySelectResult m_rdispl = " << m_rdispl[0] << " " << m_rdispl[1] << " " << m_rdispl[2] << " " << m_rdispl[3] << std::endl;
320  }
321 
322  ~MultiwaySelectResult()
323  {
324  delete[] m_sdispl;
325  delete[] m_scount;
326  delete[] m_rcount;
327  delete[] m_rdispl;
328  }
329 };
330 
331 /* ************************************************************************* */
332 // verify
333 
334 template<typename block_type, typename value_type, class Comparator>
335 void test_multiway_select(value_type* data,
336  index_type length,
337  MultiwaySelectResult* msr,
338  Comparator less,
339  SortComm* env
340  )
341 {
342  unsigned err = 0;
343 
344  index_type* lborder_global = new index_type[env->m_numpes * env->m_numpes];
345  index_type* sdispl_global = new index_type[env->m_numpes * env->m_numpes];
346  index_type* rdispl_global = new index_type[env->m_numpes * env->m_numpes];
347  index_type* scount_global = new index_type[env->m_numpes * env->m_numpes];
348  index_type* rcount_global = new index_type[env->m_numpes * env->m_numpes];
349 
350  MPI_Allgather(msr->m_lborder, env->m_numpes, MPI_INDEX_TYPE,
351  lborder_global, env->m_numpes, MPI_INDEX_TYPE,
352  env->m_comm
353  );
354 
355  CH_assert(0); // this is wrong now -- int not index_type
356  MPI_Gather(msr->m_sdispl, env->m_numpes, MPI_INDEX_TYPE,
357  sdispl_global, env->m_numpes, MPI_INDEX_TYPE,
358  0, env->m_comm
359  );
360 
361  MPI_Gather(msr->m_rdispl, env->m_numpes, MPI_INDEX_TYPE,
362  rdispl_global, env->m_numpes, MPI_INDEX_TYPE,
363  0, env->m_comm
364  );
365 
366  MPI_Gather(msr->m_scount, env->m_numpes, MPI_INDEX_TYPE,
367  scount_global, env->m_numpes, MPI_INDEX_TYPE,
368  0, env->m_comm
369  );
370 
371  MPI_Gather(msr->m_rcount, env->m_numpes, MPI_INDEX_TYPE,
372  rcount_global, env->m_numpes, MPI_INDEX_TYPE,
373  0, env->m_comm
374  );
375 
376  index_type* rlen = new index_type[env->m_numpes];
377 
378  MPI_Gather(&length, 1, MPI_INDEX_TYPE,
379  rlen, 1, MPI_INDEX_TYPE,
380  0, env->m_comm
381  );
382 
383  if (env->m_mype == 0)
384  {
385  /* ******************************************************************* */
386 
387  std::cout << "lborder = " << std::endl;
388 
389  index_type rsum = 0;
390 
391  for (int r = 0; r < env->m_numpes; ++r)
392  {
393  index_type tsum = 0;
394 
395  for (int i = 0; i < env->m_numpes; ++i)
396  {
397  std::cout << std::setw(11) << lborder_global[r * env->m_numpes + i];
398 
399  tsum += lborder_global[r * env->m_numpes + i];
400  }
401 
402  std::cout << ", sum = " << tsum
403  << " is " << (tsum == rsum ? "okay" : "WRONG")
404  << std::endl;
405 
406  err += (tsum == rsum) ? 0 : 1;
407 
408  rsum += rlen[r];
409  }
410 
411  /* ******************************************************************* */
412 
413  std::cout << "scount = " << std::endl;
414 
415  for (int r = 0; r < env->m_numpes; ++r)
416  {
417  index_type tsum = 0;
418 
419  for (int i = 0; i < env->m_numpes; ++i)
420  {
421  std::cout << std::setw(11) << scount_global[r * env->m_numpes + i];
422 
423  tsum += scount_global[r * env->m_numpes + i];
424  }
425 
426  std::cout << ", sum = " << tsum
427  << " is " << (tsum == rlen[r] ? "okay" : "WRONG")
428  << std::endl;
429 
430  err += (tsum == rlen[r]) ? 0 : 1;
431  }
432 
433 
434  /* ******************************************************************* */
435 
436  std::cout << "rcount = " << std::endl;
437 
438  for (int r = 0; r < env->m_numpes; ++r)
439  {
440  index_type tsum = 0;
441 
442  for (int i = 0; i < env->m_numpes; ++i)
443  {
444  std::cout << std::setw(11) << rcount_global[r * env->m_numpes + i];
445 
446  tsum += rcount_global[r * env->m_numpes + i];
447  }
448 
449  std::cout << ", sum = " << tsum
450  << " is " << (tsum == rlen[r] ? "okay" : "WRONG")
451  << std::endl;
452 
453  err += (tsum == rlen[r]) ? 0 : 1;
454  }
455 
456  /* ******************************************************************* */
457 
458  std::cout << "sdispl = " << std::endl;
459 
460  for (int r = 0; r < env->m_numpes; ++r)
461  {
462  err += (sdispl_global[r * env->m_numpes + 0] == 0) ? 0 : 1;
463 
464  std::cout << std::setw(11) << sdispl_global[r * env->m_numpes + 0];
465 
466  bool asc = true;
467 
468  for (int i = 1; i < env->m_numpes; ++i)
469  {
470  std::cout << std::setw(11) << sdispl_global[r * env->m_numpes + i];
471 
472  asc &= (sdispl_global[r * env->m_numpes + i]
473  >= sdispl_global[r * env->m_numpes + i - 1]
474  );
475  }
476 
477  std::cout << ", ascending: " << (asc ? "okay" : "WRONG") << std::endl;
478 
479  err += asc ? 0 : 1;
480  }
481 
482 
483  /* ******************************************************************* */
484 
485  std::cout << "rdispl = " << std::endl;
486 
487  for (int r = 0; r < env->m_numpes; ++r)
488  {
489  err += (rdispl_global[r * env->m_numpes + 0] == 0) ? 0 : 1;
490 
491  std::cout << std::setw(11) << rdispl_global[r * env->m_numpes + 0];
492 
493  bool asc = true;
494 
495  for (int i = 1; i < env->m_numpes; ++i)
496  {
497  std::cout << std::setw(11) << rdispl_global[r * env->m_numpes + i];
498 
499  asc &= (rdispl_global[r * env->m_numpes + i]
500  >= rdispl_global[r * env->m_numpes + i - 1]
501  );
502  }
503 
504  std::cout << ", ascending: " << (asc ? "okay" : "WRONG") << std::endl;
505 
506  err += asc ? 0 : 1;
507  }
508 
509  /* ******************************************************************* */
510  }
511 
512  sort_pe_out(env) << "msr local errors: " << err << std::endl;
513 
514  unsigned serr;
515 
516  MPI_Reduce(&err, &serr, 1, MPI_UNSIGNED, MPI_SUM, 0, env->m_comm);
517 
518  err = serr;
519 
520  sort_root_out(env) << "msr local error sum: " << err << std::endl;
521 
522  /* *********************************************************************** */
523 
524  value_type* blval = new value_type[env->m_numpes * env->m_numpes];
525  value_type* brval = new value_type[env->m_numpes * env->m_numpes];
526 
527  for (int r = 0; r < env->m_numpes; ++r)
528  {
529  if (lborder_global[r * env->m_numpes + env->m_mype] < length)
530  {
531  value_type val = data[lborder_global[r * env->m_numpes + env->m_mype]];
532  blval[env->m_mype * env->m_numpes + r] = val;
533  }
534  else
535  {
536  blval[env->m_mype * env->m_numpes + r] = less.max_value();
537  }
538 
539  if (lborder_global[r * env->m_numpes + env->m_mype] - 1 < 0)
540  {
541  brval[env->m_mype * env->m_numpes + r] = less.min_value();
542  }
543  else
544  {
545  value_type val = data[lborder_global[r * env->m_numpes + env->m_mype] - 1];
546  brval[env->m_mype * env->m_numpes + r] = val;
547  }
548  }
549 
550  MPI_Datatype mpi_value_type; //value_type specification for MPI
551  MPI_Type_contiguous(sizeof(value_type), MPI_BYTE, &mpi_value_type);
552  MPI_Type_commit(&mpi_value_type);
553 
554  MPI_Gather(&blval[env->m_mype * env->m_numpes], env->m_numpes, mpi_value_type,
555  blval, env->m_numpes, mpi_value_type,
556  0, env->m_comm
557  );
558  MPI_Gather(&brval[env->m_mype * env->m_numpes], env->m_numpes, mpi_value_type,
559  brval, env->m_numpes, mpi_value_type,
560  0, env->m_comm
561  );
562 
563  MPI_Type_free(const_cast<MPI_Datatype*>(&mpi_value_type));
564 
565  if (env->m_mype == 0)
566  {
567  value_type* vmin = new value_type[env->m_numpes];
568  value_type* vmax = new value_type[env->m_numpes];
569 
570  std::cout << "blval = " << std::endl;
571 
572  for (int r = 0; r < env->m_numpes; ++r)
573  {
574  vmin[r] = less.max_value();
575 
576  for (int k = 0; k < env->m_numpes; ++k)
577  {
578  std::cout << blval[k * env->m_numpes + r];
579 
580  if (less(blval[k * env->m_numpes + r], vmin[r]))
581  {
582  vmin[r] = blval[k * env->m_numpes + r];
583  }
584  }
585  std::cout << " => min[" << sort_format_pe_number(r) << "] = " << vmin[r] << std::endl;
586  }
587 
588  std::cout << "brval = " << std::endl;
589 
590  for (int r = 0; r < env->m_numpes; ++r)
591  {
592  vmax[r] = less.min_value();
593 
594  for (int k = 0; k < env->m_numpes; ++k)
595  {
596  std::cout << brval[k * env->m_numpes + r];
597 
598  if (less(vmax[r], brval[k * env->m_numpes + r]))
599  {
600  vmax[r] = brval[k * env->m_numpes + r];
601  }
602  }
603  std::cout << " => max[" << sort_format_pe_number(r - 1) << "] = " << vmax[r] << std::endl;
604  }
605 
606  for (int r = 0; r < env->m_numpes - 1; ++r)
607  {
608  std::cout << "max[" << sort_format_pe_number(r) << "] <= "
609  << "min[" << sort_format_pe_number(r + 1) << "] iff "
610  << vmax[r + 1] << " <= "
611  << vmin[r + 1] << ": "
612  << (less(vmin[r + 1], vmax[r + 1]) ? "WRONG" : "okay") << std::endl;
613 
614  err += less(vmin[r + 1], vmax[r + 1]) ? 1 : 0;
615  }
616 
617  sort_pe_out(env) << "msr global error count: " << err << std::endl;
618 
619  if (err > 0)
620  {
621  sort_pe_msg(env, "* msr is ERRONEOUS");
622  CH_assert(false);
623  }
624  else
625  {
626  sort_pe_msg(env, "* msr looks GOOD");
627  }
628 
629  delete[] vmin;
630  delete[] vmax;
631  }
632 
633  delete[] blval;
634  delete[] brval;
635  delete[] lborder_global;
636  delete[] sdispl_global;
637  delete[] rdispl_global;
638  delete[] scount_global;
639  delete[] rcount_global;
640 }
641 
642 // ugly but not used
643 #include "sort_hist.H"
644 
645 #endif // CH_MPI
646 
647 #endif // SORT_UTILS_H
IndirectIndirectInvertedLess(InIt *data, Comparator &less)
Definition: sort_utils.H:96
S * data
Definition: sort_utils.H:27
#define CH_assert(cond)
Definition: CHArray.H:37
Definition: sort_utils.H:69
Comparator & less
Definition: sort_utils.H:50
InIt * data
Definition: sort_utils.H:71
Comparator & less
Definition: sort_utils.H:28
IndirectInvertedLess(S *data, Comparator &less)
Definition: sort_utils.H:52
IndirectIndirectLess(InIt *data, Comparator &less)
Definition: sort_utils.H:74
Comparator & less
Definition: sort_utils.H:72
InIt * data
Definition: sort_utils.H:93
bool operator()(const unsigned &a, const unsigned &b)
Definition: sort_utils.H:78
Definition: sort_utils.H:91
Definition: sort_utils.H:47
Comparator & less
Definition: sort_utils.H:94
bool operator()(const index_type &a, const index_type &b)
Definition: sort_utils.H:34
bool operator()(const unsigned &a, const unsigned &b)
Definition: sort_utils.H:100
bool operator()(const index_type &a, const index_type &b)
Definition: sort_utils.H:56
Definition: sort_utils.H:25
long long index_type
Definition: sort_utils.H:22
S * data
Definition: sort_utils.H:49
IndirectLess(S *data, Comparator &less)
Definition: sort_utils.H:30
unsigned long long key_type
Definition: sort_utils.H:21