Reference documentation for deal.II version GIT d6cf33b98c 2023-09-22 19:50:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
dof_info.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2023 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.md at
12 // the top level directory of deal.II.
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/config.h>
22 
25 
28 
29 #include <array>
30 #include <memory>
31 
32 
34 
35 #ifndef DOXYGEN
36 
37 // forward declarations
38 
39 namespace internal
40 {
41  namespace MatrixFreeFunctions
42  {
43  template <int dim>
44  class HangingNodes;
45 
46  struct TaskInfo;
47 
48  template <typename Number>
49  struct ConstraintValues;
50 
51  namespace VectorDataExchange
52  {
53  class Base;
54  }
55  } // namespace MatrixFreeFunctions
56 } // namespace internal
57 
58 template <typename>
59 class AffineConstraints;
60 
62 
63 template <typename>
64 class TriaIterator;
65 
66 template <int, int, bool>
67 class DoFCellAccessor;
68 
69 namespace Utilities
70 {
71  namespace MPI
72  {
73  class Partitioner;
74  }
75 } // namespace Utilities
76 
77 #endif
78 
79 namespace internal
80 {
81  namespace MatrixFreeFunctions
82  {
87  using compressed_constraint_kind = std::uint8_t;
88 
107  struct DoFInfo
108  {
121  static const unsigned int chunk_size_zero_vector = 64;
122 
126  DoFInfo();
127 
131  DoFInfo(const DoFInfo &) = default;
132 
136  DoFInfo(DoFInfo &&) noexcept = default;
137 
141  ~DoFInfo() = default;
142 
146  DoFInfo &
147  operator=(const DoFInfo &) = default;
148 
152  DoFInfo &
153  operator=(DoFInfo &&) noexcept = default;
154 
158  void
159  clear();
160 
166  unsigned int
167  fe_index_from_degree(const unsigned int first_selected_component,
168  const unsigned int fe_degree) const;
169 
189  void
190  get_dof_indices_on_cell_batch(std::vector<unsigned int> &local_indices,
191  const unsigned int cell_batch,
192  const bool with_constraints = true) const;
193 
203  template <typename number>
204  void
206  const std::vector<types::global_dof_index> &local_indices_resolved,
207  const std::vector<types::global_dof_index> &local_indices,
208  const bool cell_has_hanging_nodes,
209  const ::AffineConstraints<number> &constraints,
210  const unsigned int cell_number,
211  ConstraintValues<double> &constraint_values,
212  bool &cell_at_boundary);
213 
218  template <int dim>
219  bool
221  const HangingNodes<dim> &hanging_nodes,
222  const std::vector<std::vector<unsigned int>> &lexicographic_mapping,
223  const unsigned int cell_number,
224  const TriaIterator<DoFCellAccessor<dim, dim, false>> &cell,
225  std::vector<types::global_dof_index> &dof_indices);
226 
234  void
235  assign_ghosts(const std::vector<unsigned int> &boundary_cells,
236  const MPI_Comm communicator_sm,
237  const bool use_vector_data_exchanger_full);
238 
245  void
246  reorder_cells(const TaskInfo &task_info,
247  const std::vector<unsigned int> &renumbering,
248  const std::vector<unsigned int> &constraint_pool_row_index,
249  const std::vector<unsigned char> &irregular_cells);
250 
255  void
257  const std::vector<unsigned char> &irregular_cells);
258 
263  template <int length>
264  void
266  const std::vector<FaceToCellTopology<length>> &faces);
267 
273  void
274  make_connectivity_graph(const TaskInfo &task_info,
275  const std::vector<unsigned int> &renumbering,
276  DynamicSparsityPattern &connectivity) const;
277 
283  void
285  const Table<2, ShapeInfo<double>> &shape_info,
286  const unsigned int n_owned_cells,
287  const unsigned int n_lanes,
288  const std::vector<FaceToCellTopology<1>> &inner_faces,
289  const std::vector<FaceToCellTopology<1>> &ghosted_faces,
290  const bool fill_cell_centric,
291  const MPI_Comm communicator_sm,
292  const bool use_vector_data_exchanger_full);
293 
299  void
301  std::array<std::vector<std::pair<unsigned int, unsigned int>>, 3>
302  &cell_indices_contiguous_sm);
303 
314  void
316  std::vector<types::global_dof_index> &renumbering);
317 
327  template <int length>
328  void
330  const TaskInfo &task_info,
331  const std::vector<FaceToCellTopology<length>> &faces);
332 
336  std::size_t
337  memory_consumption() const;
338 
343  template <typename StreamType>
344  void
345  print_memory_consumption(StreamType &out,
346  const TaskInfo &size_info) const;
347 
352  template <typename Number>
353  void
354  print(const std::vector<Number> &constraint_pool_data,
355  const std::vector<unsigned int> &constraint_pool_row_index,
356  std::ostream &out) const;
357 
369  enum class IndexStorageVariants : unsigned char
370  {
379  full,
391  interleaved,
400  contiguous,
424  interleaved_contiguous,
436  interleaved_contiguous_strided,
449  interleaved_contiguous_mixed_strides
450  };
451 
456  enum DoFAccessIndex : unsigned char
457  {
469  dof_access_cell = 2
470  };
471 
480  unsigned int vectorization_length;
481 
489  std::array<std::vector<IndexStorageVariants>, 3> index_storage_variants;
490 
497  std::vector<std::pair<unsigned int, unsigned int>> row_starts;
498 
514  std::vector<unsigned int> dof_indices;
515 
520  std::vector<std::vector<bool>> hanging_node_constraint_masks_comp;
521 
526  std::vector<compressed_constraint_kind> hanging_node_constraint_masks;
527 
537  std::vector<std::pair<unsigned short, unsigned short>>
539 
543  std::vector<unsigned int> dof_indices_interleaved;
544 
553  std::array<std::vector<unsigned int>, 3> dof_indices_contiguous;
554 
563  std::array<std::vector<std::pair<unsigned int, unsigned int>>, 3>
565 
574  std::array<std::vector<unsigned int>, 3> dof_indices_interleave_strides;
575 
585  std::array<std::vector<unsigned char>, 3> n_vectorization_lanes_filled;
586 
593  std::shared_ptr<const Utilities::MPI::Partitioner> vector_partitioner;
594 
598  std::shared_ptr<
601 
621  std::array<
622  std::shared_ptr<
624  5>
626 
631  std::vector<unsigned int> constrained_dofs;
632 
637  std::vector<unsigned int> row_starts_plain_indices;
638 
647  std::vector<unsigned int> plain_dof_indices;
648 
654 
659  unsigned int n_base_elements;
660 
665  std::vector<unsigned int> n_components;
666 
671  std::vector<unsigned int> start_components;
672 
677  std::vector<unsigned int> component_to_base_index;
678 
690  std::vector<std::vector<unsigned int>> component_dof_indices_offset;
691 
695  std::vector<unsigned int> dofs_per_cell;
696 
700  std::vector<unsigned int> dofs_per_face;
701 
706 
710  std::vector<unsigned int> cell_active_fe_index;
711 
716  unsigned int max_fe_index;
717 
723  std::vector<std::vector<unsigned int>> fe_index_conversion;
724 
730  std::vector<types::global_dof_index> ghost_dofs;
731 
737  std::vector<unsigned int> vector_zero_range_list_index;
738 
742  std::vector<std::pair<unsigned int, unsigned int>> vector_zero_range_list;
743 
749  std::vector<unsigned int> cell_loop_pre_list_index;
750 
755  std::vector<std::pair<unsigned int, unsigned int>> cell_loop_pre_list;
756 
762  std::vector<unsigned int> cell_loop_post_list_index;
763 
768  std::vector<std::pair<unsigned int, unsigned int>> cell_loop_post_list;
769  };
770 
771 
772  /*-------------------------- Inline functions ---------------------------*/
773 
774 #ifndef DOXYGEN
775 
776 
777  inline unsigned int
778  DoFInfo::fe_index_from_degree(const unsigned int first_selected_component,
779  const unsigned int fe_degree) const
780  {
781  const unsigned int n_indices = fe_index_conversion.size();
782  if (n_indices <= 1)
783  return 0;
784  for (unsigned int i = 0; i < n_indices; ++i)
785  if (fe_index_conversion[i][first_selected_component] == fe_degree)
786  return i;
788  }
789 
790 #endif // ifndef DOXYGEN
791 
792  } // end of namespace MatrixFreeFunctions
793 } // end of namespace internal
794 
796 
797 #endif
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
std::uint8_t compressed_constraint_kind
Definition: dof_info.h:87
static const unsigned int invalid_unsigned_int
Definition: types.h:213
Definition: types.h:33
unsigned int global_dof_index
Definition: types.h:82
void print(const std::vector< Number > &constraint_pool_data, const std::vector< unsigned int > &constraint_pool_row_index, std::ostream &out) const
void reorder_cells(const TaskInfo &task_info, const std::vector< unsigned int > &renumbering, const std::vector< unsigned int > &constraint_pool_row_index, const std::vector< unsigned char > &irregular_cells)
Definition: dof_info.cc:299
std::vector< std::pair< unsigned short, unsigned short > > constraint_indicator
Definition: dof_info.h:538
void compute_tight_partitioners(const Table< 2, ShapeInfo< double >> &shape_info, const unsigned int n_owned_cells, const unsigned int n_lanes, const std::vector< FaceToCellTopology< 1 >> &inner_faces, const std::vector< FaceToCellTopology< 1 >> &ghosted_faces, const bool fill_cell_centric, const MPI_Comm communicator_sm, const bool use_vector_data_exchanger_full)
Definition: dof_info.cc:896
std::vector< unsigned int > cell_loop_pre_list_index
Definition: dof_info.h:749
std::vector< std::pair< unsigned int, unsigned int > > row_starts
Definition: dof_info.h:497
std::vector< std::vector< unsigned int > > component_dof_indices_offset
Definition: dof_info.h:690
std::vector< std::vector< bool > > hanging_node_constraint_masks_comp
Definition: dof_info.h:520
void get_dof_indices_on_cell_batch(std::vector< unsigned int > &local_indices, const unsigned int cell_batch, const bool with_constraints=true) const
Definition: dof_info.cc:83
unsigned int fe_index_from_degree(const unsigned int first_selected_component, const unsigned int fe_degree) const
void compute_cell_index_compression(const std::vector< unsigned char > &irregular_cells)
Definition: dof_info.cc:531
std::vector< unsigned int > cell_loop_post_list_index
Definition: dof_info.h:762
std::vector< unsigned int > dofs_per_cell
Definition: dof_info.h:695
void compute_shared_memory_contiguous_indices(std::array< std::vector< std::pair< unsigned int, unsigned int >>, 3 > &cell_indices_contiguous_sm)
Definition: dof_info.cc:1242
void assign_ghosts(const std::vector< unsigned int > &boundary_cells, const MPI_Comm communicator_sm, const bool use_vector_data_exchanger_full)
Definition: dof_info.cc:156
void compute_dof_renumbering(std::vector< types::global_dof_index > &renumbering)
Definition: dof_info.cc:1468
static const unsigned int chunk_size_zero_vector
Definition: dof_info.h:121
std::shared_ptr< const internal::MatrixFreeFunctions::VectorDataExchange::Base > vector_exchanger
Definition: dof_info.h:600
std::vector< std::vector< unsigned int > > fe_index_conversion
Definition: dof_info.h:723
std::vector< unsigned int > dof_indices
Definition: dof_info.h:514
std::vector< std::pair< unsigned int, unsigned int > > vector_zero_range_list
Definition: dof_info.h:742
std::shared_ptr< const Utilities::MPI::Partitioner > vector_partitioner
Definition: dof_info.h:593
std::vector< compressed_constraint_kind > hanging_node_constraint_masks
Definition: dof_info.h:526
std::array< std::vector< unsigned int >, 3 > dof_indices_interleave_strides
Definition: dof_info.h:574
std::array< std::vector< std::pair< unsigned int, unsigned int > >, 3 > dof_indices_contiguous_sm
Definition: dof_info.h:564
std::vector< unsigned int > row_starts_plain_indices
Definition: dof_info.h:637
std::vector< unsigned int > dofs_per_face
Definition: dof_info.h:700
std::size_t memory_consumption() const
Definition: dof_info.cc:1526
std::vector< unsigned int > component_to_base_index
Definition: dof_info.h:677
std::vector< unsigned int > n_components
Definition: dof_info.h:665
void print_memory_consumption(StreamType &out, const TaskInfo &size_info) const
std::vector< std::pair< unsigned int, unsigned int > > cell_loop_pre_list
Definition: dof_info.h:755
std::vector< types::global_dof_index > ghost_dofs
Definition: dof_info.h:730
void read_dof_indices(const std::vector< types::global_dof_index > &local_indices_resolved, const std::vector< types::global_dof_index > &local_indices, const bool cell_has_hanging_nodes, const ::AffineConstraints< number > &constraints, const unsigned int cell_number, ConstraintValues< double > &constraint_values, bool &cell_at_boundary)
std::array< std::vector< unsigned int >, 3 > dof_indices_contiguous
Definition: dof_info.h:553
DoFInfo(const DoFInfo &)=default
std::vector< unsigned int > cell_active_fe_index
Definition: dof_info.h:710
std::array< std::shared_ptr< const internal::MatrixFreeFunctions::VectorDataExchange::Base >, 5 > vector_exchanger_face_variants
Definition: dof_info.h:625
std::vector< unsigned int > plain_dof_indices
Definition: dof_info.h:647
std::array< std::vector< unsigned char >, 3 > n_vectorization_lanes_filled
Definition: dof_info.h:585
std::vector< unsigned int > start_components
Definition: dof_info.h:671
std::vector< unsigned int > dof_indices_interleaved
Definition: dof_info.h:543
std::vector< unsigned int > constrained_dofs
Definition: dof_info.h:631
std::array< std::vector< IndexStorageVariants >, 3 > index_storage_variants
Definition: dof_info.h:489
bool process_hanging_node_constraints(const HangingNodes< dim > &hanging_nodes, const std::vector< std::vector< unsigned int >> &lexicographic_mapping, const unsigned int cell_number, const TriaIterator< DoFCellAccessor< dim, dim, false >> &cell, std::vector< types::global_dof_index > &dof_indices)
void make_connectivity_graph(const TaskInfo &task_info, const std::vector< unsigned int > &renumbering, DynamicSparsityPattern &connectivity) const
Definition: dof_info.cc:1391
std::vector< unsigned int > vector_zero_range_list_index
Definition: dof_info.h:737
void compute_face_index_compression(const std::vector< FaceToCellTopology< length >> &faces)
void compute_vector_zero_access_pattern(const TaskInfo &task_info, const std::vector< FaceToCellTopology< length >> &faces)
std::vector< std::pair< unsigned int, unsigned int > > cell_loop_post_list
Definition: dof_info.h:768
DoFInfo(DoFInfo &&) noexcept=default