AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
VectorBackend.hpp
1#pragma once
2
3#include <dune/istl/bvector.hh>
4#include <dune/functions/backends/istlvectorbackend.hh>
5
6#include <amdis/Output.hpp>
7#include <amdis/common/FakeContainer.hpp>
8#include <amdis/linearalgebra/VectorBase.hpp>
9#include <amdis/typetree/MultiIndex.hpp>
10
11namespace AMDiS
12{
13 template <class T, class = void>
15 {
16 using type = Dune::FieldVector<T,1>;
17 };
18
19 template <class T>
20 struct BlockVectorType<T, typename T::field_type>
21 {
22 using type = T;
23 };
24
25 template <class T>
27 : public VectorBase<ISTLBlockVector<T>>
28 {
29 public:
31 using BaseVector = Dune::BlockVector<typename BlockVectorType<T>::type>;
32
34 using block_type = typename BaseVector::block_type;
35
37 using value_type = T;
38
40 using field_type = typename block_type::field_type;
41
43 using size_type = typename BaseVector::size_type;
44
45 public:
47 template <class Comm,
48 Dune::disableCopyMove<ISTLBlockVector,Comm> = 0>
49 explicit ISTLBlockVector(Comm const&) {}
50
51 ISTLBlockVector() = default;
52
54 BaseVector const& vector() const
55 {
56 return vector_;
57 }
58
61 {
62 return vector_;
63 }
64
66 std::size_t size() const
67 {
68 return vector_.size();
69 }
70
72 // NOTE: this resize works recursively on the hierarchy of the vector
73 template <class Basis>
74 void init(Basis const& basis, bool clear)
75 {
76 Dune::Functions::istlVectorBackend(vector_).resize(basis);
77 if (clear)
78 vector_ = 0;
79 }
80
82 template <class MultiIndex>
84 {
85 return Dune::Functions::istlVectorBackend(vector_)[idx];
86 }
87
89 template <class MultiIndex>
90 value_type const& at(MultiIndex const& idx) const
91 {
92 return Dune::Functions::istlVectorBackend(vector_)[idx];
93 }
94
95 private:
96 BaseVector vector_;
97 };
98
99
100 namespace Recursive
101 {
102 template <class T>
103 struct ForEach<ISTLBlockVector<T>> : ForEach<VectorBase<ISTLBlockVector<T>>> {};
104
105 template <class T>
106 struct Transform<ISTLBlockVector<T>> : Transform<VectorBase<ISTLBlockVector<T>>> {};
107
108 template <class T>
109 struct InnerProduct<ISTLBlockVector<T>> : InnerProduct<VectorBase<ISTLBlockVector<T>>> {};
110
111 } // end namespace Recursive
112} // end namespace AMDiS
Definition: VectorBackend.hpp:28
value_type & at(MultiIndex const &idx)
Access the entry idx of the vector_ with write-access.
Definition: VectorBackend.hpp:83
ISTLBlockVector(Comm const &)
Constructor. Constructs new BaseVector.
Definition: VectorBackend.hpp:49
value_type const & at(MultiIndex const &idx) const
Access the entry idx of the vector_ with read-access.
Definition: VectorBackend.hpp:90
typename BaseVector::size_type size_type
The index/size - type.
Definition: VectorBackend.hpp:43
typename block_type::field_type field_type
The underlying field type.
Definition: VectorBackend.hpp:40
Dune::BlockVector< typename BlockVectorType< T >::type > BaseVector
The vector type of the underlying base vector.
Definition: VectorBackend.hpp:31
void init(Basis const &basis, bool clear)
Resize the vector_ to the size given by the sizeProvider.
Definition: VectorBackend.hpp:74
BaseVector & vector()
Return the data-vector vector.
Definition: VectorBackend.hpp:60
BaseVector const & vector() const
Return the data-vector vector.
Definition: VectorBackend.hpp:54
std::size_t size() const
Return the current size of the vector_.
Definition: VectorBackend.hpp:66
typename BaseVector::block_type block_type
The type of the elements of the DOFVector.
Definition: VectorBackend.hpp:34
T value_type
The type of the elements of the DOFVector.
Definition: VectorBackend.hpp:37
CRTP base class for flat vector backends.
Definition: VectorBase.hpp:19
constexpr bool MultiIndex
A multi-index type.
Definition: Concepts.hpp:149
Definition: VectorBackend.hpp:15
Definition: ForEach.hpp:39
General implementation of recursive inner-product.
Definition: InnerProduct.hpp:58
Definition: Transform.hpp:35