dune-istl  2.7.0
overlappingschwarz.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_OVERLAPPINGSCHWARZ_HH
4 #define DUNE_ISTL_OVERLAPPINGSCHWARZ_HH
5 #include <cassert>
6 #include <algorithm>
7 #include <functional>
8 #include <vector>
9 #include <set>
10 #include <dune/common/dynmatrix.hh>
11 #include <dune/common/sllist.hh>
12 #include <dune/common/unused.hh>
13 #include "preconditioners.hh"
14 #include "superlu.hh"
15 #include "umfpack.hh"
16 #include "bvector.hh"
17 #include "bcrsmatrix.hh"
18 #include "ilusubdomainsolver.hh"
19 #include <dune/istl/solvertype.hh>
20 
21 namespace Dune
22 {
23 
35  template<class M, class X, class TM, class TD, class TA>
36  class SeqOverlappingSchwarz;
37 
41  template<class I, class S, class D>
43  {
44  public:
46  typedef D subdomain_vector;
47 
48  typedef I InitializerList;
49  typedef typename InitializerList::value_type AtomInitializer;
50  typedef typename AtomInitializer::Matrix Matrix;
51  typedef typename Matrix::const_iterator Iter;
52  typedef typename Matrix::row_type::const_iterator CIter;
53 
54  typedef S IndexSet;
55  typedef typename IndexSet::size_type size_type;
56 
58  const IndexSet& indices,
59  const subdomain_vector& domains);
60 
61 
62  void addRowNnz(const Iter& row);
63 
64  void allocate();
65 
66  void countEntries(const Iter& row, const CIter& col) const;
67 
68  void calcColstart() const;
69 
70  void copyValue(const Iter& row, const CIter& col) const;
71 
72  void createMatrix() const;
73  private:
74  class IndexMap
75  {
76  public:
77  typedef typename S::size_type size_type;
78  typedef std::map<size_type,size_type> Map;
79  typedef typename Map::iterator iterator;
80  typedef typename Map::const_iterator const_iterator;
81 
82  IndexMap();
83 
84  void insert(size_type grow);
85 
86  const_iterator find(size_type grow) const;
87 
88  iterator find(size_type grow);
89 
90  iterator begin();
91 
92  const_iterator begin() const;
93 
94  iterator end();
95 
96  const_iterator end() const;
97 
98  private:
99  std::map<size_type,size_type> map_;
100  size_type row;
101  };
102 
103 
104  typedef typename InitializerList::iterator InitIterator;
105  typedef typename IndexSet::const_iterator IndexIteratur;
106  InitializerList* initializers;
107  const IndexSet *indices;
108  mutable std::vector<IndexMap> indexMaps;
109  const subdomain_vector& domains;
110  };
111 
116  {};
117 
122  {};
123 
129  {};
130 
135  template<class M, class X, class Y>
137 
138  // Specialization for BCRSMatrix
139  template<class K, class Al, class X, class Y>
141  {
142  typedef BCRSMatrix< K, Al> M;
143  public:
145  typedef typename std::remove_const<M>::type matrix_type;
146  typedef typename X::field_type field_type;
147  typedef typename std::remove_const<M>::type rilu_type;
149  typedef X domain_type;
151  typedef Y range_type;
152  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
153 
158  void apply (DynamicVector<field_type>& v, DynamicVector<field_type>& d)
159  {
160  assert(v.size() > 0);
161  assert(v.size() == d.size());
162  assert(A.rows() <= v.size());
163  assert(A.cols() <= v.size());
164  size_t sz = A.rows();
165  v.resize(sz);
166  d.resize(sz);
167  A.solve(v,d);
168  v.resize(v.capacity());
169  d.resize(d.capacity());
170  }
171 
179  template<class S>
180  void setSubMatrix(const M& BCRS, S& rowset)
181  {
182  size_t sz = rowset.size();
183  A.resize(sz*n,sz*n);
184  typedef typename S::const_iterator SIter;
185  size_t r = 0, c = 0;
186  for(SIter rowIdx = rowset.begin(), rowEnd=rowset.end();
187  rowIdx!= rowEnd; ++rowIdx, r++)
188  {
189  c = 0;
190  for(SIter colIdx = rowset.begin(), colEnd=rowset.end();
191  colIdx!= colEnd; ++colIdx, c++)
192  {
193  if (BCRS[*rowIdx].find(*colIdx) == BCRS[*rowIdx].end())
194  continue;
195  for (size_t i=0; i<n; i++)
196  {
197  for (size_t j=0; j<n; j++)
198  {
199  A[r*n+i][c*n+j] = Impl::asMatrix(BCRS[*rowIdx][*colIdx])[i][j];
200  }
201  }
202  }
203  }
204  }
205  private:
206  DynamicMatrix<K> A;
207  };
208 
209  template<typename T, bool tag>
211  {};
212 
213  template<typename T>
215 
216  // specialization for DynamicMatrix
217  template<class K, class Al, class X, class Y>
219  {
220  public:
222  typedef typename X::field_type field_type;
223  typedef Y range_type;
224  typedef typename range_type::block_type block_type;
226  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
234  OverlappingAssignerHelper(std::size_t maxlength, const BCRSMatrix<K, Al>& mat_, const X& b_, Y& x_);
235 
239  inline
240  void deallocate();
241 
245  inline
246  void resetIndexForNextDomain();
247 
252  inline
253  DynamicVector<field_type> & lhs();
254 
259  inline
260  DynamicVector<field_type> & rhs();
261 
266  inline
267  void relaxResult(field_type relax);
268 
273  void operator()(const size_type& domainIndex);
274 
282  inline
283  void assignResult(block_type& res);
284 
285  private:
289  const matrix_type* mat;
291  // we need a pointer, because we have to avoid deep copies
292  DynamicVector<field_type> * rhs_;
294  // we need a pointer, because we have to avoid deep copies
295  DynamicVector<field_type> * lhs_;
297  const range_type* b;
299  range_type* x;
301  std::size_t i;
303  std::size_t maxlength_;
304  };
305 
306 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
307  template<template<class> class S, typename T, typename A>
308  struct OverlappingAssignerHelper<S<BCRSMatrix<T, A>>, true>
309  {
311  typedef typename S<BCRSMatrix<T, A>>::range_type range_type;
312  typedef typename range_type::field_type field_type;
313  typedef typename range_type::block_type block_type;
314 
316 
317  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
318  static constexpr size_t m = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::cols;
326  OverlappingAssignerHelper(std::size_t maxlength, const matrix_type& mat,
327  const range_type& b, range_type& x);
333  void deallocate();
334 
335  /*
336  * @brief Resets the local index to zero.
337  */
338  void resetIndexForNextDomain();
339 
344  field_type* lhs();
345 
350  field_type* rhs();
351 
356  void relaxResult(field_type relax);
357 
362  void operator()(const size_type& domain);
363 
371  void assignResult(block_type& res);
372 
373  private:
377  const matrix_type* mat;
379  field_type* rhs_;
381  field_type* lhs_;
383  const range_type* b;
385  range_type* x;
387  std::size_t i;
389  std::size_t maxlength_;
390  };
391 
392 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
393 
394  template<class M, class X, class Y>
396  {
397  public:
398  typedef M matrix_type;
399 
400  typedef typename Y::field_type field_type;
401 
402  typedef typename Y::block_type block_type;
403 
404  typedef typename matrix_type::size_type size_type;
412  OverlappingAssignerILUBase(std::size_t maxlength, const M& mat,
413  const Y& b, X& x);
419  void deallocate();
420 
425 
430  X& lhs();
431 
436  Y& rhs();
437 
442  void relaxResult(field_type relax);
443 
448  void operator()(const size_type& domain);
449 
457  void assignResult(block_type& res);
458 
459  private:
463  const M* mat;
465  X* lhs_;
467  Y* rhs_;
469  const Y* b;
471  X* x;
473  size_type i;
474  };
475 
476  // specialization for ILU0
477  template<class M, class X, class Y>
479  : public OverlappingAssignerILUBase<M,X,Y>
480  {
481  public:
489  OverlappingAssignerHelper(std::size_t maxlength, const M& mat,
490  const Y& b, X& x)
491  : OverlappingAssignerILUBase<M,X,Y>(maxlength, mat,b,x)
492  {}
493  };
494 
495  // specialization for ILUN
496  template<class M, class X, class Y>
498  : public OverlappingAssignerILUBase<M,X,Y>
499  {
500  public:
508  OverlappingAssignerHelper(std::size_t maxlength, const M& mat,
509  const Y& b, X& x)
510  : OverlappingAssignerILUBase<M,X,Y>(maxlength, mat,b,x)
511  {}
512  };
513 
514  template<typename S, typename T>
516  {};
517 
518  template<typename S, typename T, typename A>
519  struct AdditiveAdder<S, BlockVector<T,A> >
520  {
521  typedef typename A::size_type size_type;
522  typedef typename std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::field_type field_type;
524  OverlappingAssigner<S>& assigner, const field_type& relax_);
525  void operator()(const size_type& domain);
526  void axpy();
527  static constexpr size_t n = std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::dimension;
528 
529  private:
530  BlockVector<T,A>* v;
531  BlockVector<T,A>* x;
532  OverlappingAssigner<S>* assigner;
533  field_type relax;
534  };
535 
536  template<typename S,typename T>
538  {};
539 
540  template<typename S, typename T, typename A>
542  {
543  typedef typename A::size_type size_type;
544  typedef typename std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::field_type field_type;
546  OverlappingAssigner<S>& assigner_, const field_type& relax_);
547  void operator()(const size_type& domain);
548  void axpy();
549  static constexpr size_t n = std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::dimension;
550 
551  private:
552  BlockVector<T,A>* x;
553  OverlappingAssigner<S>* assigner;
554  field_type relax;
555  };
556 
566  template<typename T, class X, class S>
568  {};
569 
570  template<class X, class S>
572  {
574  };
575 
576  template<class X, class S>
578  {
580  };
581 
582  template<class X, class S>
584  {
586  };
587 
599  template<typename T1, typename T2, bool forward>
601  {
602  typedef T1 solver_vector;
603  typedef typename solver_vector::iterator solver_iterator;
604  typedef T2 subdomain_vector;
605  typedef typename subdomain_vector::const_iterator domain_iterator;
606 
608  {
609  return sv.begin();
610  }
611 
613  {
614  return sv.end();
615  }
617  {
618  return sv.begin();
619  }
620 
622  {
623  return sv.end();
624  }
625  };
626 
627  template<typename T1, typename T2>
628  struct IteratorDirectionSelector<T1,T2,false>
629  {
630  typedef T1 solver_vector;
631  typedef typename solver_vector::reverse_iterator solver_iterator;
632  typedef T2 subdomain_vector;
633  typedef typename subdomain_vector::const_reverse_iterator domain_iterator;
634 
636  {
637  return sv.rbegin();
638  }
639 
641  {
642  return sv.rend();
643  }
645  {
646  return sv.rbegin();
647  }
648 
650  {
651  return sv.rend();
652  }
653  };
654 
663  template<class T>
665  {
666  typedef T smoother;
667  typedef typename smoother::range_type range_type;
668 
669  static void apply(smoother& sm, range_type& v, const range_type& b)
670  {
671  sm.template apply<true>(v, b);
672  }
673  };
674 
675  template<class M, class X, class TD, class TA>
677  {
680 
681  static void apply(smoother& sm, range_type& v, const range_type& b)
682  {
683  sm.template apply<true>(v, b);
684  sm.template apply<false>(v, b);
685  }
686  };
687 
688  template<class T, bool tag>
689  struct SeqOverlappingSchwarzAssemblerHelper
690  {};
691 
692  template<class T>
694 
695  template<class K, class Al, class X, class Y>
697  {
699  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
700  template<class RowToDomain, class Solvers, class SubDomains>
701  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
702  Solvers& solvers, const SubDomains& domains,
703  bool onTheFly);
704  };
705 
706  template<template<class> class S, typename T, typename A>
708  {
710  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
711  template<class RowToDomain, class Solvers, class SubDomains>
712  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
713  Solvers& solvers, const SubDomains& domains,
714  bool onTheFly);
715  };
716 
717  template<class M,class X, class Y>
719  {
720  typedef M matrix_type;
721  template<class RowToDomain, class Solvers, class SubDomains>
722  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
723  Solvers& solvers, const SubDomains& domains,
724  bool onTheFly);
725  };
726 
727  template<class M,class X, class Y>
730  {};
731 
732  template<class M,class X, class Y>
735  {};
736 
747  template<class M, class X, class TM=AdditiveSchwarzMode,
748  class TD=ILU0SubdomainSolver<M,X,X>, class TA=std::allocator<X> >
750  : public Preconditioner<X,X>
751  {
752  public:
756  typedef M matrix_type;
757 
761  typedef X domain_type;
762 
766  typedef X range_type;
767 
774  typedef TM Mode;
775 
779  typedef typename X::field_type field_type;
780 
782  typedef typename matrix_type::size_type size_type;
783 
785  typedef TA allocator;
786 
788  typedef std::set<size_type, std::less<size_type>,
789  typename TA::template rebind<size_type>::other>
791 
793  typedef std::vector<subdomain_type, typename TA::template rebind<subdomain_type>::other> subdomain_vector;
794 
796  typedef SLList<size_type, typename TA::template rebind<size_type>::other> subdomain_list;
797 
799  typedef std::vector<subdomain_list, typename TA::template rebind<subdomain_list>::other > rowtodomain_vector;
800 
802  typedef TD slu;
803 
805  typedef std::vector<slu, typename TA::template rebind<slu>::other> slu_vector;
806 
820  SeqOverlappingSchwarz(const matrix_type& mat, const subdomain_vector& subDomains,
821  field_type relaxationFactor=1, bool onTheFly_=true);
822 
834  SeqOverlappingSchwarz(const matrix_type& mat, const rowtodomain_vector& rowToDomain,
835  field_type relaxationFactor=1, bool onTheFly_=true);
836 
842  virtual void pre (X& x, X& b)
843  {
844  DUNE_UNUSED_PARAMETER(x);
845  DUNE_UNUSED_PARAMETER(b);
846  }
847 
853  virtual void apply (X& v, const X& d);
854 
860  virtual void post (X& x)
861  {
862  DUNE_UNUSED_PARAMETER(x);
863  }
864 
865  template<bool forward>
866  void apply(X& v, const X& d);
867 
870  {
872  }
873 
874  private:
875  const M& mat;
876  slu_vector solvers;
877  subdomain_vector subDomains;
878  field_type relax;
879 
880  typename M::size_type maxlength;
881 
882  bool onTheFly;
883  };
884 
885 
886 
887  template<class I, class S, class D>
889  const IndexSet& idx,
890  const subdomain_vector& domains_)
891  : initializers(&il), indices(&idx), indexMaps(il.size()), domains(domains_)
892  {}
893 
894 
895  template<class I, class S, class D>
897  {
898  typedef typename IndexSet::value_type::const_iterator iterator;
899  for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
900  (*initializers)[*domain].addRowNnz(row, domains[*domain]);
901  indexMaps[*domain].insert(row.index());
902  }
903  }
904 
905  template<class I, class S, class D>
907  {
908  for(auto&& i: *initializers)
909  i.allocateMatrixStorage();
910  for(auto&& i: *initializers)
911  i.allocateMarker();
912  }
913 
914  template<class I, class S, class D>
916  {
917  typedef typename IndexSet::value_type::const_iterator iterator;
918  for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
919  typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.index());
920  if(v!= indexMaps[*domain].end()) {
921  (*initializers)[*domain].countEntries(indexMaps[*domain].find(col.index())->second);
922  }
923  }
924  }
925 
926  template<class I, class S, class D>
928  {
929  for(auto&& i : *initializers)
930  i.calcColstart();
931  }
932 
933  template<class I, class S, class D>
935  {
936  typedef typename IndexSet::value_type::const_iterator iterator;
937  for(iterator domain=(*indices)[row.index()].begin(); domain!= (*indices)[row.index()].end(); ++domain) {
938  typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.index());
939  if(v!= indexMaps[*domain].end()) {
940  assert(indexMaps[*domain].end()!=indexMaps[*domain].find(row.index()));
941  (*initializers)[*domain].copyValue(col, indexMaps[*domain].find(row.index())->second,
942  v->second);
943  }
944  }
945  }
946 
947  template<class I, class S, class D>
949  {
950  std::vector<IndexMap>().swap(indexMaps);
951  for(auto&& i: *initializers)
952  i.createMatrix();
953  }
954 
955  template<class I, class S, class D>
957  : row(0)
958  {}
959 
960  template<class I, class S, class D>
962  {
963  assert(map_.find(grow)==map_.end());
964  map_.insert(std::make_pair(grow, row++));
965  }
966 
967  template<class I, class S, class D>
968  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
970  {
971  return map_.find(grow);
972  }
973 
974  template<class I, class S, class D>
975  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
977  {
978  return map_.find(grow);
979  }
980 
981  template<class I, class S, class D>
982  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
984  {
985  return map_.end();
986  }
987 
988  template<class I, class S, class D>
989  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
991  {
992  return map_.end();
993  }
994 
995  template<class I, class S, class D>
996  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
998  {
999  return map_.begin();
1000  }
1001 
1002  template<class I, class S, class D>
1003  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
1005  {
1006  return map_.begin();
1007  }
1008 
1009  template<class M, class X, class TM, class TD, class TA>
1011  field_type relaxationFactor, bool fly)
1012  : mat(mat_), relax(relaxationFactor), onTheFly(fly)
1013  {
1014  typedef typename rowtodomain_vector::const_iterator RowDomainIterator;
1015  typedef typename subdomain_list::const_iterator DomainIterator;
1016 #ifdef DUNE_ISTL_WITH_CHECKING
1017  assert(rowToDomain.size()==mat.N());
1018  assert(rowToDomain.size()==mat.M());
1019 
1020  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1021  assert(iter->size()>0);
1022 
1023 #endif
1024  // calculate the number of domains
1025  size_type domains=0;
1026  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1027  for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1028  domains=std::max(domains, *d);
1029  ++domains;
1030 
1031  solvers.resize(domains);
1032  subDomains.resize(domains);
1033 
1034  // initialize subdomains to row mapping from row to subdomain mapping
1035  size_type row=0;
1036  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter, ++row)
1037  for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1038  subDomains[*d].insert(row);
1039 
1040 #ifdef DUNE_ISTL_WITH_CHECKING
1041  size_type i=0;
1042  typedef typename subdomain_vector::const_iterator iterator;
1043  for(iterator iter=subDomains.begin(); iter != subDomains.end(); ++iter) {
1044  typedef typename subdomain_type::const_iterator entry_iterator;
1045  Dune::dvverb<<"domain "<<i++<<":";
1046  for(entry_iterator entry = iter->begin(); entry != iter->end(); ++entry) {
1047  Dune::dvverb<<" "<<*entry;
1048  }
1049  Dune::dvverb<<std::endl;
1050  }
1051 #endif
1053  ::assembleLocalProblems(rowToDomain, mat, solvers, subDomains, onTheFly);
1054  }
1055 
1056  template<class M, class X, class TM, class TD, class TA>
1058  const subdomain_vector& sd,
1059  field_type relaxationFactor,
1060  bool fly)
1061  : mat(mat_), solvers(sd.size()), subDomains(sd), relax(relaxationFactor),
1062  onTheFly(fly)
1063  {
1064  typedef typename subdomain_vector::const_iterator DomainIterator;
1065 
1066 #ifdef DUNE_ISTL_WITH_CHECKING
1067  size_type i=0;
1068 
1069  for(DomainIterator d=sd.begin(); d != sd.end(); ++d,++i) {
1070  //std::cout<<i<<": "<<d->size()<<std::endl;
1071  assert(d->size()>0);
1072  typedef typename DomainIterator::value_type::const_iterator entry_iterator;
1073  Dune::dvverb<<"domain "<<i<<":";
1074  for(entry_iterator entry = d->begin(); entry != d->end(); ++entry) {
1075  Dune::dvverb<<" "<<*entry;
1076  }
1077  Dune::dvverb<<std::endl;
1078  }
1079 
1080 #endif
1081 
1082  // Create a row to subdomain mapping
1083  rowtodomain_vector rowToDomain(mat.N());
1084 
1085  size_type domainId=0;
1086 
1087  for(DomainIterator domain=sd.begin(); domain != sd.end(); ++domain, ++domainId) {
1088  typedef typename subdomain_type::const_iterator iterator;
1089  for(iterator row=domain->begin(); row != domain->end(); ++row)
1090  rowToDomain[*row].push_back(domainId);
1091  }
1092 
1094  ::assembleLocalProblems(rowToDomain, mat, solvers, subDomains, onTheFly);
1095  }
1096 
1103  template<class M>
1105 
1106  template<typename T, typename A>
1108  {
1109  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
1110  static constexpr size_t m = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::cols;
1111  template<class Domain>
1112  static int size(const Domain & d)
1113  {
1114  assert(n==m);
1115  return m*d.size();
1116  }
1117  };
1118 
1119  template<class K, class Al, class X, class Y>
1120  template<class RowToDomain, class Solvers, class SubDomains>
1121  std::size_t
1122  SeqOverlappingSchwarzAssemblerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>::
1123  assembleLocalProblems(const RowToDomain& rowToDomain,
1124  const matrix_type& mat,
1125  Solvers& solvers,
1126  const SubDomains& subDomains,
1127  bool onTheFly)
1128  {
1129  DUNE_UNUSED_PARAMETER(onTheFly);
1130  DUNE_UNUSED_PARAMETER(rowToDomain);
1131  DUNE_UNUSED_PARAMETER(mat);
1132  DUNE_UNUSED_PARAMETER(solvers);
1133  typedef typename SubDomains::const_iterator DomainIterator;
1134  std::size_t maxlength = 0;
1135 
1136  assert(onTheFly);
1137 
1138  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1139  maxlength=std::max(maxlength, domain->size());
1140  maxlength*=n;
1141 
1142  return maxlength;
1143  }
1144 
1145 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1146  template<template<class> class S, typename T, typename A>
1147  template<class RowToDomain, class Solvers, class SubDomains>
1148  std::size_t SeqOverlappingSchwarzAssemblerHelper<S<BCRSMatrix<T,A>>,true>::assembleLocalProblems(const RowToDomain& rowToDomain,
1149  const matrix_type& mat,
1150  Solvers& solvers,
1151  const SubDomains& subDomains,
1152  bool onTheFly)
1153  {
1154  typedef typename S<BCRSMatrix<T,A>>::MatrixInitializer MatrixInitializer;
1155  typedef typename std::vector<MatrixInitializer>::iterator InitializerIterator;
1156  typedef typename SubDomains::const_iterator DomainIterator;
1157  typedef typename Solvers::iterator SolverIterator;
1158  std::size_t maxlength = 0;
1159 
1160  if(onTheFly) {
1161  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1162  maxlength=std::max(maxlength, domain->size());
1163  maxlength*=Impl::asMatrix(*mat[0].begin()).N();
1164  }else{
1165  // initialize the initializers
1166  DomainIterator domain=subDomains.begin();
1167 
1168  // Create the initializers list.
1169  std::vector<MatrixInitializer> initializers(subDomains.size());
1170 
1171  SolverIterator solver=solvers.begin();
1172  for(InitializerIterator initializer=initializers.begin(); initializer!=initializers.end();
1173  ++initializer, ++solver, ++domain) {
1174  solver->getInternalMatrix().N_=SeqOverlappingSchwarzDomainSize<matrix_type>::size(*domain);
1175  solver->getInternalMatrix().M_=SeqOverlappingSchwarzDomainSize<matrix_type>::size(*domain);
1176  //solver->setVerbosity(true);
1177  *initializer=MatrixInitializer(solver->getInternalMatrix());
1178  }
1179 
1180  // Set up the supermatrices according to the subdomains
1182  RowToDomain, SubDomains> Initializer;
1183 
1184  Initializer initializer(initializers, rowToDomain, subDomains);
1185  copyToColCompMatrix(initializer, mat);
1186 
1187  // Calculate the LU decompositions
1188  for(auto&& s: solvers)
1189  s.decompose();
1190  for (SolverIterator solverIt = solvers.begin(); solverIt != solvers.end(); ++solverIt)
1191  {
1192  assert(solverIt->getInternalMatrix().N() == solverIt->getInternalMatrix().M());
1193  maxlength = std::max(maxlength, solverIt->getInternalMatrix().N());
1194  }
1195  }
1196  return maxlength;
1197  }
1198 
1199 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1200 
1201  template<class M,class X,class Y>
1202  template<class RowToDomain, class Solvers, class SubDomains>
1204  const matrix_type& mat,
1205  Solvers& solvers,
1206  const SubDomains& subDomains,
1207  bool onTheFly)
1208  {
1209  DUNE_UNUSED_PARAMETER(rowToDomain);
1210  typedef typename SubDomains::const_iterator DomainIterator;
1211  typedef typename Solvers::iterator SolverIterator;
1212  std::size_t maxlength = 0;
1213 
1214  if(onTheFly) {
1215  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1216  maxlength=std::max(maxlength, domain->size());
1217  }else{
1218  // initialize the solvers of the local prolems.
1219  SolverIterator solver=solvers.begin();
1220  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end();
1221  ++domain, ++solver) {
1222  solver->setSubMatrix(mat, *domain);
1223  maxlength=std::max(maxlength, domain->size());
1224  }
1225  }
1226 
1227  return maxlength;
1228 
1229  }
1230 
1231 
1232  template<class M, class X, class TM, class TD, class TA>
1234  {
1236  }
1237 
1238  template<class M, class X, class TM, class TD, class TA>
1239  template<bool forward>
1240  void SeqOverlappingSchwarz<M,X,TM,TD,TA>::apply(X& x, const X& b)
1241  {
1242  typedef slu_vector solver_vector;
1245  domain_iterator;
1246 
1247  OverlappingAssigner<TD> assigner(maxlength, mat, b, x);
1248 
1251  X v(x); // temporary for the update
1252  v=0;
1253 
1254  typedef typename AdderSelector<TM,X,TD >::Adder Adder;
1255  Adder adder(v, x, assigner, relax);
1256 
1257  for(; domain != IteratorDirectionSelector<solver_vector,subdomain_vector,forward>::end(subDomains); ++domain) {
1258  //Copy rhs to C-array for SuperLU
1259  std::for_each(domain->begin(), domain->end(), assigner);
1260  assigner.resetIndexForNextDomain();
1261  if(onTheFly) {
1262  // Create the subdomain solver
1263  slu sdsolver;
1264  sdsolver.setSubMatrix(mat, *domain);
1265  // Apply
1266  sdsolver.apply(assigner.lhs(), assigner.rhs());
1267  }else{
1268  solver->apply(assigner.lhs(), assigner.rhs());
1269  ++solver;
1270  }
1271 
1272  //Add relaxed correction to from SuperLU to v
1273  std::for_each(domain->begin(), domain->end(), adder);
1274  assigner.resetIndexForNextDomain();
1275 
1276  }
1277 
1278  adder.axpy();
1279  assigner.deallocate();
1280  }
1281 
1282  template<class K, class Al, class X, class Y>
1283  OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1284  ::OverlappingAssignerHelper(std::size_t maxlength, const BCRSMatrix<K, Al>& mat_,
1285  const X& b_, Y& x_) :
1286  mat(&mat_),
1287  rhs_( new DynamicVector<field_type>(maxlength, 42) ),
1288  lhs_( new DynamicVector<field_type>(maxlength, -42) ),
1289  b(&b_),
1290  x(&x_),
1291  i(0),
1292  maxlength_(maxlength)
1293  {}
1294 
1295  template<class K, class Al, class X, class Y>
1296  void
1298  ::deallocate()
1299  {
1300  delete rhs_;
1301  delete lhs_;
1302  }
1303 
1304  template<class K, class Al, class X, class Y>
1305  void
1307  ::resetIndexForNextDomain()
1308  {
1309  i=0;
1310  }
1311 
1312  template<class K, class Al, class X, class Y>
1313  DynamicVector<typename X::field_type> &
1315  ::lhs()
1316  {
1317  return *lhs_;
1318  }
1319 
1320  template<class K, class Al, class X, class Y>
1321  DynamicVector<typename X::field_type> &
1323  ::rhs()
1324  {
1325  return *rhs_;
1326  }
1327 
1328  template<class K, class Al, class X, class Y>
1329  void
1331  ::relaxResult(field_type relax)
1332  {
1333  lhs() *= relax;
1334  }
1335 
1336  template<class K, class Al, class X, class Y>
1337  void
1339  ::operator()(const size_type& domainIndex)
1340  {
1341  lhs() = 0.0;
1342 #if 0
1343  //assign right hand side of current domainindex block
1344  for(size_type j=0; j<n; ++j, ++i) {
1345  assert(i<maxlength_);
1346  rhs()[i]=(*b)[domainIndex][j];
1347  }
1348 
1349  // loop over all Matrix row entries and calculate defect.
1350  typedef typename matrix_type::ConstColIterator col_iterator;
1351 
1352  // calculate defect for current row index block
1353  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1354  block_type tmp(0.0);
1355  (*col).mv((*x)[col.index()], tmp);
1356  i-=n;
1357  for(size_type j=0; j<n; ++j, ++i) {
1358  assert(i<maxlength_);
1359  rhs()[i]-=tmp[j];
1360  }
1361  }
1362 #else
1363  //assign right hand side of current domainindex block
1364  for(size_type j=0; j<n; ++j, ++i) {
1365  assert(i<maxlength_);
1366  rhs()[i]=Impl::asVector((*b)[domainIndex])[j];
1367 
1368  // loop over all Matrix row entries and calculate defect.
1369  typedef typename matrix_type::ConstColIterator col_iterator;
1370 
1371  // calculate defect for current row index block
1372  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1373  for(size_type k=0; k<n; ++k) {
1374  rhs()[i]-=Impl::asMatrix(*col)[j][k] * Impl::asVector((*x)[col.index()])[k];
1375  }
1376  }
1377  }
1378 #endif
1379  }
1380 
1381  template<class K, class Al, class X, class Y>
1382  void
1384  ::assignResult(block_type& res)
1385  {
1386  // assign the result of the local solve to the global vector
1387  for(size_type j=0; j<n; ++j, ++i) {
1388  assert(i<maxlength_);
1389  Impl::asVector(res)[j]+=lhs()[i];
1390  }
1391  }
1392 
1393 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1394 
1395  template<template<class> class S, typename T, typename A>
1397  ::OverlappingAssignerHelper(std::size_t maxlength,
1398  const BCRSMatrix<T,A>& mat_,
1399  const range_type& b_,
1400  range_type& x_)
1401  : mat(&mat_),
1402  b(&b_),
1403  x(&x_), i(0), maxlength_(maxlength)
1404  {
1405  rhs_ = new field_type[maxlength];
1406  lhs_ = new field_type[maxlength];
1407 
1408  }
1409 
1410  template<template<class> class S, typename T, typename A>
1412  {
1413  delete[] rhs_;
1414  delete[] lhs_;
1415  }
1416 
1417  template<template<class> class S, typename T, typename A>
1418  void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::operator()(const size_type& domainIndex)
1419  {
1420  //assign right hand side of current domainindex block
1421  // rhs is an array of doubles!
1422  // rhs[starti] = b[domainindex]
1423  for(size_type j=0; j<n; ++j, ++i) {
1424  assert(i<maxlength_);
1425  rhs_[i]=Impl::asVector((*b)[domainIndex])[j];
1426  }
1427 
1428 
1429  // loop over all Matrix row entries and calculate defect.
1430  typedef typename matrix_type::ConstColIterator col_iterator;
1431 
1432  // calculate defect for current row index block
1433  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1434  block_type tmp;
1435  Impl::asMatrix(*col).mv((*x)[col.index()], tmp);
1436  i-=n;
1437  for(size_type j=0; j<n; ++j, ++i) {
1438  assert(i<maxlength_);
1439  rhs_[i]-=Impl::asVector(tmp)[j];
1440  }
1441 
1442  }
1443 
1444  }
1445 
1446  template<template<class> class S, typename T, typename A>
1448  {
1449  for(size_type j=i+n; i<j; ++i) {
1450  assert(i<maxlength_);
1451  lhs_[i]*=relax;
1452  }
1453  i-=n;
1454  }
1455 
1456  template<template<class> class S, typename T, typename A>
1458  {
1459  // assign the result of the local solve to the global vector
1460  for(size_type j=0; j<n; ++j, ++i) {
1461  assert(i<maxlength_);
1462  Impl::asVector(res)[j]+=lhs_[i];
1463  }
1464  }
1465 
1466  template<template<class> class S, typename T, typename A>
1467  void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::resetIndexForNextDomain()
1468  {
1469  i=0;
1470  }
1471 
1472  template<template<class> class S, typename T, typename A>
1473  typename OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::field_type*
1475  {
1476  return lhs_;
1477  }
1478 
1479  template<template<class> class S, typename T, typename A>
1480  typename OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::field_type*
1482  {
1483  return rhs_;
1484  }
1485 
1486 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1487 
1488  template<class M, class X, class Y>
1490  const M& mat_,
1491  const Y& b_,
1492  X& x_)
1493  : mat(&mat_),
1494  b(&b_),
1495  x(&x_), i(0)
1496  {
1497  rhs_= new Y(maxlength);
1498  lhs_ = new X(maxlength);
1499  }
1500 
1501  template<class M, class X, class Y>
1503  {
1504  delete rhs_;
1505  delete lhs_;
1506  }
1507 
1508  template<class M, class X, class Y>
1510  {
1511  (*rhs_)[i]=(*b)[domainIndex];
1512 
1513  // loop over all Matrix row entries and calculate defect.
1514  typedef typename matrix_type::ConstColIterator col_iterator;
1515 
1516  // calculate defect for current row index block
1517  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1518  Impl::asMatrix(*col).mmv((*x)[col.index()], (*rhs_)[i]);
1519  }
1520  // Goto next local index
1521  ++i;
1522  }
1523 
1524  template<class M, class X, class Y>
1526  {
1527  (*lhs_)[i]*=relax;
1528  }
1529 
1530  template<class M, class X, class Y>
1532  {
1533  res+=(*lhs_)[i++];
1534  }
1535 
1536  template<class M, class X, class Y>
1538  {
1539  return *lhs_;
1540  }
1541 
1542  template<class M, class X, class Y>
1544  {
1545  return *rhs_;
1546  }
1547 
1548  template<class M, class X, class Y>
1550  {
1551  i=0;
1552  }
1553 
1554  template<typename S, typename T, typename A>
1556  BlockVector<T,A>& x_,
1557  OverlappingAssigner<S>& assigner_,
1558  const field_type& relax_)
1559  : v(&v_), x(&x_), assigner(&assigner_), relax(relax_)
1560  {}
1561 
1562  template<typename S, typename T, typename A>
1563  void AdditiveAdder<S,BlockVector<T,A> >::operator()(const size_type& domainIndex)
1564  {
1565  // add the result of the local solve to the current update
1566  assigner->assignResult((*v)[domainIndex]);
1567  }
1568 
1569 
1570  template<typename S, typename T, typename A>
1572  {
1573  // relax the update and add it to the current guess.
1574  x->axpy(relax,*v);
1575  }
1576 
1577 
1578  template<typename S, typename T, typename A>
1581  BlockVector<T,A>& x_,
1582  OverlappingAssigner<S>& assigner_, const field_type& relax_)
1583  : x(&x_), assigner(&assigner_), relax(relax_)
1584  {
1585  DUNE_UNUSED_PARAMETER(v_);
1586  }
1587 
1588 
1589  template<typename S,typename T, typename A>
1590  void MultiplicativeAdder<S,BlockVector<T,A> >::operator()(const size_type& domainIndex)
1591  {
1592  // add the result of the local solve to the current guess
1593  assigner->relaxResult(relax);
1594  assigner->assignResult((*x)[domainIndex]);
1595  }
1596 
1597 
1598  template<typename S,typename T, typename A>
1600  {
1601  // nothing to do, as the corrections already relaxed and added in operator()
1602  }
1603 
1604 
1606 }
1607 
1608 #endif
Dune::SeqOverlappingSchwarz::SeqOverlappingSchwarz
SeqOverlappingSchwarz(const matrix_type &mat, const subdomain_vector &subDomains, field_type relaxationFactor=1, bool onTheFly_=true)
Construct the overlapping Schwarz method.
Definition: overlappingschwarz.hh:1057
Dune::OverlappingAssignerILUBase::assignResult
void assignResult(block_type &res)
Assigns the block to the current local index. At the same time the local defect is calculated for the...
Definition: overlappingschwarz.hh:1531
Dune::SeqOverlappingSchwarz::allocator
TA allocator
The allocator to use.
Definition: overlappingschwarz.hh:785
Dune::OverlappingAssignerHelper
Definition: overlappingschwarz.hh:210
Dune::AdderSelector< SymmetricMultiplicativeSchwarzMode, X, S >::Adder
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:585
Dune::IteratorDirectionSelector::end
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:621
col
Col col
Definition: matrixmatrix.hh:349
Dune::SolverCategory::sequential
@ sequential
Category for sequential solvers.
Definition: solvercategory.hh:23
Dune::AdderSelector< MultiplicativeSchwarzMode, X, S >::Adder
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:579
Dune::IteratorDirectionSelector::subdomain_vector
T2 subdomain_vector
Definition: overlappingschwarz.hh:604
Dune::IteratorDirectionSelector::solver_vector
T1 solver_vector
Definition: overlappingschwarz.hh:602
Dune::SeqOverlappingSchwarzAssemblerILUBase
Definition: overlappingschwarz.hh:718
Dune::AdditiveSchwarzMode
Tag that the tells the Schwarz method to be additive.
Definition: overlappingschwarz.hh:115
Dune::OverlappingSchwarzInitializer::OverlappingSchwarzInitializer
OverlappingSchwarzInitializer(InitializerList &il, const IndexSet &indices, const subdomain_vector &domains)
Definition: overlappingschwarz.hh:888
Dune::SeqOverlappingSchwarzApplier
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:664
Dune::IteratorDirectionSelector::begin
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:616
Dune::ILU0SubdomainSolver
Exact subdomain solver using ILU(p) with appropriate p.
Definition: ilusubdomainsolver.hh:74
Dune::MultiplicativeAdder< S, BlockVector< T, A > >::field_type
std::decay_t< decltype(Impl::asVector(std::declval< T >)))>::field_type field_type
Definition: overlappingschwarz.hh:544
Dune::SeqOverlappingSchwarz
Sequential overlapping Schwarz preconditioner.
Definition: colcompmatrix.hh:150
Dune::OverlappingSchwarzInitializer::createMatrix
void createMatrix() const
Definition: overlappingschwarz.hh:948
Dune::SeqOverlappingSchwarz::pre
virtual void pre(X &x, X &b)
Prepare the preconditioner.
Definition: overlappingschwarz.hh:842
Dune::OverlappingSchwarzInitializer::Matrix
AtomInitializer::Matrix Matrix
Definition: overlappingschwarz.hh:50
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::field_type
X::field_type field_type
Definition: overlappingschwarz.hh:146
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::apply
void apply(DynamicVector< field_type > &v, DynamicVector< field_type > &d)
Apply the subdomain solver.
Definition: overlappingschwarz.hh:158
Dune::OverlappingSchwarzInitializer::AtomInitializer
InitializerList::value_type AtomInitializer
Definition: overlappingschwarz.hh:49
Dune::SeqOverlappingSchwarz::size_type
matrix_type::size_type size_type
The return type of the size method.
Definition: overlappingschwarz.hh:782
Dune::OverlappingAssignerHelper< S< BCRSMatrix< T, A > >, true >::field_type
range_type::field_type field_type
Definition: overlappingschwarz.hh:312
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::range_type
Y range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:151
Dune::SymmetricMultiplicativeSchwarzMode
Tag that tells the Schwarz method to be multiplicative and symmetric.
Definition: overlappingschwarz.hh:128
Dune::OverlappingAssignerILUBase::operator()
void operator()(const size_type &domain)
calculate one entry of the local defect.
Definition: overlappingschwarz.hh:1509
Dune::OverlappingAssignerILUBase::matrix_type
M matrix_type
Definition: overlappingschwarz.hh:398
Dune::SeqOverlappingSchwarz::slu
TD slu
The type for the subdomain solver in use.
Definition: overlappingschwarz.hh:802
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::matrix_type
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: overlappingschwarz.hh:145
Dune::MultiplicativeAdder
Definition: overlappingschwarz.hh:537
Dune::SeqOverlappingSchwarz::domain_type
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:761
Dune::SeqOverlappingSchwarzApplier< SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > >::smoother
SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > smoother
Definition: overlappingschwarz.hh:678
Dune::SeqOverlappingSchwarzAssemblerHelper< S< BCRSMatrix< T, A > >, true >::matrix_type
BCRSMatrix< T, A > matrix_type
Definition: overlappingschwarz.hh:709
Dune::SeqOverlappingSchwarz::rowtodomain_vector
std::vector< subdomain_list, typename TA::template rebind< subdomain_list >::other > rowtodomain_vector
The vector type containing the row index to subdomain mapping.
Definition: overlappingschwarz.hh:799
Dune::OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::field_type
X::field_type field_type
Definition: overlappingschwarz.hh:222
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::rilu_type
std::remove_const< M >::type rilu_type
Definition: overlappingschwarz.hh:147
Dune::AdditiveAdder
Definition: overlappingschwarz.hh:515
Dune::IteratorDirectionSelector< T1, T2, false >::solver_iterator
solver_vector::reverse_iterator solver_iterator
Definition: overlappingschwarz.hh:631
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::setSubMatrix
void setSubMatrix(const M &BCRS, S &rowset)
Set the data of the local problem.
Definition: overlappingschwarz.hh:180
Dune::IteratorDirectionSelector< T1, T2, false >::end
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:640
Dune::SeqOverlappingSchwarzApplier::range_type
smoother::range_type range_type
Definition: overlappingschwarz.hh:667
Dune::AdderSelector< AdditiveSchwarzMode, X, S >::Adder
AdditiveAdder< S, X > Adder
Definition: overlappingschwarz.hh:573
Dune::OverlappingAssignerILUBase::relaxResult
void relaxResult(field_type relax)
relax the result.
Definition: overlappingschwarz.hh:1525
Dune::SeqOverlappingSchwarz::field_type
X::field_type field_type
The field type of the preconditioner.
Definition: overlappingschwarz.hh:779
Dune::OverlappingSchwarzInitializer::size_type
IndexSet::size_type size_type
Definition: overlappingschwarz.hh:55
Dune::SeqOverlappingSchwarzApplier::smoother
T smoother
Definition: overlappingschwarz.hh:666
Dune::SeqOverlappingSchwarzAssemblerILUBase::matrix_type
M matrix_type
Definition: overlappingschwarz.hh:720
Dune::OverlappingAssignerILUBase::resetIndexForNextDomain
void resetIndexForNextDomain()
Resets the local index to zero.
Definition: overlappingschwarz.hh:1549
Dune::SeqOverlappingSchwarz::category
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: overlappingschwarz.hh:869
Dune::SeqOverlappingSchwarz::slu_vector
std::vector< slu, typename TA::template rebind< slu >::other > slu_vector
The vector type containing subdomain solvers.
Definition: overlappingschwarz.hh:805
Dune::BCRSMatrix::size_type
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:459
Dune::OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::matrix_type
BCRSMatrix< K, Al > matrix_type
Definition: overlappingschwarz.hh:221
Dune::SeqOverlappingSchwarz::subdomain_vector
std::vector< subdomain_type, typename TA::template rebind< subdomain_type >::other > subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:793
Dune::SeqOverlappingSchwarz::post
virtual void post(X &x)
Postprocess the preconditioner.
Definition: overlappingschwarz.hh:860
Dune::IteratorDirectionSelector< T1, T2, false >::end
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:649
Dune::IteratorDirectionSelector::domain_iterator
subdomain_vector::const_iterator domain_iterator
Definition: overlappingschwarz.hh:605
Dune::SeqOverlappingSchwarz::subdomain_list
SLList< size_type, typename TA::template rebind< size_type >::other > subdomain_list
The type for the row to subdomain mapping.
Definition: overlappingschwarz.hh:796
Dune::OverlappingSchwarzInitializer::copyValue
void copyValue(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:934
Dune::SeqOverlappingSchwarz::range_type
X range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:766
Dune::OverlappingAssignerILUBase::size_type
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:404
Dune::IteratorDirectionSelector< T1, T2, false >::domain_iterator
subdomain_vector::const_reverse_iterator domain_iterator
Definition: overlappingschwarz.hh:633
Dune::copyToColCompMatrix
void copyToColCompMatrix(F &initializer, const MRS &mrs)
Definition: colcompmatrix.hh:473
Dune::OverlappingAssignerILUBase::field_type
Y::field_type field_type
Definition: overlappingschwarz.hh:400
Dune::IteratorDirectionSelector< T1, T2, false >::begin
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:635
bcrsmatrix.hh
Implementation of the BCRSMatrix class.
Dune::BCRSMatrix
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:424
Dune::AdderSelector
template meta program for choosing how to add the correction.
Definition: overlappingschwarz.hh:567
Dune::OverlappingSchwarzInitializer::IndexSet
S IndexSet
Definition: overlappingschwarz.hh:54
Dune::AdditiveAdder< S, BlockVector< T, A > >::size_type
A::size_type size_type
Definition: overlappingschwarz.hh:521
Dune::SeqOverlappingSchwarzAssemblerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::matrix_type
BCRSMatrix< K, Al > matrix_type
Definition: overlappingschwarz.hh:698
Dune::SeqOverlappingSchwarz::apply
virtual void apply(X &v, const X &d)
Apply the precondtioner.
Definition: overlappingschwarz.hh:1233
Dune::OverlappingAssignerILUBase::rhs
Y & rhs()
Get the local right hand side.
Definition: overlappingschwarz.hh:1543
Dune::OverlappingSchwarzInitializer::Iter
Matrix::const_iterator Iter
Definition: overlappingschwarz.hh:51
Dune::OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::size_type
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:225
Dune::IteratorDirectionSelector::end
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:612
Dune::ILUNSubdomainSolver
Definition: ilusubdomainsolver.hh:107
Dune::IteratorDirectionSelector::begin
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:607
Dune::SeqOverlappingSchwarzAssemblerHelper
Definition: colcompmatrix.hh:153
Dune::OverlappingSchwarzInitializer::countEntries
void countEntries(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:915
Dune::IteratorDirectionSelector
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:600
Dune::SeqOverlappingSchwarzApplier< SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > >::apply
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:681
Dune::MultiplicativeAdder< S, BlockVector< T, A > >::size_type
A::size_type size_type
Definition: overlappingschwarz.hh:543
Dune
Definition: allocator.hh:7
Dune::SeqOverlappingSchwarzApplier< SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > >::range_type
smoother::range_type range_type
Definition: overlappingschwarz.hh:679
Dune::OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::block_type
range_type::block_type block_type
Definition: overlappingschwarz.hh:224
Dune::OverlappingSchwarzInitializer::allocate
void allocate()
Definition: overlappingschwarz.hh:906
Dune::OverlappingSchwarzInitializer::subdomain_vector
D subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:46
superlu.hh
Classes for using SuperLU with ISTL matrices.
Dune::OverlappingAssignerILUBase::block_type
Y::block_type block_type
Definition: overlappingschwarz.hh:402
Dune::IteratorDirectionSelector::solver_iterator
solver_vector::iterator solver_iterator
Definition: overlappingschwarz.hh:603
Dune::OverlappingAssignerILUBase
Definition: overlappingschwarz.hh:395
Dune::OverlappingSchwarzInitializer::addRowNnz
void addRowNnz(const Iter &row)
Definition: overlappingschwarz.hh:896
Dune::IteratorDirectionSelector< T1, T2, false >::begin
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:644
Dune::DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >::domain_type
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:149
Dune::OverlappingAssignerHelper< S< BCRSMatrix< T, A > >, true >::size_type
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:315
Dune::OverlappingAssignerILUBase::lhs
X & lhs()
Get the local left hand side.
Definition: overlappingschwarz.hh:1537
Dune::SeqOverlappingSchwarz::subdomain_type
std::set< size_type, std::less< size_type >, typename TA::template rebind< size_type >::other > subdomain_type
The type for the subdomain to row index mapping.
Definition: overlappingschwarz.hh:790
mat
Matrix & mat
Definition: matrixmatrix.hh:345
Dune::OverlappingAssignerILUBase::deallocate
void deallocate()
Deallocates memory of the local vector.
Definition: overlappingschwarz.hh:1502
Dune::OverlappingAssignerILUBase::OverlappingAssignerILUBase
OverlappingAssignerILUBase(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:1489
Dune::OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al >, X, Y >, false >::range_type
Y range_type
Definition: overlappingschwarz.hh:223
Dune::SeqOverlappingSchwarzDomainSize
Definition: overlappingschwarz.hh:1104
Dune::OverlappingAssignerHelper< ILU0SubdomainSolver< M, X, Y >, false >::OverlappingAssignerHelper
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:489
Dune::DynamicMatrixSubdomainSolver
Exact subdomain solver using Dune::DynamicMatrix<T>::solve.
Definition: overlappingschwarz.hh:136
umfpack.hh
Classes for using UMFPack with ISTL matrices.
Dune::Preconditioner
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
Dune::BCRSMatrix::ConstColIterator
row_type::ConstIterator ConstColIterator
Const iterator to the entries of a row.
Definition: bcrsmatrix.hh:699
solvertype.hh
Templates characterizing the type of a solver.
preconditioners.hh
Define general preconditioner interface.
Dune::SeqOverlappingSchwarzApplier::apply
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:669
Dune::BlockVector
A vector of blocks with memory management.
Definition: bvector.hh:402
Dune::MultiplicativeSchwarzMode
Tag that tells the Schwarz method to be multiplicative.
Definition: overlappingschwarz.hh:121
Dune::SolverCategory::Category
Category
Definition: solvercategory.hh:21
Dune::OverlappingAssignerHelper< S< BCRSMatrix< T, A > >, true >::range_type
S< BCRSMatrix< T, A > >::range_type range_type
Definition: overlappingschwarz.hh:311
Dune::SeqOverlappingSchwarz::matrix_type
M matrix_type
The type of the matrix to precondition.
Definition: overlappingschwarz.hh:756
Dune::OverlappingSchwarzInitializer
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:42
Dune::OverlappingAssignerHelper< S< BCRSMatrix< T, A > >, true >::block_type
range_type::block_type block_type
Definition: overlappingschwarz.hh:313
Dune::IteratorDirectionSelector< T1, T2, false >::subdomain_vector
T2 subdomain_vector
Definition: overlappingschwarz.hh:632
Dune::SeqOverlappingSchwarzAssemblerILUBase::assembleLocalProblems
static std::size_t assembleLocalProblems(const RowToDomain &rowToDomain, const matrix_type &mat, Solvers &solvers, const SubDomains &domains, bool onTheFly)
Definition: overlappingschwarz.hh:1203
bvector.hh
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Dune::OverlappingSchwarzInitializer::calcColstart
void calcColstart() const
Definition: overlappingschwarz.hh:927
Dune::SeqOverlappingSchwarzDomainSize< BCRSMatrix< T, A > >::size
static int size(const Domain &d)
Definition: overlappingschwarz.hh:1112
Dune::IteratorDirectionSelector< T1, T2, false >::solver_vector
T1 solver_vector
Definition: overlappingschwarz.hh:630
ilusubdomainsolver.hh
Various local subdomain solvers based on ILU for SeqOverlappingSchwarz.
Dune::OverlappingSchwarzInitializer::CIter
Matrix::row_type::const_iterator CIter
Definition: overlappingschwarz.hh:52
Dune::AdditiveAdder< S, BlockVector< T, A > >::field_type
std::decay_t< decltype(Impl::asVector(std::declval< T >)))>::field_type field_type
Definition: overlappingschwarz.hh:522
Dune::OverlappingSchwarzInitializer::InitializerList
I InitializerList
Definition: overlappingschwarz.hh:48
Dune::OverlappingAssignerHelper< S< BCRSMatrix< T, A > >, true >::matrix_type
BCRSMatrix< T, A > matrix_type
Definition: overlappingschwarz.hh:310
Dune::OverlappingAssignerHelper< ILUNSubdomainSolver< M, X, Y >, false >::OverlappingAssignerHelper
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:508
Dune::SeqOverlappingSchwarz::Mode
TM Mode
The mode (additive or multiplicative) of the Schwarz method.
Definition: overlappingschwarz.hh:774