dune-pdelab  2.0.0
istlvectorbackend.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_BACKEND_ISTLVECTORBACKEND_HH
4 #define DUNE_PDELAB_BACKEND_ISTLVECTORBACKEND_HH
5 
6 #include <dune/common/fvector.hh>
7 #include <dune/istl/bvector.hh>
8 #include <dune/typetree/typetree.hh>
9 
19 
20 namespace Dune {
21  namespace PDELab {
22 
23  template<typename GFS, typename C>
25  {
26 
27  public:
28  typedef typename C::field_type ElementType;
29  typedef ElementType E;
30  typedef C Container;
31  typedef GFS GridFunctionSpace;
32  typedef Container BaseT;
33  typedef typename Container::field_type field_type;
34  typedef typename Container::block_type block_type;
35  typedef typename Container::size_type size_type;
36 
37  typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
38 
41 
42 
43 #if HAVE_TEMPLATE_ALIASES
44 
45  template<typename LFSCache>
47 
48  template<typename LFSCache>
50 
51 #else
52 
53  template<typename LFSCache>
54  struct LocalView
55  : public UncachedVectorView<ISTLBlockVectorContainer,LFSCache>
56  {
57 
59  {}
60 
63  {}
64 
65  };
66 
67  template<typename LFSCache>
69  : public ConstUncachedVectorView<const ISTLBlockVectorContainer,LFSCache>
70  {
71 
73  {}
74 
77  {}
78 
79  };
80 
81 #endif // HAVE_TEMPLATE_ALIASES
82 
83 
85  : _gfs(rhs._gfs)
86  , _container(make_shared<Container>(_gfs.ordering().blockCount()))
87  {
88  istl::dispatch_vector_allocation(_gfs.ordering(),*_container,typename GFS::Ordering::ContainerAllocationTag());
89  (*_container) = rhs.base();
90  }
91 
93  : _gfs(gfs)
94  , _container(make_shared<Container>(gfs.ordering().blockCount()))
95  {
96  istl::dispatch_vector_allocation(gfs.ordering(),*_container,typename GFS::Ordering::ContainerAllocationTag());
97  }
98 
101  : _gfs(gfs)
102  {}
103 
109  ISTLBlockVectorContainer (const GFS& gfs, Container& container)
110  : _gfs(gfs)
111  , _container(stackobject_to_shared_ptr(container))
112  {
113  _container->resize(gfs.ordering().blockCount());
114  istl::dispatch_vector_allocation(gfs.ordering(),*_container,typename GFS::Ordering::ContainerAllocationTag());
115  }
116 
117  ISTLBlockVectorContainer (const GFS& gfs, const E& e)
118  : _gfs(gfs)
119  , _container(make_shared<Container>(gfs.ordering().blockCount()))
120  {
121  istl::dispatch_vector_allocation(gfs.ordering(),*_container,typename GFS::Ordering::ContainerAllocationTag());
122  (*_container)=e;
123  }
124 
125  void detach()
126  {
127  _container.reset();
128  }
129 
130  void attach(shared_ptr<Container> container)
131  {
132  _container = container;
133  }
134 
135  bool attached() const
136  {
137  return bool(_container);
138  }
139 
140  const shared_ptr<Container>& storage() const
141  {
142  return _container;
143  }
144 
145  size_type N() const
146  {
147  return _container->N();
148  }
149 
151  {
152  if (this == &r)
153  return *this;
154  if (attached())
155  {
156  (*_container) = r.base();
157  }
158  else
159  {
160  _container = make_shared<Container>(r.base());
161  }
162  return *this;
163  }
164 
166  {
167  (*_container)=e;
168  return *this;
169  }
170 
172  {
173  (*_container)*=e;
174  return *this;
175  }
176 
177 
179  {
180  (*_container)+=e;
181  return *this;
182  }
183 
185  {
186  (*_container)+= e.base();
187  return *this;
188  }
189 
191  {
192  (*_container)-= e.base();
193  return *this;
194  }
195 
196  block_type& block(std::size_t i)
197  {
198  return (*_container)[i];
199  }
200 
201  const block_type& block(std::size_t i) const
202  {
203  return (*_container)[i];
204  }
205 
207  {
208  return istl::access_vector_element(istl::container_tag(*_container),*_container,ci,ci.size()-1);
209  }
210 
211  const E& operator[](const ContainerIndex& ci) const
212  {
213  return istl::access_vector_element(istl::container_tag(*_container),*_container,ci,ci.size()-1);
214  }
215 
216  typename Dune::template FieldTraits<E>::real_type two_norm() const
217  {
218  return _container->two_norm();
219  }
220 
221  typename Dune::template FieldTraits<E>::real_type one_norm() const
222  {
223  return _container->one_norm();
224  }
225 
226  typename Dune::template FieldTraits<E>::real_type infinity_norm() const
227  {
228  return _container->infinity_norm();
229  }
230 
232  {
233  return (*_container)*y.base();
234  }
235 
237  {
238  return _container->dot(y.base());
239  }
240 
242  {
243  _container->axpy(a, y.base());
244  return *this;
245  }
246 
247  // for debugging and AMG access
249  {
250  return *_container;
251  }
252 
253  const Container& base () const
254  {
255  return *_container;
256  }
257 
258  operator Container&()
259  {
260  return *_container;
261  }
262 
263  operator const Container&() const
264  {
265  return *_container;
266  }
267 
269  {
270  return iterator(*_container,false);
271  }
272 
273 
275  {
276  return const_iterator(*_container,false);
277  }
278 
280  {
281  return iterator(*_container,true);
282  }
283 
284 
286  {
287  return const_iterator(*_container,true);
288  }
289 
290  size_t flatsize() const
291  {
292  return _container->dim();
293  }
294 
295  const GFS& gridFunctionSpace() const
296  {
297  return _gfs;
298  }
299 
300  private:
301  const GFS& _gfs;
302  shared_ptr<Container> _container;
303  };
304 
305 
306 
307 
308 
309 #ifndef DOXYGEN
310 
311  // helper struct invoking the GFS tree -> ISTL vector reduction
312  template<typename GFS, typename E>
313  struct ISTLVectorSelectorHelper
314  {
315 
316  typedef typename TypeTree::AccumulateType<
317  GFS,
318  istl::vector_creation_policy<E>
319  >::type vector_descriptor;
320 
321  typedef ISTLBlockVectorContainer<GFS,typename vector_descriptor::vector_type> Type;
322 
323  };
324 
325  template<ISTLParameters::Blocking blocking, std::size_t block_size, typename GFS, typename E>
326  struct BackendVectorSelectorHelper<ISTLVectorBackend<blocking,block_size>, GFS, E>
327  : public ISTLVectorSelectorHelper<GFS,E>
328  {};
329 
330 #endif // DOXYGEN
331 
332  } // namespace PDELab
333 } // namespace Dune
334 
335 #endif // DUNE_PDELAB_BACKEND_ISTLVECTORBACKEND_HH
E operator*(const ISTLBlockVectorContainer &y) const
Definition: istlvectorbackend.hh:231
const block_type & block(std::size_t i) const
Definition: istlvectorbackend.hh:201
const GFS & gridFunctionSpace() const
Definition: istlvectorbackend.hh:295
void detach()
Definition: istlvectorbackend.hh:125
C Container
Definition: istlvectorbackend.hh:30
size_t flatsize() const
Definition: istlvectorbackend.hh:290
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/tags.hh:27
Container & base()
Definition: istlvectorbackend.hh:248
const Container & base() const
Definition: istlvectorbackend.hh:253
Container BaseT
Definition: istlvectorbackend.hh:32
ISTLBlockVectorContainer & operator*=(const E &e)
Definition: istlvectorbackend.hh:171
iterator begin()
Definition: istlvectorbackend.hh:268
ISTLBlockVectorContainer & operator=(const ISTLBlockVectorContainer &r)
Definition: istlvectorbackend.hh:150
Container::size_type size_type
Definition: istlvectorbackend.hh:35
ConstLocalView(const ISTLBlockVectorContainer &vc)
Definition: istlvectorbackend.hh:75
E dot(const ISTLBlockVectorContainer &y) const
Definition: istlvectorbackend.hh:236
Dune::template FieldTraits< E >::real_type two_norm() const
Definition: istlvectorbackend.hh:216
GFS::Ordering::Traits::ContainerIndex ContainerIndex
Definition: istlvectorbackend.hh:37
Dune::template FieldTraits< E >::real_type infinity_norm() const
Definition: istlvectorbackend.hh:226
istl::vector_iterator< const C > const_iterator
Definition: istlvectorbackend.hh:40
Tag for requesting a vector or matrix container without a pre-attached underlying object...
Definition: backend/tags.hh:23
ConstLocalView()
Definition: istlvectorbackend.hh:72
E & operator[](const ContainerIndex &ci)
Definition: istlvectorbackend.hh:206
void attach(shared_ptr< Container > container)
Definition: istlvectorbackend.hh:130
ElementType E
Definition: istlvectorbackend.hh:29
const E & operator[](const ContainerIndex &ci) const
Definition: istlvectorbackend.hh:211
LocalView()
Definition: istlvectorbackend.hh:58
Definition: istlvectorbackend.hh:24
block_type & block(std::size_t i)
Definition: istlvectorbackend.hh:196
ISTLBlockVectorContainer & axpy(const E &a, const ISTLBlockVectorContainer &y)
Definition: istlvectorbackend.hh:241
Container::field_type field_type
Definition: istlvectorbackend.hh:33
const shared_ptr< Container > & storage() const
Definition: istlvectorbackend.hh:140
istl::vector_iterator< C > iterator
Definition: istlvectorbackend.hh:39
ISTLBlockVectorContainer(const GFS &gfs, const E &e)
Definition: istlvectorbackend.hh:117
ISTLBlockVectorContainer & operator+=(const E &e)
Definition: istlvectorbackend.hh:178
Dune::template FieldTraits< E >::real_type one_norm() const
Definition: istlvectorbackend.hh:221
iterator end()
Definition: istlvectorbackend.hh:279
size_type N() const
Definition: istlvectorbackend.hh:145
GFS GridFunctionSpace
Definition: istlvectorbackend.hh:31
Container::block_type block_type
Definition: istlvectorbackend.hh:34
ISTLBlockVectorContainer(const GFS &gfs, tags::attached_container=tags::attached_container())
Definition: istlvectorbackend.hh:92
bool attached() const
Definition: istlvectorbackend.hh:135
tags::container< T >::type container_tag(const T &)
Gets instance of container tag associated with T.
Definition: backend/istl/tags.hh:247
ISTLBlockVectorContainer(const GFS &gfs, tags::unattached_container)
Creates an ISTLBlockVectorContainer without allocating an underlying ISTL vector. ...
Definition: istlvectorbackend.hh:100
LocalView(ISTLBlockVectorContainer &vc)
Definition: istlvectorbackend.hh:61
Definition: vectoriterator.hh:182
const_iterator begin() const
Definition: istlvectorbackend.hh:274
C::field_type ElementType
Definition: istlvectorbackend.hh:28
const E & e
Definition: interpolate.hh:172
ISTLBlockVectorContainer & operator-=(const ISTLBlockVectorContainer &e)
Definition: istlvectorbackend.hh:190
ISTLBlockVectorContainer(const GFS &gfs, Container &container)
Constructs an ISTLBlockVectorContainer for an explicitly given vector object.
Definition: istlvectorbackend.hh:109
ISTLBlockVectorContainer(const ISTLBlockVectorContainer &rhs)
Definition: istlvectorbackend.hh:84
Definition: istlvectorbackend.hh:54
Various tags for influencing backend behavior.
const_iterator end() const
Definition: istlvectorbackend.hh:285