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;
Implementation of the BCRSMatrix class.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
The incomplete LU factorization kernels.
Some generic functions for pretty printing vectors and matrices.
Define general, extensible interface for operators. The available implementation wraps a matrix.
Classes providing communication interfaces for overlapping Schwarz methods.
Define general preconditioner interface.
Define base class for scalar product and norm.
Implementations of the inverse operator interface.
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