dune-istl  2.8.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 <memory>
9 #include <vector>
10 #include <set>
11 #include <dune/common/dynmatrix.hh>
12 #include <dune/common/sllist.hh>
13 
15 #include "preconditioners.hh"
16 #include "superlu.hh"
17 #include "umfpack.hh"
18 #include "bvector.hh"
19 #include "bcrsmatrix.hh"
20 #include "ilusubdomainsolver.hh"
21 #include <dune/istl/solvertype.hh>
22 
23 namespace Dune
24 {
25 
37  template<class M, class X, class TM, class TD, class TA>
38  class SeqOverlappingSchwarz;
39 
43  template<class I, class S, class D>
45  {
46  public:
48  typedef D subdomain_vector;
49 
50  typedef I InitializerList;
51  typedef typename InitializerList::value_type AtomInitializer;
52  typedef typename AtomInitializer::Matrix Matrix;
53  typedef typename Matrix::const_iterator Iter;
54  typedef typename Matrix::row_type::const_iterator CIter;
55 
56  typedef S IndexSet;
57  typedef typename IndexSet::size_type size_type;
58 
60  const IndexSet& indices,
61  const subdomain_vector& domains);
62 
63 
64  void addRowNnz(const Iter& row);
65 
66  void allocate();
67 
68  void countEntries(const Iter& row, const CIter& col) const;
69 
70  void calcColstart() const;
71 
72  void copyValue(const Iter& row, const CIter& col) const;
73 
74  void createMatrix() const;
75  private:
76  class IndexMap
77  {
78  public:
79  typedef typename S::size_type size_type;
80  typedef std::map<size_type,size_type> Map;
81  typedef typename Map::iterator iterator;
82  typedef typename Map::const_iterator const_iterator;
83 
84  IndexMap();
85 
86  void insert(size_type grow);
87 
88  const_iterator find(size_type grow) const;
89 
90  iterator find(size_type grow);
91 
92  iterator begin();
93 
94  const_iterator begin() const;
95 
96  iterator end();
97 
98  const_iterator end() const;
99 
100  private:
101  std::map<size_type,size_type> map_;
102  size_type row;
103  };
104 
105 
106  typedef typename InitializerList::iterator InitIterator;
107  typedef typename IndexSet::const_iterator IndexIteratur;
108  InitializerList* initializers;
109  const IndexSet *indices;
110  mutable std::vector<IndexMap> indexMaps;
111  const subdomain_vector& domains;
112  };
113 
118  {};
119 
124  {};
125 
131  {};
132 
137  template<class M, class X, class Y>
139 
140  // Specialization for BCRSMatrix
141  template<class K, class Al, class X, class Y>
143  {
144  typedef BCRSMatrix< K, Al> M;
145  public:
147  typedef typename std::remove_const<M>::type matrix_type;
148  typedef typename X::field_type field_type;
149  typedef typename std::remove_const<M>::type rilu_type;
151  typedef X domain_type;
153  typedef Y range_type;
154  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
155 
160  void apply (DynamicVector<field_type>& v, DynamicVector<field_type>& d)
161  {
162  assert(v.size() > 0);
163  assert(v.size() == d.size());
164  assert(A.rows() <= v.size());
165  assert(A.cols() <= v.size());
166  size_t sz = A.rows();
167  v.resize(sz);
168  d.resize(sz);
169  A.solve(v,d);
170  v.resize(v.capacity());
171  d.resize(d.capacity());
172  }
173 
181  template<class S>
182  void setSubMatrix(const M& BCRS, S& rowset)
183  {
184  size_t sz = rowset.size();
185  A.resize(sz*n,sz*n);
186  typedef typename S::const_iterator SIter;
187  size_t r = 0, c = 0;
188  for(SIter rowIdx = rowset.begin(), rowEnd=rowset.end();
189  rowIdx!= rowEnd; ++rowIdx, r++)
190  {
191  c = 0;
192  for(SIter colIdx = rowset.begin(), colEnd=rowset.end();
193  colIdx!= colEnd; ++colIdx, c++)
194  {
195  if (BCRS[*rowIdx].find(*colIdx) == BCRS[*rowIdx].end())
196  continue;
197  for (size_t i=0; i<n; i++)
198  {
199  for (size_t j=0; j<n; j++)
200  {
201  A[r*n+i][c*n+j] = Impl::asMatrix(BCRS[*rowIdx][*colIdx])[i][j];
202  }
203  }
204  }
205  }
206  }
207  private:
208  DynamicMatrix<K> A;
209  };
210 
211  template<typename T, bool tag>
213  {};
214 
215  template<typename T>
217 
218  // specialization for DynamicMatrix
219  template<class K, class Al, class X, class Y>
221  {
222  public:
224  typedef typename X::field_type field_type;
225  typedef Y range_type;
226  typedef typename range_type::block_type block_type;
228  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
236  OverlappingAssignerHelper(std::size_t maxlength, const BCRSMatrix<K, Al>& mat_, const X& b_, Y& x_);
237 
241  inline
242  void deallocate();
243 
247  inline
248  void resetIndexForNextDomain();
249 
254  inline
255  DynamicVector<field_type> & lhs();
256 
261  inline
262  DynamicVector<field_type> & rhs();
263 
268  inline
269  void relaxResult(field_type relax);
270 
275  void operator()(const size_type& domainIndex);
276 
284  inline
285  void assignResult(block_type& res);
286 
287  private:
291  const matrix_type* mat;
293  // we need a pointer, because we have to avoid deep copies
294  DynamicVector<field_type> * rhs_;
296  // we need a pointer, because we have to avoid deep copies
297  DynamicVector<field_type> * lhs_;
299  const range_type* b;
301  range_type* x;
303  std::size_t i;
305  std::size_t maxlength_;
306  };
307 
308 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
309  template<template<class> class S, typename T, typename A>
310  struct OverlappingAssignerHelper<S<BCRSMatrix<T, A>>, true>
311  {
313  typedef typename S<BCRSMatrix<T, A>>::range_type range_type;
314  typedef typename range_type::field_type field_type;
315  typedef typename range_type::block_type block_type;
316 
318 
319  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
320  static constexpr size_t m = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::cols;
328  OverlappingAssignerHelper(std::size_t maxlength, const matrix_type& mat,
329  const range_type& b, range_type& x);
335  void deallocate();
336 
337  /*
338  * @brief Resets the local index to zero.
339  */
340  void resetIndexForNextDomain();
341 
346  field_type* lhs();
347 
352  field_type* rhs();
353 
358  void relaxResult(field_type relax);
359 
364  void operator()(const size_type& domain);
365 
373  void assignResult(block_type& res);
374 
375  private:
379  const matrix_type* mat;
381  field_type* rhs_;
383  field_type* lhs_;
385  const range_type* b;
387  range_type* x;
389  std::size_t i;
391  std::size_t maxlength_;
392  };
393 
394 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
395 
396  template<class M, class X, class Y>
398  {
399  public:
400  typedef M matrix_type;
401 
402  typedef typename Y::field_type field_type;
403 
404  typedef typename Y::block_type block_type;
405 
406  typedef typename matrix_type::size_type size_type;
414  OverlappingAssignerILUBase(std::size_t maxlength, const M& mat,
415  const Y& b, X& x);
421  void deallocate();
422 
427 
432  X& lhs();
433 
438  Y& rhs();
439 
444  void relaxResult(field_type relax);
445 
450  void operator()(const size_type& domain);
451 
459  void assignResult(block_type& res);
460 
461  private:
465  const M* mat;
467  X* lhs_;
469  Y* rhs_;
471  const Y* b;
473  X* x;
475  size_type i;
476  };
477 
478  // specialization for ILU0
479  template<class M, class X, class Y>
481  : public OverlappingAssignerILUBase<M,X,Y>
482  {
483  public:
491  OverlappingAssignerHelper(std::size_t maxlength, const M& mat,
492  const Y& b, X& x)
493  : OverlappingAssignerILUBase<M,X,Y>(maxlength, mat,b,x)
494  {}
495  };
496 
497  // specialization for ILUN
498  template<class M, class X, class Y>
500  : public OverlappingAssignerILUBase<M,X,Y>
501  {
502  public:
510  OverlappingAssignerHelper(std::size_t maxlength, const M& mat,
511  const Y& b, X& x)
512  : OverlappingAssignerILUBase<M,X,Y>(maxlength, mat,b,x)
513  {}
514  };
515 
516  template<typename S, typename T>
518  {};
519 
520  template<typename S, typename T, typename A>
521  struct AdditiveAdder<S, BlockVector<T,A> >
522  {
523  typedef typename A::size_type size_type;
524  typedef typename std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::field_type field_type;
526  OverlappingAssigner<S>& assigner, const field_type& relax_);
527  void operator()(const size_type& domain);
528  void axpy();
529  static constexpr size_t n = std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::dimension;
530 
531  private:
532  BlockVector<T,A>* v;
533  BlockVector<T,A>* x;
534  OverlappingAssigner<S>* assigner;
535  field_type relax;
536  };
537 
538  template<typename S,typename T>
540  {};
541 
542  template<typename S, typename T, typename A>
544  {
545  typedef typename A::size_type size_type;
546  typedef typename std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::field_type field_type;
548  OverlappingAssigner<S>& assigner_, const field_type& relax_);
549  void operator()(const size_type& domain);
550  void axpy();
551  static constexpr size_t n = std::decay_t<decltype(Impl::asVector(std::declval<T>()))>::dimension;
552 
553  private:
554  BlockVector<T,A>* x;
555  OverlappingAssigner<S>* assigner;
556  field_type relax;
557  };
558 
568  template<typename T, class X, class S>
570  {};
571 
572  template<class X, class S>
574  {
576  };
577 
578  template<class X, class S>
580  {
582  };
583 
584  template<class X, class S>
586  {
588  };
589 
601  template<typename T1, typename T2, bool forward>
603  {
604  typedef T1 solver_vector;
605  typedef typename solver_vector::iterator solver_iterator;
606  typedef T2 subdomain_vector;
607  typedef typename subdomain_vector::const_iterator domain_iterator;
608 
610  {
611  return sv.begin();
612  }
613 
615  {
616  return sv.end();
617  }
619  {
620  return sv.begin();
621  }
622 
624  {
625  return sv.end();
626  }
627  };
628 
629  template<typename T1, typename T2>
630  struct IteratorDirectionSelector<T1,T2,false>
631  {
632  typedef T1 solver_vector;
633  typedef typename solver_vector::reverse_iterator solver_iterator;
634  typedef T2 subdomain_vector;
635  typedef typename subdomain_vector::const_reverse_iterator domain_iterator;
636 
638  {
639  return sv.rbegin();
640  }
641 
643  {
644  return sv.rend();
645  }
647  {
648  return sv.rbegin();
649  }
650 
652  {
653  return sv.rend();
654  }
655  };
656 
665  template<class T>
667  {
668  typedef T smoother;
669  typedef typename smoother::range_type range_type;
670 
671  static void apply(smoother& sm, range_type& v, const range_type& b)
672  {
673  sm.template apply<true>(v, b);
674  }
675  };
676 
677  template<class M, class X, class TD, class TA>
679  {
682 
683  static void apply(smoother& sm, range_type& v, const range_type& b)
684  {
685  sm.template apply<true>(v, b);
686  sm.template apply<false>(v, b);
687  }
688  };
689 
690  template<class T, bool tag>
692  {};
693 
694  template<class T>
696 
697  template<class K, class Al, class X, class Y>
699  {
701  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<K>()))>::rows;
702  template<class RowToDomain, class Solvers, class SubDomains>
703  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
704  Solvers& solvers, const SubDomains& domains,
705  bool onTheFly);
706  };
707 
708  template<template<class> class S, typename T, typename A>
710  {
712  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
713  template<class RowToDomain, class Solvers, class SubDomains>
714  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
715  Solvers& solvers, const SubDomains& domains,
716  bool onTheFly);
717  };
718 
719  template<class M,class X, class Y>
721  {
722  typedef M matrix_type;
723  template<class RowToDomain, class Solvers, class SubDomains>
724  static std::size_t assembleLocalProblems(const RowToDomain& rowToDomain, const matrix_type& mat,
725  Solvers& solvers, const SubDomains& domains,
726  bool onTheFly);
727  };
728 
729  template<class M,class X, class Y>
732  {};
733 
734  template<class M,class X, class Y>
737  {};
738 
749  template<class M, class X, class TM=AdditiveSchwarzMode,
750  class TD=ILU0SubdomainSolver<M,X,X>, class TA=std::allocator<X> >
752  : public Preconditioner<X,X>
753  {
754  public:
758  typedef M matrix_type;
759 
763  typedef X domain_type;
764 
768  typedef X range_type;
769 
776  typedef TM Mode;
777 
781  typedef typename X::field_type field_type;
782 
784  typedef typename matrix_type::size_type size_type;
785 
787  typedef TA allocator;
788 
790  typedef std::set<size_type, std::less<size_type>,
791  typename std::allocator_traits<TA>::template rebind_alloc<size_type> >
793 
795  typedef std::vector<subdomain_type, typename std::allocator_traits<TA>::template rebind_alloc<subdomain_type> > subdomain_vector;
796 
798  typedef SLList<size_type, typename std::allocator_traits<TA>::template rebind_alloc<size_type> > subdomain_list;
799 
801  typedef std::vector<subdomain_list, typename std::allocator_traits<TA>::template rebind_alloc<subdomain_list> > rowtodomain_vector;
802 
804  typedef TD slu;
805 
807  typedef std::vector<slu, typename std::allocator_traits<TA>::template rebind_alloc<slu> > slu_vector;
808 
822  SeqOverlappingSchwarz(const matrix_type& mat, const subdomain_vector& subDomains,
823  field_type relaxationFactor=1, bool onTheFly_=true);
824 
836  SeqOverlappingSchwarz(const matrix_type& mat, const rowtodomain_vector& rowToDomain,
837  field_type relaxationFactor=1, bool onTheFly_=true);
838 
844  virtual void pre ([[maybe_unused]] X& x, [[maybe_unused]] X& b)
845  {}
846 
852  virtual void apply (X& v, const X& d);
853 
859  virtual void post ([[maybe_unused]] X& x)
860  {}
861 
862  template<bool forward>
863  void apply(X& v, const X& d);
864 
867  {
869  }
870 
871  private:
872  const M& mat;
873  slu_vector solvers;
874  subdomain_vector subDomains;
875  field_type relax;
876 
877  typename M::size_type maxlength;
878 
879  bool onTheFly;
880  };
881 
882 
883 
884  template<class I, class S, class D>
886  const IndexSet& idx,
887  const subdomain_vector& domains_)
888  : initializers(&il), indices(&idx), indexMaps(il.size()), domains(domains_)
889  {}
890 
891 
892  template<class I, class S, class D>
894  {
895  typedef typename IndexSet::value_type::const_iterator iterator;
896  for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
897  (*initializers)[*domain].addRowNnz(row, domains[*domain]);
898  indexMaps[*domain].insert(row.index());
899  }
900  }
901 
902  template<class I, class S, class D>
904  {
905  for(auto&& i: *initializers)
906  i.allocateMatrixStorage();
907  for(auto&& i: *initializers)
908  i.allocateMarker();
909  }
910 
911  template<class I, class S, class D>
913  {
914  typedef typename IndexSet::value_type::const_iterator iterator;
915  for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
916  typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.index());
917  if(v!= indexMaps[*domain].end()) {
918  (*initializers)[*domain].countEntries(indexMaps[*domain].find(col.index())->second);
919  }
920  }
921  }
922 
923  template<class I, class S, class D>
925  {
926  for(auto&& i : *initializers)
927  i.calcColstart();
928  }
929 
930  template<class I, class S, class D>
932  {
933  typedef typename IndexSet::value_type::const_iterator iterator;
934  for(iterator domain=(*indices)[row.index()].begin(); domain!= (*indices)[row.index()].end(); ++domain) {
935  typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.index());
936  if(v!= indexMaps[*domain].end()) {
937  assert(indexMaps[*domain].end()!=indexMaps[*domain].find(row.index()));
938  (*initializers)[*domain].copyValue(col, indexMaps[*domain].find(row.index())->second,
939  v->second);
940  }
941  }
942  }
943 
944  template<class I, class S, class D>
946  {
947  std::vector<IndexMap>().swap(indexMaps);
948  for(auto&& i: *initializers)
949  i.createMatrix();
950  }
951 
952  template<class I, class S, class D>
954  : row(0)
955  {}
956 
957  template<class I, class S, class D>
959  {
960  assert(map_.find(grow)==map_.end());
961  map_.insert(std::make_pair(grow, row++));
962  }
963 
964  template<class I, class S, class D>
965  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
967  {
968  return map_.find(grow);
969  }
970 
971  template<class I, class S, class D>
972  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
974  {
975  return map_.find(grow);
976  }
977 
978  template<class I, class S, class D>
979  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
981  {
982  return map_.end();
983  }
984 
985  template<class I, class S, class D>
986  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
988  {
989  return map_.end();
990  }
991 
992  template<class I, class S, class D>
993  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
995  {
996  return map_.begin();
997  }
998 
999  template<class I, class S, class D>
1000  typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
1002  {
1003  return map_.begin();
1004  }
1005 
1006  template<class M, class X, class TM, class TD, class TA>
1008  field_type relaxationFactor, bool fly)
1009  : mat(mat_), relax(relaxationFactor), onTheFly(fly)
1010  {
1011  typedef typename rowtodomain_vector::const_iterator RowDomainIterator;
1012  typedef typename subdomain_list::const_iterator DomainIterator;
1013 #ifdef DUNE_ISTL_WITH_CHECKING
1014  assert(rowToDomain.size()==mat.N());
1015  assert(rowToDomain.size()==mat.M());
1016 
1017  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1018  assert(iter->size()>0);
1019 
1020 #endif
1021  // calculate the number of domains
1022  size_type domains=0;
1023  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1024  for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1025  domains=std::max(domains, *d);
1026  ++domains;
1027 
1028  solvers.resize(domains);
1029  subDomains.resize(domains);
1030 
1031  // initialize subdomains to row mapping from row to subdomain mapping
1032  size_type row=0;
1033  for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter, ++row)
1034  for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1035  subDomains[*d].insert(row);
1036 
1037 #ifdef DUNE_ISTL_WITH_CHECKING
1038  size_type i=0;
1039  typedef typename subdomain_vector::const_iterator iterator;
1040  for(iterator iter=subDomains.begin(); iter != subDomains.end(); ++iter) {
1041  typedef typename subdomain_type::const_iterator entry_iterator;
1042  Dune::dvverb<<"domain "<<i++<<":";
1043  for(entry_iterator entry = iter->begin(); entry != iter->end(); ++entry) {
1044  Dune::dvverb<<" "<<*entry;
1045  }
1046  Dune::dvverb<<std::endl;
1047  }
1048 #endif
1050  ::assembleLocalProblems(rowToDomain, mat, solvers, subDomains, onTheFly);
1051  }
1052 
1053  template<class M, class X, class TM, class TD, class TA>
1055  const subdomain_vector& sd,
1056  field_type relaxationFactor,
1057  bool fly)
1058  : mat(mat_), solvers(sd.size()), subDomains(sd), relax(relaxationFactor),
1059  onTheFly(fly)
1060  {
1061  typedef typename subdomain_vector::const_iterator DomainIterator;
1062 
1063 #ifdef DUNE_ISTL_WITH_CHECKING
1064  size_type i=0;
1065 
1066  for(DomainIterator d=sd.begin(); d != sd.end(); ++d,++i) {
1067  //std::cout<<i<<": "<<d->size()<<std::endl;
1068  assert(d->size()>0);
1069  typedef typename DomainIterator::value_type::const_iterator entry_iterator;
1070  Dune::dvverb<<"domain "<<i<<":";
1071  for(entry_iterator entry = d->begin(); entry != d->end(); ++entry) {
1072  Dune::dvverb<<" "<<*entry;
1073  }
1074  Dune::dvverb<<std::endl;
1075  }
1076 
1077 #endif
1078 
1079  // Create a row to subdomain mapping
1080  rowtodomain_vector rowToDomain(mat.N());
1081 
1082  size_type domainId=0;
1083 
1084  for(DomainIterator domain=sd.begin(); domain != sd.end(); ++domain, ++domainId) {
1085  typedef typename subdomain_type::const_iterator iterator;
1086  for(iterator row=domain->begin(); row != domain->end(); ++row)
1087  rowToDomain[*row].push_back(domainId);
1088  }
1089 
1091  ::assembleLocalProblems(rowToDomain, mat, solvers, subDomains, onTheFly);
1092  }
1093 
1100  template<class M>
1102 
1103  template<typename T, typename A>
1105  {
1106  static constexpr size_t n = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::rows;
1107  static constexpr size_t m = std::decay_t<decltype(Impl::asMatrix(std::declval<T>()))>::cols;
1108  template<class Domain>
1109  static int size(const Domain & d)
1110  {
1111  assert(n==m);
1112  return m*d.size();
1113  }
1114  };
1115 
1116  template<class K, class Al, class X, class Y>
1117  template<class RowToDomain, class Solvers, class SubDomains>
1118  std::size_t
1120  assembleLocalProblems([[maybe_unused]] const RowToDomain& rowToDomain,
1121  [[maybe_unused]] const matrix_type& mat,
1122  [[maybe_unused]] Solvers& solvers,
1123  const SubDomains& subDomains,
1124  [[maybe_unused]] bool onTheFly)
1125  {
1126  typedef typename SubDomains::const_iterator DomainIterator;
1127  std::size_t maxlength = 0;
1128 
1129  assert(onTheFly);
1130 
1131  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1132  maxlength=std::max(maxlength, domain->size());
1133  maxlength*=n;
1134 
1135  return maxlength;
1136  }
1137 
1138 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1139  template<template<class> class S, typename T, typename A>
1140  template<class RowToDomain, class Solvers, class SubDomains>
1141  std::size_t SeqOverlappingSchwarzAssemblerHelper<S<BCRSMatrix<T,A>>,true>::assembleLocalProblems(const RowToDomain& rowToDomain,
1142  const matrix_type& mat,
1143  Solvers& solvers,
1144  const SubDomains& subDomains,
1145  bool onTheFly)
1146  {
1147  typedef typename S<BCRSMatrix<T,A>>::MatrixInitializer MatrixInitializer;
1148  typedef typename std::vector<MatrixInitializer>::iterator InitializerIterator;
1149  typedef typename SubDomains::const_iterator DomainIterator;
1150  typedef typename Solvers::iterator SolverIterator;
1151  std::size_t maxlength = 0;
1152 
1153  if(onTheFly) {
1154  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1155  maxlength=std::max(maxlength, domain->size());
1156  maxlength*=Impl::asMatrix(*mat[0].begin()).N();
1157  }else{
1158  // initialize the initializers
1159  DomainIterator domain=subDomains.begin();
1160 
1161  // Create the initializers list.
1162  std::vector<MatrixInitializer> initializers(subDomains.size());
1163 
1164  SolverIterator solver=solvers.begin();
1165  for(InitializerIterator initializer=initializers.begin(); initializer!=initializers.end();
1166  ++initializer, ++solver, ++domain) {
1167  solver->getInternalMatrix().N_=SeqOverlappingSchwarzDomainSize<matrix_type>::size(*domain);
1168  solver->getInternalMatrix().M_=SeqOverlappingSchwarzDomainSize<matrix_type>::size(*domain);
1169  //solver->setVerbosity(true);
1170  *initializer=MatrixInitializer(solver->getInternalMatrix());
1171  }
1172 
1173  // Set up the supermatrices according to the subdomains
1175  RowToDomain, SubDomains> Initializer;
1176 
1177  Initializer initializer(initializers, rowToDomain, subDomains);
1178  copyToBCCSMatrix(initializer, mat);
1179 
1180  // Calculate the LU decompositions
1181  for(auto&& s: solvers)
1182  s.decompose();
1183  for (SolverIterator solverIt = solvers.begin(); solverIt != solvers.end(); ++solverIt)
1184  {
1185  assert(solverIt->getInternalMatrix().N() == solverIt->getInternalMatrix().M());
1186  maxlength = std::max(maxlength, solverIt->getInternalMatrix().N());
1187  }
1188  }
1189  return maxlength;
1190  }
1191 
1192 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1193 
1194  template<class M,class X,class Y>
1195  template<class RowToDomain, class Solvers, class SubDomains>
1196  std::size_t SeqOverlappingSchwarzAssemblerILUBase<M,X,Y>::assembleLocalProblems([[maybe_unused]] const RowToDomain& rowToDomain,
1197  const matrix_type& mat,
1198  Solvers& solvers,
1199  const SubDomains& subDomains,
1200  bool onTheFly)
1201  {
1202  typedef typename SubDomains::const_iterator DomainIterator;
1203  typedef typename Solvers::iterator SolverIterator;
1204  std::size_t maxlength = 0;
1205 
1206  if(onTheFly) {
1207  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1208  maxlength=std::max(maxlength, domain->size());
1209  }else{
1210  // initialize the solvers of the local prolems.
1211  SolverIterator solver=solvers.begin();
1212  for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end();
1213  ++domain, ++solver) {
1214  solver->setSubMatrix(mat, *domain);
1215  maxlength=std::max(maxlength, domain->size());
1216  }
1217  }
1218 
1219  return maxlength;
1220 
1221  }
1222 
1223 
1224  template<class M, class X, class TM, class TD, class TA>
1226  {
1228  }
1229 
1230  template<class M, class X, class TM, class TD, class TA>
1231  template<bool forward>
1232  void SeqOverlappingSchwarz<M,X,TM,TD,TA>::apply(X& x, const X& b)
1233  {
1234  typedef slu_vector solver_vector;
1237  domain_iterator;
1238 
1239  OverlappingAssigner<TD> assigner(maxlength, mat, b, x);
1240 
1243  X v(x); // temporary for the update
1244  v=0;
1245 
1246  typedef typename AdderSelector<TM,X,TD >::Adder Adder;
1247  Adder adder(v, x, assigner, relax);
1248 
1249  for(; domain != IteratorDirectionSelector<solver_vector,subdomain_vector,forward>::end(subDomains); ++domain) {
1250  //Copy rhs to C-array for SuperLU
1251  std::for_each(domain->begin(), domain->end(), assigner);
1252  assigner.resetIndexForNextDomain();
1253  if(onTheFly) {
1254  // Create the subdomain solver
1255  slu sdsolver;
1256  sdsolver.setSubMatrix(mat, *domain);
1257  // Apply
1258  sdsolver.apply(assigner.lhs(), assigner.rhs());
1259  }else{
1260  solver->apply(assigner.lhs(), assigner.rhs());
1261  ++solver;
1262  }
1263 
1264  //Add relaxed correction to from SuperLU to v
1265  std::for_each(domain->begin(), domain->end(), adder);
1266  assigner.resetIndexForNextDomain();
1267 
1268  }
1269 
1270  adder.axpy();
1271  assigner.deallocate();
1272  }
1273 
1274  template<class K, class Al, class X, class Y>
1275  OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1276  ::OverlappingAssignerHelper(std::size_t maxlength, const BCRSMatrix<K, Al>& mat_,
1277  const X& b_, Y& x_) :
1278  mat(&mat_),
1279  rhs_( new DynamicVector<field_type>(maxlength, 42) ),
1280  lhs_( new DynamicVector<field_type>(maxlength, -42) ),
1281  b(&b_),
1282  x(&x_),
1283  i(0),
1284  maxlength_(maxlength)
1285  {}
1286 
1287  template<class K, class Al, class X, class Y>
1288  void
1290  ::deallocate()
1291  {
1292  delete rhs_;
1293  delete lhs_;
1294  }
1295 
1296  template<class K, class Al, class X, class Y>
1297  void
1299  ::resetIndexForNextDomain()
1300  {
1301  i=0;
1302  }
1303 
1304  template<class K, class Al, class X, class Y>
1305  DynamicVector<typename X::field_type> &
1307  ::lhs()
1308  {
1309  return *lhs_;
1310  }
1311 
1312  template<class K, class Al, class X, class Y>
1313  DynamicVector<typename X::field_type> &
1315  ::rhs()
1316  {
1317  return *rhs_;
1318  }
1319 
1320  template<class K, class Al, class X, class Y>
1321  void
1323  ::relaxResult(field_type relax)
1324  {
1325  lhs() *= relax;
1326  }
1327 
1328  template<class K, class Al, class X, class Y>
1329  void
1331  ::operator()(const size_type& domainIndex)
1332  {
1333  lhs() = 0.0;
1334 #if 0
1335  //assign right hand side of current domainindex block
1336  for(size_type j=0; j<n; ++j, ++i) {
1337  assert(i<maxlength_);
1338  rhs()[i]=(*b)[domainIndex][j];
1339  }
1340 
1341  // loop over all Matrix row entries and calculate defect.
1342  typedef typename matrix_type::ConstColIterator col_iterator;
1343 
1344  // calculate defect for current row index block
1345  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1346  block_type tmp(0.0);
1347  (*col).mv((*x)[col.index()], tmp);
1348  i-=n;
1349  for(size_type j=0; j<n; ++j, ++i) {
1350  assert(i<maxlength_);
1351  rhs()[i]-=tmp[j];
1352  }
1353  }
1354 #else
1355  //assign right hand side of current domainindex block
1356  for(size_type j=0; j<n; ++j, ++i) {
1357  assert(i<maxlength_);
1358  rhs()[i]=Impl::asVector((*b)[domainIndex])[j];
1359 
1360  // loop over all Matrix row entries and calculate defect.
1361  typedef typename matrix_type::ConstColIterator col_iterator;
1362 
1363  // calculate defect for current row index block
1364  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1365  for(size_type k=0; k<n; ++k) {
1366  rhs()[i]-=Impl::asMatrix(*col)[j][k] * Impl::asVector((*x)[col.index()])[k];
1367  }
1368  }
1369  }
1370 #endif
1371  }
1372 
1373  template<class K, class Al, class X, class Y>
1374  void
1376  ::assignResult(block_type& res)
1377  {
1378  // assign the result of the local solve to the global vector
1379  for(size_type j=0; j<n; ++j, ++i) {
1380  assert(i<maxlength_);
1381  Impl::asVector(res)[j]+=lhs()[i];
1382  }
1383  }
1384 
1385 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1386 
1387  template<template<class> class S, typename T, typename A>
1389  ::OverlappingAssignerHelper(std::size_t maxlength,
1390  const BCRSMatrix<T,A>& mat_,
1391  const range_type& b_,
1392  range_type& x_)
1393  : mat(&mat_),
1394  b(&b_),
1395  x(&x_), i(0), maxlength_(maxlength)
1396  {
1397  rhs_ = new field_type[maxlength];
1398  lhs_ = new field_type[maxlength];
1399 
1400  }
1401 
1402  template<template<class> class S, typename T, typename A>
1404  {
1405  delete[] rhs_;
1406  delete[] lhs_;
1407  }
1408 
1409  template<template<class> class S, typename T, typename A>
1410  void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::operator()(const size_type& domainIndex)
1411  {
1412  //assign right hand side of current domainindex block
1413  // rhs is an array of doubles!
1414  // rhs[starti] = b[domainindex]
1415  for(size_type j=0; j<n; ++j, ++i) {
1416  assert(i<maxlength_);
1417  rhs_[i]=Impl::asVector((*b)[domainIndex])[j];
1418  }
1419 
1420 
1421  // loop over all Matrix row entries and calculate defect.
1422  typedef typename matrix_type::ConstColIterator col_iterator;
1423 
1424  // calculate defect for current row index block
1425  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1426  block_type tmp;
1427  Impl::asMatrix(*col).mv((*x)[col.index()], tmp);
1428  i-=n;
1429  for(size_type j=0; j<n; ++j, ++i) {
1430  assert(i<maxlength_);
1431  rhs_[i]-=Impl::asVector(tmp)[j];
1432  }
1433 
1434  }
1435 
1436  }
1437 
1438  template<template<class> class S, typename T, typename A>
1440  {
1441  for(size_type j=i+n; i<j; ++i) {
1442  assert(i<maxlength_);
1443  lhs_[i]*=relax;
1444  }
1445  i-=n;
1446  }
1447 
1448  template<template<class> class S, typename T, typename A>
1450  {
1451  // assign the result of the local solve to the global vector
1452  for(size_type j=0; j<n; ++j, ++i) {
1453  assert(i<maxlength_);
1454  Impl::asVector(res)[j]+=lhs_[i];
1455  }
1456  }
1457 
1458  template<template<class> class S, typename T, typename A>
1459  void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::resetIndexForNextDomain()
1460  {
1461  i=0;
1462  }
1463 
1464  template<template<class> class S, typename T, typename A>
1465  typename OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::field_type*
1467  {
1468  return lhs_;
1469  }
1470 
1471  template<template<class> class S, typename T, typename A>
1472  typename OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::field_type*
1474  {
1475  return rhs_;
1476  }
1477 
1478 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1479 
1480  template<class M, class X, class Y>
1482  const M& mat_,
1483  const Y& b_,
1484  X& x_)
1485  : mat(&mat_),
1486  b(&b_),
1487  x(&x_), i(0)
1488  {
1489  rhs_= new Y(maxlength);
1490  lhs_ = new X(maxlength);
1491  }
1492 
1493  template<class M, class X, class Y>
1495  {
1496  delete rhs_;
1497  delete lhs_;
1498  }
1499 
1500  template<class M, class X, class Y>
1502  {
1503  (*rhs_)[i]=(*b)[domainIndex];
1504 
1505  // loop over all Matrix row entries and calculate defect.
1506  typedef typename matrix_type::ConstColIterator col_iterator;
1507 
1508  // calculate defect for current row index block
1509  for(col_iterator col=(*mat)[domainIndex].begin(); col!=(*mat)[domainIndex].end(); ++col) {
1510  Impl::asMatrix(*col).mmv((*x)[col.index()], (*rhs_)[i]);
1511  }
1512  // Goto next local index
1513  ++i;
1514  }
1515 
1516  template<class M, class X, class Y>
1518  {
1519  (*lhs_)[i]*=relax;
1520  }
1521 
1522  template<class M, class X, class Y>
1524  {
1525  res+=(*lhs_)[i++];
1526  }
1527 
1528  template<class M, class X, class Y>
1530  {
1531  return *lhs_;
1532  }
1533 
1534  template<class M, class X, class Y>
1536  {
1537  return *rhs_;
1538  }
1539 
1540  template<class M, class X, class Y>
1542  {
1543  i=0;
1544  }
1545 
1546  template<typename S, typename T, typename A>
1548  BlockVector<T,A>& x_,
1549  OverlappingAssigner<S>& assigner_,
1550  const field_type& relax_)
1551  : v(&v_), x(&x_), assigner(&assigner_), relax(relax_)
1552  {}
1553 
1554  template<typename S, typename T, typename A>
1555  void AdditiveAdder<S,BlockVector<T,A> >::operator()(const size_type& domainIndex)
1556  {
1557  // add the result of the local solve to the current update
1558  assigner->assignResult((*v)[domainIndex]);
1559  }
1560 
1561 
1562  template<typename S, typename T, typename A>
1564  {
1565  // relax the update and add it to the current guess.
1566  x->axpy(relax,*v);
1567  }
1568 
1569 
1570  template<typename S, typename T, typename A>
1573  BlockVector<T,A>& x_,
1574  OverlappingAssigner<S>& assigner_, const field_type& relax_)
1575  : x(&x_), assigner(&assigner_), relax(relax_)
1576  {}
1577 
1578 
1579  template<typename S,typename T, typename A>
1580  void MultiplicativeAdder<S,BlockVector<T,A> >::operator()(const size_type& domainIndex)
1581  {
1582  // add the result of the local solve to the current guess
1583  assigner->relaxResult(relax);
1584  assigner->assignResult((*x)[domainIndex]);
1585  }
1586 
1587 
1588  template<typename S,typename T, typename A>
1590  {
1591  // nothing to do, as the corrections already relaxed and added in operator()
1592  }
1593 
1594 
1596 }
1597 
1598 #endif
Implementation of the BCRSMatrix class.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Various local subdomain solvers based on ILU for SeqOverlappingSchwarz.
Define general preconditioner interface.
Templates characterizing the type of a solver.
Classes for using SuperLU with ISTL matrices.
Classes for using UMFPack with ISTL matrices.
Col col
Definition: matrixmatrix.hh:349
Matrix & mat
Definition: matrixmatrix.hh:345
void addRowNnz(const Iter &row)
Definition: overlappingschwarz.hh:893
X & lhs()
Get the local left hand side.
Definition: overlappingschwarz.hh:1529
void calcColstart() const
Definition: overlappingschwarz.hh:924
Y & rhs()
Get the local right hand side.
Definition: overlappingschwarz.hh:1535
void resetIndexForNextDomain()
Resets the local index to zero.
Definition: overlappingschwarz.hh:1541
void copyValue(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:931
void createMatrix() const
Definition: overlappingschwarz.hh:945
OverlappingSchwarzInitializer(InitializerList &il, const IndexSet &indices, const subdomain_vector &domains)
Definition: overlappingschwarz.hh:885
virtual void apply(X &v, const X &d)
Apply the precondtioner.
Definition: overlappingschwarz.hh:1225
OverlappingAssignerILUBase(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:1481
void operator()(const size_type &domain)
calculate one entry of the local defect.
Definition: overlappingschwarz.hh:1501
SeqOverlappingSchwarz(const matrix_type &mat, const subdomain_vector &subDomains, field_type relaxationFactor=1, bool onTheFly_=true)
Construct the overlapping Schwarz method.
Definition: overlappingschwarz.hh:1054
void allocate()
Definition: overlappingschwarz.hh:903
void deallocate()
Deallocates memory of the local vector.
Definition: overlappingschwarz.hh:1494
void countEntries(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:912
void relaxResult(field_type relax)
relax the result.
Definition: overlappingschwarz.hh:1517
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:1523
SeqOverlappingSchwarz(const matrix_type &mat, const rowtodomain_vector &rowToDomain, field_type relaxationFactor=1, bool onTheFly_=true)
Definition: overlappingschwarz.hh:1007
Definition: allocator.hh:9
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:45
Matrix::row_type::const_iterator CIter
Definition: overlappingschwarz.hh:54
S IndexSet
Definition: overlappingschwarz.hh:56
Matrix::const_iterator Iter
Definition: overlappingschwarz.hh:53
IndexSet::size_type size_type
Definition: overlappingschwarz.hh:57
I InitializerList
Definition: overlappingschwarz.hh:50
AtomInitializer::Matrix Matrix
Definition: overlappingschwarz.hh:52
InitializerList::value_type AtomInitializer
Definition: overlappingschwarz.hh:51
D subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:48
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:464
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:498
row_type::ConstIterator ConstColIterator
Const iterator to the entries of a row.
Definition: bcrsmatrix.hh:739
A vector of blocks with memory management.
Definition: bvector.hh:393
Exact subdomain solver using ILU(p) with appropriate p.
Definition: ilusubdomainsolver.hh:76
Definition: ilusubdomainsolver.hh:109
Sequential overlapping Schwarz preconditioner.
Definition: overlappingschwarz.hh:753
X::field_type field_type
The field type of the preconditioner.
Definition: overlappingschwarz.hh:781
void apply(X &v, const X &d)
Apply one step of the preconditioner to the system A(v)=d.
SLList< size_type, typename std::allocator_traits< TA >::template rebind_alloc< size_type > > subdomain_list
The type for the row to subdomain mapping.
Definition: overlappingschwarz.hh:798
std::vector< slu, typename std::allocator_traits< TA >::template rebind_alloc< slu > > slu_vector
The vector type containing subdomain solvers.
Definition: overlappingschwarz.hh:807
M matrix_type
The type of the matrix to precondition.
Definition: overlappingschwarz.hh:758
TM Mode
The mode (additive or multiplicative) of the Schwarz method.
Definition: overlappingschwarz.hh:776
X range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:768
std::set< size_type, std::less< size_type >, typename std::allocator_traits< TA >::template rebind_alloc< size_type > > subdomain_type
The type for the subdomain to row index mapping.
Definition: overlappingschwarz.hh:792
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:763
TD slu
The type for the subdomain solver in use.
Definition: overlappingschwarz.hh:804
virtual void post([[maybe_unused]] X &x)
Postprocess the preconditioner.
Definition: overlappingschwarz.hh:859
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: overlappingschwarz.hh:866
virtual void pre([[maybe_unused]] X &x, [[maybe_unused]] X &b)
Prepare the preconditioner.
Definition: overlappingschwarz.hh:844
TA allocator
The allocator to use.
Definition: overlappingschwarz.hh:787
std::vector< subdomain_type, typename std::allocator_traits< TA >::template rebind_alloc< subdomain_type > > subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:795
std::vector< subdomain_list, typename std::allocator_traits< TA >::template rebind_alloc< subdomain_list > > rowtodomain_vector
The vector type containing the row index to subdomain mapping.
Definition: overlappingschwarz.hh:801
matrix_type::size_type size_type
The return type of the size method.
Definition: overlappingschwarz.hh:784
Definition: overlappingschwarz.hh:692
Tag that the tells the Schwarz method to be additive.
Definition: overlappingschwarz.hh:118
Tag that tells the Schwarz method to be multiplicative.
Definition: overlappingschwarz.hh:124
Tag that tells the Schwarz method to be multiplicative and symmetric.
Definition: overlappingschwarz.hh:131
Exact subdomain solver using Dune::DynamicMatrix<T>::solve.
Definition: overlappingschwarz.hh:138
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: overlappingschwarz.hh:147
X::field_type field_type
Definition: overlappingschwarz.hh:148
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:151
Y range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:153
void setSubMatrix(const M &BCRS, S &rowset)
Set the data of the local problem.
Definition: overlappingschwarz.hh:182
void apply(DynamicVector< field_type > &v, DynamicVector< field_type > &d)
Apply the subdomain solver.
Definition: overlappingschwarz.hh:160
std::remove_const< M >::type rilu_type
Definition: overlappingschwarz.hh:149
Definition: overlappingschwarz.hh:213
S< BCRSMatrix< T, A > >::range_type range_type
Definition: overlappingschwarz.hh:313
range_type::block_type block_type
Definition: overlappingschwarz.hh:315
range_type::field_type field_type
Definition: overlappingschwarz.hh:314
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:317
BCRSMatrix< T, A > matrix_type
Definition: overlappingschwarz.hh:312
Definition: overlappingschwarz.hh:398
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:406
Y::field_type field_type
Definition: overlappingschwarz.hh:402
Y::block_type block_type
Definition: overlappingschwarz.hh:404
M matrix_type
Definition: overlappingschwarz.hh:400
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:491
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:510
Definition: overlappingschwarz.hh:518
std::decay_t< decltype(Impl::asVector(std::declval< T >)))>::field_type field_type
Definition: overlappingschwarz.hh:524
A::size_type size_type
Definition: overlappingschwarz.hh:523
Definition: overlappingschwarz.hh:540
std::decay_t< decltype(Impl::asVector(std::declval< T >)))>::field_type field_type
Definition: overlappingschwarz.hh:546
A::size_type size_type
Definition: overlappingschwarz.hh:545
template meta program for choosing how to add the correction.
Definition: overlappingschwarz.hh:570
AdditiveAdder< S, X > Adder
Definition: overlappingschwarz.hh:575
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:581
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:587
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:603
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:609
solver_vector::iterator solver_iterator
Definition: overlappingschwarz.hh:605
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:623
subdomain_vector::const_iterator domain_iterator
Definition: overlappingschwarz.hh:607
T1 solver_vector
Definition: overlappingschwarz.hh:604
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:618
T2 subdomain_vector
Definition: overlappingschwarz.hh:606
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:614
T2 subdomain_vector
Definition: overlappingschwarz.hh:634
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:642
solver_vector::reverse_iterator solver_iterator
Definition: overlappingschwarz.hh:633
subdomain_vector::const_reverse_iterator domain_iterator
Definition: overlappingschwarz.hh:635
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:637
T1 solver_vector
Definition: overlappingschwarz.hh:632
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:646
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:651
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:667
smoother::range_type range_type
Definition: overlappingschwarz.hh:669
T smoother
Definition: overlappingschwarz.hh:668
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:671
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:683
SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > smoother
Definition: overlappingschwarz.hh:680
static std::size_t assembleLocalProblems(const RowToDomain &rowToDomain, const matrix_type &mat, Solvers &solvers, const SubDomains &domains, bool onTheFly)
BCRSMatrix< T, A > matrix_type
Definition: overlappingschwarz.hh:711
Definition: overlappingschwarz.hh:721
M matrix_type
Definition: overlappingschwarz.hh:722
static std::size_t assembleLocalProblems(const RowToDomain &rowToDomain, const matrix_type &mat, Solvers &solvers, const SubDomains &domains, bool onTheFly)
Definition: overlappingschwarz.hh:1101
static int size(const Domain &d)
Definition: overlappingschwarz.hh:1109
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
Category
Definition: solvercategory.hh:21
@ sequential
Category for sequential solvers.
Definition: solvercategory.hh:23