16 #ifndef dealii_mg_constrained_dofs_h 17 #define dealii_mg_constrained_dofs_h 19 #include <deal.II/base/config.h> 21 #include <deal.II/base/subscriptor.h> 23 #include <deal.II/lac/affine_constraints.h> 25 #include <deal.II/multigrid/mg_tools.h> 30 DEAL_II_NAMESPACE_OPEN
32 template <
int dim,
int spacedim>
45 using size_dof = std::vector<std::set<types::global_dof_index>>::size_type;
65 template <
int dim,
int spacedim>
79 template <
int dim,
int spacedim>
80 DEAL_II_DEPRECATED
void 96 template <
int dim,
int spacedim>
100 const std::set<types::boundary_id> &boundary_ids,
186 template <
int dim,
int spacedim>
199 for (
unsigned int l = 0; l < nlevels; ++l)
209 for (; cell != endc; ++cell)
212 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
213 if (cell->has_periodic_neighbor(f))
215 if (cell->is_locally_owned_on_level())
218 cell->periodic_neighbor(f)->level_subdomain_id() !=
221 "Periodic neighbor of a locally owned cell must either be owned or ghost."));
225 else if (cell->periodic_neighbor(f)->level_subdomain_id() ==
228 Assert(cell->is_locally_owned_on_level() ==
false,
233 const unsigned int dofs_per_face =
234 cell->face(f)->
get_fe(0).dofs_per_face;
235 std::vector<types::global_dof_index> dofs_1(dofs_per_face);
236 std::vector<types::global_dof_index> dofs_2(dofs_per_face);
238 cell->periodic_neighbor(f)
239 ->face(cell->periodic_neighbor_face_no(f))
240 ->get_mg_dof_indices(l, dofs_1, 0);
241 cell->face(f)->get_mg_dof_indices(l, dofs_2, 0);
247 for (
unsigned int i = 0; i < dofs_per_face; ++i)
270 template <
int dim,
int spacedim>
292 template <
int dim,
int spacedim>
296 const std::set<types::boundary_id> &boundary_ids,
343 const unsigned int level,
347 const IndexSet &interface_dofs_on_level =
393 DEAL_II_NAMESPACE_CLOSE
bool is_boundary_index(const unsigned int level, const types::global_dof_index index) const
bool is_interface_matrix_entry(const unsigned int level, const types::global_dof_index i, const types::global_dof_index j) const
bool at_refinement_edge(const unsigned int level, const types::global_dof_index index) const
const AffineConstraints< double > & get_level_constraint_matrix(const unsigned int level) const
const IndexSet & get_boundary_indices(const unsigned int level) const
#define AssertIndexRange(index, range)
std::vector< IndexSet > boundary_indices
const Triangulation< dim, spacedim > & get_triangulation() const
unsigned long long int global_dof_index
static::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
types::global_dof_index n_dofs() const
cell_iterator end() const
void make_zero_boundary_constraints(const DoFHandler< dim, spacedim > &dof, const std::set< types::boundary_id > &boundary_ids, const ComponentMask &component_mask=ComponentMask())
std::vector< IndexSet > refinement_edge_indices
const types::subdomain_id artificial_subdomain_id
cell_iterator begin(const unsigned int level=0) const
const IndexSet & get_refinement_edge_indices(unsigned int level) const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
bool have_boundary_indices() const
bool is_element(const size_type index) const
typename ActiveSelector::cell_iterator cell_iterator
std::vector< AffineConstraints< double > > level_constraints
static::ExceptionBase & ExcInternalError()
void initialize(const DoFHandler< dim, spacedim > &dof)