3 #ifndef DUNE_ISTL_OVERLAPPINGSCHWARZ_HH 4 #define DUNE_ISTL_OVERLAPPINGSCHWARZ_HH 10 #include <dune/common/dynmatrix.hh> 11 #include <dune/common/sllist.hh> 12 #include <dune/common/unused.hh> 35 template<
class M,
class X,
class TM,
class TD,
class TA>
36 class SeqOverlappingSchwarz;
41 template<
class I,
class S,
class D>
50 typedef typename AtomInitializer::Matrix
Matrix;
51 typedef typename Matrix::const_iterator
Iter;
58 const IndexSet& indices,
59 const subdomain_vector& domains);
78 typedef std::map<size_type,size_type> Map;
79 typedef typename Map::iterator iterator;
80 typedef typename Map::const_iterator const_iterator;
84 void insert(size_type grow);
86 const_iterator find(size_type grow)
const;
88 iterator find(size_type grow);
92 const_iterator begin()
const;
96 const_iterator end()
const;
99 std::map<size_type,size_type> map_;
104 typedef typename InitializerList::iterator InitIterator;
105 typedef typename IndexSet::const_iterator IndexIteratur;
106 InitializerList* initializers;
107 const IndexSet *indices;
108 mutable std::vector<IndexMap> indexMaps;
109 const subdomain_vector& domains;
135 template<
class M,
class X,
class Y>
139 template<
class K,
int n,
class Al,
class X,
class Y>
157 void apply (DynamicVector<field_type>& v, DynamicVector<field_type>& d)
159 assert(v.size() > 0);
160 assert(v.size() == d.size());
161 assert(A.rows() <= v.size());
162 assert(A.cols() <= v.size());
163 size_t sz = A.rows();
167 v.resize(v.capacity());
168 d.resize(d.capacity());
181 size_t sz = rowset.size();
183 typedef typename S::const_iterator SIter;
185 for(SIter rowIdx = rowset.begin(), rowEnd=rowset.end();
186 rowIdx!= rowEnd; ++rowIdx, r++)
189 for(SIter colIdx = rowset.begin(), colEnd=rowset.end();
190 colIdx!= colEnd; ++colIdx, c++)
192 if (BCRS[*rowIdx].find(*colIdx) == BCRS[*rowIdx].
end())
194 for (
size_t i=0; i<n; i++)
196 for (
size_t j=0; j<n; j++)
198 A[r*n+i][c*n+j] = BCRS[*rowIdx][*colIdx][i][j];
208 template<
typename T,
bool tag>
216 template<
class K,
int n,
class Al,
class X,
class Y>
245 void resetIndexForNextDomain();
252 DynamicVector<K> & lhs();
259 DynamicVector<K> & rhs();
266 void relaxResult(field_type relax);
272 void operator()(
const size_type& domainIndex);
282 void assignResult(block_type& res);
288 const matrix_type*
mat;
291 DynamicVector<field_type> * rhs_;
294 DynamicVector<field_type> * lhs_;
302 std::size_t maxlength_;
305 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 306 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
324 const range_type& b, range_type& x);
335 void resetIndexForNextDomain();
353 void relaxResult(field_type relax);
359 void operator()(
const size_type& domain);
368 void assignResult(block_type& res);
374 const matrix_type*
mat;
386 std::size_t maxlength_;
389 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 391 template<
class M,
class X,
class Y>
421 void resetIndexForNextDomain();
439 void relaxResult(field_type relax);
445 void operator()(
const size_type& domain);
454 void assignResult(block_type& res);
474 template<
class M,
class X,
class Y>
493 template<
class M,
class X,
class Y>
511 template<
typename S,
typename T>
515 template<
typename S,
typename T,
typename A,
int n>
521 void operator()(
const size_type& domain);
531 template<
typename S,
typename T>
535 template<
typename S,
typename T,
typename A,
int n>
541 void operator()(
const size_type& domain);
559 template<
typename T,
class X,
class S>
563 template<
class X,
class S>
569 template<
class X,
class S>
575 template<
class X,
class S>
592 template<
typename T1,
typename T2,
bool forward>
600 static solver_iterator
begin(solver_vector& sv)
605 static solver_iterator
end(solver_vector& sv)
609 static domain_iterator
begin(
const subdomain_vector& sv)
614 static domain_iterator
end(
const subdomain_vector& sv)
620 template<
typename T1,
typename T2>
628 static solver_iterator
begin(solver_vector& sv)
633 static solver_iterator
end(solver_vector& sv)
637 static domain_iterator
begin(
const subdomain_vector& sv)
642 static domain_iterator
end(
const subdomain_vector& sv)
662 static void apply(smoother& sm, range_type& v,
const range_type& b)
664 sm.template apply<true>(v, b);
668 template<
class M,
class X,
class TD,
class TA>
674 static void apply(smoother& sm, range_type& v,
const range_type& b)
676 sm.template apply<true>(v, b);
677 sm.template apply<false>(v, b);
681 template<
class T,
bool tag>
688 template<
class K,
int n,
class Al,
class X,
class Y>
692 template<
class RowToDomain,
class Solvers,
class SubDomains>
693 static std::size_t assembleLocalProblems(
const RowToDomain& rowToDomain,
const matrix_type&
mat,
694 Solvers& solvers,
const SubDomains& domains,
698 template<
template<
class>
class S,
typename T,
typename A,
int m,
int n>
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,
708 template<
class M,
class X,
class Y>
712 template<
class RowToDomain,
class Solvers,
class SubDomains>
713 static std::size_t assembleLocalProblems(
const RowToDomain& rowToDomain,
const matrix_type&
mat,
714 Solvers& solvers,
const SubDomains& domains,
718 template<
class M,
class X,
class Y>
723 template<
class M,
class X,
class Y>
779 typedef std::set<size_type, std::less<size_type>,
780 typename TA::template rebind<size_type>::other>
784 typedef std::vector<subdomain_type, typename TA::template rebind<subdomain_type>::other>
subdomain_vector;
787 typedef SLList<size_type, typename TA::template rebind<size_type>::other>
subdomain_list;
790 typedef std::vector<subdomain_list, typename TA::template rebind<subdomain_list>::other >
rowtodomain_vector;
796 typedef std::vector<slu, typename TA::template rebind<slu>::other>
slu_vector;
817 field_type relaxationFactor=1,
bool onTheFly_=
true);
831 field_type relaxationFactor=1,
bool onTheFly_=
true);
838 virtual void pre (X& x, X& b)
840 DUNE_UNUSED_PARAMETER(x);
841 DUNE_UNUSED_PARAMETER(b);
849 virtual void apply (X& v,
const X& d);
858 DUNE_UNUSED_PARAMETER(x);
861 template<
bool forward>
862 void apply(X& v,
const X& d);
867 subdomain_vector subDomains;
870 typename M::size_type maxlength;
877 template<
class I,
class S,
class D>
881 : initializers(&il), indices(&idx), indexMaps(il.size()), domains(domains_)
885 template<
class I,
class S,
class D>
888 typedef typename IndexSet::value_type::const_iterator iterator;
889 for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
890 (*initializers)[*domain].addRowNnz(row, domains[*domain]);
891 indexMaps[*domain].insert(row.index());
895 template<
class I,
class S,
class D>
898 std::for_each(initializers->begin(), initializers->end(),
899 std::mem_fun_ref(&AtomInitializer::allocateMatrixStorage));
900 std::for_each(initializers->begin(), initializers->end(),
901 std::mem_fun_ref(&AtomInitializer::allocateMarker));
904 template<
class I,
class S,
class D>
907 typedef typename IndexSet::value_type::const_iterator iterator;
908 for(iterator domain=(*indices)[row.index()].begin(); domain != (*indices)[row.index()].end(); ++domain) {
909 typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.
index());
910 if(v!= indexMaps[*domain].end()) {
911 (*initializers)[*domain].countEntries(indexMaps[*domain].find(col.
index())->second);
916 template<
class I,
class S,
class D>
919 std::for_each(initializers->begin(), initializers->end(),
920 std::mem_fun_ref(&AtomInitializer::calcColstart));
923 template<
class I,
class S,
class D>
926 typedef typename IndexSet::value_type::const_iterator iterator;
927 for(iterator domain=(*indices)[row.index()].begin(); domain!= (*indices)[row.index()].end(); ++domain) {
928 typename std::map<size_type,size_type>::const_iterator v = indexMaps[*domain].find(col.
index());
929 if(v!= indexMaps[*domain].end()) {
930 assert(indexMaps[*domain].end()!=indexMaps[*domain].find(row.index()));
931 (*initializers)[*domain].copyValue(col, indexMaps[*domain].find(row.index())->second,
937 template<
class I,
class S,
class D>
940 std::vector<IndexMap>().swap(indexMaps);
941 std::for_each(initializers->begin(), initializers->end(),
942 std::mem_fun_ref(&AtomInitializer::createMatrix));
945 template<
class I,
class S,
class D>
950 template<
class I,
class S,
class D>
953 assert(map_.find(grow)==map_.end());
954 map_.insert(std::make_pair(grow, row++));
957 template<
class I,
class S,
class D>
958 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
961 return map_.find(grow);
964 template<
class I,
class S,
class D>
965 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
968 return map_.find(grow);
971 template<
class I,
class S,
class D>
972 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
978 template<
class I,
class S,
class D>
979 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
985 template<
class I,
class S,
class D>
986 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::const_iterator
992 template<
class I,
class S,
class D>
993 typename OverlappingSchwarzInitializer<I,S,D>::IndexMap::iterator
999 template<
class M,
class X,
class TM,
class TD,
class TA>
1002 :
mat(mat_), relax(relaxationFactor), onTheFly(fly)
1004 typedef typename rowtodomain_vector::const_iterator RowDomainIterator;
1005 typedef typename subdomain_list::const_iterator DomainIterator;
1006 #ifdef DUNE_ISTL_WITH_CHECKING 1007 assert(rowToDomain.size()==
mat.N());
1008 assert(rowToDomain.size()==
mat.M());
1010 for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1011 assert(iter->size()>0);
1016 for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter)
1017 for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1018 domains=std::max(domains, *d);
1021 solvers.resize(domains);
1022 subDomains.resize(domains);
1026 for(RowDomainIterator iter=rowToDomain.begin(); iter != rowToDomain.end(); ++iter, ++row)
1027 for(DomainIterator d=iter->begin(); d != iter->end(); ++d)
1028 subDomains[*d].insert(row);
1030 #ifdef DUNE_ISTL_WITH_CHECKING 1032 typedef typename subdomain_vector::const_iterator iterator;
1033 for(iterator iter=subDomains.begin(); iter != subDomains.end(); ++iter) {
1034 typedef typename subdomain_type::const_iterator entry_iterator;
1035 Dune::dvverb<<
"domain "<<i++<<
":";
1036 for(entry_iterator entry = iter->begin(); entry != iter->end(); ++entry) {
1037 Dune::dvverb<<
" "<<*entry;
1039 Dune::dvverb<<std::endl;
1046 template<
class M,
class X,
class TM,
class TD,
class TA>
1051 :
mat(mat_), solvers(sd.size()), subDomains(sd), relax(relaxationFactor),
1054 typedef typename subdomain_vector::const_iterator DomainIterator;
1056 #ifdef DUNE_ISTL_WITH_CHECKING 1059 for(DomainIterator d=sd.begin(); d != sd.end(); ++d,++i) {
1061 assert(d->size()>0);
1062 typedef typename DomainIterator::value_type::const_iterator entry_iterator;
1063 Dune::dvverb<<
"domain "<<i<<
":";
1064 for(entry_iterator entry = d->begin(); entry != d->end(); ++entry) {
1065 Dune::dvverb<<
" "<<*entry;
1067 Dune::dvverb<<std::endl;
1077 for(DomainIterator domain=sd.begin(); domain != sd.end(); ++domain, ++domainId) {
1078 typedef typename subdomain_type::const_iterator iterator;
1079 for(iterator row=domain->begin(); row != domain->end(); ++row)
1080 rowToDomain[*row].push_back(domainId);
1096 template<
typename T,
typename A,
int n,
int m>
1099 template<
class Domain>
1107 template<
class K,
int n,
class Al,
class X,
class Y>
1108 template<
class RowToDomain,
class Solvers,
class SubDomains>
1111 assembleLocalProblems(
const RowToDomain& rowToDomain,
1114 const SubDomains& subDomains,
1117 DUNE_UNUSED_PARAMETER(onTheFly);
1118 DUNE_UNUSED_PARAMETER(rowToDomain);
1119 DUNE_UNUSED_PARAMETER(mat);
1120 DUNE_UNUSED_PARAMETER(solvers);
1121 typedef typename SubDomains::const_iterator DomainIterator;
1122 std::size_t maxlength = 0;
1126 for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1127 maxlength=std::max(maxlength, domain->size());
1133 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 1134 template<
template<
class>
class S,
typename T,
typename A,
int m,
int n>
1135 template<
class RowToDomain,
class Solvers,
class SubDomains>
1139 const SubDomains& subDomains,
1142 typedef typename S<BCRSMatrix<FieldMatrix<T,m,n>,A> >::MatrixInitializer MatrixInitializer;
1143 typedef typename std::vector<MatrixInitializer>::iterator InitializerIterator;
1144 typedef typename SubDomains::const_iterator DomainIterator;
1145 typedef typename Solvers::iterator SolverIterator;
1146 std::size_t maxlength = 0;
1149 for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1150 maxlength=std::max(maxlength, domain->size());
1151 maxlength*=mat[0].
begin()->N();
1154 DomainIterator domain=subDomains.begin();
1157 std::vector<MatrixInitializer> initializers(subDomains.size());
1159 SolverIterator solver=solvers.begin();
1160 for(InitializerIterator initializer=initializers.begin(); initializer!=initializers.end();
1161 ++initializer, ++solver, ++domain) {
1165 *initializer=MatrixInitializer(solver->getInternalMatrix());
1170 RowToDomain, SubDomains> Initializer;
1172 Initializer initializer(initializers, rowToDomain, subDomains);
1177 for (SolverIterator solverIt = solvers.begin(); solverIt != solvers.end(); ++solverIt)
1179 assert(solverIt->getInternalMatrix().N() == solverIt->getInternalMatrix().M());
1180 maxlength = std::max(maxlength, solverIt->getInternalMatrix().N());
1186 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 1188 template<
class M,
class X,
class Y>
1189 template<
class RowToDomain,
class Solvers,
class SubDomains>
1193 const SubDomains& subDomains,
1196 DUNE_UNUSED_PARAMETER(rowToDomain);
1197 typedef typename SubDomains::const_iterator DomainIterator;
1198 typedef typename Solvers::iterator SolverIterator;
1199 std::size_t maxlength = 0;
1202 for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end(); ++domain)
1203 maxlength=std::max(maxlength, domain->size());
1206 SolverIterator solver=solvers.begin();
1207 for(DomainIterator domain=subDomains.begin(); domain!=subDomains.end();
1208 ++domain, ++solver) {
1209 solver->setSubMatrix(mat, *domain);
1210 maxlength=std::max(maxlength, domain->size());
1219 template<
class M,
class X,
class TM,
class TD,
class TA>
1225 template<
class M,
class X,
class TM,
class TD,
class TA>
1226 template<
bool forward>
1242 Adder adder(v, x, assigner, relax);
1246 std::for_each(domain->begin(), domain->end(), assigner);
1247 assigner.resetIndexForNextDomain();
1251 sdsolver.setSubMatrix(
mat, *domain);
1253 sdsolver.apply(assigner.lhs(), assigner.rhs());
1255 solver->apply(assigner.lhs(), assigner.rhs());
1260 std::for_each(domain->begin(), domain->end(), adder);
1261 assigner.resetIndexForNextDomain();
1266 assigner.deallocate();
1269 template<
class K,
int n,
class Al,
class X,
class Y>
1272 const X& b_, Y& x_) :
1274 rhs_( new DynamicVector<
field_type>(maxlength, 42) ),
1275 lhs_( new DynamicVector<
field_type>(maxlength, -42) ),
1279 maxlength_(maxlength)
1282 template<
class K,
int n,
class Al,
class X,
class Y>
1284 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1291 template<
class K,
int n,
class Al,
class X,
class Y>
1293 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1294 ::resetIndexForNextDomain()
1299 template<
class K,
int n,
class Al,
class X,
class Y>
1301 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1307 template<
class K,
int n,
class Al,
class X,
class Y>
1309 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1315 template<
class K,
int n,
class Al,
class X,
class Y>
1317 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1323 template<
class K,
int n,
class Al,
class X,
class Y>
1325 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1332 assert(i<maxlength_);
1333 rhs()[i]=(*b)[domainIndex][j];
1340 for(col_iterator
col=(*
mat)[domainIndex].begin();
col!=(*mat)[domainIndex].end(); ++
col) {
1342 (*col).mv((*x)[
col.index()], tmp);
1345 assert(i<maxlength_);
1352 assert(i<maxlength_);
1353 rhs()[i]=(*b)[domainIndex][j];
1359 for(col_iterator
col=(*
mat)[domainIndex].begin();
col!=(*mat)[domainIndex].end(); ++
col) {
1361 rhs()[i]-=(*col)[j][k] * (*x)[
col.index()][k];
1368 template<
class K,
int n,
class Al,
class X,
class Y>
1370 OverlappingAssignerHelper< DynamicMatrixSubdomainSolver< BCRSMatrix< FieldMatrix<K,n,n>, Al>, X, Y >,
false>
1375 assert(i<maxlength_);
1380 #if HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 1382 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1390 x(&x_), i(0), maxlength_(maxlength)
1397 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1398 void OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::deallocate()
1404 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1405 void OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::operator()(
const size_type& domainIndex)
1411 assert(i<maxlength_);
1412 rhs_[i]=(*b)[domainIndex][j];
1420 for(col_iterator
col=(*
mat)[domainIndex].begin();
col!=(*mat)[domainIndex].end(); ++
col) {
1422 (*col).mv((*x)[
col.index()], tmp);
1425 assert(i<maxlength_);
1433 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1434 void OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::relaxResult(
field_type relax)
1437 assert(i<maxlength_);
1443 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1444 void OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::assignResult(
block_type& res)
1448 assert(i<maxlength_);
1453 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1454 void OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::resetIndexForNextDomain()
1459 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1460 typename OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>
::field_type*
1461 OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::lhs()
1466 template<
template<
class>
class S,
int n,
int m,
typename T,
typename A>
1467 typename OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>
::field_type*
1468 OverlappingAssignerHelper<S<BCRSMatrix<FieldMatrix<T,n,m>,A> >,
true>::rhs()
1473 #endif // HAVE_SUPERLU || HAVE_SUITESPARSE_UMFPACK 1475 template<
class M,
class X,
class Y>
1484 rhs_=
new Y(maxlength);
1485 lhs_ =
new X(maxlength);
1488 template<
class M,
class X,
class Y>
1495 template<
class M,
class X,
class Y>
1498 (*rhs_)[i]=(*b)[domainIndex];
1501 typedef typename matrix_type::ConstColIterator col_iterator;
1504 for(col_iterator
col=(*mat)[domainIndex].begin();
col!=(*mat)[domainIndex].end(); ++
col) {
1505 (*col).mmv((*x)[
col.index()], (*rhs_)[i]);
1511 template<
class M,
class X,
class Y>
1517 template<
class M,
class X,
class Y>
1523 template<
class M,
class X,
class Y>
1529 template<
class M,
class X,
class Y>
1535 template<
class M,
class X,
class Y>
1541 template<
typename S,
typename T,
typename A,
int n>
1546 : v(&v_), x(&x_), assigner(&assigner_), relax(relax_)
1549 template<
typename S,
typename T,
typename A,
int n>
1553 assigner->assignResult((*v)[domainIndex]);
1557 template<
typename S,
typename T,
typename A,
int n>
1558 void AdditiveAdder<S,BlockVector<FieldVector<T,n>,A> >::axpy()
1565 template<
typename S,
typename T,
typename A,
int n>
1570 : x(&x_), assigner(&assigner_), relax(relax_)
1572 DUNE_UNUSED_PARAMETER(v_);
1576 template<
typename S,
typename T,
typename A,
int n>
1577 void MultiplicativeAdder<S,BlockVector<FieldVector<T,n>,A> >::operator()(
const size_type& domainIndex)
1580 assigner->relaxResult(relax);
1581 assigner->assignResult((*x)[domainIndex]);
1585 template<
typename S,
typename T,
typename A,
int n>
1586 void MultiplicativeAdder<S,BlockVector<FieldVector<T,n>,A> >::axpy()
Tag that tells the Schwarz method to be multiplicative.
Definition: overlappingschwarz.hh:121
std::vector< subdomain_list, typename TA::template rebind< subdomain_list >::other > rowtodomain_vector
The vector type containing the row index to subdomain mapping.
Definition: overlappingschwarz.hh:790
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:578
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:42
X & lhs()
Get the local left hand side.
Definition: overlappingschwarz.hh:1524
void allocate()
Definition: overlappingschwarz.hh:896
solver_vector::reverse_iterator solver_iterator
Definition: overlappingschwarz.hh:624
Col col
Definition: matrixmatrix.hh:349
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:486
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:224
X::field_type field_type
The field type of the preconditioner.
Definition: overlappingschwarz.hh:770
solver_vector::iterator solver_iterator
Definition: overlappingschwarz.hh:596
void resetIndexForNextDomain()
Resets the local index to zero.
Definition: overlappingschwarz.hh:1536
Y::block_type block_type
Definition: overlappingschwarz.hh:399
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:149
range_type::block_type block_type
Definition: overlappingschwarz.hh:312
void copyToColCompMatrix(F &initializer, const MRS &mrs)
Definition: colcompmatrix.hh:428
std::vector< subdomain_type, typename TA::template rebind< subdomain_type >::other > subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:784
A::size_type size_type
The type for the index access and the size.
Definition: bcrsmatrix.hh:457
A vector of blocks with memory management.
Definition: bvector.hh:312
Exact subdomain solver using Dune::DynamicMatrix<T>::solve.
Definition: overlappingschwarz.hh:136
virtual void pre(X &x, X &b)
Prepare the preconditioner.
Definition: overlappingschwarz.hh:838
Definition: overlappingschwarz.hh:1094
Tag that the tells the schwarz method to be additive.
Definition: overlappingschwarz.hh:115
void setSubMatrix(const M &BCRS, S &rowset)
Set the data of the local problem.
Definition: overlappingschwarz.hh:179
Definition: overlappingschwarz.hh:709
void apply(DynamicVector< field_type > &v, DynamicVector< field_type > &d)
Apply the subdomain solver.
Definition: overlappingschwarz.hh:157
BCRSMatrix< FieldMatrix< T, m, n >, A > matrix_type
Definition: overlappingschwarz.hh:701
Y & rhs()
Get the local right hand side.
Definition: overlappingschwarz.hh:1530
void calcColstart() const
Definition: overlappingschwarz.hh:917
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:633
Iterator implementation class.
Definition: basearray.hh:84
Definition: overlappingschwarz.hh:512
A::size_type size_type
Definition: overlappingschwarz.hh:538
Iterator end()
Get iterator to one beyond last row.
Definition: bcrsmatrix.hh:640
Exact subdomain solver using ILU(p) with appropriate p.
Definition: ilusubdomainsolver.hh:74
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:637
Definition: overlappingschwarz.hh:392
Iterator begin()
Get iterator to first row.
Definition: bcrsmatrix.hh:634
T2 subdomain_vector
Definition: overlappingschwarz.hh:625
Definition: overlappingschwarz.hh:209
subdomain_vector::const_reverse_iterator domain_iterator
Definition: overlappingschwarz.hh:626
std::vector< slu, typename TA::template rebind< slu >::other > slu_vector
The vector type containing subdomain solvers.
Definition: overlappingschwarz.hh:796
Y range_type
Definition: overlappingschwarz.hh:222
Templates characterizing the type of a solver.
void deallocate()
Deallocates memory of the local vector.
Definition: overlappingschwarz.hh:1489
void copyValue(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:924
BCRSMatrix< FieldMatrix< K, n, n >, Al > matrix_type
Definition: overlappingschwarz.hh:691
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:642
Matrix::const_iterator Iter
Definition: overlappingschwarz.hh:51
std::remove_const< M >::type rilu_type
Definition: overlappingschwarz.hh:147
Definition: colcompmatrix.hh:160
Category for sequential solvers.
Definition: solvercategory.hh:21
Define general preconditioner interface.
Definition: overlappingschwarz.hh:532
BCRSMatrix< FieldMatrix< K, n, n >, Al > matrix_type
Definition: overlappingschwarz.hh:220
T2 subdomain_vector
Definition: overlappingschwarz.hh:597
subdomain_vector::const_iterator domain_iterator
Definition: overlappingschwarz.hh:598
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:1518
T1 solver_vector
Definition: overlappingschwarz.hh:623
SeqOverlappingSchwarz< M, X, SymmetricMultiplicativeSchwarzMode, TD, TA > smoother
Definition: overlappingschwarz.hh:671
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:662
smoother::range_type range_type
Definition: overlappingschwarz.hh:672
T1 solver_vector
Definition: overlappingschwarz.hh:595
static int size(const Domain &d)
Definition: overlappingschwarz.hh:1100
Various local subdomain solvers based on ILU for SeqOverlappingSchwarz.
template meta program for choosing how to add the correction.
Definition: overlappingschwarz.hh:560
Definition: matrixutils.hh:25
SeqOverlappingSchwarz(const matrix_type &mat, const subdomain_vector &subDomains, field_type relaxationFactor=1, bool onTheFly_=true)
Construct the overlapping Schwarz method.
Definition: overlappingschwarz.hh:1047
size_type index() const
return index
Definition: basearray.hh:109
M matrix_type
Definition: overlappingschwarz.hh:395
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:401
static solver_iterator end(solver_vector &sv)
Definition: overlappingschwarz.hh:605
OverlappingAssignerILUBase(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:1476
A sparse block matrix with compressed row storage.
Definition: bcrsmatrix.hh:422
M::field_type field_type
Definition: overlappingschwarz.hh:397
K field_type
Definition: overlappingschwarz.hh:221
InitializerList::value_type AtomInitializer
Definition: overlappingschwarz.hh:49
Matrix & mat
Definition: matrixmatrix.hh:345
Definition: basearray.hh:19
TA allocator
The allocator to use.
Definition: overlappingschwarz.hh:776
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:752
matrix_type::size_type size_type
Definition: overlappingschwarz.hh:314
std::set< size_type, std::less< size_type >, typename TA::template rebind< size_type >::other > subdomain_type
The type for the subdomain to row index mapping.
Definition: overlappingschwarz.hh:781
X range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:757
This file implements a vector space as a tensor product of a given vector space. The number of compon...
M matrix_type
Definition: overlappingschwarz.hh:711
Helper template meta program for application of overlapping schwarz.
Definition: overlappingschwarz.hh:593
AdditiveAdder< S, X > Adder
Definition: overlappingschwarz.hh:566
range_type::block_type block_type
Definition: overlappingschwarz.hh:223
IndexSet::size_type size_type
Definition: overlappingschwarz.hh:55
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:600
Implementation of the BCRSMatrix class.
range_type::field_type field_type
Definition: overlappingschwarz.hh:311
Tag that tells the Schwarz method to be multiplicative and symmetric.
Definition: overlappingschwarz.hh:128
static domain_iterator begin(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:609
T smoother
Definition: overlappingschwarz.hh:659
void relaxResult(field_type relax)
relax the result.
Definition: overlappingschwarz.hh:1512
I InitializerList
Definition: overlappingschwarz.hh:48
virtual void post(X &x)
Postprocess the preconditioner.
Definition: overlappingschwarz.hh:856
S IndexSet
Definition: overlappingschwarz.hh:54
smoother::range_type range_type
Definition: overlappingschwarz.hh:660
OverlappingSchwarzInitializer(InitializerList &il, const IndexSet &indices, const subdomain_vector &domains)
Definition: overlappingschwarz.hh:878
S< BCRSMatrix< FieldMatrix< T, n, m >, A > >::range_type range_type
Definition: overlappingschwarz.hh:310
Matrix::row_type::const_iterator CIter
Definition: overlappingschwarz.hh:52
Classes for using UMFPack with ISTL matrices.
TM Mode
The mode (additive or multiplicative) of the Schwarz method.
Definition: overlappingschwarz.hh:765
MultiplicativeAdder< S, X > Adder
Definition: overlappingschwarz.hh:572
AtomInitializer::Matrix Matrix
Definition: overlappingschwarz.hh:50
iterator class for sequential access
Definition: basearray.hh:563
TD slu
The type for the subdomain solver in use.
Definition: overlappingschwarz.hh:793
D subdomain_vector
The vector type containing the subdomain to row index mapping.
Definition: overlappingschwarz.hh:46
static solver_iterator begin(solver_vector &sv)
Definition: overlappingschwarz.hh:628
Sequential overlapping Schwarz preconditioner.
Definition: colcompmatrix.hh:157
void operator()(const size_type &domain)
calculate one entry of the local defect.
Definition: overlappingschwarz.hh:1496
Definition: ilusubdomainsolver.hh:107
A::size_type size_type
Definition: overlappingschwarz.hh:518
Helper template meta program for application of overlapping schwarz.
Definition: overlappingschwarz.hh:657
K field_type
Definition: overlappingschwarz.hh:146
Y range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:151
OverlappingAssignerHelper(std::size_t maxlength, const M &mat, const Y &b, X &x)
Constructor.
Definition: overlappingschwarz.hh:505
M matrix_type
The type of the matrix to precondition.
Definition: overlappingschwarz.hh:747
std::remove_const< M >::type matrix_type
The matrix type the preconditioner is for.
Definition: overlappingschwarz.hh:145
static void apply(smoother &sm, range_type &v, const range_type &b)
Definition: overlappingschwarz.hh:674
void createMatrix() const
Definition: overlappingschwarz.hh:938
static domain_iterator end(const subdomain_vector &sv)
Definition: overlappingschwarz.hh:614
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
Classes for using SuperLU with ISTL matrices.
void addRowNnz(const Iter &row)
Definition: overlappingschwarz.hh:886
void countEntries(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:905
SLList< size_type, typename TA::template rebind< size_type >::other > subdomain_list
The type for the row to subdomain mapping.
Definition: overlappingschwarz.hh:787
matrix_type::size_type size_type
The return type of the size method.
Definition: overlappingschwarz.hh:773