dune-common  2.2.1
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 
363  const const_iterator end=remoteIndices.end();
364 
365  int rank;
366 
367  MPI_Comm_rank(remoteIndices.communicator(), &rank);
368 
369  // Allocate memory for the type construction.
370  for(const_iterator process=remoteIndices.begin(); process != end; ++process){
371  // Messure the number of indices send to the remote process first
372  int size=0;
373  typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
374  const RemoteIterator remoteEnd = send ? process->second.first->end() :
375  process->second.second->end();
376  RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
377 
378  while(remote!=remoteEnd){
379  if( send ? destFlags.contains(remote->attribute()) :
380  sourceFlags.contains(remote->attribute())){
381 
382  // do we send the index?
383  if( send ? sourceFlags.contains(remote->localIndexPair().local().attribute()) :
384  destFlags.contains(remote->localIndexPair().local().attribute()))
385  ++size;
386  }
387  ++remote;
388  }
389  interfaceInformation.reserve(process->first, size);
390  }
391 
392  // compare the local and remote indices and set up the types
393 
394  for(const_iterator process=remoteIndices.begin(); process != end; ++process){
395  typedef typename RemoteIndices::RemoteIndexList::const_iterator RemoteIterator;
396  const RemoteIterator remoteEnd = send ? process->second.first->end() :
397  process->second.second->end();
398  RemoteIterator remote = send ? process->second.first->begin() : process->second.second->begin();
399 
400  while(remote!=remoteEnd){
401  if( send ? destFlags.contains(remote->attribute()) :
402  sourceFlags.contains(remote->attribute())){
403  // do we send the index?
404  if( send ? sourceFlags.contains(remote->localIndexPair().local().attribute()) :
405  destFlags.contains(remote->localIndexPair().local().attribute()))
406  interfaceInformation.add(process->first,remote->localIndexPair().local().local());
407  }
408  ++remote;
409  }
410  }
411  }
412 
413  inline MPI_Comm Interface::communicator() const
414  {
415  return communicator_;
416 
417  }
418 
419 
420  inline const std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces() const
421  {
422  return interfaces_;
423  }
424 
425  inline std::map<int,std::pair<InterfaceInformation,InterfaceInformation> >& Interface::interfaces()
426  {
427  return interfaces_;
428  }
429 
430  inline void Interface::print() const
431  {
432  typedef InformationMap::const_iterator const_iterator;
433  const const_iterator end=interfaces_.end();
434  int rank;
435  MPI_Comm_rank(communicator(), &rank);
436 
437  for(const_iterator infoPair=interfaces_.begin(); infoPair!=end; ++infoPair){
438  {
439  std::cout<<rank<<": send for process "<<infoPair->first<<": ";
440  const InterfaceInformation& info(infoPair->second.first);
441  for(size_t i=0; i < info.size(); i++)
442  std::cout<<info[i]<<" ";
443  std::cout<<std::endl;
444  }{
445 
446  std::cout<<rank<<": receive for process "<<infoPair->first<<": ";
447  const InterfaceInformation& info(infoPair->second.second);
448  for(size_t i=0; i < info.size(); i++)
449  std::cout<<info[i]<<" ";
450  std::cout<<std::endl;
451  }
452 
453  }
454  }
455 
456  template<typename R, typename T1, typename T2>
457  inline void Interface::build(const R& remoteIndices, const T1& sourceFlags,
458  const T2& destFlags)
459  {
460  communicator_=remoteIndices.communicator();
461 
462  assert(interfaces_.empty());
463 
464  // Build the send interface
465  InformationBuilder<true> sendInformation(interfaces_);
466  this->template buildInterface<R,T1,T2,InformationBuilder<true>,true>(remoteIndices, sourceFlags,
467  destFlags, sendInformation);
468 
469  // Build the receive interface
470  InformationBuilder<false> recvInformation(interfaces_);
471  this->template buildInterface<R,T1,T2,InformationBuilder<false>,false>(remoteIndices,sourceFlags,
472  destFlags, recvInformation);
473  strip();
474  }
475  inline void Interface::strip()
476  {
477  typedef InformationMap::iterator const_iterator;
478  for(const_iterator interfacePair = interfaces_.begin(); interfacePair != interfaces_.end();)
479  if(interfacePair->second.first.size()==0 && interfacePair->second.second.size()==0){
480  interfacePair->second.first.free();
481  interfacePair->second.second.free();
482  const_iterator toerase=interfacePair++;
483  interfaces_.erase(toerase);
484  }else
485  ++interfacePair;
486  }
487 
488  inline void Interface::free()
489  {
490  typedef InformationMap::iterator iterator;
491  typedef InformationMap::const_iterator const_iterator;
492  const const_iterator end = interfaces_.end();
493  for(iterator interfacePair = interfaces_.begin(); interfacePair != end; ++interfacePair){
494  interfacePair->second.first.free();
495  interfacePair->second.second.free();
496  }
497  interfaces_.clear();
498  }
499 
501  {
502  free();
503  }
505 }
506 
507 namespace std
508 {
509  inline ostream& operator<<(ostream& os, const Dune::Interface& interface)
510  {
511  typedef Dune::Interface::InformationMap InfoMap;
512  typedef InfoMap::const_iterator Iter;
513  for(Iter i=interface.interfaces().begin(), end = interface.interfaces().end();
514  i!=end; ++i)
515  {
516  os<<i->first<<": [ source=[";
517  for(std::size_t j=0; j < i->second.first.size(); ++j)
518  os<<i->second.first[j]<<" ";
519  os<<"] size="<<i->second.first.size()<<", target=[";
520  for(std::size_t j=0; j < i->second.second.size(); ++j)
521  os<<i->second.second[j]<<" ";
522  os<<"] size="<<i->second.second.size()<<"\n";
523  }
524  return os;
525  }
526 }// end namespace std
527 #endif // HAVE_MPI
528 
529 #endif