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