dune-pdelab  2.5-dev
simple/vector.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_SIMPLE_VECTOR_HH
4 #define DUNE_PDELAB_BACKEND_SIMPLE_VECTOR_HH
5 
6 #include <algorithm>
7 #include <functional>
8 #include <memory>
9 #include <numeric>
10 
11 #include <dune/common/fvector.hh>
12 #include <dune/istl/bvector.hh>
13 
20 
21 namespace Dune {
22  namespace PDELab {
23 
24  namespace Simple {
25 
26  namespace {
27 
28  // For some reason std::bind cannot directly deduce the correct version
29  // of Dune::fvmeta::abs2, so we package it into a functor to help it along.
30  template<typename K>
31  struct abs2
32  {
33  typename FieldTraits<K>::real_type operator()(const K& k) const
34  {
35  return Dune::fvmeta::abs2(k);
36  }
37  };
38 
39  }
40 
41  template<typename GFS, typename C>
43  : public Backend::impl::Wrapper<C>
44  {
45 
46  friend Backend::impl::Wrapper<C>;
47 
48  public:
49  typedef C Container;
50  typedef typename Container::value_type ElementType;
51  typedef ElementType E;
52 
53  // for ISTL solver compatibility
55 
56  typedef GFS GridFunctionSpace;
57  typedef typename Container::size_type size_type;
58 
59  typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
60 
61  typedef typename Container::iterator iterator;
62  typedef typename Container::const_iterator const_iterator;
63 
64  template<typename LFSCache>
66 
67  template<typename LFSCache>
69 
70 
72  : _gfs(rhs._gfs)
73  , _container(std::make_shared<Container>(*(rhs._container)))
74  {}
75 
77  : _gfs(gfs)
78  , _container(std::make_shared<Container>(gfs.ordering().blockCount()))
79  {}
80 
83  : _gfs(gfs)
84  {}
85 
91  VectorContainer (const GFS& gfs, Container& container)
92  : _gfs(gfs)
93  , _container(stackobject_to_shared_ptr(container))
94  {
95  _container->resize(gfs.ordering().blockCount());
96  }
97 
98  VectorContainer (const GFS& gfs, const E& e)
99  : _gfs(gfs)
100  , _container(std::make_shared<Container>(gfs.ordering().blockCount(),e))
101  {}
102 
103  void detach()
104  {
105  _container.reset();
106  }
107 
108  void attach(std::shared_ptr<Container> container)
109  {
110  _container = container;
111  }
112 
113  bool attached() const
114  {
115  return bool(_container);
116  }
117 
118  const std::shared_ptr<Container>& storage() const
119  {
120  return _container;
121  }
122 
123  size_type N() const
124  {
125  return _container->size();
126  }
127 
129  {
130  if (this == &r)
131  return *this;
132  if (attached())
133  {
134  (*_container) = (*r._container);
135  }
136  else
137  {
138  _container = std::make_shared<Container>(*(r._container));
139  }
140  return *this;
141  }
142 
144  {
145  std::fill(_container->begin(),_container->end(),e);
146  return *this;
147  }
148 
150  {
151  std::transform(_container->begin(),_container->end(),_container->begin(),
152  std::bind(std::multiplies<E>(),e,std::placeholders::_1));
153  return *this;
154  }
155 
156 
158  {
159  std::transform(_container->begin(),_container->end(),_container->begin(),
160  std::bind(std::plus<E>(),e,std::placeholders::_1));
161  return *this;
162  }
163 
165  {
166  std::transform(_container->begin(),_container->end(),y._container->begin(),
167  _container->begin(),std::plus<E>());
168  return *this;
169  }
170 
172  {
173  std::transform(_container->begin(),_container->end(),y._container->begin(),
174  _container->begin(),std::minus<E>());
175  return *this;
176  }
177 
179  {
180  return (*_container)[ci[0]];
181  }
182 
183  const E& operator[](const ContainerIndex& ci) const
184  {
185  return (*_container)[ci[0]];
186  }
187 
188  typename Dune::template FieldTraits<E>::real_type two_norm() const
189  {
190  using namespace std::placeholders;
191  typedef typename Dune::template FieldTraits<E>::real_type Real;
192  return std::sqrt(std::accumulate(_container->begin(),_container->end(),Real(0),std::bind(std::plus<Real>(),_1,std::bind(abs2<E>(),_2))));
193  }
194 
195  typename Dune::template FieldTraits<E>::real_type one_norm() const
196  {
197  typedef typename Dune::template FieldTraits<E>::real_type Real;
198  return std::accumulate(_container->begin(),_container->end(),Real(0),
199  [](const auto& n, const auto& e) {
200  using std::abs;
201  return n + abs(e);
202  });
203  }
204 
205  typename Dune::template FieldTraits<E>::real_type infinity_norm() const
206  {
207  if (_container->size() == 0)
208  return 0;
209  using std::abs;
210  typedef typename Dune::template FieldTraits<E>::real_type Real;
211  return abs(*std::max_element(_container->begin(),_container->end(),
212  [](const auto& a, const auto& b) {
213  using std::abs;
214  return abs(a) < abs(b);
215  }));
216  }
217 
218  E operator*(const VectorContainer& y) const
219  {
220  return std::inner_product(_container->begin(),_container->end(),y._container->begin(),E(0));
221  }
222 
223  E dot(const VectorContainer& y) const
224  {
225  return std::inner_product(_container->begin(),_container->end(),y._container->begin(),E(0),std::plus<E>(),Dune::dot<E,E>);
226  }
227 
228  VectorContainer& axpy(const E& a, const VectorContainer& y)
229  {
230  using namespace std::placeholders;
231  std::transform(_container->begin(),_container->end(),y._container->begin(),
232  _container->begin(),std::bind(std::plus<E>(),_1,std::bind(std::multiplies<E>(),a,_2)));
233  return *this;
234  }
235 
236  // for debugging and AMG access
238  {
239  return *_container;
240  }
241 
242  const Container& base () const
243  {
244  return *_container;
245  }
246 
247  private:
248 
249  Container& native ()
250  {
251  return *_container;
252  }
253 
254  const Container& native () const
255  {
256  return *_container;
257  }
258 
259  public:
260 
262  {
263  return _container->begin();
264  }
265 
267  {
268  return _container->begin();
269  }
270 
272  {
273  return _container->end();
274  }
275 
277  {
278  return _container->end();
279  }
280 
281  size_t flatsize() const
282  {
283  return _container->size();
284  }
285 
286  const GFS& gridFunctionSpace() const
287  {
288  return _gfs;
289  }
290 
291  private:
292  const GFS& _gfs;
293  std::shared_ptr<Container> _container;
294  };
295 
296 
297  } // namespace Simple
298 
299 
300 #ifndef DOXYGEN
301 
302  template<typename GFS, typename E>
303  struct SimpleVectorSelectorHelper
304  {
305 
306  using vector_type = typename GFS::Traits::Backend::template vector_type<E>;
307 
308  using Type = Simple::VectorContainer<GFS,vector_type>;
309 
310  };
311 
312  namespace Backend {
313  namespace impl {
314 
315  template<template<typename> class Container, typename GFS, typename E>
316  struct BackendVectorSelectorHelper<Simple::VectorBackend<Container>, GFS, E>
317  : public SimpleVectorSelectorHelper<GFS,E>
318  {};
319 
320  } // namespace impl
321  } // namespace Backend
322 
323 #endif // DOXYGEN
324 
325  } // namespace PDELab
326 } // namespace Dune
327 
328 #endif // DUNE_PDELAB_BACKEND_SIMPLE_VECTOR_HH
Dune::PDELab::Simple::VectorContainer::Container
C Container
Definition: simple/vector.hh:49
Dune::PDELab::Simple::VectorContainer::detach
void detach()
Definition: simple/vector.hh:103
Dune::PDELab::Simple::VectorContainer
Definition: simple/vector.hh:42
Dune::PDELab::Simple::VectorContainer::attached
bool attached() const
Definition: simple/vector.hh:113
Dune::PDELab::Simple::VectorContainer::operator*=
VectorContainer & operator*=(const E &e)
Definition: simple/vector.hh:149
Dune::PDELab::Simple::VectorContainer::VectorContainer
VectorContainer(const GFS &gfs, const E &e)
Definition: simple/vector.hh:98
Dune::PDELab::Simple::VectorContainer::field_type
ElementType field_type
Definition: simple/vector.hh:54
Dune::PDELab::Simple::VectorContainer::begin
iterator begin()
Definition: simple/vector.hh:261
lfsindexcache.hh
Dune::PDELab::Simple::VectorContainer::one_norm
Dune::template FieldTraits< E >::real_type one_norm() const
Definition: simple/vector.hh:195
Dune::PDELab::Simple::VectorContainer::N
size_type N() const
Definition: simple/vector.hh:123
gridfunctionspace.hh
Dune::PDELab::Backend::attached_container
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/common/tags.hh:27
Dune::PDELab::Simple::VectorContainer::begin
const_iterator begin() const
Definition: simple/vector.hh:266
Dune::PDELab::Simple::VectorContainer::base
Container & base()
Definition: simple/vector.hh:237
Dune::PDELab::Simple::VectorContainer::ContainerIndex
GFS::Ordering::Traits::ContainerIndex ContainerIndex
Definition: simple/vector.hh:59
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::Simple::VectorContainer::operator=
VectorContainer & operator=(const VectorContainer &r)
Definition: simple/vector.hh:128
Dune::PDELab::Simple::VectorContainer::operator+=
VectorContainer & operator+=(const VectorContainer &y)
Definition: simple/vector.hh:164
Dune::PDELab::Simple::VectorContainer::end
iterator end()
Definition: simple/vector.hh:271
Dune::PDELab::Simple::VectorContainer::size_type
Container::size_type size_type
Definition: simple/vector.hh:57
Dune::PDELab::Simple::VectorContainer::operator[]
const E & operator[](const ContainerIndex &ci) const
Definition: simple/vector.hh:183
Dune::PDELab::Simple::VectorContainer::VectorContainer
VectorContainer(const GFS &gfs, Backend::attached_container=Backend::attached_container())
Definition: simple/vector.hh:76
Dune::PDELab::Simple::VectorContainer::infinity_norm
Dune::template FieldTraits< E >::real_type infinity_norm() const
Definition: simple/vector.hh:205
Dune::PDELab::Simple::VectorContainer::E
ElementType E
Definition: simple/vector.hh:51
Dune::PDELab::UncachedVectorView
Definition: uncachedvectorview.hh:134
Dune::PDELab::Simple::VectorContainer::const_iterator
Container::const_iterator const_iterator
Definition: simple/vector.hh:62
Dune::PDELab::Simple::VectorContainer::operator+=
VectorContainer & operator+=(const E &e)
Definition: simple/vector.hh:157
Dune::PDELab::Simple::VectorContainer::gridFunctionSpace
const GFS & gridFunctionSpace() const
Definition: simple/vector.hh:286
Dune::PDELab::Simple::VectorContainer::two_norm
Dune::template FieldTraits< E >::real_type two_norm() const
Definition: simple/vector.hh:188
Dune::PDELab::Simple::VectorContainer::VectorContainer
VectorContainer(const GFS &gfs, Container &container)
Constructs an VectorContainer for an explicitly given vector object.
Definition: simple/vector.hh:91
Dune::PDELab::Simple::VectorContainer::axpy
VectorContainer & axpy(const E &a, const VectorContainer &y)
Definition: simple/vector.hh:228
Dune::PDELab::Simple::VectorContainer::VectorContainer
VectorContainer(const VectorContainer &rhs)
Definition: simple/vector.hh:71
descriptors.hh
Dune::PDELab::Simple::VectorContainer::ElementType
Container::value_type ElementType
Definition: simple/vector.hh:50
Dune::PDELab::Simple::VectorContainer::operator-=
VectorContainer & operator-=(const VectorContainer &y)
Definition: simple/vector.hh:171
e
const Entity & e
Definition: localfunctionspace.hh:120
tags.hh
Various tags for influencing backend behavior.
Dune::PDELab::Simple::VectorContainer::VectorContainer
VectorContainer(const GFS &gfs, Backend::unattached_container)
Creates a VectorContainer without allocating storage.
Definition: simple/vector.hh:82
Dune::PDELab::Simple::VectorContainer::flatsize
size_t flatsize() const
Definition: simple/vector.hh:281
Dune::PDELab::Simple::VectorContainer::operator=
VectorContainer & operator=(const E &e)
Definition: simple/vector.hh:143
Dune::PDELab::Simple::VectorContainer::operator*
E operator*(const VectorContainer &y) const
Definition: simple/vector.hh:218
Dune::PDELab::Simple::VectorContainer::attach
void attach(std::shared_ptr< Container > container)
Definition: simple/vector.hh:108
Dune::PDELab::Simple::VectorContainer::GridFunctionSpace
GFS GridFunctionSpace
Definition: simple/vector.hh:56
Dune::PDELab::ConstUncachedVectorView
Definition: uncachedvectorview.hh:14
Dune::PDELab::Simple::VectorContainer::operator[]
E & operator[](const ContainerIndex &ci)
Definition: simple/vector.hh:178
Dune::PDELab::Backend::unattached_container
Tag for requesting a vector or matrix container without a pre-attached underlying object.
Definition: backend/common/tags.hh:23
Dune::PDELab::Simple::VectorContainer::iterator
Container::iterator iterator
Definition: simple/vector.hh:61
interface.hh
uncachedvectorview.hh
Dune::PDELab::Simple::VectorContainer::end
const_iterator end() const
Definition: simple/vector.hh:276
Dune::PDELab::Simple::VectorContainer::dot
E dot(const VectorContainer &y) const
Definition: simple/vector.hh:223
Dune::PDELab::Simple::VectorContainer::storage
const std::shared_ptr< Container > & storage() const
Definition: simple/vector.hh:118
Dune::PDELab::Simple::VectorContainer::base
const Container & base() const
Definition: simple/vector.hh:242