DUNE PDELab (2.8)

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>
12#include <dune/common/sllist.hh>
13
14#include <dune/istl/bccsmatrixinitializer.hh>
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
59 OverlappingSchwarzInitializer(InitializerList& il,
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>
142 class DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, 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
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>
212 class OverlappingAssignerHelper
213 {};
214
215 template<typename T>
216 using OverlappingAssigner = OverlappingAssignerHelper<T, Dune::StoresColumnCompressed<T>::value>;
217
218 // specialization for DynamicMatrix
219 template<class K, class Al, class X, class Y>
220 class OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix<K, Al>, X, Y >,false>
221 {
222 public:
223 typedef BCRSMatrix< K, Al> matrix_type;
224 typedef typename X::field_type field_type;
225 typedef Y range_type;
226 typedef typename range_type::block_type block_type;
227 typedef typename matrix_type::size_type size_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 {
312 typedef BCRSMatrix<T, A> matrix_type;
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
317 typedef typename matrix_type::size_type size_type;
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>
397 class OverlappingAssignerILUBase
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
426 void resetIndexForNextDomain();
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>
480 class OverlappingAssignerHelper<ILU0SubdomainSolver<M,X,Y>, false>
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>
499 class OverlappingAssignerHelper<ILUNSubdomainSolver<M,X,Y>,false>
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>
517 struct AdditiveAdder
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;
525 AdditiveAdder(BlockVector<T,A>& v, BlockVector<T,A>& x,
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>
539 struct MultiplicativeAdder
540 {};
541
542 template<typename S, typename T, typename A>
543 struct MultiplicativeAdder<S, BlockVector<T,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;
547 MultiplicativeAdder(BlockVector<T,A>& v, BlockVector<T,A>& x,
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 {
575 typedef AdditiveAdder<S,X> Adder;
576 };
577
578 template<class X, class S>
579 struct AdderSelector<MultiplicativeSchwarzMode,X,S>
580 {
581 typedef MultiplicativeAdder<S,X> Adder;
582 };
583
584 template<class X, class S>
585 struct AdderSelector<SymmetricMultiplicativeSchwarzMode,X,S>
586 {
587 typedef MultiplicativeAdder<S,X> Adder;
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
609 static solver_iterator begin(solver_vector& sv)
610 {
611 return sv.begin();
612 }
613
614 static solver_iterator end(solver_vector& sv)
615 {
616 return sv.end();
617 }
618 static domain_iterator begin(const subdomain_vector& sv)
619 {
620 return sv.begin();
621 }
622
623 static domain_iterator end(const subdomain_vector& sv)
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
637 static solver_iterator begin(solver_vector& sv)
638 {
639 return sv.rbegin();
640 }
641
642 static solver_iterator end(solver_vector& sv)
643 {
644 return sv.rend();
645 }
646 static domain_iterator begin(const subdomain_vector& sv)
647 {
648 return sv.rbegin();
649 }
650
651 static domain_iterator end(const subdomain_vector& sv)
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 {
681 typedef typename smoother::range_type range_type;
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>
691 struct SeqOverlappingSchwarzAssemblerHelper
692 {};
693
694 template<class T>
695 using SeqOverlappingSchwarzAssembler = SeqOverlappingSchwarzAssemblerHelper<T,Dune::StoresColumnCompressed<T>::value>;
696
697 template<class K, class Al, class X, class Y>
698 struct SeqOverlappingSchwarzAssemblerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
699 {
700 typedef BCRSMatrix< K, Al> matrix_type;
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>
709 struct SeqOverlappingSchwarzAssemblerHelper<S<BCRSMatrix<T,A>>,true>
710 {
711 typedef BCRSMatrix<T,A> matrix_type;
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>
720 struct SeqOverlappingSchwarzAssemblerILUBase
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>
730 struct SeqOverlappingSchwarzAssemblerHelper<ILU0SubdomainSolver<M,X,Y>,false>
731 : public SeqOverlappingSchwarzAssemblerILUBase<M,X,Y>
732 {};
733
734 template<class M,class X, class Y>
735 struct SeqOverlappingSchwarzAssemblerHelper<ILUNSubdomainSolver<M,X,Y>,false>
736 : public SeqOverlappingSchwarzAssemblerILUBase<M,X,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
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>
885 OverlappingSchwarzInitializer<I,S,D>::OverlappingSchwarzInitializer(InitializerList& il,
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>
893 void OverlappingSchwarzInitializer<I,S,D>::addRowNnz(const Iter& row)
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>
903 void OverlappingSchwarzInitializer<I,S,D>::allocate()
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>
912 void OverlappingSchwarzInitializer<I,S,D>::countEntries(const Iter& row, const CIter& col) const
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>
924 void OverlappingSchwarzInitializer<I,S,D>::calcColstart() const
925 {
926 for(auto&& i : *initializers)
927 i.calcColstart();
928 }
929
930 template<class I, class S, class D>
931 void OverlappingSchwarzInitializer<I,S,D>::copyValue(const Iter& row, const CIter& col) const
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>
945 void OverlappingSchwarzInitializer<I,S,D>::createMatrix() const
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>
953 OverlappingSchwarzInitializer<I,S,D>::IndexMap::IndexMap()
954 : row(0)
955 {}
956
957 template<class I, class S, class D>
958 void OverlappingSchwarzInitializer<I,S,D>::IndexMap::insert(size_type grow)
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
966 OverlappingSchwarzInitializer<I,S,D>::IndexMap::find(size_type grow) const
967 {
968 return map_.find(grow);
969 }
970
971 template<class I, class S, class D>
972 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
973 OverlappingSchwarzInitializer<I,S,D>::IndexMap::find(size_type grow)
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
980 OverlappingSchwarzInitializer<I,S,D>::IndexMap::end() const
981 {
982 return map_.end();
983 }
984
985 template<class I, class S, class D>
986 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
987 OverlappingSchwarzInitializer<I,S,D>::IndexMap::end()
988 {
989 return map_.end();
990 }
991
992 template<class I, class S, class D>
993 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
994 OverlappingSchwarzInitializer<I,S,D>::IndexMap::begin() const
995 {
996 return map_.begin();
997 }
998
999 template<class I, class S, class D>
1000 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
1001 OverlappingSchwarzInitializer<I,S,D>::IndexMap::begin()
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
1049 maxlength = SeqOverlappingSchwarzAssembler<slu>
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
1090 maxlength = SeqOverlappingSchwarzAssembler<slu>
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
1119 SeqOverlappingSchwarzAssemblerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>::
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
1174 typedef OverlappingSchwarzInitializer<std::vector<MatrixInitializer>,
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;
1235 typedef typename IteratorDirectionSelector<solver_vector,subdomain_vector,forward>::solver_iterator iterator;
1236 typedef typename IteratorDirectionSelector<solver_vector,subdomain_vector,forward>::domain_iterator
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
1289 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1290 ::deallocate()
1291 {
1292 delete rhs_;
1293 delete lhs_;
1294 }
1295
1296 template<class K, class Al, class X, class Y>
1297 void
1298 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1299 ::resetIndexForNextDomain()
1300 {
1301 i=0;
1302 }
1303
1304 template<class K, class Al, class X, class Y>
1306 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1307 ::lhs()
1308 {
1309 return *lhs_;
1310 }
1311
1312 template<class K, class Al, class X, class Y>
1314 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1315 ::rhs()
1316 {
1317 return *rhs_;
1318 }
1319
1320 template<class K, class Al, class X, class Y>
1321 void
1322 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
1323 ::relaxResult(field_type relax)
1324 {
1325 lhs() *= relax;
1326 }
1327
1328 template<class K, class Al, class X, class Y>
1329 void
1330 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
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
1375 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< K, Al>, X, Y >,false>
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>
1388 OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>
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>
1403 void OverlappingAssignerHelper<S<BCRSMatrix<T,A> >,true>::deallocate()
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>
1439 void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::relaxResult(field_type relax)
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>
1449 void OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::assignResult(block_type& res)
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*
1466 OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::lhs()
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*
1473 OverlappingAssignerHelper<S<BCRSMatrix<T,A>>,true>::rhs()
1474 {
1475 return rhs_;
1476 }
1477
1478#endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK
1479
1480 template<class M, class X, class Y>
1481 OverlappingAssignerILUBase<M,X,Y>::OverlappingAssignerILUBase(std::size_t maxlength,
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>
1494 void OverlappingAssignerILUBase<M,X,Y>::deallocate()
1495 {
1496 delete rhs_;
1497 delete lhs_;
1498 }
1499
1500 template<class M, class X, class Y>
1501 void OverlappingAssignerILUBase<M,X,Y>::operator()(const size_type& domainIndex)
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>
1517 void OverlappingAssignerILUBase<M,X,Y>::relaxResult(field_type relax)
1518 {
1519 (*lhs_)[i]*=relax;
1520 }
1521
1522 template<class M, class X, class Y>
1523 void OverlappingAssignerILUBase<M,X,Y>::assignResult(block_type& res)
1524 {
1525 res+=(*lhs_)[i++];
1526 }
1527
1528 template<class M, class X, class Y>
1529 X& OverlappingAssignerILUBase<M,X,Y>::lhs()
1530 {
1531 return *lhs_;
1532 }
1533
1534 template<class M, class X, class Y>
1535 Y& OverlappingAssignerILUBase<M,X,Y>::rhs()
1536 {
1537 return *rhs_;
1538 }
1539
1540 template<class M, class X, class Y>
1541 void OverlappingAssignerILUBase<M,X,Y>::resetIndexForNextDomain()
1542 {
1543 i=0;
1544 }
1545
1546 template<typename S, typename T, typename A>
1547 AdditiveAdder<S,BlockVector<T,A> >::AdditiveAdder(BlockVector<T,A>& v_,
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>
1563 void AdditiveAdder<S,BlockVector<T,A> >::axpy()
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>
1571 MultiplicativeAdder<S,BlockVector<T,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>
1589 void MultiplicativeAdder<S,BlockVector<T,A> >::axpy()
1590 {
1591 // nothing to do, as the corrections already relaxed and added in operator()
1592 }
1593
1594
1596}
1597
1598#endif
This file implements a vector space as a tensor product of a given vector space. The number of compon...
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:464
row_type::ConstIterator ConstColIterator
Const iterator to the entries of a row.
Definition: bcrsmatrix.hh:739
Iterator begin()
begin iterator
Definition: densevector.hh:348
Iterator end()
end iterator
Definition: densevector.hh:354
Iterator find(size_type i)
return iterator to given element or end()
Definition: densevector.hh:374
Exact subdomain solver using Dune::DynamicMatrix<T>::solve.
Definition: overlappingschwarz.hh:138
size_type capacity() const
Number of elements for which memory has been allocated.
Definition: dynvector.hh:135
Index Set Interface base class.
Definition: indexidset.hh:76
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:45
D subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:48
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
A constant iterator for the SLList.
Definition: sllist.hh:369
A single linked list.
Definition: sllist.hh:42
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
This file implements a dense matrix with dynamic numbers of rows and columns.
virtual void apply(X &v, const X &d)
Apply the precondtioner.
Definition: overlappingschwarz.hh:1225
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
SeqOverlappingSchwarz(const matrix_type &mat, const rowtodomain_vector &rowToDomain, field_type relaxationFactor=1, bool onTheFly_=true)
Definition: overlappingschwarz.hh:1007
auto max(ADLTag< 0 >, const V &v1, const V &v2)
implements binary Simd::max()
Definition: defaults.hh:79
DVVerbType dvverb(std::cout)
stream for very verbose output.
Definition: stdstreams.hh:93
Various local subdomain solvers based on ILU for SeqOverlappingSchwarz.
Dune namespace.
Definition: alignedallocator.hh:11
Define general preconditioner interface.
Implements a singly linked list together with the necessary iterators.
Templates characterizing the type of a solver.
template meta program for choosing how to add the correction.
Definition: overlappingschwarz.hh:570
Tag that the tells the Schwarz method to be additive.
Definition: overlappingschwarz.hh:118
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:603
Tag that tells the Schwarz method to be multiplicative.
Definition: overlappingschwarz.hh:124
Helper template meta program for application of overlapping Schwarz.
Definition: overlappingschwarz.hh:667
Definition: overlappingschwarz.hh:1101
Category
Definition: solvercategory.hh:21
@ sequential
Category for sequential solvers.
Definition: solvercategory.hh:23
Tag that tells the Schwarz method to be multiplicative and symmetric.
Definition: overlappingschwarz.hh:131
Classes for using SuperLU with ISTL matrices.
Classes for using UMFPack with ISTL matrices.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Dec 21, 23:30, 2024)