dune-istl  2.8.0
matrixredistribute.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_ISTL_MATRIXREDISTRIBUTE_HH
4 #define DUNE_ISTL_MATRIXREDISTRIBUTE_HH
5 #include <memory>
6 #include "repartition.hh"
7 #include <dune/common/exceptions.hh>
8 #include <dune/common/parallel/indexset.hh>
10 #include <dune/istl/paamg/pinfo.hh>
16 namespace Dune
17 {
18  template<typename T>
20  {
21  bool isSetup() const
22  {
23  return false;
24  }
25  template<class D>
26  void redistribute([[maybe_unused]] const D& from, [[maybe_unused]] D& to) const
27  {}
28 
29  template<class D>
30  void redistributeBackward([[maybe_unused]] D& from, [[maybe_unused]]const D& to) const
31  {}
32 
33  void resetSetup()
34  {}
35 
36  void setNoRows([[maybe_unused]] std::size_t size)
37  {}
38 
39  void setNoCopyRows([[maybe_unused]] std::size_t size)
40  {}
41 
42  void setNoBackwardsCopyRows([[maybe_unused]] std::size_t size)
43  {}
44 
45  std::size_t getRowSize([[maybe_unused]] std::size_t index) const
46  {
47  return -1;
48  }
49 
50  std::size_t getCopyRowSize([[maybe_unused]] std::size_t index) const
51  {
52  return -1;
53  }
54 
55  std::size_t getBackwardsCopyRowSize([[maybe_unused]] std::size_t index) const
56  {
57  return -1;
58  }
59 
60  };
61 
62 #if HAVE_MPI
63  template<typename T, typename T1>
65  {
66  public:
68 
70  : interface(), setup_(false)
71  {}
72 
74  {
75  return interface;
76  }
77  template<typename IS>
78  void checkInterface(const IS& source,
79  const IS& target, MPI_Comm comm)
80  {
81  auto ri = std::make_unique<RemoteIndices<IS> >(source, target, comm);
82  ri->template rebuild<true>();
83  Interface inf;
85  int rank;
86  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
87  inf.free();
88  inf.build(*ri, flags, flags);
89 
90 
91 #ifdef DEBUG_REPART
92  if(inf!=interface) {
93 
94  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
95  if(rank==0)
96  std::cout<<"Interfaces do not match!"<<std::endl;
97  std::cout<<rank<<": redist interface new :"<<inf<<std::endl;
98  std::cout<<rank<<": redist interface :"<<interface<<std::endl;
99 
100  throw "autsch!";
101  }
102 #endif
103  }
104  void setSetup()
105  {
106  setup_=true;
107  interface.strip();
108  }
109 
110  void resetSetup()
111  {
112  setup_=false;
113  }
114 
115  template<class GatherScatter, class D>
116  void redistribute(const D& from, D& to) const
117  {
118  BufferedCommunicator communicator;
119  communicator.template build<D>(from,to, interface);
120  communicator.template forward<GatherScatter>(from, to);
121  communicator.free();
122  }
123  template<class GatherScatter, class D>
124  void redistributeBackward(D& from, const D& to) const
125  {
126 
127  BufferedCommunicator communicator;
128  communicator.template build<D>(from,to, interface);
129  communicator.template backward<GatherScatter>(from, to);
130  communicator.free();
131  }
132 
133  template<class D>
134  void redistribute(const D& from, D& to) const
135  {
136  redistribute<CopyGatherScatter<D> >(from,to);
137  }
138  template<class D>
139  void redistributeBackward(D& from, const D& to) const
140  {
141  redistributeBackward<CopyGatherScatter<D> >(from,to);
142  }
143  bool isSetup() const
144  {
145  return setup_;
146  }
147 
148  void reserve(std::size_t size)
149  {}
150 
151  std::size_t& getRowSize(std::size_t index)
152  {
153  return rowSize[index];
154  }
155 
156  std::size_t getRowSize(std::size_t index) const
157  {
158  return rowSize[index];
159  }
160 
161  std::size_t& getCopyRowSize(std::size_t index)
162  {
163  return copyrowSize[index];
164  }
165 
166  std::size_t getCopyRowSize(std::size_t index) const
167  {
168  return copyrowSize[index];
169  }
170 
171  std::size_t& getBackwardsCopyRowSize(std::size_t index)
172  {
173  return backwardscopyrowSize[index];
174  }
175 
176  std::size_t getBackwardsCopyRowSize(std::size_t index) const
177  {
178  return backwardscopyrowSize[index];
179  }
180 
181  void setNoRows(std::size_t rows)
182  {
183  rowSize.resize(rows, 0);
184  }
185 
186  void setNoCopyRows(std::size_t rows)
187  {
188  copyrowSize.resize(rows, 0);
189  }
190 
191  void setNoBackwardsCopyRows(std::size_t rows)
192  {
193  backwardscopyrowSize.resize(rows, 0);
194  }
195 
196  private:
197  std::vector<std::size_t> rowSize;
198  std::vector<std::size_t> copyrowSize;
199  std::vector<std::size_t> backwardscopyrowSize;
200  RedistributeInterface interface;
201  bool setup_;
202  };
203 
212  template<class M, class RI>
214  {
215  // Make the default communication policy work.
216  typedef typename M::size_type value_type;
217  typedef typename M::size_type size_type;
218 
224  CommMatrixRowSize(const M& m_, RI& rowsize_)
225  : matrix(m_), rowsize(rowsize_)
226  {}
227  const M& matrix;
228  RI& rowsize;
229 
230  };
231 
232 
241  template<class M, class I>
243  {
244  typedef typename M::size_type size_type;
245 
252  CommMatrixSparsityPattern(const M& m_, const Dune::GlobalLookupIndexSet<I>& idxset_, const I& aggidxset_)
253  : matrix(m_), idxset(idxset_), aggidxset(aggidxset_), rowsize()
254  {}
255 
263  CommMatrixSparsityPattern(const M& m_, const Dune::GlobalLookupIndexSet<I>& idxset_, const I& aggidxset_,
264  const std::vector<typename M::size_type>& rowsize_)
265  : matrix(m_), idxset(idxset_), aggidxset(aggidxset_), sparsity(aggidxset_.size()), rowsize(&rowsize_)
266  {}
267 
275  {
276  // insert diagonal to overlap rows
277  typedef typename Dune::GlobalLookupIndexSet<I>::const_iterator IIter;
278  typedef typename Dune::OwnerOverlapCopyCommunication<int>::OwnerSet OwnerSet;
279  std::size_t nnz=0;
280 #ifdef DEBUG_REPART
281  int rank;
282 
283  MPI_Comm_rank(MPI_COMM_WORLD, &rank);
284 #endif
285  for(IIter i= aggidxset.begin(), end=aggidxset.end(); i!=end; ++i) {
286  if(!OwnerSet::contains(i->local().attribute())) {
287 #ifdef DEBUG_REPART
288  std::cout<<rank<<" Inserting diagonal for"<<i->local()<<std::endl;
289 #endif
290  sparsity[i->local()].insert(i->local());
291  }
292 
293  nnz+=sparsity[i->local()].size();
294  }
295  assert( aggidxset.size()==sparsity.size());
296 
297  if(nnz>0) {
298  m.setSize(aggidxset.size(), aggidxset.size(), nnz);
299  m.setBuildMode(M::row_wise);
300  typename M::CreateIterator citer=m.createbegin();
301 #ifdef DEBUG_REPART
302  std::size_t idx=0;
303  bool correct=true;
304  Dune::GlobalLookupIndexSet<I> global(aggidxset);
305 #endif
306  typedef typename std::vector<std::set<size_type> >::const_iterator Iter;
307  for(Iter i=sparsity.begin(), end=sparsity.end(); i!=end; ++i, ++citer)
308  {
309  typedef typename std::set<size_type>::const_iterator SIter;
310  for(SIter si=i->begin(), send=i->end(); si!=send; ++si)
311  citer.insert(*si);
312 #ifdef DEBUG_REPART
313  if(i->find(idx)==i->end()) {
314  const typename I::IndexPair* gi=global.pair(idx);
315  assert(gi);
316  std::cout<<rank<<": row "<<idx<<" is missing a diagonal entry! global="<<gi->global()<<" attr="<<gi->local().attribute()<<" "<<
317  OwnerSet::contains(gi->local().attribute())<<
318  " row size="<<i->size()<<std::endl;
319  correct=false;
320  }
321  ++idx;
322 #endif
323  }
324 #ifdef DEBUG_REPART
325  if(!correct)
326  throw "bla";
327 #endif
328  }
329  }
330 
338  void completeSparsityPattern(std::vector<std::set<size_type> > add_sparsity)
339  {
340  for (unsigned int i = 0; i != sparsity.size(); ++i) {
341  if (add_sparsity[i].size() != 0) {
342  typedef std::set<size_type> Set;
343  Set tmp_set;
344  std::insert_iterator<Set> tmp_insert (tmp_set, tmp_set.begin());
345  std::set_union(add_sparsity[i].begin(), add_sparsity[i].end(),
346  sparsity[i].begin(), sparsity[i].end(), tmp_insert);
347  sparsity[i].swap(tmp_set);
348  }
349  }
350  }
351 
352  const M& matrix;
353  typedef Dune::GlobalLookupIndexSet<I> LookupIndexSet;
354  const Dune::GlobalLookupIndexSet<I>& idxset;
355  const I& aggidxset;
356  std::vector<std::set<size_type> > sparsity;
357  const std::vector<size_type>* rowsize;
358  };
359 
360  template<class M, class I>
361  struct CommPolicy<CommMatrixSparsityPattern<M,I> >
362  {
364 
369  typedef typename I::GlobalIndex IndexedType;
370 
372  typedef VariableSize IndexedTypeFlag;
373 
374  static typename M::size_type getSize(const Type& t, std::size_t i)
375  {
376  if(!t.rowsize)
377  return t.matrix[i].size();
378  else
379  {
380  assert((*t.rowsize)[i]>0);
381  return (*t.rowsize)[i];
382  }
383  }
384  };
385 
392  template<class M, class I>
394  {
403  CommMatrixRow(M& m_, const Dune::GlobalLookupIndexSet<I>& idxset_, const I& aggidxset_)
404  : matrix(m_), idxset(idxset_), aggidxset(aggidxset_), rowsize()
405  {}
406 
410  CommMatrixRow(M& m_, const Dune::GlobalLookupIndexSet<I>& idxset_, const I& aggidxset_,
411  std::vector<size_t>& rowsize_)
412  : matrix(m_), idxset(idxset_), aggidxset(aggidxset_), rowsize(&rowsize_)
413  {}
420  {
421  typedef typename Dune::GlobalLookupIndexSet<I>::const_iterator Iter;
422  typedef typename Dune::OwnerOverlapCopyCommunication<int>::OwnerSet OwnerSet;
423 
424  for(Iter i= aggidxset.begin(), end=aggidxset.end(); i!=end; ++i)
425  if(!OwnerSet::contains(i->local().attribute())) {
426  // Set to Dirchlet
427  typedef typename M::ColIterator CIter;
428  for(CIter c=matrix[i->local()].begin(), cend= matrix[i->local()].end();
429  c!= cend; ++c)
430  {
431  *c=0;
432  if(c.index()==i->local()) {
433  auto setDiagonal = [](auto&& scalarOrMatrix, const auto& value) {
434  auto&& matrixView = Dune::Impl::asMatrix(scalarOrMatrix);
435  for (auto rowIt = matrixView.begin(); rowIt != matrixView.end(); ++rowIt)
436  (*rowIt)[rowIt.index()] = value;
437  };
438  setDiagonal(*c, 1);
439  }
440  }
441  }
442  }
444  M& matrix;
446  const Dune::GlobalLookupIndexSet<I>& idxset;
448  const I& aggidxset;
450  std::vector<size_t>* rowsize; // row sizes differ from sender side in overlap!
451  };
452 
453  template<class M, class I>
454  struct CommPolicy<CommMatrixRow<M,I> >
455  {
457 
462  typedef std::pair<typename I::GlobalIndex,typename M::block_type> IndexedType;
463 
465  typedef VariableSize IndexedTypeFlag;
466 
467  static std::size_t getSize(const Type& t, std::size_t i)
468  {
469  if(!t.rowsize)
470  return t.matrix[i].size();
471  else
472  {
473  assert((*t.rowsize)[i]>0);
474  return (*t.rowsize)[i];
475  }
476  }
477  };
478 
479  template<class M, class I, class RI>
481  {
483 
484  static const typename M::size_type gather(const Container& cont, std::size_t i)
485  {
486  return cont.matrix[i].size();
487  }
488  static void scatter(Container& cont, const typename M::size_type& rowsize, std::size_t i)
489  {
490  assert(rowsize);
491  cont.rowsize.getRowSize(i)=rowsize;
492  }
493 
494  };
495 
496  template<class M, class I, class RI>
498  {
500 
501  static const typename M::size_type gather(const Container& cont, std::size_t i)
502  {
503  return cont.matrix[i].size();
504  }
505  static void scatter(Container& cont, const typename M::size_type& rowsize, std::size_t i)
506  {
507  assert(rowsize);
508  if (rowsize > cont.rowsize.getCopyRowSize(i))
509  cont.rowsize.getCopyRowSize(i)=rowsize;
510  }
511 
512  };
513 
514  template<class M, class I>
516  {
517  typedef typename I::GlobalIndex GlobalIndex;
519  typedef typename M::ConstColIterator ColIter;
520 
521  static ColIter col;
523 
524  static const GlobalIndex& gather(const Container& cont, std::size_t i, std::size_t j)
525  {
526  if(j==0)
527  col=cont.matrix[i].begin();
528  else if (col!=cont.matrix[i].end())
529  ++col;
530 
531  //copy communication: different row sizes for copy rows with the same global index
532  //are possible. If all values of current matrix row are sent, send
533  //std::numeric_limits<GlobalIndex>::max()
534  //and receiver will ignore it.
535  if (col==cont.matrix[i].end()) {
536  numlimits = std::numeric_limits<GlobalIndex>::max();
537  return numlimits;
538  }
539  else {
540  const typename I::IndexPair* index=cont.idxset.pair(col.index());
541  assert(index);
542  // Only send index if col is no ghost
543  if ( index->local().attribute() != 2)
544  return index->global();
545  else {
546  numlimits = std::numeric_limits<GlobalIndex>::max();
547  return numlimits;
548  }
549  }
550  }
551  static void scatter(Container& cont, const GlobalIndex& gi, std::size_t i, [[maybe_unused]] std::size_t j)
552  {
553  try{
554  if (gi != std::numeric_limits<GlobalIndex>::max()) {
555  const typename I::IndexPair& ip=cont.aggidxset.at(gi);
556  assert(ip.global()==gi);
557  std::size_t column = ip.local();
558  cont.sparsity[i].insert(column);
559 
560  typedef typename Dune::OwnerOverlapCopyCommunication<int>::OwnerSet OwnerSet;
561  if(!OwnerSet::contains(ip.local().attribute()))
562  // preserve symmetry for overlap
563  cont.sparsity[column].insert(i);
564  }
565  }
566  catch(const Dune::RangeError&) {
567  // Entry not present in the new index set. Ignore!
568 #ifdef DEBUG_REPART
569  typedef typename Container::LookupIndexSet GlobalLookup;
570  typedef typename GlobalLookup::IndexPair IndexPair;
571  typedef typename Dune::OwnerOverlapCopyCommunication<int>::OwnerSet OwnerSet;
572 
573  GlobalLookup lookup(cont.aggidxset);
574  const IndexPair* pi=lookup.pair(i);
575  assert(pi);
576  if(OwnerSet::contains(pi->local().attribute())) {
577  int rank;
578  MPI_Comm_rank(MPI_COMM_WORLD,&rank);
579  std::cout<<rank<<cont.aggidxset<<std::endl;
580  std::cout<<rank<<": row "<<i<<" (global="<<gi <<") not in index set for owner index "<<pi->global()<<std::endl;
581  throw;
582  }
583 #endif
584  }
585  }
586 
587  };
588  template<class M, class I>
590 
591  template<class M, class I>
593 
594 
595  template<class M, class I>
597  {
598  typedef typename I::GlobalIndex GlobalIndex;
600  typedef typename M::ConstColIterator ColIter;
601  typedef typename std::pair<GlobalIndex,typename M::block_type> Data;
602  static ColIter col;
603  static Data datastore;
605 
606  static const Data& gather(const Container& cont, std::size_t i, std::size_t j)
607  {
608  if(j==0)
609  col=cont.matrix[i].begin();
610  else if (col!=cont.matrix[i].end())
611  ++col;
612  // copy communication: different row sizes for copy rows with the same global index
613  // are possible. If all values of current matrix row are sent, send
614  // std::numeric_limits<GlobalIndex>::max()
615  // and receiver will ignore it.
616  if (col==cont.matrix[i].end()) {
617  numlimits = std::numeric_limits<GlobalIndex>::max();
619  return datastore;
620  }
621  else {
622  // convert local column index to global index
623  const typename I::IndexPair* index=cont.idxset.pair(col.index());
624  assert(index);
625  // Store the data to prevent reference to temporary
626  // Only send index if col is no ghost
627  if ( index->local().attribute() != 2)
628  datastore = Data(index->global(),*col);
629  else {
630  numlimits = std::numeric_limits<GlobalIndex>::max();
632  }
633  return datastore;
634  }
635  }
636  static void scatter(Container& cont, const Data& data, std::size_t i, [[maybe_unused]] std::size_t j)
637  {
638  try{
639  if (data.first != std::numeric_limits<GlobalIndex>::max()) {
640  typename M::size_type column=cont.aggidxset.at(data.first).local();
641  cont.matrix[i][column]=data.second;
642  }
643  }
644  catch(const Dune::RangeError&) {
645  // This an overlap row and might therefore lack some entries!
646  }
647 
648  }
649  };
650 
651  template<class M, class I>
653 
654  template<class M, class I>
656 
657  template<class M, class I>
659 
660  template<typename M, typename C>
661  void redistributeSparsityPattern(M& origMatrix, M& newMatrix, C& origComm, C& newComm,
663  {
664  typename C::CopySet copyflags;
665  typename C::OwnerSet ownerflags;
666  typedef typename C::ParallelIndexSet IndexSet;
667  typedef RedistributeInformation<C> RI;
668  std::vector<typename M::size_type> rowsize(newComm.indexSet().size(), 0);
669  std::vector<typename M::size_type> copyrowsize(newComm.indexSet().size(), 0);
670  std::vector<typename M::size_type> backwardscopyrowsize(origComm.indexSet().size(), 0);
671 
672  // get owner rowsizes
673  CommMatrixRowSize<M,RI> commRowSize(origMatrix, ri);
674  ri.template redistribute<MatrixRowSizeGatherScatter<M,IndexSet,RI> >(commRowSize,commRowSize);
675 
676  origComm.buildGlobalLookup();
677 
678  for (std::size_t i=0; i < newComm.indexSet().size(); i++) {
679  rowsize[i] = ri.getRowSize(i);
680  }
681  // get sparsity pattern from owner rows
683  origsp(origMatrix, origComm.globalLookup(), newComm.indexSet());
685  newsp(origMatrix, origComm.globalLookup(), newComm.indexSet(), rowsize);
686 
687  ri.template redistribute<MatrixSparsityPatternGatherScatter<M,IndexSet> >(origsp,newsp);
688 
689  // build copy to owner interface to get missing matrix values for novlp case
691  RemoteIndices<IndexSet> *ris = new RemoteIndices<IndexSet>(origComm.indexSet(),
692  newComm.indexSet(),
693  origComm.communicator());
694  ris->template rebuild<true>();
695 
696  ri.getInterface().free();
697  ri.getInterface().build(*ris,copyflags,ownerflags);
698 
699  // get copy rowsizes
700  CommMatrixRowSize<M,RI> commRowSize_copy(origMatrix, ri);
701  ri.template redistribute<MatrixCopyRowSizeGatherScatter<M,IndexSet,RI> >(commRowSize_copy,
702  commRowSize_copy);
703 
704  for (std::size_t i=0; i < newComm.indexSet().size(); i++) {
705  copyrowsize[i] = ri.getCopyRowSize(i);
706  }
707  //get copy rowsizes for sender
708  ri.redistributeBackward(backwardscopyrowsize,copyrowsize);
709  for (std::size_t i=0; i < origComm.indexSet().size(); i++) {
710  ri.getBackwardsCopyRowSize(i) = backwardscopyrowsize[i];
711  }
712 
713  // get sparsity pattern from copy rows
714  CommMatrixSparsityPattern<M,IndexSet> origsp_copy(origMatrix,
715  origComm.globalLookup(),
716  newComm.indexSet(),
717  backwardscopyrowsize);
718  CommMatrixSparsityPattern<M,IndexSet> newsp_copy(origMatrix, origComm.globalLookup(),
719  newComm.indexSet(), copyrowsize);
720  ri.template redistribute<MatrixSparsityPatternGatherScatter<M,IndexSet> >(origsp_copy,
721  newsp_copy);
722 
723  newsp.completeSparsityPattern(newsp_copy.sparsity);
724  newsp.storeSparsityPattern(newMatrix);
725  }
726  else
727  newsp.storeSparsityPattern(newMatrix);
728 
729 #ifdef DUNE_ISTL_WITH_CHECKING
730  // Check for symmetry
731  int ret=0;
732  typedef typename M::ConstRowIterator RIter;
733  for(RIter row=newMatrix.begin(), rend=newMatrix.end(); row != rend; ++row) {
734  typedef typename M::ConstColIterator CIter;
735  for(CIter col=row->begin(), cend=row->end(); col!=cend; ++col)
736  {
737  try{
738  newMatrix[col.index()][row.index()];
739  }catch(const Dune::ISTLError&) {
740  std::cerr<<newComm.communicator().rank()<<": entry ("
741  <<col.index()<<","<<row.index()<<") missing! for symmetry!"<<std::endl;
742  ret=1;
743 
744  }
745 
746  }
747  }
748 
749  if(ret)
750  DUNE_THROW(ISTLError, "Matrix not symmetric!");
751 #endif
752  }
753 
754  template<typename M, typename C>
755  void redistributeMatrixEntries(M& origMatrix, M& newMatrix, C& origComm, C& newComm,
757  {
758  typedef typename C::ParallelIndexSet IndexSet;
759  typename C::OwnerSet ownerflags;
760  std::vector<typename M::size_type> rowsize(newComm.indexSet().size(), 0);
761  std::vector<typename M::size_type> copyrowsize(newComm.indexSet().size(), 0);
762  std::vector<typename M::size_type> backwardscopyrowsize(origComm.indexSet().size(), 0);
763 
764  for (std::size_t i=0; i < newComm.indexSet().size(); i++) {
765  rowsize[i] = ri.getRowSize(i);
767  copyrowsize[i] = ri.getCopyRowSize(i);
768  }
769  }
770 
771  for (std::size_t i=0; i < origComm.indexSet().size(); i++)
773  backwardscopyrowsize[i] = ri.getBackwardsCopyRowSize(i);
774 
775 
777  // fill sparsity pattern from copy rows
778  CommMatrixRow<M,IndexSet> origrow_copy(origMatrix, origComm.globalLookup(),
779  newComm.indexSet(), backwardscopyrowsize);
780  CommMatrixRow<M,IndexSet> newrow_copy(newMatrix, origComm.globalLookup(),
781  newComm.indexSet(),copyrowsize);
782  ri.template redistribute<MatrixRowGatherScatter<M,IndexSet> >(origrow_copy,
783  newrow_copy);
784  ri.getInterface().free();
785  RemoteIndices<IndexSet> *ris = new RemoteIndices<IndexSet>(origComm.indexSet(),
786  newComm.indexSet(),
787  origComm.communicator());
788  ris->template rebuild<true>();
789  ri.getInterface().build(*ris,ownerflags,ownerflags);
790  }
791 
793  origrow(origMatrix, origComm.globalLookup(), newComm.indexSet());
795  newrow(newMatrix, origComm.globalLookup(), newComm.indexSet(),rowsize);
796  ri.template redistribute<MatrixRowGatherScatter<M,IndexSet> >(origrow,newrow);
797  if (SolverCategory::category(origComm) != static_cast<int>(SolverCategory::nonoverlapping))
798  newrow.setOverlapRowsToDirichlet();
799  }
800 
817  template<typename M, typename C>
818  void redistributeMatrix(M& origMatrix, M& newMatrix, C& origComm, C& newComm,
820  {
821  ri.setNoRows(newComm.indexSet().size());
822  ri.setNoCopyRows(newComm.indexSet().size());
823  ri.setNoBackwardsCopyRows(origComm.indexSet().size());
824  redistributeSparsityPattern(origMatrix, newMatrix, origComm, newComm, ri);
825  redistributeMatrixEntries(origMatrix, newMatrix, origComm, newComm, ri);
826  }
827 #endif
828 
829 template<typename M>
830  void redistributeMatrixEntries(M& origMatrix, M& newMatrix,
834  {
835  DUNE_THROW(InvalidStateException, "Trying to redistribute in sequential program!");
836  }
837  template<typename M>
838  void redistributeMatrix(M& origMatrix, M& newMatrix,
842  {
843  DUNE_THROW(InvalidStateException, "Trying to redistribute in sequential program!");
844  }
845 }
846 #endif
Classes providing communication interfaces for overlapping Schwarz methods.
Functionality for redistributing a parallel index set using graph partitioning.
Col col
Definition: matrixmatrix.hh:349
Definition: allocator.hh:9
void redistributeMatrixEntries(M &origMatrix, M &newMatrix, C &origComm, C &newComm, RedistributeInformation< C > &ri)
Definition: matrixredistribute.hh:755
void redistributeSparsityPattern(M &origMatrix, M &newMatrix, C &origComm, C &newComm, RedistributeInformation< C > &ri)
Definition: matrixredistribute.hh:661
void redistributeMatrix(M &origMatrix, M &newMatrix, C &origComm, C &newComm, RedistributeInformation< C > &ri)
Redistribute a matrix according to given domain decompositions.
Definition: matrixredistribute.hh:818
derive error class from the base class in common
Definition: istlexception.hh:17
Definition: matrixredistribute.hh:20
std::size_t getBackwardsCopyRowSize([[maybe_unused]] std::size_t index) const
Definition: matrixredistribute.hh:55
std::size_t getRowSize([[maybe_unused]] std::size_t index) const
Definition: matrixredistribute.hh:45
void setNoCopyRows([[maybe_unused]] std::size_t size)
Definition: matrixredistribute.hh:39
void redistributeBackward([[maybe_unused]] D &from, [[maybe_unused]]const D &to) const
Definition: matrixredistribute.hh:30
void setNoRows([[maybe_unused]] std::size_t size)
Definition: matrixredistribute.hh:36
std::size_t getCopyRowSize([[maybe_unused]] std::size_t index) const
Definition: matrixredistribute.hh:50
void resetSetup()
Definition: matrixredistribute.hh:33
bool isSetup() const
Definition: matrixredistribute.hh:21
void setNoBackwardsCopyRows([[maybe_unused]] std::size_t size)
Definition: matrixredistribute.hh:42
void redistribute([[maybe_unused]] const D &from, [[maybe_unused]] D &to) const
Definition: matrixredistribute.hh:26
std::size_t getRowSize(std::size_t index) const
Definition: matrixredistribute.hh:156
void redistribute(const D &from, D &to) const
Definition: matrixredistribute.hh:134
void setNoBackwardsCopyRows(std::size_t rows)
Definition: matrixredistribute.hh:191
std::size_t & getBackwardsCopyRowSize(std::size_t index)
Definition: matrixredistribute.hh:171
std::size_t getCopyRowSize(std::size_t index) const
Definition: matrixredistribute.hh:166
void setNoRows(std::size_t rows)
Definition: matrixredistribute.hh:181
RedistributeInterface & getInterface()
Definition: matrixredistribute.hh:73
void reserve(std::size_t size)
Definition: matrixredistribute.hh:148
OwnerOverlapCopyCommunication< T, T1 > Comm
Definition: matrixredistribute.hh:67
void setNoCopyRows(std::size_t rows)
Definition: matrixredistribute.hh:186
void redistributeBackward(D &from, const D &to) const
Definition: matrixredistribute.hh:139
std::size_t getBackwardsCopyRowSize(std::size_t index) const
Definition: matrixredistribute.hh:176
void redistribute(const D &from, D &to) const
Definition: matrixredistribute.hh:116
void redistributeBackward(D &from, const D &to) const
Definition: matrixredistribute.hh:124
std::size_t & getCopyRowSize(std::size_t index)
Definition: matrixredistribute.hh:161
void checkInterface(const IS &source, const IS &target, MPI_Comm comm)
Definition: matrixredistribute.hh:78
bool isSetup() const
Definition: matrixredistribute.hh:143
std::size_t & getRowSize(std::size_t index)
Definition: matrixredistribute.hh:151
Utility class to communicate and set the row sizes of a redistributed matrix.
Definition: matrixredistribute.hh:214
M::size_type size_type
Definition: matrixredistribute.hh:217
M::size_type value_type
Definition: matrixredistribute.hh:216
RI & rowsize
Definition: matrixredistribute.hh:228
const M & matrix
Definition: matrixredistribute.hh:227
CommMatrixRowSize(const M &m_, RI &rowsize_)
Constructor.
Definition: matrixredistribute.hh:224
Utility class to communicate and build the sparsity pattern of a redistributed matrix.
Definition: matrixredistribute.hh:243
M::size_type size_type
Definition: matrixredistribute.hh:244
const Dune::GlobalLookupIndexSet< I > & idxset
Definition: matrixredistribute.hh:354
void storeSparsityPattern(M &m)
Creates and stores the sparsity pattern of the redistributed matrix.
Definition: matrixredistribute.hh:274
const I & aggidxset
Definition: matrixredistribute.hh:355
const std::vector< size_type > * rowsize
Definition: matrixredistribute.hh:357
void completeSparsityPattern(std::vector< std::set< size_type > > add_sparsity)
Completes the sparsity pattern of the redistributed matrix with data from copy rows for the novlp cas...
Definition: matrixredistribute.hh:338
CommMatrixSparsityPattern(const M &m_, const Dune::GlobalLookupIndexSet< I > &idxset_, const I &aggidxset_)
Constructor for the original side.
Definition: matrixredistribute.hh:252
const M & matrix
Definition: matrixredistribute.hh:352
CommMatrixSparsityPattern(const M &m_, const Dune::GlobalLookupIndexSet< I > &idxset_, const I &aggidxset_, const std::vector< typename M::size_type > &rowsize_)
Constructor for the redistruted side.
Definition: matrixredistribute.hh:263
std::vector< std::set< size_type > > sparsity
Definition: matrixredistribute.hh:356
Dune::GlobalLookupIndexSet< I > LookupIndexSet
Definition: matrixredistribute.hh:353
static M::size_type getSize(const Type &t, std::size_t i)
Definition: matrixredistribute.hh:374
CommMatrixSparsityPattern< M, I > Type
Definition: matrixredistribute.hh:363
I::GlobalIndex IndexedType
The indexed type we send. This is the global index indentitfying the column.
Definition: matrixredistribute.hh:369
VariableSize IndexedTypeFlag
Each row varies in size.
Definition: matrixredistribute.hh:372
Utility class for comunicating the matrix entries.
Definition: matrixredistribute.hh:394
std::vector< size_t > * rowsize
row size information for the receiving side.
Definition: matrixredistribute.hh:450
M & matrix
The matrix to communicate the values of.
Definition: matrixredistribute.hh:444
CommMatrixRow(M &m_, const Dune::GlobalLookupIndexSet< I > &idxset_, const I &aggidxset_, std::vector< size_t > &rowsize_)
Constructor.
Definition: matrixredistribute.hh:410
const Dune::GlobalLookupIndexSet< I > & idxset
Index set for the original matrix.
Definition: matrixredistribute.hh:446
void setOverlapRowsToDirichlet()
Sets the non-owner rows correctly as Dirichlet boundaries.
Definition: matrixredistribute.hh:419
const I & aggidxset
Index set for the redistributed matrix.
Definition: matrixredistribute.hh:448
CommMatrixRow(M &m_, const Dune::GlobalLookupIndexSet< I > &idxset_, const I &aggidxset_)
Constructor.
Definition: matrixredistribute.hh:403
std::pair< typename I::GlobalIndex, typename M::block_type > IndexedType
The indexed type we send. This is the pair of global index indentitfying the column and the value its...
Definition: matrixredistribute.hh:462
CommMatrixRow< M, I > Type
Definition: matrixredistribute.hh:456
static std::size_t getSize(const Type &t, std::size_t i)
Definition: matrixredistribute.hh:467
VariableSize IndexedTypeFlag
Each row varies in size.
Definition: matrixredistribute.hh:465
Definition: matrixredistribute.hh:481
static void scatter(Container &cont, const typename M::size_type &rowsize, std::size_t i)
Definition: matrixredistribute.hh:488
static const M::size_type gather(const Container &cont, std::size_t i)
Definition: matrixredistribute.hh:484
CommMatrixRowSize< M, RI > Container
Definition: matrixredistribute.hh:482
Definition: matrixredistribute.hh:498
static const M::size_type gather(const Container &cont, std::size_t i)
Definition: matrixredistribute.hh:501
static void scatter(Container &cont, const typename M::size_type &rowsize, std::size_t i)
Definition: matrixredistribute.hh:505
CommMatrixRowSize< M, RI > Container
Definition: matrixredistribute.hh:499
Definition: matrixredistribute.hh:516
M::ConstColIterator ColIter
Definition: matrixredistribute.hh:519
CommMatrixSparsityPattern< M, I > Container
Definition: matrixredistribute.hh:518
static GlobalIndex numlimits
Definition: matrixredistribute.hh:522
static ColIter col
Definition: matrixredistribute.hh:521
I::GlobalIndex GlobalIndex
Definition: matrixredistribute.hh:517
static const GlobalIndex & gather(const Container &cont, std::size_t i, std::size_t j)
Definition: matrixredistribute.hh:524
static void scatter(Container &cont, const GlobalIndex &gi, std::size_t i, [[maybe_unused]] std::size_t j)
Definition: matrixredistribute.hh:551
Definition: matrixredistribute.hh:597
I::GlobalIndex GlobalIndex
Definition: matrixredistribute.hh:598
static Data datastore
Definition: matrixredistribute.hh:603
static GlobalIndex numlimits
Definition: matrixredistribute.hh:604
static const Data & gather(const Container &cont, std::size_t i, std::size_t j)
Definition: matrixredistribute.hh:606
M::ConstColIterator ColIter
Definition: matrixredistribute.hh:600
std::pair< GlobalIndex, typename M::block_type > Data
Definition: matrixredistribute.hh:601
static void scatter(Container &cont, const Data &data, std::size_t i, [[maybe_unused]] std::size_t j)
Definition: matrixredistribute.hh:636
static ColIter col
Definition: matrixredistribute.hh:602
CommMatrixRow< M, I > Container
Definition: matrixredistribute.hh:599
A class setting up standard communication for a two-valued attribute set with owner/overlap/copy sema...
Definition: owneroverlapcopy.hh:172
EnumItem< AttributeSet, OwnerOverlapCopyAttributeSet::owner > OwnerSet
Definition: owneroverlapcopy.hh:192
Definition: pinfo.hh:26
Definition: repartition.hh:258
@ nonoverlapping
Category for non-overlapping solvers.
Definition: solvercategory.hh:25
static Category category(const OP &op, decltype(op.category()) *=nullptr)
Helperfunction to extract the solver category either from an enum, or from the newly introduced virtu...
Definition: solvercategory.hh:32