dune-functions  2.7.1
flatvectorview.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH
4 #define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH
5 
6 
7 #include <dune/common/concept.hh>
8 
10 
11 
12 
13 
14 namespace Dune {
15 namespace Functions {
16 namespace Impl {
17 
18 
19 
20 template<class V>
21 struct FlatVectorBackend
22 {
23 
24  template<class VV, class Index,
25  typename std::enable_if< models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0>
26  static auto getEntry(VV&& v, const Index& i)
27  ->decltype(v[i])
28  {
29  return v[i];
30  }
31 
32  template<class VV, class Index,
33  typename std::enable_if< not models<Concept::HasIndexAccess, VV, Index>(), int>::type = 0>
34  static auto getEntry(VV&& v, const Index& i)
35  ->decltype(v)
36  {
37  return std::forward<VV>(v);
38  }
39 
40  template<class VV,
41  typename std::enable_if< models<Concept::HasSizeMethod, VV>(), int>::type = 0>
42  static auto size(VV&& v)
43  ->decltype(v.size())
44  {
45  return v.size();
46  }
47 
48  template<class VV,
49  typename std::enable_if< not models<Concept::HasSizeMethod, VV>(), int>::type = 0>
50  static std::size_t size(VV&& v)
51  {
52  return 1;
53  }
54 
55 };
56 
57 
58 
59 
60 
61 template<class K, int n, int m>
62 struct FlatVectorBackend<typename Dune::FieldMatrix<K, n, m> >
63 {
64 
65  template<class VV, class Index>
66  static auto getEntry(VV&& v, const Index& i) -> decltype(v[i/m][i%m])
67  {
68  return v[i/m][i%m];
69  }
70 
71  template<class VV>
72  static int size(VV&& v)
73  {
74  return n*m;
75  }
76 };
77 
78 
79 
80 
81 template<class T>
82 class FlatVectorView
83 {
84  using Backend = FlatVectorBackend<std::decay_t<T>>;
85 public:
86  FlatVectorView(T& t) :
87  t_(&t)
88  {}
89 
90  auto size() const
91  {
92  return Backend::size(*t_);
93  }
94 
95  template<class Index>
96  decltype(auto) operator[](const Index& i) const
97  {
98  return Backend::getEntry(*t_, i);
99  }
100 
101  template<class Index>
102  decltype(auto) operator[](const Index& i)
103  {
104  return Backend::getEntry(*t_, i);
105  }
106 
107 private:
108  T* t_;
109 };
110 
111 
112 template<class T>
113 class FlatVectorView<T&&>
114 {
115  using Backend = FlatVectorBackend<std::decay_t<T>>;
116 public:
117  FlatVectorView(T&& t) :
118  t_(std::move(t))
119  {}
120 
121  auto size() const
122  {
123  return Backend::size(t_);
124  }
125 
126  template<class Index>
127  decltype(auto) operator[](const Index& i) const
128  {
129  return Backend::getEntry(t_, i);
130  }
131 
132  template<class Index>
133  decltype(auto) operator[](const Index& i)
134  {
135  return Backend::getEntry(t_, i);
136  }
137 
138 private:
139  T t_;
140 };
141 
142 } // namespace Impl
143 
144 
145 
158 template<class T>
159 auto flatVectorView(T& t)
160 {
161  return Impl::FlatVectorView<T>(t);
162 }
163 
176 template<class T>
177 auto flatVectorView(const T& t)
178 {
179  return Impl::FlatVectorView<const T>(t);
180 }
181 
194 template<class T>
195 auto flatVectorView(T&& t)
196 {
197  return Impl::FlatVectorView<T&&>(std::move(t));
198 }
199 
200 
201 } // namespace Dune::Functions
202 } // namespace Dune
203 
204 
205 #endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_FLATVECTORVIEW_HH
Definition: polynomial.hh:10
auto flatVectorView(T &t)
Create flat vector view of passed mutable container.
Definition: flatvectorview.hh:159