Reference documentation for deal.II version 8.5.1
All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
dof_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2016 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii__matrix_free_dof_info_h
18 #define dealii__matrix_free_dof_info_h
19 
20 
21 #include <deal.II/base/exceptions.h>
22 #include <deal.II/base/vectorization.h>
23 #include <deal.II/base/partitioner.h>
24 #include <deal.II/lac/constraint_matrix.h>
25 #include <deal.II/lac/dynamic_sparsity_pattern.h>
26 #include <deal.II/dofs/dof_handler.h>
27 #include <deal.II/matrix_free/helper_functions.h>
28 
29 #include <deal.II/base/std_cxx11/array.h>
30 
31 #include <memory>
32 
33 
34 DEAL_II_NAMESPACE_OPEN
35 
36 namespace internal
37 {
38  namespace MatrixFreeFunctions
39  {
40  template <typename Number> struct ConstraintValues;
41 
60  struct DoFInfo
61  {
65  DoFInfo ();
66 
70  DoFInfo (const DoFInfo &dof_info);
71 
75  void clear ();
76 
77 
81  const unsigned int *begin_indices (const unsigned int row) const;
82 
87  const unsigned int *end_indices (const unsigned int row) const;
88 
92  unsigned int row_length_indices (const unsigned int row) const;
93 
98  const std::pair<unsigned short,unsigned short> *
99  begin_indicators (const unsigned int row) const;
100 
105  const std::pair<unsigned short,unsigned short> *
106  end_indicators (const unsigned int row) const;
107 
112  unsigned int row_length_indicators (const unsigned int row) const;
113 
118  const unsigned int *begin_indices_plain (const unsigned int row) const;
119 
124  const unsigned int *end_indices_plain (const unsigned int row) const;
125 
132  unsigned int fe_index_from_degree (const unsigned int fe_degree) const;
133 
134 
141  unsigned int
142  fe_index_from_dofs_per_cell (const unsigned int dofs_per_cell) const;
143 
152  void read_dof_indices (const std::vector<types::global_dof_index> &local_indices,
153  const std::vector<unsigned int> &lexicographic_inv,
154  const ConstraintMatrix &constraints,
155  const unsigned int cell_number,
156  ConstraintValues<double> &constraint_values,
157  bool &cell_at_boundary);
158 
166  void assign_ghosts(const std::vector<unsigned int> &boundary_cells);
167 
174  void compute_renumber_serial (const std::vector<unsigned int> &boundary_cells,
175  const SizeInfo &size_info,
176  std::vector<unsigned int> &renumbering);
177 
184  void compute_renumber_hp_serial (SizeInfo &size_info,
185  std::vector<unsigned int> &renumbering,
186  std::vector<unsigned int> &irregular_cells);
187 
193  void compute_renumber_parallel (const std::vector<unsigned int> &boundary_cells,
194  SizeInfo &size_info,
195  std::vector<unsigned int> &renumbering);
196 
203  void reorder_cells (const SizeInfo &size_info,
204  const std::vector<unsigned int> &renumbering,
205  const std::vector<unsigned int> &constraint_pool_row_index,
206  const std::vector<unsigned int> &irregular_cells,
207  const unsigned int vectorization_length);
208 
215  void guess_block_size (const SizeInfo &size_info,
216  TaskInfo &task_info);
217 
231  void
232  make_thread_graph_partition_color (SizeInfo &size_info,
233  TaskInfo &task_info,
234  std::vector<unsigned int> &renumbering,
235  std::vector<unsigned int> &irregular_cells,
236  const bool hp_bool);
237 
252  void
253  make_thread_graph_partition_partition (SizeInfo &size_info,
254  TaskInfo &task_info,
255  std::vector<unsigned int> &renumbering,
256  std::vector<unsigned int> &irregular_cells,
257  const bool hp_bool);
258 
265  void
266  make_connectivity_graph (const SizeInfo &size_info,
267  const TaskInfo &task_info,
268  const std::vector<unsigned int> &renumbering,
269  const std::vector<unsigned int> &irregular_cells,
270  const bool do_blocking,
271  DynamicSparsityPattern &connectivity) const;
272 
276  void renumber_dofs (std::vector<types::global_dof_index> &renumbering);
277 
281  std::size_t memory_consumption() const;
282 
287  template <typename StreamType>
288  void print_memory_consumption(StreamType &out,
289  const SizeInfo &size_info) const;
290 
295  template <typename Number>
296  void print (const std::vector<Number> &constraint_pool_data,
297  const std::vector<unsigned int> &constraint_pool_row_index,
298  std::ostream &out) const;
299 
310  std::vector<std_cxx11::array<unsigned int, 3> > row_starts;
311 
327  std::vector<unsigned int> dof_indices;
328 
338  std::vector<std::pair<unsigned short,unsigned short> > constraint_indicator;
339 
346  std_cxx11::shared_ptr<const Utilities::MPI::Partitioner> vector_partitioner;
347 
352  std::vector<unsigned int> constrained_dofs;
353 
358  std::vector<unsigned int> row_starts_plain_indices;
359 
368  std::vector<unsigned int> plain_dof_indices;
369 
375  unsigned int dimension;
376 
381  unsigned int n_components;
382 
386  std::vector<unsigned int> dofs_per_cell;
387 
391  std::vector<unsigned int> dofs_per_face;
392 
396  bool store_plain_indices;
397 
401  std::vector<unsigned int> cell_active_fe_index;
402 
407  unsigned int max_fe_index;
408 
414  std::vector<std::pair<unsigned int,unsigned int> > fe_index_conversion;
415 
421  std::vector<types::global_dof_index> ghost_dofs;
422  };
423 
424 
425  /*----------------------- Inline functions ----------------------------------*/
426 
427 #ifndef DOXYGEN
428 
429  inline
430  const unsigned int *
431  DoFInfo::begin_indices (const unsigned int row) const
432  {
433  AssertIndexRange (row, row_starts.size()-1);
434  const unsigned int index = row_starts[row][0];
435  AssertIndexRange(index, dof_indices.size()+1);
436  return dof_indices.empty() ? 0 : &dof_indices[0] + index;
437  }
438 
439 
440 
441  inline
442  const unsigned int *
443  DoFInfo::end_indices (const unsigned int row) const
444  {
445  AssertIndexRange (row, row_starts.size()-1);
446  const unsigned int index = row_starts[row+1][0];
447  AssertIndexRange(index, dof_indices.size()+1);
448  return dof_indices.empty() ? 0 : &dof_indices[0] + index;
449  }
450 
451 
452 
453  inline
454  unsigned int
455  DoFInfo::row_length_indices (const unsigned int row) const
456  {
457  AssertIndexRange (row, row_starts.size()-1);
458  return (row_starts[row+1][0] - row_starts[row][0]);
459  }
460 
461 
462 
463  inline
464  const std::pair<unsigned short,unsigned short> *
465  DoFInfo::begin_indicators (const unsigned int row) const
466  {
467  AssertIndexRange (row, row_starts.size()-1);
468  const unsigned int index = row_starts[row][1];
469  AssertIndexRange (index, constraint_indicator.size()+1);
470  return constraint_indicator.empty() ? 0 : &constraint_indicator[0] + index;
471  }
472 
473 
474 
475  inline
476  const std::pair<unsigned short,unsigned short> *
477  DoFInfo::end_indicators (const unsigned int row) const
478  {
479  AssertIndexRange (row, row_starts.size()-1);
480  const unsigned int index = row_starts[row+1][1];
481  AssertIndexRange (index, constraint_indicator.size()+1);
482  return constraint_indicator.empty() ? 0 : &constraint_indicator[0] + index;
483  }
484 
485 
486 
487  inline
488  unsigned int
489  DoFInfo::row_length_indicators (const unsigned int row) const
490  {
491  AssertIndexRange (row, row_starts.size()-1);
492  return (row_starts[row+1][1] - row_starts[row][1]);
493  }
494 
495 
496 
497  inline
498  const unsigned int *
499  DoFInfo::begin_indices_plain (const unsigned int row) const
500  {
501  // if we have no constraints, should take the data from dof_indices
502  if (row_length_indicators(row) == 0)
503  {
504  Assert (row_starts_plain_indices[row]==numbers::invalid_unsigned_int,
505  ExcInternalError());
506  return begin_indices(row);
507  }
508  else
509  {
510  AssertDimension (row_starts.size(), row_starts_plain_indices.size());
511  const unsigned int index = row_starts_plain_indices[row];
512  AssertIndexRange(index, plain_dof_indices.size()+1);
513  return plain_dof_indices.empty() ? 0 : &plain_dof_indices[0] + index;
514  }
515  }
516 
517 
518 
519  inline
520  const unsigned int *
521  DoFInfo::end_indices_plain (const unsigned int row) const
522  {
523  return begin_indices_plain(row) +
524  dofs_per_cell[(cell_active_fe_index.size()==0)?
525  0:cell_active_fe_index[row]];
526  }
527 
528 
529 
530  inline
531  unsigned int
532  DoFInfo::fe_index_from_degree (const unsigned int fe_degree) const
533  {
534  const unsigned int n_indices = fe_index_conversion.size();
535  for (unsigned int i=0; i<n_indices; ++i)
536  if (fe_index_conversion[i].first == fe_degree)
537  return i;
538  return n_indices;
539  }
540 
541 
542 
543  inline
544  unsigned int
545  DoFInfo::fe_index_from_dofs_per_cell (const unsigned int dofs_per_cell) const
546  {
547  for (unsigned int i=0; i<fe_index_conversion.size(); ++i)
548  if (fe_index_conversion[i].second == dofs_per_cell)
549  return i;
550  return 0;
551  }
552 
553  } // end of namespace MatrixFreeFunctions
554 } // end of namespace internal
555 
556 #endif // ifndef DOXYGEN
557 
558 DEAL_II_NAMESPACE_CLOSE
559 
560 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:170
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1146
#define AssertIndexRange(index, range)
Definition: exceptions.h:1170
#define Assert(cond, exc)
Definition: exceptions.h:313
std_cxx11::enable_if< std_cxx11::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static ::ExceptionBase & ExcInternalError()