3 #ifndef DUNE_NOVLPSCHWARZ_HH
4 #define DUNE_NOVLPSCHWARZ_HH
13 #include <dune/common/timer.hh>
58 template<
class M,
class X,
class Y,
class C>
73 typedef typename C::PIS
PIS;
74 typedef typename C::RI
RI;
75 typedef typename RI::RemoteIndexList
RIL;
80 typedef std::multimap<int,int>
MM;
81 typedef std::multimap<int,std::pair<int,RILIterator> >
RIMap;
97 : _A_(A), communication(com), buildcomm(true)
101 virtual void apply (
const X& x, Y& y)
const
105 communication.addOwnerCopyToOwnerCopy(y,y);
116 communication.addOwnerCopyToOwnerCopy(y,y);
129 const PIS& pis=communication.indexSet();
130 const RI& ri = communication.remoteIndices();
135 if (buildcomm ==
true){
138 if (mask.size()!=
static_cast<typename std::vector<double>::size_type
>(x.size())){
139 mask.resize(x.size());
140 for (
typename std::vector<double>::size_type i=0; i<mask.size(); i++)
142 for (
typename PIS::const_iterator i=pis.begin(); i!=pis.end(); ++i)
144 mask[i->local().local()] = 0;
146 mask[i->local().local()] = 2;
149 for (MM::iterator iter = bordercontribution.begin();
150 iter != bordercontribution.end(); ++iter)
151 bordercontribution.erase(iter);
152 std::map<int,int> owner;
157 for (
RowIterator i = _A_.begin(); i != _A_.end(); ++i)
158 if (mask[i.index()] == 0)
159 for (
RIIterator remote = ri.begin(); remote != ri.end(); ++remote){
160 RIL& ril = *(remote->second.first);
161 for (
RILIterator rindex = ril.begin(); rindex != ril.end(); ++rindex)
163 if (rindex->localIndexPair().local().local() == i.index()){
165 (std::make_pair(i.index(),
166 std::pair<int,RILIterator>(remote->first, rindex)));
168 owner.insert(std::make_pair(i.index(),remote->first));
173 for (
RowIterator i = _A_.begin(); i != _A_.end(); ++i){
174 if (mask[i.index()] == 0){
175 std::map<int,int>::iterator it = owner.find(i.index());
177 std::pair<RIMapit, RIMapit> foundiit = rimap.equal_range(i.index());
178 for (
ColIterator j = _A_[i.index()].begin(); j != _A_[i.index()].end(); ++j){
179 if (mask[j.index()] == 0){
181 for (
RIMapit foundi = foundiit.first; foundi != foundiit.second; ++foundi){
182 std::pair<RIMapit, RIMapit> foundjit = rimap.equal_range(j.index());
183 for (
RIMapit foundj = foundjit.first; foundj != foundjit.second; ++foundj)
184 if (foundj->second.first == foundi->second.first)
186 || foundj->second.first == iowner
187 || foundj->second.first < communication.communicator().rank()){
202 bordercontribution.insert(std::pair<int,int>(i.index(),j.index()));
211 for (
RowIterator i = _A_.begin(); i != _A_.end(); ++i){
212 if (mask[i.index()] == 0){
214 for (
ColIterator j = _A_[i.index()].begin(); j != _A_[i.index()].end(); ++j){
215 if (mask[j.index()] == 1)
216 (*j).usmv(alpha,x[j.index()],y[i.index()]);
217 else if (mask[j.index()] == 0){
218 std::pair<MM::iterator, MM::iterator> itp =
219 bordercontribution.equal_range(i.index());
220 for (MM::iterator it = itp.first; it != itp.second; ++it)
221 if ((*it).second == (
int)j.index())
222 (*j).usmv(alpha,x[j.index()],y[i.index()]);
226 else if (mask[i.index()] == 1){
227 for (
ColIterator j = _A_[i.index()].begin(); j != _A_[i.index()].end(); ++j)
228 if (mask[j.index()] != 2)
229 (*j).usmv(alpha,x[j.index()],y[i.index()]);
237 mutable bool buildcomm;
238 mutable std::vector<double> mask;
239 mutable std::multimap<int,int> bordercontribution;
253 template<
class X,
class C>
282 communication.dot(x,y,result);
289 virtual double norm (
const X& x)
291 return communication.norm(x);
298 communication.copyOwnerToAll(x,x);
305 template<
class X,
class C>
338 template<
class C,
class P>
363 : preconditioner(prec), communication(c)
373 preconditioner.pre(x,b);
386 preconditioner.apply(v,d);
387 communication.addOwnerCopyToOwnerCopy(v,v);
397 preconditioner.post(x);