3#ifndef DUNE_ISTL_NOVLPSCHWARZ_HH
4#define DUNE_ISTL_NOVLPSCHWARZ_HH
13#include <dune/common/timer.hh>
15#include <dune/common/hybridutilities.hh>
60 template<
class M,
class X,
class Y,
class C>
75 typedef typename C::PIS
PIS;
76 typedef typename C::RI
RI;
77 typedef typename RI::RemoteIndexList
RIL;
82 typedef std::multimap<int,int>
MM;
83 typedef std::multimap<int,std::pair<int,RILIterator> >
RIMap;
94 : _A_(stackobject_to_shared_ptr(A)), communication(com), buildcomm(true)
98 : _A_(A), communication(com), buildcomm(true)
102 virtual void apply (
const X& x, Y& y)
const
106 communication.addOwnerCopyToOwnerCopy(y,y);
117 communication.addOwnerCopyToOwnerCopy(y,y);
130 const PIS& pis=communication.indexSet();
131 const RI& ri = communication.remoteIndices();
136 if (buildcomm ==
true) {
139 if (mask.size()!=
static_cast<typename std::vector<double>::size_type
>(x.size())) {
140 mask.resize(x.size());
141 for (
typename std::vector<double>::size_type i=0; i<mask.size(); i++)
143 for (
typename PIS::const_iterator i=pis.begin(); i!=pis.end(); ++i)
145 mask[i->local().local()] = 0;
147 mask[i->local().local()] = 2;
150 for (MM::iterator iter = bordercontribution.begin();
151 iter != bordercontribution.end(); ++iter)
152 bordercontribution.erase(iter);
153 std::map<int,int> owner;
158 for (
RowIterator i = _A_->begin(); i != _A_->end(); ++i)
159 if (mask[i.index()] == 0)
160 for (
RIIterator remote = ri.begin(); remote != ri.end(); ++remote) {
161 RIL& ril = *(remote->second.first);
162 for (
RILIterator rindex = ril.begin(); rindex != ril.end(); ++rindex)
164 if (rindex->localIndexPair().local().local() == i.index()) {
166 (std::make_pair(i.index(),
167 std::pair<int,RILIterator>(remote->first, rindex)));
169 owner.insert(std::make_pair(i.index(),remote->first));
174 for (
RowIterator i = _A_->begin(); i != _A_->end(); ++i) {
175 if (mask[i.index()] == 0) {
176 std::map<int,int>::iterator it = owner.find(i.index());
178 std::pair<RIMapit, RIMapit> foundiit = rimap.equal_range(i.index());
179 for (
ColIterator j = (*_A_)[i.index()].begin(); j != (*_A_)[i.index()].end(); ++j) {
180 if (mask[j.index()] == 0) {
182 for (
RIMapit foundi = foundiit.first; foundi != foundiit.second; ++foundi) {
183 std::pair<RIMapit, RIMapit> foundjit = rimap.equal_range(j.index());
184 for (
RIMapit foundj = foundjit.first; foundj != foundjit.second; ++foundj)
185 if (foundj->second.first == foundi->second.first)
187 || foundj->second.first == iowner
188 || foundj->second.first < communication.communicator().rank()) {
203 bordercontribution.insert(std::pair<int,int>(i.index(),j.index()));
212 for (
RowIterator i = _A_->begin(); i != _A_->end(); ++i) {
213 if (mask[i.index()] == 0) {
215 for (
ColIterator j = (*_A_)[i.index()].begin(); j != (*_A_)[i.index()].end(); ++j) {
216 if (mask[j.index()] == 1)
217 Impl::asMatrix(*j).usmv(alpha,x[j.index()],y[i.index()]);
218 else if (mask[j.index()] == 0) {
219 std::pair<MM::iterator, MM::iterator> itp =
220 bordercontribution.equal_range(i.index());
221 for (MM::iterator it = itp.first; it != itp.second; ++it)
222 if ((*it).second == (
int)j.index())
223 Impl::asMatrix(*j).usmv(alpha,x[j.index()],y[i.index()]);
227 else if (mask[i.index()] == 1) {
228 for (
ColIterator j = (*_A_)[i.index()].begin(); j != (*_A_)[i.index()].end(); ++j)
229 if (mask[j.index()] != 2)
230 Impl::asMatrix(*j).usmv(alpha,x[j.index()],y[i.index()]);
244 return communication;
247 std::shared_ptr<const matrix_type> _A_;
249 mutable bool buildcomm;
250 mutable std::vector<double> mask;
251 mutable std::multimap<int,int> bordercontribution;
275 template<
class C,
class P>
277 :
public Preconditioner<typename P::domain_type,typename P::range_type> {
279 using X =
typename P::domain_type;
280 using Y =
typename P::range_type;
304 : _preconditioner(stackobject_to_shared_ptr(p)), _communication(c)
315 : _preconditioner(p), _communication(c)
325 _preconditioner->pre(x,b);
338 _preconditioner->apply(v,d);
339 _communication.addOwnerCopyToOwnerCopy(v,v);
342 template<
bool forward>
345 _preconditioner->template apply<forward>(v,d);
346 _communication.addOwnerCopyToOwnerCopy(v,v);
356 _preconditioner->post(x);
367 std::shared_ptr<P> _preconditioner;
Implementations of the inverse operator interface.
Define general preconditioner interface.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Classes providing communication interfaces for overlapping Schwarz methods.
Implementation of the BCRSMatrix class.
Some generic functions for pretty printing vectors and matrices.
Define base class for scalar product and norm.
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
Define general, extensible interface for operators. The available implementation wraps a matrix.
The incomplete LU factorization kernels.
Definition: allocator.hh:9
A nonoverlapping operator with communication object.
Definition: novlpschwarz.hh:62
C::PIS PIS
Definition: novlpschwarz.hh:75
C communication_type
The type of the communication object.
Definition: novlpschwarz.hh:73
std::multimap< int, std::pair< int, RILIterator > > RIMap
Definition: novlpschwarz.hh:83
C::RI RI
Definition: novlpschwarz.hh:76
void novlp_op_apply(const X &x, Y &y, field_type alpha) const
Definition: novlpschwarz.hh:127
NonoverlappingSchwarzOperator(std::shared_ptr< const matrix_type > A, const communication_type &com)
Definition: novlpschwarz.hh:97
M::ConstColIterator ColIterator
Definition: novlpschwarz.hh:80
virtual void apply(const X &x, Y &y) const
apply operator to x:
Definition: novlpschwarz.hh:102
X domain_type
The type of the domain.
Definition: novlpschwarz.hh:67
virtual SolverCategory::Category category() const
Category of the linear operator (see SolverCategory::Category)
Definition: novlpschwarz.hh:236
RIMap::iterator RIMapit
Definition: novlpschwarz.hh:84
virtual const matrix_type & getmat() const
get matrix via *
Definition: novlpschwarz.hh:122
RIL::const_iterator RILIterator
Definition: novlpschwarz.hh:79
std::multimap< int, int > MM
Definition: novlpschwarz.hh:82
Y range_type
The type of the range.
Definition: novlpschwarz.hh:69
M matrix_type
The type of the matrix we operate on.
Definition: novlpschwarz.hh:65
M::ConstRowIterator RowIterator
Definition: novlpschwarz.hh:81
const communication_type & getCommunication() const
Get the object responsible for communication.
Definition: novlpschwarz.hh:242
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const
apply operator to x, scale and add:
Definition: novlpschwarz.hh:110
RI::RemoteIndexList RIL
Definition: novlpschwarz.hh:77
X::field_type field_type
The field type of the range.
Definition: novlpschwarz.hh:71
NonoverlappingSchwarzOperator(const matrix_type &A, const communication_type &com)
constructor: just store a reference to a matrix.
Definition: novlpschwarz.hh:93
RI::const_iterator RIIterator
Definition: novlpschwarz.hh:78
Traits class for generically constructing non default constructable types.
Definition: construction.hh:37
Nonoverlapping parallel preconditioner.
Definition: novlpschwarz.hh:277
NonoverlappingBlockPreconditioner(P &p, const communication_type &c)
Constructor.
Definition: novlpschwarz.hh:303
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: novlpschwarz.hh:360
virtual void apply(domain_type &v, const range_type &d)
Apply the preconditioner.
Definition: novlpschwarz.hh:333
P::range_type range_type
The range type of the preconditioner.
Definition: novlpschwarz.hh:285
NonoverlappingBlockPreconditioner(const std::shared_ptr< P > &p, const communication_type &c)
Constructor.
Definition: novlpschwarz.hh:314
virtual void post(domain_type &x)
Clean up.
Definition: novlpschwarz.hh:354
C communication_type
The type of the communication object.
Definition: novlpschwarz.hh:287
virtual void pre(domain_type &x, range_type &b)
Prepare the preconditioner.
Definition: novlpschwarz.hh:323
void apply(X &v, const Y &d)
Apply one step of the preconditioner to the system A(v)=d.
Definition: novlpschwarz.hh:343
P::domain_type domain_type
The domain type of the preconditioner.
Definition: novlpschwarz.hh:283
A linear operator exporting itself in matrix form.
Definition: operators.hh:107
@ owner
Definition: owneroverlapcopy.hh:59
@ copy
Definition: owneroverlapcopy.hh:59
@ overlap
Definition: owneroverlapcopy.hh:59
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:30
Category
Definition: solvercategory.hh:21
@ nonoverlapping
Category for non-overlapping solvers.
Definition: solvercategory.hh:25