dune-istl  2.2.0
transfer.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set ts=8 sw=2 et sts=2:
3 #ifndef DUNE_AMGTRANSFER_HH
4 #define DUNE_AMGTRANSFER_HH
5 
6 #include<dune/istl/bvector.hh>
11 #include<dune/common/exceptions.hh>
12 
13 namespace Dune
14 {
15  namespace Amg
16  {
17 
28  template<class V1, class V2, class T>
29  class Transfer
30  {
31 
32  public:
33  typedef V1 Vertex;
34  typedef V2 Vector;
35 
36  template<typename T1, typename R>
37  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
38  Vector& fineRedist,T1 damp, R& redistributor=R());
39 
40  template<typename T1, typename R>
41  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
42  T1 damp);
43 
44  static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
45  T& comm);
46  };
47 
48  template<class V,class V1>
50  {
51  public:
52  typedef V Vertex;
53  typedef V1 Vector;
55  template<typename T1>
56  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
57  Vector& fineRedist, T1 damp,
59  const Redist& redist=Redist());
60  template<typename T1>
61  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
62  T1 damp,
64 
65 
66  static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
67  const SequentialInformation& comm);
68  };
69 
70 #if HAVE_MPI
71 
72  template<class V,class V1, class T1, class T2>
74  {
75  public:
76  typedef V Vertex;
77  typedef V1 Vector;
79  template<typename T3>
80  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
81  Vector& fineRedist, T3 damp, OwnerOverlapCopyCommunication<T1,T2>& comm,
82  const Redist& redist);
83  template<typename T3>
84  static void prolongate(const AggregatesMap<Vertex>& aggregates, Vector& coarse, Vector& fine,
86 
87  static void restrict(const AggregatesMap<Vertex>& aggregates, Vector& coarse, const Vector& fine,
89  };
90 
91 #endif
92 
93  template<class V, class V1>
94  template<typename T>
95  inline void
97  Vector& coarse, Vector& fine, Vector& fineRedist,
98  T damp,
99  const SequentialInformation& comm,
100  const Redist& redist)
101  {
102  prolongate(aggregates, coarse, fine, damp);
103  }
104  template<class V, class V1>
105  template<typename T>
106  inline void
108  Vector& coarse, Vector& fine,
109  T damp,
110  const SequentialInformation& comm)
111  {
112  typedef typename Vector::iterator Iterator;
113 
114  Iterator end = coarse.end();
115  Iterator begin= coarse.begin();
116  for(;begin!=end;++begin)
117  *begin*=damp;
118  end=fine.end();
119  begin=fine.begin();
120 
121  for(Iterator block=begin; block != end; ++block){
122  std::ptrdiff_t index=block-begin;
123  const Vertex& vertex = aggregates[index];
124  if(vertex != AggregatesMap<Vertex>::ISOLATED)
125  *block += coarse[aggregates[index]];
126  }
127  }
128 
129  template<class V, class V1>
130  inline void
132  Vector& coarse,
133  const Vector& fine,
134  const SequentialInformation& comm)
135  {
136  // Set coarse vector to zero
137  coarse=0;
138 
139  typedef typename Vector::const_iterator Iterator;
140  Iterator end = fine.end();
141  Iterator begin=fine.begin();
142 
143  for(Iterator block=begin; block != end; ++block){
144  const Vertex& vertex = aggregates[block-begin];
145  if(vertex != AggregatesMap<Vertex>::ISOLATED)
146  coarse[vertex] += *block;
147  }
148  }
149 
150 #if HAVE_MPI
151  template<class V, class V1, class T1, class T2>
152  template<typename T3>
154  Vector& coarse, Vector& fine,
155  Vector& fineRedist, T3 damp,
157  const Redist& redist)
158  {
159  if(fineRedist.size()>0)
160  // we operated on the coarse level
161  Transfer<V,V1,SequentialInformation>::prolongate(aggregates, coarse, fineRedist, damp);
162 
163  // TODO This could be accomplished with one communication, too!
164  redist.redistributeBackward(fine, fineRedist);
165  comm.copyOwnerToAll(fine,fine);
166  }
167 
168  template<class V, class V1, class T1, class T2>
169  template<typename T3>
171  Vector& coarse, Vector& fine,
172  T3 damp,
174  {
175  Transfer<V,V1,SequentialInformation>::prolongate(aggregates, coarse, fine, damp);
176  }
177  template<class V, class V1, class T1, class T2>
179  Vector& coarse, const Vector& fine,
181  {
183  // We need this here to avoid it in the smoothers on the coarse level.
184  // There (in the preconditioner d is const.
185  comm.project(coarse);
186  }
187 #endif
188 
189  }// namspace Amg
190  } // namspace Dune
191 #endif