3 #ifndef DUNE_PDELAB_BACKEND_SIMPLE_VECTOR_HH
4 #define DUNE_PDELAB_BACKEND_SIMPLE_VECTOR_HH
11 #include <dune/common/fvector.hh>
12 #include <dune/istl/bvector.hh>
33 typename FieldTraits<K>::real_type operator()(
const K& k)
const
35 return Dune::fvmeta::abs2(k);
41 template<
typename GFS,
typename C>
43 :
public Backend::impl::Wrapper<C>
46 friend Backend::impl::Wrapper<C>;
64 template<
typename LFSCache>
67 template<
typename LFSCache>
73 , _container(std::make_shared<
Container>(*(rhs._container)))
78 , _container(std::make_shared<
Container>(gfs.ordering().blockCount()))
93 , _container(stackobject_to_shared_ptr(container))
95 _container->resize(gfs.ordering().blockCount());
100 , _container(std::make_shared<
Container>(gfs.ordering().blockCount(),
e))
108 void attach(std::shared_ptr<Container> container)
110 _container = container;
115 return bool(_container);
118 const std::shared_ptr<Container>&
storage()
const
125 return _container->size();
134 (*_container) = (*r._container);
138 _container = std::make_shared<Container>(*(r._container));
145 std::fill(_container->begin(),_container->end(),
e);
151 std::transform(_container->begin(),_container->end(),_container->begin(),
152 std::bind(std::multiplies<E>(),
e,std::placeholders::_1));
159 std::transform(_container->begin(),_container->end(),_container->begin(),
160 std::bind(std::plus<E>(),
e,std::placeholders::_1));
166 std::transform(_container->begin(),_container->end(),y._container->begin(),
167 _container->begin(),std::plus<E>());
173 std::transform(_container->begin(),_container->end(),y._container->begin(),
174 _container->begin(),std::minus<E>());
180 return (*_container)[ci[0]];
185 return (*_container)[ci[0]];
188 typename Dune::template FieldTraits<E>::real_type
two_norm()
const
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))));
195 typename Dune::template FieldTraits<E>::real_type
one_norm()
const
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) {
207 if (_container->size() == 0)
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) {
214 return abs(a) < abs(b);
220 return std::inner_product(_container->begin(),_container->end(),y._container->begin(),
E(0));
225 return std::inner_product(_container->begin(),_container->end(),y._container->begin(),
E(0),std::plus<E>(),Dune::dot<E,E>);
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)));
263 return _container->begin();
268 return _container->begin();
273 return _container->end();
278 return _container->end();
283 return _container->size();
293 std::shared_ptr<Container> _container;
302 template<
typename GFS,
typename E>
303 struct SimpleVectorSelectorHelper
306 using vector_type =
typename GFS::Traits::Backend::template vector_type<E>;
308 using Type = Simple::VectorContainer<GFS,vector_type>;
315 template<
template<
typename>
class Container,
typename GFS,
typename E>
316 struct BackendVectorSelectorHelper<Simple::VectorBackend<Container>, GFS,
E>
317 :
public SimpleVectorSelectorHelper<GFS,E>
328 #endif // DUNE_PDELAB_BACKEND_SIMPLE_VECTOR_HH