2 #ifndef DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH 3 #define DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH 5 #include <dune/common/typetraits.hh> 17 template<
typename GFSV,
typename GFSU,
typename C,
typename Stats>
19 :
public Backend::impl::Wrapper<C>
22 friend Backend::impl::Wrapper<C>;
27 typedef ElementType
E;
36 typedef typename GFSV::Ordering::Traits::ContainerIndex
RowIndex;
37 typedef typename GFSU::Ordering::Traits::ContainerIndex
ColIndex;
39 typedef typename ISTL::build_pattern_type<C,GFSV,GFSU,typename GFSV::Ordering::ContainerAllocationTag>::type
Pattern;
49 typedef typename std::conditional<
51 std::vector<PatternStatistics>,
53 >::type StatisticsReturnType;
57 template<
typename RowCache,
typename ColCache>
60 template<
typename RowCache,
typename ColCache>
63 template<
typename RowCache,
typename ColCache>
66 template<
typename RowCache,
typename ColCache>
71 : _container(
std::make_shared<Container>())
73 _stats = go.matrixBackend().buildPattern(go,*
this);
87 : _container(
Dune::stackobject_to_shared_ptr(container))
89 _stats = go.matrixBackend().buildPattern(go,*
this);
94 : _container(
std::make_shared<Container>())
96 _stats = go.matrixBackend().buildPattern(go,*
this);
106 : _container(
std::make_shared<Container>())
110 : _container(
std::make_shared<Container>(*(rhs._container)))
120 (*_container) = (*(rhs._container));
124 _container = std::make_shared<Container>(*(rhs._container));
142 DUNE_THROW(InvalidStateException,
"no pattern statistics available");
146 const std::vector<PatternStatistics>&
patternStatistics(std::true_type multiple)
const 149 DUNE_THROW(InvalidStateException,
"no pattern statistics available");
163 void attach(std::shared_ptr<Container> container)
165 _container = container;
170 return bool(_container);
173 const std::shared_ptr<Container>&
storage()
const 180 return _container->N();
185 return _container->M();
202 return ISTL::access_matrix_element(
ISTL::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
205 const E&
operator()(
const RowIndex& ri,
const ColIndex& ci)
const 207 return ISTL::access_matrix_element(
ISTL::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
212 const Container& native()
const 224 template<
typename RowCache,
typename ColCache>
227 return &((*this)(row_cache.containerIndex(0),col_cache.containerIndex(0)));
230 template<
typename RowCache,
typename ColCache>
231 const value_type*
data(
const RowCache& row_cache,
const ColCache& col_cache)
const 233 return &((*this)(row_cache.containerIndex(0),col_cache.containerIndex(0)));
242 void clear_row(
const RowIndex& ri,
const E& diagonal_entry)
245 ISTL::write_matrix_element_if_exists(diagonal_entry,
ISTL::container_tag(*_container),*_container,ri,ri,ri.size()-1,ri.size()-1);
251 ISTL::write_matrix_element_if_exists_to_block(diagonal_entry,
ISTL::container_tag(*_container),*_container,ri,ri,ri.size()-1,ri.size()-1);
256 std::shared_ptr<Container> _container;
257 std::vector<PatternStatistics> _stats;
266 #endif // DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH Definition: uncachedmatrixview.hh:12
BCRSMatrix(const GO &go)
Definition: bcrsmatrix.hh:70
Definition: bcrsmatrix.hh:18
tags::container< T >::type container_tag(const T &)
Gets instance of container tag associated with T.
Definition: backend/istl/tags.hh:246
const Entity & e
Definition: localfunctionspace.hh:120
Stats PatternStatistics
Definition: bcrsmatrix.hh:41
BCRSMatrix(const GO &go, const E &e)
Definition: bcrsmatrix.hh:93
const value_type * data(const RowCache &row_cache, const ColCache &col_cache) const
Definition: bcrsmatrix.hh:231
BCRSMatrix(Backend::attached_container)
Creates an BCRSMatrix with an empty underlying ISTL matrix.
Definition: bcrsmatrix.hh:105
void clear_row_block(const RowIndex &ri, const E &diagonal_entry)
Definition: bcrsmatrix.hh:248
Definition: aliasedmatrixview.hh:12
BCRSMatrix & operator=(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:113
GFSV::Ordering::Traits::ContainerIndex RowIndex
Definition: bcrsmatrix.hh:36
C::size_type size_type
Definition: bcrsmatrix.hh:31
E & operator()(const RowIndex &ri, const ColIndex &ci)
Definition: bcrsmatrix.hh:200
C::block_type block_type
Definition: bcrsmatrix.hh:30
GFSU::Ordering::Traits::ContainerIndex ColIndex
Definition: bcrsmatrix.hh:37
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/common/tags.hh:27
void clear_row(const RowIndex &ri, const E &diagonal_entry)
Definition: bcrsmatrix.hh:242
void flush()
Definition: bcrsmatrix.hh:236
BCRSMatrix & operator*=(const E &e)
Definition: bcrsmatrix.hh:194
void detach()
Definition: bcrsmatrix.hh:157
Definition: uncachedmatrixview.hh:165
bool attached() const
Definition: bcrsmatrix.hh:168
E value_type
Definition: bcrsmatrix.hh:43
BCRSMatrix(const GO &go, Container &container)
Construct matrix container using an externally given matrix as storage.
Definition: bcrsmatrix.hh:86
size_type M() const
Definition: bcrsmatrix.hh:183
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
GFSV TestGridFunctionSpace
Definition: bcrsmatrix.hh:34
C::field_type field_type
Definition: bcrsmatrix.hh:29
void finalize()
Definition: bcrsmatrix.hh:239
const E & operator()(const RowIndex &ri, const ColIndex &ci) const
Definition: bcrsmatrix.hh:205
ElementType E
Definition: bcrsmatrix.hh:27
const std::shared_ptr< Container > & storage() const
Definition: bcrsmatrix.hh:173
C::field_type ElementType
Definition: bcrsmatrix.hh:26
size_type N() const
Definition: bcrsmatrix.hh:178
BCRSMatrix(Backend::unattached_container=Backend::unattached_container())
Creates an BCRSMatrix without allocating an underlying ISTL matrix.
Definition: bcrsmatrix.hh:101
ISTL::build_pattern_type< C, GFSV, GFSU, typename GFSV::Ordering::ContainerAllocationTag >::type Pattern
Definition: bcrsmatrix.hh:39
GFSU TrialGridFunctionSpace
Definition: bcrsmatrix.hh:33
value_type * data(const RowCache &row_cache, const ColCache &col_cache)
Definition: bcrsmatrix.hh:225
Tag for requesting a vector or matrix container without a pre-attached underlying object...
Definition: backend/common/tags.hh:23
C Container
Definition: bcrsmatrix.hh:28
Definition: aliasedmatrixview.hh:166
Various tags for influencing backend behavior.
void attach(std::shared_ptr< Container > container)
Definition: bcrsmatrix.hh:163
BCRSMatrix(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:109
const StatisticsReturnType & patternStatistics() const
Returns pattern statistics for all contained BCRSMatrix objects.
Definition: bcrsmatrix.hh:130