dune-pdelab  2.5-dev
bcrsmatrixbackend.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil -*-
2 #ifndef DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIXBACKEND_HH
3 #define DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIXBACKEND_HH
4 
5 // this is here for backwards compatibility and deprecation warnings, remove after 2.5.0
6 #include "ensureistlinclude.hh"
7 
11 
12 namespace Dune {
13  namespace PDELab {
14  namespace ISTL {
15 
16  // ********************************************************************************
17  // infrastructure for deducing the pattern type from row and column orderings
18  // ********************************************************************************
19 
20  namespace {
21 
22  template<typename M, typename RowOrdering, typename ColOrdering, bool pattern>
23  struct _build_bcrs_pattern_type
24  {
25  // we use void as a marker for a nonexisting subpattern
26  typedef void type;
27  };
28 
29  // central TMP that is invoked recursively while traversing the ordering hierarchy
30  // and builds up the possibly nested pattern.
31  template<typename M, typename RowOrdering, typename ColOrdering>
32  struct _build_bcrs_pattern_type<M,RowOrdering,ColOrdering,true>
33  {
34 
35  // descend into blocks
36  typedef typename _build_bcrs_pattern_type<
37  typename M::block_type,
38  RowOrdering,
39  ColOrdering,
40  requires_pattern<
41  typename M::block_type
42  >::value
43  >::type BlockOrdering;
44 
45  // Leafs -> BCRSPattern, interior nodes -> NestedPattern
46  typedef typename std::conditional<
48  BCRSPattern<
49  RowOrdering,
50  ColOrdering
51  >,
52  NestedPattern<
53  RowOrdering,
54  ColOrdering,
55  BlockOrdering
56  >
57  >::type type;
58 
59  };
60 
61  // Wrapper TMP for constructing OrderingBase types from function spaces and for
62  // shielding user from recursive implementation
63  template<typename M, typename GFSV, typename GFSU, typename Tag>
64  struct build_bcrs_pattern_type
65  {
66 
67  typedef OrderingBase<
68  typename GFSV::Ordering::Traits::DOFIndex,
69  typename GFSV::Ordering::Traits::ContainerIndex
70  > RowOrdering;
71 
72  typedef OrderingBase<
73  typename GFSU::Ordering::Traits::DOFIndex,
74  typename GFSU::Ordering::Traits::ContainerIndex
75  > ColOrdering;
76 
78  };
79 
80  // Specialization for forcibly flat backends
81  template<typename M, typename GFSV, typename GFSU>
82  struct build_bcrs_pattern_type<M,GFSV,GFSU,FlatContainerAllocationTag>
83  {
84  typedef BCRSPattern<typename GFSV::Ordering, typename GFSU::Ordering> type;
85  };
86 
87 
88  // leaf BCRSMatrix
89  template<typename OrderingV, typename OrderingU, typename Pattern, typename Container, typename StatsVector>
90  typename std::enable_if<
92  >::type
93  allocate_bcrs_matrix(const OrderingV& ordering_v,
94  const OrderingU& ordering_u,
95  Pattern& p,
96  Container& c,
97  StatsVector& stats)
98  {
99  c.setSize(ordering_v.blockCount(),ordering_u.blockCount(),0);
100  c.setBuildMode(Container::random);
101 
102  std::vector<typename Pattern::size_type> row_sizes(p.sizes());
103 
104  typename Pattern::size_type nnz = 0;
105  typename Pattern::size_type longest_row = 0;
106 
107  for (typename Pattern::size_type i = 0; i < c.N(); ++i)
108  {
109  nnz += row_sizes[i];
110  longest_row = std::max(longest_row,row_sizes[i]);
111  c.setrowsize(i,row_sizes[i]);
112  }
113  c.endrowsizes();
114 
115  stats.push_back(typename StatsVector::value_type(nnz,longest_row,p.overflowCount(),p.entriesPerRow(),ordering_v.blockCount()));
116 
117  for (typename Pattern::size_type i = 0; i < c.N(); ++i)
118  c.setIndices(i,p.begin(i),p.end(i));
119 
120  // free temporary index storage in pattern before allocating data array in matrix
121  p.clear();
122  // allocate data array
123  c.endindices();
124  }
125 
126 
127  // ********************************************************************************
128  // nested matrix allocation
129  // In contrast to the older implementation, this code still uses BCRSMatrix for nested matrices,
130  // but we do not attempt to keep the pattern of those interior matrices sparse and always allocate all
131  // blocks. That greatly simplifies the code, and as those interior matrices really shouldn't be very
132  // large, any performance impact is minimal.
133  // The code also collects statistics about the pattern of the leaf BCRSMatrices and collates those stats
134  // in a row-major ordering.
135  // ********************************************************************************
136 
137  // interior BCRSMatrix
138  template<typename OrderingV, typename OrderingU, typename Pattern, typename Container, typename StatsVector>
139  typename std::enable_if<
142  >::type
143  allocate_bcrs_matrix(const OrderingV& ordering_v,
144  const OrderingU& ordering_u,
145  Pattern& p,
146  Container& c,
147  StatsVector& stats)
148  {
149  c.setSize(ordering_v.blockCount(),ordering_u.blockCount(),ordering_v.blockCount()*ordering_u.blockCount());
150  c.setBuildMode(Container::random);
151 
152  for (std::size_t i = 0; i < c.N(); ++i)
153  c.setrowsize(i,ordering_u.blockCount());
154  c.endrowsizes();
155 
156  for (std::size_t i = 0; i < c.N(); ++i)
157  for (std::size_t j = 0; j < c.M(); ++j)
158  c.addindex(i,j);
159  c.endindices();
160 
161  for (std::size_t i = 0; i < c.N(); ++i)
162  for (std::size_t j = 0; j < c.M(); ++j)
163  {
164  allocate_bcrs_matrix(ordering_v.childOrdering(i),
165  ordering_u.childOrdering(j),
166  p.subPattern(i,j),
167  c[i][j],
168  stats);
169  }
170  }
171 
172  } // anonymous namespace
173 
174 
175 
177 
189  template<typename EntriesPerRow = std::size_t>
191  {
192 
194  typedef std::size_t size_type;
195 
198 
200  template<typename Matrix, typename GFSV, typename GFSU>
201  using Pattern = typename build_bcrs_pattern_type<
202  typename Matrix::Container,
203  GFSV,
204  GFSU,
205  typename GFSV::Ordering::ContainerAllocationTag
206  >::type;
207 
208  template<typename VV, typename VU, typename E>
210  {
211  typedef BCRSMatrix<
212  typename VV::GridFunctionSpace,
213  typename VU::GridFunctionSpace,
214  typename build_matrix_type<
215  E,
216  typename VV::Container,
217  typename VU::Container
218  >::type,
219  Statistics
220  > type;
221  };
222 
224 
227  template<typename GridOperator, typename Matrix>
228  std::vector<Statistics> buildPattern(const GridOperator& grid_operator, Matrix& matrix) const
229  {
230  Pattern<
231  Matrix,
234  > pattern(grid_operator.testGridFunctionSpace().ordering(),grid_operator.trialGridFunctionSpace().ordering(),_entries_per_row);
235  grid_operator.fill_pattern(pattern);
236  std::vector<Statistics> stats;
237  allocate_bcrs_matrix(grid_operator.testGridFunctionSpace().ordering(),
238  grid_operator.trialGridFunctionSpace().ordering(),
239  pattern,
240  Backend::native(matrix),
241  stats
242  );
243  return stats;
244  }
245 
247 
252  BCRSMatrixBackend(const EntriesPerRow& entries_per_row)
253  : _entries_per_row(entries_per_row)
254  {}
255 
256  private:
257 
258  EntriesPerRow _entries_per_row;
259 
260  };
261 
262  } // namespace ISTL
263  } // namespace PDELab
264 } // namespace Dune
265 
266 #endif // DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIXBACKEND_HH
PatternStatistics< size_type > Statistics
The type of the object holding the statistics generated during pattern construction.
Definition: bcrsmatrixbackend.hh:197
Definition: bcrsmatrix.hh:20
Definition: bcrsmatrixbackend.hh:209
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
BCRSMatrixBackend(const EntriesPerRow &entries_per_row)
Constructs a BCRSMatrixBackend.
Definition: bcrsmatrixbackend.hh:252
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:120
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:191
Standard grid operator implementation.
Definition: gridoperator.hh:57
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:126
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:127
std::size_t size_type
The size type of the BCRSMatrix.
Definition: bcrsmatrixbackend.hh:194
std::vector< Statistics > buildPattern(const GridOperator &grid_operator, Matrix &matrix) const
Builds the matrix pattern associated with grid_operator and initializes matrix with it...
Definition: bcrsmatrixbackend.hh:228
const P & p
Definition: constraints.hh:147
Backend using (possibly nested) ISTL BCRSMatrices.
Definition: bcrsmatrixbackend.hh:190
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > &>::type native(T &t)
Definition: backend/interface.hh:192
typename build_bcrs_pattern_type< typename Matrix::Container, GFSV, GFSU, typename GFSV::Ordering::ContainerAllocationTag >::type Pattern
The type of the pattern object passed to the GridOperator for pattern construction.
Definition: bcrsmatrixbackend.hh:206
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
BCRSMatrix< typename VV::GridFunctionSpace, typename VU::GridFunctionSpace, typename build_matrix_type< E, typename VV::Container, typename VU::Container >::type, Statistics > type
Definition: bcrsmatrixbackend.hh:220
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
Statistics about the pattern of a BCRSMatrix.
Definition: patternstatistics.hh:16