dune-pdelab  2.5-dev
uncachedmatrixview.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_BACKEND_COMMON_UNCACHEDMATRIXVIEW_HH
3 #define DUNE_PDELAB_BACKEND_COMMON_UNCACHEDMATRIXVIEW_HH
4 
5 #include <type_traits>
6 
7 namespace Dune {
8  namespace PDELab {
9 
10 
11  template<typename M_, typename RowCache, typename ColCache>
13  {
14 
15  public:
16 
17  typedef typename std::remove_const<M_>::type Container;
18 
19  static_assert(
20  (std::is_same<
21  typename RowCache::LocalFunctionSpace::Traits::GridFunctionSpace,
22  typename Container::TestGridFunctionSpace
23  >::value),
24  "The RowCache passed to LocalView must belong to the underlying GFSV"
25  );
26 
27  static_assert(
28  (std::is_same<
29  typename ColCache::LocalFunctionSpace::Traits::GridFunctionSpace,
30  typename Container::TrialGridFunctionSpace
31  >::value),
32  "The ColCache passed to LocalView must belong to the underlying GFSU"
33  );
34 
35  public:
36 
37  typedef typename Container::field_type E;
38  typedef typename Container::size_type size_type;
39 
40  typedef E ElementType;
41 
42  typedef RowCache RowIndexCache;
43  typedef ColCache ColIndexCache;
44 
45  typedef typename RowCache::LocalFunctionSpace LFSV;
46  typedef typename ColCache::LocalFunctionSpace LFSU;
47 
48  typedef typename LFSV::Traits::DOFIndex RowDOFIndex;
49  typedef typename LFSV::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex RowContainerIndex;
50 
51  typedef typename LFSU::Traits::DOFIndex ColDOFIndex;
52  typedef typename LFSU::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex ColContainerIndex;
53 
55  : _container(nullptr)
56  , _row_cache(nullptr)
57  , _col_cache(nullptr)
58  {}
59 
62  , _row_cache(nullptr)
63  , _col_cache(nullptr)
64  {}
65 
67  {
68  assert(_row_cache);
69  return *_row_cache;
70  }
71 
73  {
74  assert(_col_cache);
75  return *_col_cache;
76  }
77 
78  void attach(M_& container)
79  {
81  }
82 
83  void detach()
84  {
85  _container = nullptr;
86  }
87 
88  void bind(const RowCache& row_cache, const ColCache& col_cache)
89  {
90  _row_cache = &row_cache;
91  _col_cache = &col_cache;
92  }
93 
94  void unbind()
95  {}
96 
97  size_type N() const
98  {
99  return rowIndexCache().size();
100  }
101 
102  size_type M() const
103  {
104  return colIndexCache().size();
105  }
106 
107  template<typename LC>
108  void read(LC& local_container) const
109  {
110  for (size_type i = 0; i < N(); ++i)
111  for (size_type j = 0; j < M(); ++j)
112  local_container.getEntry(i,j) = container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j));
113  }
114 
115 
116 
118  {
119  return container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j));
120  }
121 
122  // disable this function if DOFIndex and ContainerIndex have the same type - required for interoperability
123  // with function spaces based on dune-functions bases
124  template<typename RDI, typename CDI>
125  std::enable_if_t<
126  (std::is_same<RDI,RowDOFIndex>{} and std::is_same<CDI,ColDOFIndex>{} and not
127  (std::is_same<RDI,RowContainerIndex>{} and std::is_same<CDI,ColContainerIndex>{})),
128  const ElementType&
129  >
130  operator()(const RDI& i, const CDI& j) const
131  {
132  return container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j));
133  }
134 
136  {
137  return container()(i,j);
138  }
139 
141  {
142  return container()(i,colIndexCache().containerIndex(j));
143  }
144 
146  {
147  return container()(rowIndexCache().containerIndex(i),j);
148  }
149 
150  const Container& container() const
151  {
152  return *_container;
153  }
154 
155  protected:
156 
158  const RowCache* _row_cache;
159  const ColCache* _col_cache;
160 
161  };
162 
163 
164  template<typename M_, typename RowCache, typename ColCache>
166  : public ConstUncachedMatrixView<M_,RowCache,ColCache>
167  {
168 
170 
171  public:
172 
173  typedef M_ Container;
174  typedef typename Container::ElementType ElementType;
175  typedef typename Container::size_type size_type;
176 
177  typedef RowCache RowIndexCache;
178  typedef ColCache ColIndexCache;
179 
180  typedef typename RowCache::LocalFunctionSpace LFSV;
181  typedef typename ColCache::LocalFunctionSpace LFSU;
182 
183  typedef typename LFSV::Traits::DOFIndex RowDOFIndex;
184  typedef typename LFSV::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex RowContainerIndex;
185 
186  typedef typename LFSU::Traits::DOFIndex ColDOFIndex;
187  typedef typename LFSU::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex ColContainerIndex;
188 
189  using BaseT::rowIndexCache;
190  using BaseT::colIndexCache;
191  using BaseT::N;
192  using BaseT::M;
193 
194  // Explicitly pull in operator() from the base class to work around a problem
195  // with clang not finding the const overloads of the operator from the base class.
196  using BaseT::operator();
197 
199  {}
200 
202  : BaseT(container)
203  {}
204 
205  void commit()
206  {}
207 
208  template<typename LC>
209  void write(const LC& local_container)
210  {
211  for (size_type i = 0; i < N(); ++i)
212  for (size_type j = 0; j < M(); ++j)
213  container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j)) = local_container.getEntry(i,j);
214  }
215 
216  template<typename LC>
217  void add(const LC& local_container)
218  {
219  for (size_type i = 0; i < N(); ++i)
220  for (size_type j = 0; j < M(); ++j)
221  container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j)) += local_container.getEntry(i,j);
222  }
223 
224 
225 
227  {
228  return container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j));
229  }
230 
231  // disable this function if DOFIndex and ContainerIndex have the same type - required for interoperability
232  // with function spaces based on dune-functions bases
233  template<typename RDI, typename CDI>
234  std::enable_if_t<
235  (std::is_same<RDI,RowDOFIndex>{} and std::is_same<CDI,ColDOFIndex>{} and not
236  (std::is_same<RDI,RowContainerIndex>{} and std::is_same<CDI,ColContainerIndex>{})),
237  const ElementType&
238  >
239  operator()(const RDI& i, const CDI& j)
240  {
241  return container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j));
242  }
243 
245  {
246  return container()(i,j);
247  }
248 
250  {
251  return container()(i,colIndexCache().containerIndex(j));
252  }
253 
255  {
256  return container()(rowIndexCache().containerIndex(i),j);
257  }
258 
259  void add(size_type i, size_type j, const ElementType& v)
260  {
261  container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j)) += v;
262  }
263 
264  // disable this function if DOFIndex and ContainerIndex have the same type - required for interoperability
265  // with function spaces based on dune-functions bases
266  template<typename RDI, typename CDI>
267  std::enable_if_t<
268  (std::is_same<RDI,RowDOFIndex>{} and std::is_same<CDI,ColDOFIndex>{} and not
269  (std::is_same<RDI,RowContainerIndex>{} and std::is_same<CDI,ColContainerIndex>{}))
270  >
271  add(const RDI& i, const CDI& j, const ElementType& v)
272  {
273  container()(rowIndexCache().containerIndex(i),colIndexCache().containerIndex(j)) += v;
274  }
275 
276  void add(const RowContainerIndex& i, const ColContainerIndex& j, const ElementType& v)
277  {
278  container()(i,j) += v;
279  }
280 
281  void add(const RowContainerIndex& i, size_type j, const ElementType& v)
282  {
283  container()(i,colIndexCache().containerIndex(j)) += v;
284  }
285 
286  void add(size_type i, const ColContainerIndex& j, const ElementType& v)
287  {
288  container()(rowIndexCache().containerIndex(i),j) += v;
289  }
290 
292  {
293  return *(this->_container);
294  }
295 
296  };
297 
298 
299  } // namespace PDELab
300 } // namespace Dune
301 
302 #endif // DUNE_PDELAB_BACKEND_COMMON_UNCACHEDMATRIXVIEW_HH
Dune::PDELab::UncachedMatrixView::LFSV
RowCache::LocalFunctionSpace LFSV
Definition: uncachedmatrixview.hh:180
Dune::PDELab::ConstUncachedMatrixView
Definition: uncachedmatrixview.hh:12
Dune::PDELab::ConstUncachedMatrixView::ConstUncachedMatrixView
ConstUncachedMatrixView(M_ &container)
Definition: uncachedmatrixview.hh:60
Dune::PDELab::ConstUncachedMatrixView::E
Container::field_type E
Definition: uncachedmatrixview.hh:25
Dune::PDELab::ConstUncachedMatrixView::detach
void detach()
Definition: uncachedmatrixview.hh:83
Dune::PDELab::ConstUncachedMatrixView::rowIndexCache
const RowIndexCache & rowIndexCache() const
Definition: uncachedmatrixview.hh:66
Dune::PDELab::ConstUncachedMatrixView::attach
void attach(M_ &container)
Definition: uncachedmatrixview.hh:78
Dune::PDELab::ConstUncachedMatrixView::ElementType
E ElementType
Definition: uncachedmatrixview.hh:40
Dune::PDELab::ConstUncachedMatrixView::ColContainerIndex
LFSU::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex ColContainerIndex
Definition: uncachedmatrixview.hh:52
Dune::PDELab::ConstUncachedMatrixView::_container
M_ * _container
Definition: uncachedmatrixview.hh:157
Dune::PDELab::UncachedMatrixView::add
void add(const RowContainerIndex &i, size_type j, const ElementType &v)
Definition: uncachedmatrixview.hh:281
Dune::PDELab::UncachedMatrixView::UncachedMatrixView
UncachedMatrixView()
Definition: uncachedmatrixview.hh:198
Dune::PDELab::ConstUncachedMatrixView::_col_cache
const ColCache * _col_cache
Definition: uncachedmatrixview.hh:159
Dune
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::PDELab::ConstUncachedMatrixView::ConstUncachedMatrixView
ConstUncachedMatrixView()
Definition: uncachedmatrixview.hh:54
Dune::PDELab::ConstUncachedMatrixView::Container
std::remove_const< M_ >::type Container
Definition: uncachedmatrixview.hh:17
Dune::PDELab::UncachedMatrixView::ElementType
Container::ElementType ElementType
Definition: uncachedmatrixview.hh:174
Dune::PDELab::ConstUncachedMatrixView::N
size_type N() const
Definition: uncachedmatrixview.hh:97
Dune::PDELab::UncachedMatrixView::RowDOFIndex
LFSV::Traits::DOFIndex RowDOFIndex
Definition: uncachedmatrixview.hh:183
Dune::PDELab::ConstUncachedMatrixView::container
const Container & container() const
Definition: uncachedmatrixview.hh:150
Dune::PDELab::UncachedMatrixView::RowContainerIndex
LFSV::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex RowContainerIndex
Definition: uncachedmatrixview.hh:184
Dune::PDELab::UncachedMatrixView::ColDOFIndex
LFSU::Traits::DOFIndex ColDOFIndex
Definition: uncachedmatrixview.hh:186
Dune::PDELab::UncachedMatrixView
Definition: uncachedmatrixview.hh:165
Dune::PDELab::ConstUncachedMatrixView::operator()
std::enable_if_t<(std::is_same< RDI, RowDOFIndex >{} and std::is_same< CDI, ColDOFIndex >{} and not(std::is_same< RDI, RowContainerIndex >{} and std::is_same< CDI, ColContainerIndex >{})), const ElementType & > operator()(const RDI &i, const CDI &j) const
Definition: uncachedmatrixview.hh:130
Dune::PDELab::UncachedMatrixView::ColContainerIndex
LFSU::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex ColContainerIndex
Definition: uncachedmatrixview.hh:187
Dune::PDELab::UncachedMatrixView::add
void add(size_type i, const ColContainerIndex &j, const ElementType &v)
Definition: uncachedmatrixview.hh:286
Dune::PDELab::UncachedMatrixView::add
void add(const LC &local_container)
Definition: uncachedmatrixview.hh:217
Dune::PDELab::ConstUncachedMatrixView::read
void read(LC &local_container) const
Definition: uncachedmatrixview.hh:108
Dune::PDELab::ConstUncachedMatrixView::operator()
const ElementType & operator()(size_type i, const ColContainerIndex &j) const
Definition: uncachedmatrixview.hh:145
Dune::PDELab::ConstUncachedMatrixView::LFSU
ColCache::LocalFunctionSpace LFSU
Definition: uncachedmatrixview.hh:46
Dune::PDELab::UncachedMatrixView::operator()
ElementType & operator()(const RowContainerIndex &i, size_type j)
Definition: uncachedmatrixview.hh:249
Dune::PDELab::UncachedMatrixView::add
std::enable_if_t<(std::is_same< RDI, RowDOFIndex >{} and std::is_same< CDI, ColDOFIndex >{} and not(std::is_same< RDI, RowContainerIndex >{} and std::is_same< CDI, ColContainerIndex >{})) > add(const RDI &i, const CDI &j, const ElementType &v)
Definition: uncachedmatrixview.hh:271
Dune::PDELab::ConstUncachedMatrixView::ColIndexCache
ColCache ColIndexCache
Definition: uncachedmatrixview.hh:43
Dune::PDELab::UncachedMatrixView::RowIndexCache
RowCache RowIndexCache
Definition: uncachedmatrixview.hh:177
Dune::PDELab::ConstUncachedMatrixView::unbind
void unbind()
Definition: uncachedmatrixview.hh:94
Dune::PDELab::UncachedMatrixView::operator()
ElementType & operator()(const RowContainerIndex &i, const ColContainerIndex &j)
Definition: uncachedmatrixview.hh:244
Dune::PDELab::UncachedMatrixView::Container
M_ Container
Definition: uncachedmatrixview.hh:173
Dune::PDELab::UncachedMatrixView::add
void add(size_type i, size_type j, const ElementType &v)
Definition: uncachedmatrixview.hh:259
Dune::PDELab::ConstUncachedMatrixView::RowIndexCache
RowCache RowIndexCache
Definition: uncachedmatrixview.hh:42
Dune::PDELab::ConstUncachedMatrixView::ColDOFIndex
LFSU::Traits::DOFIndex ColDOFIndex
Definition: uncachedmatrixview.hh:51
Dune::PDELab::ConstUncachedMatrixView::operator()
const ElementType & operator()(const RowContainerIndex &i, size_type j) const
Definition: uncachedmatrixview.hh:140
Dune::PDELab::ConstUncachedMatrixView::colIndexCache
const ColIndexCache & colIndexCache() const
Definition: uncachedmatrixview.hh:72
Dune::PDELab::ConstUncachedMatrixView::_row_cache
const RowCache * _row_cache
Definition: uncachedmatrixview.hh:158
Dune::PDELab::ConstUncachedMatrixView::bind
void bind(const RowCache &row_cache, const ColCache &col_cache)
Definition: uncachedmatrixview.hh:88
Dune::PDELab::UncachedMatrixView::write
void write(const LC &local_container)
Definition: uncachedmatrixview.hh:209
Dune::PDELab::UncachedMatrixView::operator()
ElementType & operator()(size_type i, size_type j)
Definition: uncachedmatrixview.hh:226
Dune::PDELab::UncachedMatrixView::ColIndexCache
ColCache ColIndexCache
Definition: uncachedmatrixview.hh:178
Dune::PDELab::ConstUncachedMatrixView::LFSV
RowCache::LocalFunctionSpace LFSV
Definition: uncachedmatrixview.hh:45
Dune::PDELab::ConstUncachedMatrixView::RowDOFIndex
LFSV::Traits::DOFIndex RowDOFIndex
Definition: uncachedmatrixview.hh:48
Dune::PDELab::UncachedMatrixView::LFSU
ColCache::LocalFunctionSpace LFSU
Definition: uncachedmatrixview.hh:181
Dune::PDELab::UncachedMatrixView::operator()
ElementType & operator()(size_type i, const ColContainerIndex &j)
Definition: uncachedmatrixview.hh:254
Dune::PDELab::UncachedMatrixView::add
void add(const RowContainerIndex &i, const ColContainerIndex &j, const ElementType &v)
Definition: uncachedmatrixview.hh:276
Dune::PDELab::ConstUncachedMatrixView::RowContainerIndex
LFSV::Traits::GridFunctionSpace::Ordering::Traits::ContainerIndex RowContainerIndex
Definition: uncachedmatrixview.hh:49
Dune::PDELab::ConstUncachedMatrixView::operator()
const ElementType & operator()(const RowContainerIndex &i, const ColContainerIndex &j) const
Definition: uncachedmatrixview.hh:135
value
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Dune::PDELab::UncachedMatrixView::operator()
std::enable_if_t<(std::is_same< RDI, RowDOFIndex >{} and std::is_same< CDI, ColDOFIndex >{} and not(std::is_same< RDI, RowContainerIndex >{} and std::is_same< CDI, ColContainerIndex >{})), const ElementType & > operator()(const RDI &i, const CDI &j)
Definition: uncachedmatrixview.hh:239
Dune::PDELab::ConstUncachedMatrixView::operator()
const ElementType & operator()(size_type i, size_type j) const
Definition: uncachedmatrixview.hh:117
Dune::PDELab::UncachedMatrixView::commit
void commit()
Definition: uncachedmatrixview.hh:205
Dune::PDELab::UncachedMatrixView::container
Container & container()
Definition: uncachedmatrixview.hh:291
Dune::PDELab::UncachedMatrixView::size_type
Container::size_type size_type
Definition: uncachedmatrixview.hh:175
Dune::PDELab::UncachedMatrixView::UncachedMatrixView
UncachedMatrixView(Container &container)
Definition: uncachedmatrixview.hh:201
Dune::PDELab::ConstUncachedMatrixView::M
size_type M() const
Definition: uncachedmatrixview.hh:102
Dune::PDELab::ConstUncachedMatrixView::size_type
Container::size_type size_type
Definition: uncachedmatrixview.hh:38