2 #ifndef DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
3 #define DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
5 #include <dune/common/typetraits.hh>
16 template<
typename GFSV,
typename GFSU,
typename C,
typename Stats>
18 :
public Backend::impl::Wrapper<C>
21 friend Backend::impl::Wrapper<C>;
26 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;
47 typedef typename conditional<
49 std::vector<PatternStatistics>,
51 >::type StatisticsReturnType;
55 template<
typename RowCache,
typename ColCache>
58 template<
typename RowCache,
typename ColCache>
63 : _container(
std::make_shared<Container>())
65 _stats = go.matrixBackend().buildPattern(go,*
this);
79 : _container(
Dune::stackobject_to_shared_ptr(container))
81 _stats = go.matrixBackend().buildPattern(go,*
this);
86 : _container(
std::make_shared<Container>())
88 _stats = go.matrixBackend().buildPattern(go,*
this);
98 : _container(
std::make_shared<Container>())
102 : _container(
std::make_shared<Container>(*(rhs._container)))
112 (*_container) = (*(rhs._container));
116 _container = std::make_shared<Container>(*(rhs._container));
134 DUNE_THROW(InvalidStateException,
"no pattern statistics available");
138 const std::vector<PatternStatistics>&
patternStatistics(true_type multiple)
const
141 DUNE_THROW(InvalidStateException,
"no pattern statistics available");
155 void attach(std::shared_ptr<Container> container)
157 _container = container;
162 return bool(_container);
165 const std::shared_ptr<Container>&
storage()
const
172 return _container->N();
177 return _container->M();
194 return istl::access_matrix_element(
istl::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
197 const E&
operator()(
const RowIndex& ri,
const ColIndex& ci)
const
199 return istl::access_matrix_element(
istl::container_tag(*_container),*_container,ri,ci,ri.size()-1,ci.size()-1);
203 DUNE_DEPRECATED_MSG(
"base() is deprecated and will be removed after PDELab 2.4. Use Backend::native() instead.")
210 DUNE_DEPRECATED_MSG(
"base() is deprecated and will be removed after PDELab 2.4. Use Backend::native() instead.")
218 const Container& native()
const
236 void clear_row(
const RowIndex& ri,
const E& diagonal_entry)
239 istl::write_matrix_element_if_exists(diagonal_entry,
istl::container_tag(*_container),*_container,ri,ri,ri.size()-1,ri.size()-1);
244 std::shared_ptr<Container> _container;
245 std::vector<PatternStatistics> _stats;
254 #endif // DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIX_HH
C::block_type block_type
Definition: bcrsmatrix.hh:30
GFSV TestGridFunctionSpace
Definition: bcrsmatrix.hh:34
C::size_type size_type
Definition: bcrsmatrix.hh:31
const std::shared_ptr< Container > & storage() const
Definition: bcrsmatrix.hh:165
const Container & base() const
Definition: bcrsmatrix.hh:204
Definition: bcrsmatrix.hh:17
BCRSMatrix(const GO &go)
Definition: bcrsmatrix.hh:62
void clear_row(const RowIndex &ri, const E &diagonal_entry)
Definition: bcrsmatrix.hh:236
tags::container< T >::type container_tag(const T &)
Gets instance of container tag associated with T.
Definition: backend/istl/tags.hh:246
Stats PatternStatistics
Definition: bcrsmatrix.hh:41
GFSV::Ordering::Traits::ContainerIndex RowIndex
Definition: bcrsmatrix.hh:36
Various tags for influencing backend behavior.
C Container
Definition: bcrsmatrix.hh:27
void detach()
Definition: bcrsmatrix.hh:149
const E & e
Definition: interpolate.hh:172
Definition: uncachedmatrixview.hh:158
ElementType E
Definition: bcrsmatrix.hh:26
bool attached() const
Definition: bcrsmatrix.hh:160
BCRSMatrix(Backend::attached_container)
Creates an BCRSMatrix with an empty underlying ISTL matrix.
Definition: bcrsmatrix.hh:97
istl::build_pattern_type< C, GFSV, GFSU, typename GFSV::Ordering::ContainerAllocationTag >::type Pattern
Definition: bcrsmatrix.hh:39
C::field_type field_type
Definition: bcrsmatrix.hh:29
size_type M() const
Definition: bcrsmatrix.hh:175
void attach(std::shared_ptr< Container > container)
Definition: bcrsmatrix.hh:155
const StatisticsReturnType & patternStatistics() const
Returns pattern statistics for all contained BCRSMatrix objects.
Definition: bcrsmatrix.hh:122
C::field_type ElementType
Definition: bcrsmatrix.hh:25
void finalize()
Definition: bcrsmatrix.hh:233
Tag for requesting a vector or matrix container without a pre-attached underlying object...
Definition: backend/common/tags.hh:23
BCRSMatrix(const GO &go, const E &e)
Definition: bcrsmatrix.hh:85
Tag for requesting a vector or matrix container with a pre-attached underlying object.
Definition: backend/common/tags.hh:27
BCRSMatrix(Backend::unattached_container=Backend::unattached_container())
Creates an BCRSMatrix without allocating an underlying ISTL matrix.
Definition: bcrsmatrix.hh:93
const E & operator()(const RowIndex &ri, const ColIndex &ci) const
Definition: bcrsmatrix.hh:197
Definition: uncachedmatrixview.hh:13
BCRSMatrix & operator=(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:105
Definition: adaptivity.hh:27
BCRSMatrix & operator*=(const E &e)
Definition: bcrsmatrix.hh:186
GFSU TrialGridFunctionSpace
Definition: bcrsmatrix.hh:33
BCRSMatrix(const GO &go, Container &container)
Construct matrix container using an externally given matrix as storage.
Definition: bcrsmatrix.hh:78
GFSU::Ordering::Traits::ContainerIndex ColIndex
Definition: bcrsmatrix.hh:37
C BaseT
Definition: bcrsmatrix.hh:28
E & operator()(const RowIndex &ri, const ColIndex &ci)
Definition: bcrsmatrix.hh:192
void flush()
Definition: bcrsmatrix.hh:230
BCRSMatrix(const BCRSMatrix &rhs)
Definition: bcrsmatrix.hh:101
size_type N() const
Definition: bcrsmatrix.hh:170