Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+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\}}\)
mg_constrained_dofs.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2010 - 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 #ifndef dealii_mg_constrained_dofs_h
17 #define dealii_mg_constrained_dofs_h
18 
19 #include <deal.II/base/config.h>
20 
23 
25 
26 #include <set>
27 #include <vector>
28 
30 
31 // Forward declaration
32 #ifndef DOXYGEN
33 template <int dim, int spacedim>
34 DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
35 class DoFHandler;
36 #endif
37 
38 
46 {
47 public:
48  using size_dof = std::vector<std::set<types::global_dof_index>>::size_type;
74  template <int dim, int spacedim>
75  void
77  const MGLevelObject<IndexSet> &level_relevant_dofs =
79 
90  template <int dim, int spacedim>
91  void
93  const DoFHandler<dim, spacedim> &dof,
94  const std::set<types::boundary_id> &boundary_ids,
95  const ComponentMask &component_mask = {});
96 
103  template <int dim, int spacedim>
104  void
106  const unsigned int level,
107  const IndexSet &boundary_indices);
108 
123  void
124  add_user_constraints(const unsigned int level,
125  const AffineConstraints<double> &constraints_on_level);
126 
141  template <int dim, int spacedim>
142  void
144  const types::boundary_id bid,
145  const unsigned int first_vector_component);
146 
150  void
152 
156  void
157  clear();
158 
162  bool
163  is_boundary_index(const unsigned int level,
164  const types::global_dof_index index) const;
165 
169  bool
170  at_refinement_edge(const unsigned int level,
171  const types::global_dof_index index) const;
172 
173 
180  bool
181  is_interface_matrix_entry(const unsigned int level,
182  const types::global_dof_index i,
183  const types::global_dof_index j) const;
184 
191  const IndexSet &
192  get_boundary_indices(const unsigned int level) const;
193 
194 
199  const IndexSet &
200  get_refinement_edge_indices(unsigned int level) const;
201 
202 
206  bool
207  have_boundary_indices() const;
208 
214  get_level_constraints(const unsigned int level) const;
215 
224  get_user_constraint_matrix(const unsigned int level) const;
225 
238  template <typename Number>
239  void
241  const unsigned int level,
242  const bool add_boundary_indices,
243  const bool add_refinement_edge_indices,
244  const bool add_level_constraints,
245  const bool add_user_constraints) const;
246 
247 private:
251  std::vector<IndexSet> boundary_indices;
252 
257  std::vector<IndexSet> refinement_edge_indices;
258 
263  std::vector<AffineConstraints<double>> level_constraints;
264 
268  std::vector<AffineConstraints<double>> user_constraints;
269 };
270 
271 
272 
273 inline bool
275  const types::global_dof_index index) const
276 {
277  if (boundary_indices.empty())
278  return false;
279 
281  return boundary_indices[level].is_element(index);
282 }
283 
284 
285 
286 inline bool
288  const types::global_dof_index index) const
289 {
291 
292  return refinement_edge_indices[level].is_element(index);
293 }
294 
295 
296 
297 inline bool
299  const unsigned int level,
300  const types::global_dof_index i,
301  const types::global_dof_index j) const
302 {
303  const IndexSet &interface_dofs_on_level =
304  this->get_refinement_edge_indices(level);
305 
306  return interface_dofs_on_level.is_element(i) // at_refinement_edge(i)
307  && !interface_dofs_on_level.is_element(j) // !at_refinement_edge(j)
308  && !this->is_boundary_index(level, i) // !on_boundary(i)
309  && !this->is_boundary_index(level, j); // !on_boundary(j)
310 }
311 
312 
313 
314 inline const IndexSet &
316 {
318  return boundary_indices[level];
319 }
320 
321 
322 
323 inline const IndexSet &
325 {
328 }
329 
330 
331 
332 inline bool
334 {
335  return boundary_indices.size() != 0;
336 }
337 
338 
339 
340 inline const AffineConstraints<double> &
342 {
344  return level_constraints[level];
345 }
346 
347 
348 
349 inline const AffineConstraints<double> &
351 {
353  return user_constraints[level];
354 }
355 
356 
357 
359 
360 #endif
bool is_element(const size_type index) const
Definition: index_set.h:1879
std::vector< IndexSet > refinement_edge_indices
std::vector< AffineConstraints< double > > user_constraints
void make_no_normal_flux_constraints(const DoFHandler< dim, spacedim > &dof, const types::boundary_id bid, const unsigned int first_vector_component)
bool is_interface_matrix_entry(const unsigned int level, const types::global_dof_index i, const types::global_dof_index j) const
void make_zero_boundary_constraints(const DoFHandler< dim, spacedim > &dof, const std::set< types::boundary_id > &boundary_ids, const ComponentMask &component_mask={})
std::vector< IndexSet > boundary_indices
void add_boundary_indices(const DoFHandler< dim, spacedim > &dof, const unsigned int level, const IndexSet &boundary_indices)
bool at_refinement_edge(const unsigned int level, const types::global_dof_index index) const
void initialize(const DoFHandler< dim, spacedim > &dof, const MGLevelObject< IndexSet > &level_relevant_dofs=MGLevelObject< IndexSet >())
bool is_boundary_index(const unsigned int level, const types::global_dof_index index) const
void merge_constraints(AffineConstraints< Number > &constraints, const unsigned int level, const bool add_boundary_indices, const bool add_refinement_edge_indices, const bool add_level_constraints, const bool add_user_constraints) const
std::vector< AffineConstraints< double > > level_constraints
bool have_boundary_indices() const
void add_user_constraints(const unsigned int level, const AffineConstraints< double > &constraints_on_level)
const IndexSet & get_refinement_edge_indices(unsigned int level) const
const AffineConstraints< double > & get_user_constraint_matrix(const unsigned int level) const
std::vector< std::set< types::global_dof_index > >::size_type size_dof
const AffineConstraints< double > & get_level_constraints(const unsigned int level) const
const IndexSet & get_boundary_indices(const unsigned int level) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_CXX20_REQUIRES(condition)
Definition: config.h:166
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
unsigned int level
Definition: grid_out.cc:4617
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
types::global_dof_index size_type
Definition: cuda_kernels.h:45
unsigned int global_dof_index
Definition: types.h:82
unsigned int boundary_id
Definition: types.h:145