1 #ifndef DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
2 #define DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
6 #include <dune/common/shared_ptr.hh>
10 #include <Eigen/Dense>
17 template<
typename GFS,
typename ET>
21 typedef Eigen::Matrix<ET, Eigen::Dynamic, 1> Container;
22 typedef ET ElementType;
26 typedef ElementType field_type;
28 typedef GFS GridFunctionSpace;
29 typedef std::size_t size_type;
31 typedef typename GFS::Ordering::Traits::ContainerIndex ContainerIndex;
33 typedef ElementType* iterator;
34 typedef const ElementType* const_iterator;
39 #if HAVE_TEMPLATE_ALIASES
41 template<
typename LFSCache>
42 using LocalView = UncachedVectorView<VectorContainer,LFSCache>;
44 template<
typename LFSCache>
45 using ConstLocalView = ConstUncachedVectorView<const VectorContainer,LFSCache>;
49 template<
typename LFSCache>
51 :
public UncachedVectorView<VectorContainer,LFSCache>
57 LocalView(VectorContainer& vc)
58 : UncachedVectorView<VectorContainer,LFSCache>(vc)
63 template<
typename LFSCache>
65 :
public ConstUncachedVectorView<const VectorContainer,LFSCache>
71 ConstLocalView(
const VectorContainer& vc)
72 : ConstUncachedVectorView<const VectorContainer,LFSCache>(vc)
77 #endif // HAVE_TEMPLATE_ALIASES
79 VectorContainer(
const VectorContainer& rhs)
81 , _container(make_shared<Container>(*(rhs._container)))
84 VectorContainer (
const GFS& gfs, tags::attached_container = tags::attached_container())
86 , _container(make_shared<Container>(gfs.ordering().blockCount()))
90 VectorContainer(
const GFS& gfs, tags::unattached_container)
99 VectorContainer (
const GFS& gfs, Container& container)
101 , _container(stackobject_to_shared_ptr(container))
103 _container->resize(gfs.ordering().blockCount());
106 VectorContainer (
const GFS& gfs,
const E&
e)
108 , _container(make_shared<Container>(Container::Constant(gfs.ordering().blockCount(),e)))
116 void attach(shared_ptr<Container> container)
118 _container = container;
121 bool attached()
const
123 return bool(_container);
126 const shared_ptr<Container>& storage()
const
133 return _container->size();
136 VectorContainer& operator=(
const VectorContainer& r)
142 (*_container) = (*r._container);
146 _container = make_shared<Container>(*(r._container));
151 VectorContainer& operator=(
const E&
e)
153 (*_container) = Container::Constant(N(),e);
157 VectorContainer& operator*=(
const E& e)
164 VectorContainer& operator+=(
const E& e)
166 (*_container) += Container::Constant(N(),e);
170 VectorContainer& operator+=(
const VectorContainer& y)
172 (*_container) += (*y._container);
176 VectorContainer& operator-= (
const VectorContainer& y)
178 (*_container) -= (*y._container);
182 E& operator[](
const ContainerIndex& ci)
184 return (*_container)(ci[0]);
187 const E& operator[](
const ContainerIndex& ci)
const
189 return (*_container)(ci[0]);
192 typename Dune::template FieldTraits<E>::real_type two_norm()
const
194 return _container->norm();
197 typename Dune::template FieldTraits<E>::real_type one_norm()
const
199 return _container->template lpNorm<1>();
202 typename Dune::template FieldTraits<E>::real_type infinity_norm()
const
204 return _container->template lpNorm<Eigen::Infinity>();
208 E operator*(
const VectorContainer& y)
const
210 return _container->transpose() * (*y._container);
214 E dot(
const VectorContainer& y)
const
216 return _container->dot(*y._container);
220 VectorContainer& axpy(
const E& a,
const VectorContainer& y)
222 (*_container) += a * (*y._container);
232 const Container& base ()
const
239 return _container->data();
242 const_iterator begin()
const
244 return _container->data();
249 return _container->data() + N();
252 const_iterator end()
const
254 return _container->data() + N();
257 size_t flatsize()
const
259 return _container->size();
262 const GFS& gridFunctionSpace()
const
269 shared_ptr<Container> _container;
278 template<
typename GFS,
typename E>
279 struct EigenVectorSelectorHelper
281 using Type = EIGEN::VectorContainer<GFS, E>;
284 template<
typename GFS,
typename E>
285 struct BackendVectorSelectorHelper<EigenVectorBackend, GFS, E>
286 :
public EigenVectorSelectorHelper<GFS,E>
296 #endif // DUNE_PDELAB_BACKEND_EIGEN_VECTOR_HH
const E & e
Definition: interpolate.hh:172
Various tags for influencing backend behavior.