dune-common  2.2.0
interface.hh
Go to the documentation of this file.
1 // $Id$
2 #ifndef DUNE_INTERFACE_HH
3 #define DUNE_INTERFACE_HH
4 
5 #if HAVE_MPI
6 
7 #include"remoteindices.hh"
9 
10 namespace Dune
11 {
32  {
33  public:
35  {};
36 
38  {}
39 
40  protected:
45  {}
46 
84  template<class R, class T1, class T2, class Op, bool send>
85  void buildInterface (const R& remoteIndices,
86  const T1& sourceFlags, const T2& destFlags,
87  Op& functor) const;
88  };
89 
98  {
99 
100  public:
101 
105  size_t size() const
106  {
107  return size_;
108  }
113  std::size_t& operator[](size_t i)
114  {
115  assert(i<size_);
116  return indices_[i];
117  }
122  std::size_t operator[](size_t i) const
123  {
124  assert(i<size_);
125  return indices_[i];
126  }
131  void reserve(size_t size)
132  {
133  indices_ = new std::size_t[size];
134  maxSize_ = size;
135 
136  }
140  void free()
141  {
142  if(indices_)
143  delete[] indices_;
144  maxSize_ = 0;
145  size_=0;
146  indices_=0;
147  }
151  void add(std::size_t index)
152  {
153  assert(size_<maxSize_);
154  indices_[size_++]=index;
155  }
156 
158  : size_(0), maxSize_(0), indices_(0)
159  {}
160 
162  {
163  }
164 
165  bool operator!=(const InterfaceInformation& o) const
166  {
167  return !operator==(o);
168  }
169 
170  bool operator==(const InterfaceInformation& o) const
171  {
172  if(size_!=o.size_)
173  return false;
174  for(std::size_t i=0; i< size_; ++i)
175  if(indices_[i]!=o.indices_[i])
176  return false;
177  return true;
178  }
179 
180  private:
184  size_t size_;
188  size_t maxSize_;
192  std::size_t* indices_;
193  };
194 
207  {
208 
209  public:
211 
212  typedef std::map<int,std::pair<Information,Information> > InformationMap;
213 
230  template<typename R, typename T1, typename T2>
231  void build(const R& remoteIndices, const T1& sourceFlags,
232  const T2& destFlags);
233 
237  void free();
238 
242  MPI_Comm communicator() const;
243 
252  const InformationMap& interfaces() const;
253 
254  Interface(MPI_Comm comm)
255  : communicator_(comm), interfaces_()
256  {}
257 
259  : communicator_(MPI_COMM_NULL), interfaces_()
260  {}
261 
265  void print() const;
266 
267  bool operator!=(const Interface& o) const
268  {
269  return ! operator==(o);
270  }
271 
272  bool operator==(const Interface& o) const
273  {
275  return false;
276  if(interfaces_.size()!=o.interfaces_.size())
277  return false;
278  typedef InformationMap::const_iterator MIter;
279 
280  for(MIter m=interfaces_.begin(), om=o.interfaces_.begin();
281  m!=interfaces_.end(); ++m, ++om)
282  {
283  if(om->first!=m->first)
284  return false;
285  if(om->second.first!=om->second.first)
286  return false;
287  if(om->second.second!=om->second.second)
288  return false;
289  }
290  return true;
291  }
292 
296  virtual ~Interface();
297 
298  void strip();
299  protected:
300 
310 
312  MPI_Comm communicator_;
313 
314  private:
322  InformationMap interfaces_;
323 
324  template<bool send>
325  class InformationBuilder
326  {
327  public:
328  InformationBuilder(InformationMap& interfaces)
329  : interfaces_(interfaces)
330  {}
331 
332  void reserve(int proc, int size)
333  {
334  if(send)
335  interfaces_[proc].first.reserve(size);
336  else
337  interfaces_[proc].second.reserve(size);
338  }
339  void add(int proc, std::size_t local)
340  {
341  if(send){
342  interfaces_[proc].first.add(local);
343  }else{
344  interfaces_[proc].second.add(local);
345  }
346  }
347 
348  private:
349  InformationMap& interfaces_;
350  };
351  };
352 
353  template<class R, class T1, class T2, class Op, bool send>
354  void InterfaceBuilder::buildInterface(const R& remoteIndices, const T1& sourceFlags, const T2& destFlags, Op& interfaceInformation) const
355  {
356 
357  if(!remoteIndices.isSynced())
358  DUNE_THROW(RemotexIndicesStateError,"RemoteIndices is not in sync with the index set. Call RemoteIndices::rebuild first!");
359  // Allocate the memory for the data type construction.
360  typedef R RemoteIndices;
361  typedef typename RemoteIndices::RemoteIndexMap::const_iterator const_iterator;
362  typedef typename RemoteIndices::ParallelIndexSet::const_iterator LocalIterator;
363 
364  const const_iterator end=remoteIndices.end();
365 
366  int rank;
367 
368  MPI_Comm_rank(remoteIndices.communicator(), &rank);
369 
370  // Allocate memory for the type construction.
371  for(const_iterator process=remoteIndices.begin(); process != end; ++process){
372  // Messure the number of indices send to the remote process first
373  int size=0;
374  typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
375  const RemoteIterator remoteEnd = send ? process->second.first->end() :
376  process->second.second->end();
377  RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
378 
379  while(remote!=remoteEnd){
380  if( send ? destFlags.contains(remote->attribute()) :
381  sourceFlags.contains(remote->attribute())){
382 
383  // do we send the index?
384  if( send ? sourceFlags.contains(remote->localIndexPair().local().attribute()) :
385  destFlags.contains(remote->localIndexPair().local().attribute()))
386  ++size;
387  }
388  ++remote;
389  }
390  interfaceInformation.reserve(process->first, size);
391  }
392 
393  // compare the local and remote indices and set up the types
394 
395  for(const_iterator process=remoteIndices.begin(); process != end; ++process){
396  typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
397  const RemoteIterator remoteEnd = send ? process->second.first->end() :
398  process->second.second->end();
399  RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
400 
401  while(remote!=remoteEnd){
402  if( send ? destFlags.contains(remote->attribute()) :
403  sourceFlags.contains(remote->attribute())){
404  // do we send the index?
405  if( send ? sourceFlags.contains(remote->localIndexPair().local().attribute()) :
406  destFlags.contains(remote->localIndexPair().local().attribute()))
407  interfaceInformation.add(process->first,remote->localIndexPair().local().local());
408  }
409  ++remote;
410  }
411  }
412  }
413 
414  inline MPI_Comm Interface::communicator() const
415  {
416  return communicator_;
417 
418  }
419 
420 
421  inline const std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces() const
422  {
423  return interfaces_;
424  }
425 
426  inline std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces()
427  {
428  return interfaces_;
429  }
430 
431  inline void Interface::print() const
432  {
433  typedef InformationMap::const_iterator const_iterator;
434  const const_iterator end=interfaces_.end();
435  int rank;
436  MPI_Comm_rank(communicator(), &rank);
437 
438  for(const_iterator infoPair=interfaces_.begin(); infoPair!=end; ++infoPair){
439  {
440  std::cout<<rank<<": send for process "<<infoPair->first<<": ";
441  const InterfaceInformation& info(infoPair->second.first);
442  for(size_t i=0; i < info.size(); i++)
443  std::cout<<info[i]<<" ";
444  std::cout<<std::endl;
445  }{
446 
447  std::cout<<rank<<": receive for process "<<infoPair->first<<": ";
448  const InterfaceInformation& info(infoPair->second.second);
449  for(size_t i=0; i < info.size(); i++)
450  std::cout<<info[i]<<" ";
451  std::cout<<std::endl;
452  }
453 
454  }
455  }
456 
457  template<typename R, typename T1, typename T2>
458  inline void Interface::build(const R& remoteIndices, const T1& sourceFlags,
459  const T2& destFlags)
460  {
461  communicator_=remoteIndices.communicator();
462 
463  assert(interfaces_.empty());
464 
465  // Build the send interface
466  InformationBuilder<true> sendInformation(interfaces_);
467  this->template buildInterface<R,T1,T2,InformationBuilder<true>,true>(remoteIndices, sourceFlags,
468  destFlags, sendInformation);
469 
470  // Build the receive interface
471  InformationBuilder<false> recvInformation(interfaces_);
472  this->template buildInterface<R,T1,T2,InformationBuilder<false>,false>(remoteIndices,sourceFlags,
473  destFlags, recvInformation);
474  strip();
475  }
476  inline void Interface::strip()
477  {
478  typedef InformationMap::iterator const_iterator;
479  for(const_iterator interfacePair = interfaces_.begin(); interfacePair != interfaces_.end();)
480  if(interfacePair->second.first.size()==0 && interfacePair->second.second.size()==0){
481  interfacePair->second.first.free();
482  interfacePair->second.second.free();
483  const_iterator toerase=interfacePair++;
484  interfaces_.erase(toerase);
485  }else
486  ++interfacePair;
487  }
488 
489  inline void Interface::free()
490  {
491  typedef InformationMap::iterator iterator;
492  typedef InformationMap::const_iterator const_iterator;
493  const const_iterator end = interfaces_.end();
494  for(iterator interfacePair = interfaces_.begin(); interfacePair != end; ++interfacePair){
495  interfacePair->second.first.free();
496  interfacePair->second.second.free();
497  }
498  interfaces_.clear();
499  }
500 
502  {
503  free();
504  }
506 }
507 
508 namespace std
509 {
510  inline ostream& operator<<(ostream& os, const Dune::Interface& interface)
511  {
512  typedef Dune::Interface::InformationMap InfoMap;
513  typedef InfoMap::const_iterator Iter;
514  for(Iter i=interface.interfaces().begin(), end = interface.interfaces().end();
515  i!=end; ++i)
516  {
517  os<<i->first<<": [ source=[";
518  for(std::size_t j=0; j < i->second.first.size(); ++j)
519  os<<i->second.first[j]<<" ";
520  os<<"] size="<<i->second.first.size()<<", target=[";
521  for(std::size_t j=0; j < i->second.second.size(); ++j)
522  os<<i->second.second[j]<<" ";
523  os<<"] size="<<i->second.second.size()<<"\n";
524  }
525  return os;
526  }
527 }// end namespace std
528 #endif // HAVE_MPI
529 
530 #endif