3 #ifndef DUNE_PDELAB_BACKEND_ISTLVECTORBACKEND_HH
4 #define DUNE_PDELAB_BACKEND_ISTLVECTORBACKEND_HH
6 #include <dune/common/fvector.hh>
7 #include <dune/istl/bvector.hh>
8 #include <dune/typetree/typetree.hh>
23 template<
typename GFS,
typename C>
43 #if HAVE_TEMPLATE_ALIASES
45 template<
typename LFSCache>
48 template<
typename LFSCache>
53 template<
typename LFSCache>
67 template<
typename LFSCache>
81 #endif // HAVE_TEMPLATE_ALIASES
86 , _container(make_shared<
Container>(_gfs.ordering().blockCount()))
88 istl::dispatch_vector_allocation(_gfs.ordering(),*_container,
typename GFS::Ordering::ContainerAllocationTag());
89 (*_container) = rhs.
base();
94 , _container(make_shared<
Container>(gfs.ordering().blockCount()))
96 istl::dispatch_vector_allocation(gfs.ordering(),*_container,
typename GFS::Ordering::ContainerAllocationTag());
111 , _container(stackobject_to_shared_ptr(container))
113 _container->resize(gfs.ordering().blockCount());
114 istl::dispatch_vector_allocation(gfs.ordering(),*_container,
typename GFS::Ordering::ContainerAllocationTag());
119 , _container(make_shared<
Container>(gfs.ordering().blockCount()))
121 istl::dispatch_vector_allocation(gfs.ordering(),*_container,
typename GFS::Ordering::ContainerAllocationTag());
130 void attach(shared_ptr<Container> container)
132 _container = container;
137 return bool(_container);
147 return _container->N();
156 (*_container) = r.
base();
160 _container = make_shared<Container>(r.
base());
186 (*_container)+= e.
base();
192 (*_container)-= e.
base();
198 return (*_container)[i];
203 return (*_container)[i];
208 return istl::access_vector_element(
istl::container_tag(*_container),*_container,ci,ci.size()-1);
213 return istl::access_vector_element(
istl::container_tag(*_container),*_container,ci,ci.size()-1);
216 typename Dune::template FieldTraits<E>::real_type
two_norm()
const
218 return _container->two_norm();
221 typename Dune::template FieldTraits<E>::real_type
one_norm()
const
223 return _container->one_norm();
228 return _container->infinity_norm();
233 return (*_container)*y.
base();
238 return _container->dot(y.
base());
243 _container->axpy(a, y.
base());
292 return _container->dim();
302 shared_ptr<Container> _container;
312 template<
typename GFS,
typename E>
313 struct ISTLVectorSelectorHelper
316 typedef typename TypeTree::AccumulateType<
318 istl::vector_creation_policy<E>
319 >::type vector_descriptor;
321 typedef ISTLBlockVectorContainer<GFS,typename vector_descriptor::vector_type> Type;
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>
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
LFSCache LFSCache
Definition: uncachedvectorview.hh:19
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
Definition: istlvectorbackend.hh:68
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