dune-pdelab  2.0.0
eigen/vector.hh
Go to the documentation of this file.
1 #ifndef DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
2 #define DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
3 
4 #if HAVE_EIGEN
5 
6 #include <dune/common/shared_ptr.hh>
10 #include <Eigen/Dense>
11 
12 namespace Dune {
13  namespace PDELab {
14 
15  namespace EIGEN {
16 
17  template<typename GFS, typename ET>
18  class VectorContainer
19  {
20  public:
21  typedef Eigen::Matrix<ET, Eigen::Dynamic, 1> Container;
22  typedef ET ElementType;
23  typedef ET E;
24 
25  // for ISTL solver compatibility
26  typedef ElementType field_type;
27 
28  typedef GFS GridFunctionSpace;
29  typedef std::size_t size_type;
30 
31  typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
32 
33  typedef ElementType* iterator;
34  typedef const ElementType* const_iterator;
35  // #warning iterators does not work well with Eigen
36  // typedef typename Container::iterator iterator;
37  // typedef typename Container::const_iterator const_iterator;
38 
39 #if HAVE_TEMPLATE_ALIASES
40 
41  template<typename LFSCache>
42  using LocalView = UncachedVectorView<VectorContainer,LFSCache>;
43 
44  template<typename LFSCache>
45  using ConstLocalView = ConstUncachedVectorView<const VectorContainer,LFSCache>;
46 
47 #else
48 
49  template<typename LFSCache>
50  struct LocalView
51  : public UncachedVectorView<VectorContainer,LFSCache>
52  {
53 
54  LocalView()
55  {}
56 
57  LocalView(VectorContainer& vc)
58  : UncachedVectorView<VectorContainer,LFSCache>(vc)
59  {}
60 
61  };
62 
63  template<typename LFSCache>
64  struct ConstLocalView
65  : public ConstUncachedVectorView<const VectorContainer,LFSCache>
66  {
67 
68  ConstLocalView()
69  {}
70 
71  ConstLocalView(const VectorContainer& vc)
72  : ConstUncachedVectorView<const VectorContainer,LFSCache>(vc)
73  {}
74 
75  };
76 
77 #endif // HAVE_TEMPLATE_ALIASES
78 
79  VectorContainer(const VectorContainer& rhs)
80  : _gfs(rhs._gfs)
81  , _container(make_shared<Container>(*(rhs._container)))
82  {}
83 
84  VectorContainer (const GFS& gfs, tags::attached_container = tags::attached_container())
85  : _gfs(gfs)
86  , _container(make_shared<Container>(gfs.ordering().blockCount()))
87  {}
88 
90  VectorContainer(const GFS& gfs, tags::unattached_container)
91  : _gfs(gfs)
92  {}
93 
99  VectorContainer (const GFS& gfs, Container& container)
100  : _gfs(gfs)
101  , _container(stackobject_to_shared_ptr(container))
102  {
103  _container->resize(gfs.ordering().blockCount());
104  }
105 
106  VectorContainer (const GFS& gfs, const E& e)
107  : _gfs(gfs)
108  , _container(make_shared<Container>(Container::Constant(gfs.ordering().blockCount(),e)))
109  {}
110 
111  void detach()
112  {
113  _container.reset();
114  }
115 
116  void attach(shared_ptr<Container> container)
117  {
118  _container = container;
119  }
120 
121  bool attached() const
122  {
123  return bool(_container);
124  }
125 
126  const shared_ptr<Container>& storage() const
127  {
128  return _container;
129  }
130 
131  size_type N() const
132  {
133  return _container->size();
134  }
135 
136  VectorContainer& operator=(const VectorContainer& r)
137  {
138  if (this == &r)
139  return *this;
140  if (attached())
141  {
142  (*_container) = (*r._container);
143  }
144  else
145  {
146  _container = make_shared<Container>(*(r._container));
147  }
148  return *this;
149  }
150 
151  VectorContainer& operator=(const E& e)
152  {
153  (*_container) = Container::Constant(N(),e);
154  return *this;
155  }
156 
157  VectorContainer& operator*=(const E& e)
158  {
159  (*_container) *= e;
160  return *this;
161  }
162 
163 
164  VectorContainer& operator+=(const E& e)
165  {
166  (*_container) += Container::Constant(N(),e);
167  return *this;
168  }
169 
170  VectorContainer& operator+=(const VectorContainer& y)
171  {
172  (*_container) += (*y._container);
173  return *this;
174  }
175 
176  VectorContainer& operator-= (const VectorContainer& y)
177  {
178  (*_container) -= (*y._container);
179  return *this;
180  }
181 
182  E& operator[](const ContainerIndex& ci)
183  {
184  return (*_container)(ci[0]);
185  }
186 
187  const E& operator[](const ContainerIndex& ci) const
188  {
189  return (*_container)(ci[0]);
190  }
191 
192  typename Dune::template FieldTraits<E>::real_type two_norm() const
193  {
194  return _container->norm();
195  }
196 
197  typename Dune::template FieldTraits<E>::real_type one_norm() const
198  {
199  return _container->template lpNorm<1>();
200  }
201 
202  typename Dune::template FieldTraits<E>::real_type infinity_norm() const
203  {
204  return _container->template lpNorm<Eigen::Infinity>();
205  }
206 
208  E operator*(const VectorContainer& y) const
209  {
210  return _container->transpose() * (*y._container);
211  }
212 
214  E dot(const VectorContainer& y) const
215  {
216  return _container->dot(*y._container);
217  }
218 
220  VectorContainer& axpy(const E& a, const VectorContainer& y)
221  {
222  (*_container) += a * (*y._container);
223  return *this;
224  }
225 
226  // for debugging and AMG access
227  Container& base ()
228  {
229  return *_container;
230  }
231 
232  const Container& base () const
233  {
234  return *_container;
235  }
236 
237  iterator begin()
238  {
239  return _container->data();
240  }
241 
242  const_iterator begin() const
243  {
244  return _container->data();
245  }
246 
247  iterator end()
248  {
249  return _container->data() + N();
250  }
251 
252  const_iterator end() const
253  {
254  return _container->data() + N();
255  }
256 
257  size_t flatsize() const
258  {
259  return _container->size();
260  }
261 
262  const GFS& gridFunctionSpace() const
263  {
264  return _gfs;
265  }
266 
267  private:
268  const GFS& _gfs;
269  shared_ptr<Container> _container;
270 
271  };
272 
273  } // end namespace EIGEN
274 
275 
276 #ifndef DOXYGEN
277 
278  template<typename GFS, typename E>
279  struct EigenVectorSelectorHelper
280  {
281  using Type = EIGEN::VectorContainer<GFS, E>;
282  };
283 
284  template<typename GFS, typename E>
285  struct BackendVectorSelectorHelper<EigenVectorBackend, GFS, E>
286  : public EigenVectorSelectorHelper<GFS,E>
287  {};
288 
289 #endif // DOXYGEN
290 
291  } // namespace PDELab
292 } // namespace Dune
293 
294 #endif
295 
296 #endif // DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
297 
298 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
299 // vi: set et ts=4 sw=2 sts=2:
const E & e
Definition: interpolate.hh:172
Various tags for influencing backend behavior.