dune-istl 2.8.0
Loading...
Searching...
No Matches
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"
22
23namespace 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:
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>
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:
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:
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;
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;
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
823 field_type relaxationFactor=1, bool onTheFly_=true);
824
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>
1572 ::MultiplicativeAdder([[maybe_unused]] BlockVector<T,A>& v_,
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
Templates characterizing the type of a solver.
Define general preconditioner interface.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Implementation of the BCRSMatrix class.
Various local subdomain solvers based on ILU for SeqOverlappingSchwarz.
Classes for using UMFPack with ISTL matrices.
Classes for using SuperLU 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
iterator end()
Definition: overlappingschwarz.hh:987
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
iterator begin()
Definition: overlappingschwarz.hh:1001
OverlappingSchwarzInitializer(InitializerList &il, const IndexSet &indices, const subdomain_vector &domains)
Definition: overlappingschwarz.hh:885
IndexMap()
Definition: overlappingschwarz.hh:953
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
const_iterator find(size_type grow) const
Definition: overlappingschwarz.hh:966
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 insert(size_type grow)
Definition: overlappingschwarz.hh:958
void countEntries(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:912
static std::size_t assembleLocalProblems(const RowToDomain &rowToDomain, const matrix_type &mat, Solvers &solvers, const SubDomains &domains, bool onTheFly)
Definition: overlappingschwarz.hh:1196
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 SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: overlappingschwarz.hh:866
virtual void pre(X &x, X &b)
Prepare the preconditioner.
Definition: overlappingschwarz.hh:844
virtual void post(X &x)
Postprocess the preconditioner.
Definition: overlappingschwarz.hh:859
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
A::size_type size_type
Definition: overlappingschwarz.hh:545
std::decay_t< decltype(Impl::asVector(std::declval< T >()))>::field_type field_type
Definition: overlappingschwarz.hh:546
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
BCRSMatrix< T, A > matrix_type
Definition: overlappingschwarz.hh:711
Definition: overlappingschwarz.hh:721
M matrix_type
Definition: overlappingschwarz.hh:722
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