dune-pdelab  2.0.0
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 
8 
9 namespace Dune {
10  namespace PDELab {
11  namespace istl {
12 
13  // ********************************************************************************
14  // infrastructure for deducing the pattern type from row and column orderings
15  // ********************************************************************************
16 
17  namespace {
18 
19  template<typename M, typename RowOrdering, typename ColOrdering, bool pattern>
20  struct _build_bcrs_pattern_type
21  {
22  // we use void as a marker for a nonexisting subpattern
23  typedef void type;
24  };
25 
26  // central TMP that is invoked recursively while traversing the ordering hierarchy
27  // and builds up the possibly nested pattern.
28  template<typename M, typename RowOrdering, typename ColOrdering>
29  struct _build_bcrs_pattern_type<M,RowOrdering,ColOrdering,true>
30  {
31 
32  // descend into blocks
33  typedef typename _build_bcrs_pattern_type<
34  typename M::block_type,
35  RowOrdering,
36  ColOrdering,
37  requires_pattern<
38  typename M::block_type
39  >::value
40  >::type BlockOrdering;
41 
42  // Leafs -> BCRSPattern, interior nodes -> NestedPattern
43  typedef typename conditional<
46  RowOrdering,
47  ColOrdering
48  >,
50  RowOrdering,
51  ColOrdering,
52  BlockOrdering
53  >
54  >::type type;
55 
56  };
57 
58  // Wrapper TMP for constructing OrderingBase types from function spaces and for
59  // shielding user from recursive implementation
60  template<typename M, typename GFSV, typename GFSU, typename Tag>
61  struct build_bcrs_pattern_type
62  {
63 
64  typedef OrderingBase<
65  typename GFSV::Ordering::Traits::DOFIndex,
66  typename GFSV::Ordering::Traits::ContainerIndex
67  > RowOrdering;
68 
69  typedef OrderingBase<
70  typename GFSU::Ordering::Traits::DOFIndex,
71  typename GFSU::Ordering::Traits::ContainerIndex
72  > ColOrdering;
73 
75  };
76 
77  // Specialization for forcibly flat backends
78  template<typename M, typename GFSV, typename GFSU>
79  struct build_bcrs_pattern_type<M,GFSV,GFSU,FlatContainerAllocationTag>
80  {
82  };
83 
84 
85 
86  // ********************************************************************************
87  // nested matrix allocation
88  // In contrast to the older implementation, this code still uses BCRSMatrix for nested matrices,
89  // but we do not attempt to keep the pattern of those interior matrices sparse and always allocate all
90  // blocks. That greatly simplifies the code, and as those interior matrices really shouldn't be very
91  // large, any performance impact is minimal.
92  // The code also collects statistics about the pattern of the leaf BCRSMatrices and collates those stats
93  // in a row-major ordering.
94  // ********************************************************************************
95 
96  // interior BCRSMatrix
97  template<typename OrderingV, typename OrderingU, typename Pattern, typename Container, typename StatsVector>
98  typename enable_if<
101  >::type
102  allocate_bcrs_matrix(const OrderingV& ordering_v,
103  const OrderingU& ordering_u,
104  Pattern& p,
105  Container& c,
106  StatsVector& stats)
107  {
108  c.setSize(ordering_v.blockCount(),ordering_u.blockCount(),ordering_v.blockCount()*ordering_u.blockCount());
109  c.setBuildMode(Container::random);
110 
111  for (std::size_t i = 0; i < c.N(); ++i)
112  c.setrowsize(i,ordering_u.blockCount());
113  c.endrowsizes();
114 
115  for (std::size_t i = 0; i < c.N(); ++i)
116  for (std::size_t j = 0; j < c.M(); ++j)
117  c.addindex(i,j);
118  c.endindices();
119 
120  for (std::size_t i = 0; i < c.N(); ++i)
121  for (std::size_t j = 0; j < c.M(); ++j)
122  {
123  allocate_matrix(ordering_v.childOrdering(i),
124  ordering_u.childOrdering(j),
125  p.subPattern(i,j),
126  c[i][j],
127  stats);
128  }
129  }
130 
131  // leaf BCRSMatrix
132  template<typename OrderingV, typename OrderingU, typename Pattern, typename Container, typename StatsVector>
133  typename enable_if<
135  >::type
136  allocate_bcrs_matrix(const OrderingV& ordering_v,
137  const OrderingU& ordering_u,
138  Pattern& p,
139  Container& c,
140  StatsVector& stats)
141  {
142  c.setSize(ordering_v.blockCount(),ordering_u.blockCount(),0);
143  c.setBuildMode(Container::random);
144 
145  std::vector<typename Pattern::size_type> row_sizes(p.sizes());
146 
147  typename Pattern::size_type nnz = 0;
148  typename Pattern::size_type longest_row = 0;
149 
150  for (typename Pattern::size_type i = 0; i < c.N(); ++i)
151  {
152  nnz += row_sizes[i];
153  longest_row = std::max(longest_row,row_sizes[i]);
154  c.setrowsize(i,row_sizes[i]);
155  }
156  c.endrowsizes();
157 
158  stats.push_back(typename StatsVector::value_type(nnz,longest_row,p.overflowCount(),p.entriesPerRow(),ordering_v.blockCount()));
159 
160  for (typename Pattern::size_type i = 0; i < c.N(); ++i)
161  c.setIndices(i,p.begin(i),p.end(i));
162 
163  // free temporary index storage in pattern before allocating data array in matrix
164  p.clear();
165  // allocate data array
166  c.endindices();
167  }
168 
169  } // anonymous namespace
170 
171 
172 
174 
186  template<typename EntriesPerRow = std::size_t>
188  {
189 
191  typedef std::size_t size_type;
192 
195 
196 #if HAVE_TEMPLATE_ALIASES || DOXYGEN
197 
199  template<typename Matrix, typename GFSV, typename GFSU>
200  using Pattern = typename build_bcrs_pattern_type<
201  typename Matrix::Container,
202  GFSV,
203  GFSU,
204  typename GFSV::Ordering::ContainerAllocationTag
205  >::type;
206 
207 #else // HAVE_TEMPLATE_ALIASES
208 
209  template<typename Matrix, typename GFSV, typename GFSU>
210  struct Pattern
211  : public build_bcrs_pattern_type<typename Matrix::Container,
212  GFSV,
213  GFSU,
214  typename GFSV::Ordering::ContainerAllocationTag
215  >::type
216  {
217 
218  typedef OrderingBase<
219  typename GFSV::Ordering::Traits::DOFIndex,
220  typename GFSV::Ordering::Traits::ContainerIndex
221  > RowOrdering;
222 
223  typedef OrderingBase<
224  typename GFSU::Ordering::Traits::DOFIndex,
225  typename GFSU::Ordering::Traits::ContainerIndex
226  > ColOrdering;
227 
228  typedef typename build_bcrs_pattern_type<
229  typename Matrix::Container,
230  GFSV,
231  GFSU,
232  typename GFSV::Ordering::ContainerAllocationTag
233  >::type BaseT;
234 
235  template<typename EntriesPerRow_>
236  Pattern(const RowOrdering& row_ordering, const ColOrdering& col_ordering, const EntriesPerRow_& entries_per_row)
237  : BaseT(row_ordering,col_ordering,entries_per_row)
238  {}
239 
240  };
241 
242 #endif // HAVE_TEMPLATE_ALIASES
243 
244  template<typename VV, typename VU, typename E>
246  {
247  typedef ISTLMatrixContainer<
248  typename VV::GridFunctionSpace,
249  typename VU::GridFunctionSpace,
250  typename build_matrix_type<
251  E,
252  typename VV::Container,
253  typename VU::Container
254  >::type,
255  Statistics
256  > type;
257  };
258 
260 
263  template<typename GridOperator, typename Matrix>
264  std::vector<Statistics> buildPattern(const GridOperator& grid_operator, Matrix& matrix) const
265  {
266  Pattern<
267  Matrix,
270  > pattern(grid_operator.testGridFunctionSpace().ordering(),grid_operator.trialGridFunctionSpace().ordering(),_entries_per_row);
271  grid_operator.fill_pattern(pattern);
272  std::vector<Statistics> stats;
273  allocate_bcrs_matrix(grid_operator.testGridFunctionSpace().ordering(),
274  grid_operator.trialGridFunctionSpace().ordering(),
275  pattern,
276  istl::raw(matrix),
277  stats
278  );
279 #if HAVE_RVALUE_REFERENCES
280  return std::move(stats);
281 #else
282  return stats;
283 #endif
284  }
285 
287 
292  BCRSMatrixBackend(const EntriesPerRow& entries_per_row)
293  : _entries_per_row(entries_per_row)
294  {}
295 
296  private:
297 
298  EntriesPerRow _entries_per_row;
299 
300  };
301 
302  } // namespace istl
303  } // namespace PDELab
304 } // namespace Dune
305 
306 #endif // DUNE_PDELAB_BACKEND_ISTL_BCRSMATRIXBACKEND_HH
ISTLMatrixContainer< typename VV::GridFunctionSpace, typename VU::GridFunctionSpace, typename build_matrix_type< E, typename VV::Container, typename VU::Container >::type, Statistics > type
Definition: bcrsmatrixbackend.hh:256
GFSU TrialGridFunctionSpace
The trial grid function space.
Definition: gridoperatorutilities.hh:37
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: gridoperator.hh:158
Pattern builder for generic BCRS-like sparse matrices.
Definition: bcrspattern.hh:43
Pattern builder for nested hierarchies of generic BCRS-like sparse matrices.
Definition: bcrspattern.hh:328
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: gridoperator.hh:87
Definition: orderingbase.hh:22
Backend using (possibly nested) ISTL BCRSMatrices.
Definition: bcrsmatrixbackend.hh:187
Definition: bcrsmatrixbackend.hh:245
Standard grid operator implementation.
Definition: gridoperator.hh:34
V & raw(V &v)
Returns the raw ISTL object associated with v, or v itself it is already an ISTL object.
Definition: backend/istl/utility.hh:26
PatternStatistics< size_type > Statistics
The type of the object holding the statistics generated during pattern construction.
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:264
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: gridoperator.hh:93
Statistics about the pattern of a BCRSMatrix.
Definition: patternstatistics.hh:13
GFSV TestGridFunctionSpace
The test grid function space.
Definition: gridoperatorutilities.hh:40
R p(int i, D x)
Definition: qkdg.hh:62
Definition: istlmatrixbackend.hh:15
std::size_t size_type
The size type of the BCRSMatrix.
Definition: bcrsmatrixbackend.hh:191
BCRSMatrixBackend(const EntriesPerRow &entries_per_row)
Constructs a BCRSMatrixBackend.
Definition: bcrsmatrixbackend.hh:292
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:205