Sierra Toolkit  Version of the Day
ParallelIndex.hpp
1 /*------------------------------------------------------------------------*/
2 /* Copyright 2010 Sandia Corporation. */
3 /* Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive */
4 /* license for use of this work by or on behalf of the U.S. Government. */
5 /* Export of this program may require a license from the */
6 /* United States Government. */
7 /*------------------------------------------------------------------------*/
8 
9 #ifndef stk_util_parallel_ParallelIndex_hpp
10 #define stk_util_parallel_ParallelIndex_hpp
11 
12 #include <utility>
13 #include <vector>
14 #include <algorithm>
15 
16 #include <stdint.h>
17 #include <stk_util/parallel/Parallel.hpp>
18 #include <stk_util/parallel/ParallelComm.hpp>
19 #include <stk_util/util/PairIter.hpp>
20 
21 namespace stk_classic {
22 namespace util {
23 
24 template <class K, class P>
25 struct ParallelIndexDecomp
26 {
27  typedef K Key ;
28  typedef P Proc ;
29 
30  Proc operator()( const Proc proc_size, const Key key ) const {
31  return ( key >> 8 ) % proc_size ;
32  };
33 };
34 
35 
45 template <class K = uint64_t, class P = unsigned, class D = ParallelIndexDecomp<K, P> >
47 public:
48  typedef K Key ;
49  typedef P Proc ;
50  typedef D Decomp ;
51  typedef std::pair<Key, Proc> KeyProc ;
52 
53 #ifndef DOXYGEN_COMPILE
54  struct LessKeyProc {
55  bool operator()( const KeyProc & lhs, const KeyProc & rhs ) const {
56  return lhs < rhs ;
57  }
58 
59  bool operator()( const KeyProc & lhs, const Key rhs ) const {
60  return lhs.first < rhs ;
61  }
62  };
63 
64  struct EqualKeyProc {
65  bool operator()( const KeyProc & lhs, const KeyProc & rhs ) const {
66  return lhs == rhs ;
67  }
68  };
69 #endif /* DOXYGEN_COMPILE */
70 
71 
74  const std::vector<Key> & local )
75  : m_comm( comm ),
76  m_key_proc(),
77  m_decomp()
78  {
79  const unsigned p_size = parallel_machine_size( comm );
80 
81  CommAll all( comm );
82 
83  pack_map( all, local );
84 
85  all.allocate_buffers( p_size / 4, false );
86 
87  pack_map( all, local );
88 
89  all.communicate();
90 
91  unpack_map( all, m_key_proc );
92 
93  sort_unique( m_key_proc );
94  }
95 
96  ~ParallelIndex()
97  {}
98 
99 
103  void query(std::vector<KeyProc> & global ) const
104  {
105  const unsigned p_size = parallel_machine_size( m_comm );
106 
107  CommAll all( m_comm );
108 
109  pack_query( all, m_key_proc );
110 
111  all.allocate_buffers( p_size / 4, false );
112 
113  pack_query( all, m_key_proc );
114 
115  all.communicate();
116 
117  unpack_query( all, global );
118 
119  sort_unique( global );
120  }
121 
122 
127  void query(const std::vector<Key> & local,
128  std::vector<KeyProc> & global ) const
129  {
130  const unsigned p_size = parallel_machine_size( m_comm );
131 
132  std::vector<KeyProc> tmp ;
133 
134  {
135  CommAll all( m_comm );
136 
137  pack_map( all, local );
138 
139  all.allocate_buffers( p_size / 4, false );
140 
141  pack_map( all, local );
142 
143  all.communicate();
144 
145  unpack_map( all, tmp ); // { ( key, querying_processor ) }
146 
147  sort_unique( tmp );
148  }
149 
150  {
151  CommAll all( m_comm );
152 
153  pack_query( all, m_key_proc, tmp );
154 
155  all.allocate_buffers( p_size / 4, false );
156 
157  pack_query( all, m_key_proc, tmp );
158 
159  all.communicate();
160 
161  unpack_query( all, global );
162 
163  sort_unique( global );
164  }
165  }
166 
167 private:
168  void sort_unique( std::vector<KeyProc> & key_proc ) const
169  {
170  typename std::vector<KeyProc>::iterator i = key_proc.begin();
171  typename std::vector<KeyProc>::iterator j = key_proc.end();
172 
173  std::sort( i, j, LessKeyProc() );
174  i = std::unique( i, j, EqualKeyProc() );
175  key_proc.erase( i, j );
176  }
177 
178  void pack_map( CommAll & all, const std::vector<Key> & local ) const
179  {
180  const unsigned p_size = all.parallel_size();
181 
182  typename std::vector<Key>::const_iterator i ;
183 
184  for ( i = local.begin() ; i != local.end() ; ++i ) {
185  const Key value = *i ;
186  const unsigned proc = m_decomp( p_size, value );
187  CommBuffer & buf = all.send_buffer(proc);
188  buf.pack<Key>( value );
189  }
190  }
191 
192  void unpack_map( CommAll & all, std::vector< KeyProc > & key_proc ) const
193  {
194  const unsigned p_size = all.parallel_size();
195 
196  unsigned count = 0 ;
197  for ( unsigned p = 0 ; p < p_size ; ++p ) {
198  count += all.recv_buffer( p ).capacity() / sizeof(Key);
199  }
200 
201  key_proc.clear();
202  key_proc.reserve( count );
203 
204  KeyProc value ;
205 
206  for ( unsigned p = 0 ; p < p_size ; ++p ) {
207  CommBuffer & buf = all.recv_buffer( p );
208  value.second = p ;
209  while ( buf.remaining() ) {
210  buf.unpack<Key>( value.first );
211  key_proc.push_back( value );
212  }
213  }
214  }
215 
216  void pack_query( CommAll & all, const std::vector< KeyProc > & key_proc ) const
217  {
218  KeyProc value ;
219 
220  typename std::vector< KeyProc >::const_iterator i ;
221 
222  for ( i = key_proc.begin() ; i != key_proc.end() ; ) {
223  value.first = i->first ;
224 
225  const typename std::vector< KeyProc >::const_iterator i_beg = i ;
226 
227  for ( ; i != key_proc.end() && value.first == i->first ; ++i )
228  ;
229 
230  const typename std::vector< KeyProc >::const_iterator i_end = i ;
231 
232  for ( i = i_beg ; i != i_end ; ++i ) {
233  CommBuffer & buf = all.send_buffer( i->second );
234 
235  typename std::vector< KeyProc >::const_iterator j ;
236 
237  for ( j = i_beg ; j != i_end ; ++j ) {
238  if ( j != i ) {
239  value.second = j->second ;
240  buf.pack<KeyProc>( value );
241  }
242  }
243  }
244  }
245  }
246 
247  void pack_query( CommAll & all,
248  const std::vector< KeyProc > & key_proc_map,
249  const std::vector< KeyProc > & query ) const
250  {
251  KeyProc value ;
252 
253  for (typename std::vector< KeyProc >::const_iterator i = query.begin() ; i != query.end() ; ) {
254  value.first = i->first ;
255 
256  typename std::vector< KeyProc >::const_iterator key_begin = std::lower_bound( key_proc_map.begin(), key_proc_map.end(), value.first, LessKeyProc() );
257 
258  typename std::vector< KeyProc >::const_iterator key_end = key_begin;
259  while ( key_end != key_proc_map.end() && key_end->first == value.first)
260  ++key_end;
261 
262  for ( ; i != query.end() && value.first == i->first ; ++i ) {
263  CommBuffer & buf = all.send_buffer( i->second );
264 
265  for ( typename std::vector< KeyProc >::const_iterator j = key_begin ; j != key_end ; ++j ) {
266  value.second = j->second ;
267  buf.pack<KeyProc>( value );
268  }
269  }
270  }
271  }
272 
273  void unpack_query( CommAll & all, std::vector< KeyProc > & key_proc ) const
274  {
275  const unsigned p_size = all.parallel_size();
276 
277  KeyProc entry ;
278 
279  for ( unsigned p = 0 ; p < p_size ; ++p ) {
280  CommBuffer & buf = all.recv_buffer( p );
281  while ( buf.remaining() ) {
282  buf.unpack<KeyProc>( entry );
283  key_proc.push_back( entry );
284  }
285  }
286  }
287 
288 private:
289  ParallelIndex();
290  ParallelIndex( const ParallelIndex & );
291  ParallelIndex & operator = ( const ParallelIndex & );
292 
293 private:
294  ParallelMachine m_comm ;
295  std::vector<KeyProc> m_key_proc ;
296  Decomp m_decomp;
297 };
298 
299 } // namespace util
300 } // namespace stk_classic
301 
302 //----------------------------------------------------------------------
303 
304 #endif
305 
ParallelIndex(ParallelMachine comm, const std::vector< Key > &local)
Construct with locally-submitted keys.
unsigned parallel_machine_size(ParallelMachine parallel_machine)
Member function parallel_machine_size ...
Definition: Parallel.cpp:18
Sierra Toolkit.
MPI_Comm ParallelMachine
Definition: Parallel.hpp:32
void query(const std::vector< Key > &local, std::vector< KeyProc > &global) const
Query which processors submitted the given keys. The local processor is in the output if it submitted...
Parallel cross-reference index for a collection of &#39;Key&#39; keys.
void query(std::vector< KeyProc > &global) const
Query which other processors submitted the same keys that the local processor submitted.
eastl::iterator_traits< InputIterator >::difference_type count(InputIterator first, InputIterator last, const T &value)