AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
VectorBackend.hpp
1#pragma once
2
3#include <Eigen/Dense>
4
5#include <dune/common/ftraits.hh>
6
7#include <amdis/Output.hpp>
8#include <amdis/algorithm/ForEach.hpp>
9#include <amdis/algorithm/Transform.hpp>
10#include <amdis/common/FakeContainer.hpp>
11#include <amdis/linearalgebra/VectorBase.hpp>
12#include <amdis/typetree/MultiIndex.hpp>
13
14namespace AMDiS
15{
16 class DefaultIndexDistribution;
17
19 template <class T>
21 : public VectorBase<EigenVector<T>>
22 {
23 public:
25 using BaseVector = Eigen::Matrix<T, Eigen::Dynamic, 1>;
26
28 using value_type = typename BaseVector::Scalar;
29
31 using size_type = typename BaseVector::Index;
32
33 private:
35 using block_type = value_type;
36
38 using field_type = typename Dune::FieldTraits<value_type>::field_type;
39
40 public:
43
45 BaseVector const& vector() const
46 {
47 return vector_;
48 }
49
52 {
53 return vector_;
54 }
55
58 {
59 return vector_.size();
60 }
61
62
64 template <class Basis>
65 void init(Basis const& basis, bool clear)
66 {
67 vector_.resize(basis.dimension());
68 if (clear)
69 vector_.setZero();
70 }
71
72
74 template <class MultiIndex>
76 {
77 const size_type i = flatMultiIndex(idx);
78 test_exit_dbg(i < size(), "Index {} out of range [0,{})", i, size());
79 return vector_.coeffRef(i);
80 }
81
83 template <class MultiIndex>
84 value_type const& at(MultiIndex const& idx) const
85 {
86 const size_type i = flatMultiIndex(idx);
87 test_exit_dbg(i < size(), "Index {} out of range [0,{})", i, size());
88 return vector_.coeff(i);
89 }
90
91 private:
93 BaseVector vector_;
94 };
95
96
97 namespace Recursive
98 {
99 template <class T>
100 struct ForEach<EigenVector<T>> : ForEach<VectorBase<EigenVector<T>>> {};
101
102 template <class T>
103 struct Transform<EigenVector<T>> : Transform<VectorBase<EigenVector<T>>> {};
104
105 template <class T>
106 struct InnerProduct<EigenVector<T>> : InnerProduct<VectorBase<EigenVector<T>>> {};
107
108 } // end namespace Recursive
109} // end namespace AMDiS
Dummy implementation for sequential index "distribution".
Definition: IndexDistribution.hpp:7
The basic container that stores a base vector and a corresponding basis.
Definition: VectorBackend.hpp:22
value_type & at(MultiIndex const &idx)
Access the entry i of the vector with write-access.
Definition: VectorBackend.hpp:75
value_type const & at(MultiIndex const &idx) const
Access the entry i of the vector with read-access.
Definition: VectorBackend.hpp:84
typename BaseVector::Index size_type
The index/size - type.
Definition: VectorBackend.hpp:31
typename BaseVector::Scalar value_type
The type of the elements of the DOFVector.
Definition: VectorBackend.hpp:28
size_type size() const
Return the current size of the vector_.
Definition: VectorBackend.hpp:57
EigenVector(DefaultIndexDistribution const &)
Constructor. Constructs new BaseVector.
Definition: VectorBackend.hpp:42
void init(Basis const &basis, bool clear)
Resize the vector_ to the size s.
Definition: VectorBackend.hpp:65
Eigen::Matrix< T, Eigen::Dynamic, 1 > BaseVector
The type of the base vector.
Definition: VectorBackend.hpp:25
BaseVector & vector()
Return the data-vector vector_.
Definition: VectorBackend.hpp:51
BaseVector const & vector() const
Return the data-vector vector_.
Definition: VectorBackend.hpp:45
CRTP base class for flat vector backends.
Definition: VectorBase.hpp:19
constexpr bool MultiIndex
A multi-index type.
Definition: Concepts.hpp:149
Definition: ForEach.hpp:39
General implementation of recursive inner-product.
Definition: InnerProduct.hpp:58
Definition: Transform.hpp:35