dune-pdelab  2.5-dev
bcrsmatrix.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
3 #define DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
4 
5 // this is here for backwards compatibility and deprecation warnings, remove after 2.5.0
6 #include "ensureistlinclude.hh"
7 
8 #include <dune/common/typetraits.hh>
13 
14 namespace Dune {
15  namespace PDELab {
16 
17  namespace ISTL {
18 
19  template<typename GFSV, typename GFSU, typename C, typename Stats>
20  class BCRSMatrix
21  : public Backend::impl::Wrapper<C>
22  {
23 
24  friend Backend::impl::Wrapper<C>;
25 
26  public:
27 
28  typedef typename C::field_type ElementType;
29  typedef ElementType E;
30  typedef C Container;
31  typedef typename C::field_type field_type;
32  typedef typename C::block_type block_type;
33  typedef typename C::size_type size_type;
34 
35  typedef GFSU TrialGridFunctionSpace;
36  typedef GFSV TestGridFunctionSpace;
37 
38  typedef typename GFSV::Ordering::Traits::ContainerIndex RowIndex;
39  typedef typename GFSU::Ordering::Traits::ContainerIndex ColIndex;
40 
41  typedef typename ISTL::build_pattern_type<C,GFSV,GFSU,typename GFSV::Ordering::ContainerAllocationTag>::type Pattern;
42 
43  typedef Stats PatternStatistics;
44 
45 #ifndef DOXYGEN
46 
47  // some trickery to avoid exposing average users to the fact that there might
48  // be multiple statistics objects
49  typedef typename std::conditional<
50  (C::blocklevel > 2),
51  std::vector<PatternStatistics>,
52  PatternStatistics
53  >::type StatisticsReturnType;
54 
55 #endif // DOXYGEN
56 
57  template<typename RowCache, typename ColCache>
59 
60  template<typename RowCache, typename ColCache>
62 
63  template<typename GO>
64  explicit BCRSMatrix (const GO& go)
65  : _container(std::make_shared<Container>())
66  {
67  _stats = go.matrixBackend().buildPattern(go,*this);
68  }
69 
79  template<typename GO>
80  BCRSMatrix (const GO& go, Container& container)
81  : _container(Dune::stackobject_to_shared_ptr(container))
82  {
83  _stats = go.matrixBackend().buildPattern(go,*this);
84  }
85 
86  template<typename GO>
87  BCRSMatrix (const GO& go, const E& e)
88  : _container(std::make_shared<Container>())
89  {
90  _stats = go.matrixBackend().buildPattern(go,*this);
91  (*_container) = e;
92  }
93 
96  {}
97 
100  : _container(std::make_shared<Container>())
101  {}
102 
103  BCRSMatrix(const BCRSMatrix& rhs)
104  : _container(std::make_shared<Container>(*(rhs._container)))
105  {}
106 
108  {
109  if (this == &rhs)
110  return *this;
111  _stats.clear();
112  if (attached())
113  {
114  (*_container) = (*(rhs._container));
115  }
116  else
117  {
118  _container = std::make_shared<Container>(*(rhs._container));
119  }
120  return *this;
121  }
122 
124  const StatisticsReturnType& patternStatistics() const
125  {
126  return patternStatistics(std::integral_constant<bool,(C::blocklevel > 2)>());
127  }
128 
129 #ifndef DOXYGEN
130 
131  private:
132 
133  const PatternStatistics& patternStatistics(std::false_type multiple) const
134  {
135  if (_stats.empty())
136  DUNE_THROW(InvalidStateException,"no pattern statistics available");
137  return _stats[0];
138  }
139 
140  const std::vector<PatternStatistics>& patternStatistics(std::true_type multiple) const
141  {
142  if (_stats.empty())
143  DUNE_THROW(InvalidStateException,"no pattern statistics available");
144  return _stats;
145  }
146 
147  public:
148 
149 #endif
150 
151  void detach()
152  {
153  _container.reset();
154  _stats.clear();
155  }
156 
157  void attach(std::shared_ptr<Container> container)
158  {
159  _container = container;
160  }
161 
162  bool attached() const
163  {
164  return bool(_container);
165  }
166 
167  const std::shared_ptr<Container>& storage() const
168  {
169  return _container;
170  }
171 
172  size_type N() const
173  {
174  return _container->N();
175  }
176 
177  size_type M() const
178  {
179  return _container->M();
180  }
181 
183  {
184  (*_container) = e;
185  return *this;
186  }
187 
189  {
190  (*_container) *= e;
191  return *this;
192  }
193 
194  E& operator()(const RowIndex& ri, const ColIndex& ci)
195  {
196  return ISTL::access_matrix_element(ISTL::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
197  }
198 
199  const E& operator()(const RowIndex& ri, const ColIndex& ci) const
200  {
201  return ISTL::access_matrix_element(ISTL::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
202  }
203 
204  private:
205 
206  const Container& native() const
207  {
208  return *_container;
209  }
210 
211  Container& native()
212  {
213  return *_container;
214  }
215 
216  public:
217 
218  void flush()
219  {}
220 
221  void finalize()
222  {}
223 
224  void clear_row(const RowIndex& ri, const E& diagonal_entry)
225  {
226  ISTL::clear_matrix_row(ISTL::container_tag(*_container),*_container,ri,ri.size()-1);
227  ISTL::write_matrix_element_if_exists(diagonal_entry,ISTL::container_tag(*_container),*_container,ri,ri,ri.size()-1,ri.size()-1);
228  }
229 
230  private:
231 
232  std::shared_ptr<Container> _container;
233  std::vector<PatternStatistics> _stats;
234 
235  };
236 
237  } // namespace ISTL
238 
239  } // namespace PDELab
240 } // namespace Dune
241 
242 #endif // DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
size_type N() const
Definition: bcrsmatrix.hh:172
tags::container< T >::type container_tag(const T &)
Gets instance of container tag associated with T.
Definition: backend/istl/tags.hh:249
E & operator()(const RowIndex &ri, const ColIndex &ci)
Definition: bcrsmatrix.hh:194
GFSV TestGridFunctionSpace
Definition: bcrsmatrix.hh:36
Tag for requesting a vector or matrix container without a pre-attached underlying object...
Definition: backend/common/tags.hh:25
void attach(std::shared_ptr< Container > container)
Definition: bcrsmatrix.hh:157
const E & operator()(const RowIndex &ri, const ColIndex &ci) const
Definition: bcrsmatrix.hh:199
BCRSMatrix(Backend::unattached_container=Backend::unattached_container())
Creates an BCRSMatrix without allocating an underlying ISTL matrix.
Definition: bcrsmatrix.hh:95
C::block_type block_type
Definition: bcrsmatrix.hh:32
Definition: uncachedmatrixview.hh:12
GFSU TrialGridFunctionSpace
Definition: bcrsmatrix.hh:35
BCRSMatrix(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:103
C::size_type size_type
Definition: bcrsmatrix.hh:33
Various tags for influencing backend behavior.
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void clear_row(const RowIndex &ri, const E &diagonal_entry)
Definition: bcrsmatrix.hh:224
size_type M() const
Definition: bcrsmatrix.hh:177
BCRSMatrix & operator*=(const E &e)
Definition: bcrsmatrix.hh:188
Stats PatternStatistics
Definition: bcrsmatrix.hh:43
const Entity & e
Definition: localfunctionspace.hh:111
const StatisticsReturnType & patternStatistics() const
Returns pattern statistics for all contained BCRSMatrix objects.
Definition: bcrsmatrix.hh:124
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/common/tags.hh:29
C Container
Definition: bcrsmatrix.hh:30
BCRSMatrix(const GO &go, const E &e)
Definition: bcrsmatrix.hh:87
BCRSMatrix(Backend::attached_container)
Creates an BCRSMatrix with an empty underlying ISTL matrix.
Definition: bcrsmatrix.hh:99
Definition: uncachedmatrixview.hh:157
const std::shared_ptr< Container > & storage() const
Definition: bcrsmatrix.hh:167
C::field_type ElementType
Definition: bcrsmatrix.hh:28
STL namespace.
ISTL::build_pattern_type< C, GFSV, GFSU, typename GFSV::Ordering::ContainerAllocationTag >::type Pattern
Definition: bcrsmatrix.hh:41
C::field_type field_type
Definition: bcrsmatrix.hh:31
BCRSMatrix(const GO &go)
Definition: bcrsmatrix.hh:64
Definition: bcrsmatrix.hh:20
void finalize()
Definition: bcrsmatrix.hh:221
ElementType E
Definition: bcrsmatrix.hh:29
void flush()
Definition: bcrsmatrix.hh:218
BCRSMatrix(const GO &go, Container &container)
Construct matrix container using an externally given matrix as storage.
Definition: bcrsmatrix.hh:80
GFSV::Ordering::Traits::ContainerIndex RowIndex
Definition: bcrsmatrix.hh:38
BCRSMatrix & operator=(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:107
GFSU::Ordering::Traits::ContainerIndex ColIndex
Definition: bcrsmatrix.hh:39
bool attached() const
Definition: bcrsmatrix.hh:162
void detach()
Definition: bcrsmatrix.hh:151