36#include <deal.II/multigrid/mg_transfer_component.templates.h>
79 template <
int dim,
typename number,
int spacedim>
81 reinit_vector_by_components(
82 const ::DoFHandler<dim, spacedim> & mg_dof,
84 const std::vector<bool> & sel,
85 const std::vector<unsigned int> & target_comp,
86 std::vector<std::vector<types::global_dof_index>> &ndofs)
88 std::vector<bool> selected = sel;
89 std::vector<unsigned int> target_component = target_comp;
90 const unsigned int ncomp = mg_dof.get_fe(0).n_components();
101 if (target_component.size() == 0)
103 target_component.resize(ncomp);
104 for (
unsigned int i = 0; i < ncomp; ++i)
105 target_component[i] = i;
110 if (selected.size() == 0)
112 selected.resize(target_component.size());
113 std::fill_n(selected.begin(), ncomp,
false);
114 for (
const unsigned int component : target_component)
115 selected[component] =
true;
118 Assert(selected.size() == target_component.size(),
122 const unsigned int n_selected =
123 std::accumulate(selected.begin(), selected.end(), 0u);
125 if (ndofs.size() == 0)
127 std::vector<std::vector<types::global_dof_index>> new_dofs(
128 mg_dof.get_triangulation().n_levels(),
129 std::vector<types::global_dof_index>(target_component.size()));
134 for (
unsigned int level = v.min_level();
level <= v.max_level(); ++
level)
136 v[
level].reinit(n_selected, 0);
138 for (
unsigned int i = 0;
146 v[
level].collect_sizes();
178 template <
int dim,
typename number,
int spacedim>
180 reinit_vector_by_components(
181 const ::DoFHandler<dim, spacedim> & mg_dof,
184 const std::vector<unsigned int> & target_component,
185 std::vector<std::vector<types::global_dof_index>> &ndofs)
188 ExcMessage(
"The component mask does not have the correct size."));
190 unsigned int selected_block = 0;
191 for (
unsigned int i = 0; i < target_component.size(); ++i)
192 if (component_mask[i])
193 selected_block = target_component[i];
195 if (ndofs.size() == 0)
197 std::vector<std::vector<types::global_dof_index>> new_dofs(
198 mg_dof.get_triangulation().n_levels(),
199 std::vector<types::global_dof_index>(target_component.size()));
204 for (
unsigned int level = v.min_level();
level <= v.max_level(); ++
level)
206 v[
level].reinit(ndofs[
level][selected_block]);
212template <
typename number>
213template <
int dim,
class InVector,
int spacedim>
218 const InVector & src)
const
223 ExcMatricesNotBuilt());
225 reinit_vector_by_components(
226 mg_dof_handler, dst, mg_component_mask, mg_target_component, sizes);
244 using IT = std::vector<
245 std::pair<types::global_dof_index, unsigned int>>::const_iterator;
246 for (IT i = copy_to_and_from_indices[
level].
begin();
247 i != copy_to_and_from_indices[
level].end();
249 dst[
level](i->second) = src(i->first);
255template <
int dim,
int spacedim>
265template <
int dim,
int spacedim>
318 const unsigned int n_components =
325 ExcMessage(
"Component mask has wrong size."));
328 sizes.resize(n_levels);
329 for (
unsigned int l = 0;
l < n_levels; ++
l)
330 sizes[
l].resize(n_components);
344 level_component_start = k;
385 for (
unsigned int i = 0; i < n_levels - 1; ++i)
394 std::vector<types::global_dof_index> dof_indices_parent(dofs_per_cell);
395 std::vector<types::global_dof_index> dof_indices_child(dofs_per_cell);
415 for (
unsigned int i = 0; i < n_components; ++i)
416 for (
unsigned int j = 0; j < n_components; ++j)
430 if (cell->has_children())
432 cell->get_mg_dof_indices(dof_indices_parent);
434 for (
unsigned int child = 0; child < cell->n_children(); ++child)
441 child, cell->refinement_case());
443 cell->child(child)->get_mg_dof_indices(dof_indices_child);
448 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
449 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
450 if (prolongation(i, j) != 0)
452 const unsigned int icomp =
454 const unsigned int jcomp =
458 dof_indices_child[i], dof_indices_parent[j]);
470 if (cell->has_children())
472 cell->get_mg_dof_indices(dof_indices_parent);
474 for (
unsigned int child = 0; child < cell->n_children(); ++child)
481 child, cell->refinement_case());
483 cell->child(child)->get_mg_dof_indices(dof_indices_child);
487 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
488 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
489 if (prolongation(i, j) != 0)
491 const unsigned int icomp =
493 const unsigned int jcomp =
497 dof_indices_child[i],
498 dof_indices_parent[j],
511 std::vector<std::vector<types::global_dof_index>> dofs_per_component(
513 std::vector<types::global_dof_index>(n_components));
523 for (
unsigned int iblock = 0; iblock < n_components; ++iblock)
524 for (
unsigned int jblock = 0; jblock < n_components; ++jblock)
525 if (iblock == jblock)
533 ->block(iblock, jblock)
536 ->block(iblock, jblock)
545 dofs_per_component[
level]);
550 std::set<types::global_dof_index>::const_iterator
553 const bool is_boundary_index =
556 if (is_boundary_index)
559 ->block(iblock, jblock)
560 .set(i, column_number, 0);
571template <
typename number>
572template <
int dim,
int spacedim>
578 unsigned int mg_select,
579 const std::vector<unsigned int> & t_component,
580 const std::vector<unsigned int> & mg_t_component,
581 const std::vector<std::set<types::global_dof_index>> &bdry_indices)
583 build(mg_dof, select, mg_select, t_component, mg_t_component, bdry_indices);
588template <
typename number>
589template <
int dim,
int spacedim>
594 unsigned int mg_select,
595 const std::vector<unsigned int> & t_component,
596 const std::vector<unsigned int> & mg_t_component,
597 const std::vector<std::set<types::global_dof_index>> &bdry_indices)
602 target_component = t_component;
603 mg_target_component = mg_t_component;
604 boundary_indices = bdry_indices;
606 selected_component = select;
607 mg_selected_component = mg_select;
610 std::vector<bool> tmp(ncomp,
false);
611 for (
unsigned int c = 0; c < ncomp; ++c)
612 if (t_component[c] == selected_component)
618 std::vector<bool> tmp(ncomp,
false);
619 for (
unsigned int c = 0; c < ncomp; ++c)
620 if (mg_t_component[c] == mg_selected_component)
629 for (
unsigned int i = 0; i < target_component.size(); ++i)
631 if (target_component[i] == select)
633 selected_component = i;
638 for (
unsigned int i = 0; i < mg_target_component.size(); ++i)
640 if (mg_target_component[i] == mg_select)
642 mg_selected_component = i;
652 interface_dofs[
l].clear();
653 interface_dofs[
l].set_size(mg_dof.
n_dofs(
l));
659 std::vector<types::global_dof_index> temp_copy_indices;
660 std::vector<types::global_dof_index> global_dof_indices(fe.
n_dofs_per_cell());
661 std::vector<types::global_dof_index> level_dof_indices(fe.
n_dofs_per_cell());
665 copy_to_and_from_indices[
level].clear();
671 temp_copy_indices.resize(0);
677 for (; level_cell != level_end; ++level_cell)
683 level_cell->get_dof_indices(global_dof_indices);
684 level_cell->get_mg_dof_indices(level_dof_indices);
688 const unsigned int component =
690 if (component_mask[component] &&
691 !interface_dofs[
level].is_element(level_dof_indices[i]))
694 mg_component_start[
level][mg_target_component[component]];
696 component_start[target_component[component]];
697 temp_copy_indices[level_dof_indices[i] - level_start] =
698 global_dof_indices[i] - global_start;
706 std::count_if(temp_copy_indices.begin(),
707 temp_copy_indices.end(),
709 return index != numbers::invalid_dof_index;
711 copy_to_and_from_indices[
level].resize(n_active_dofs);
715 copy_to_and_from_indices[
level][counter++] =
716 std::pair<types::global_dof_index, unsigned int>(
717 temp_copy_indices[i], i);
725#include "mg_transfer_component.inst"
size_type local_to_global(const unsigned int block, const size_type index) const
bool represents_n_components(const unsigned int n) const
cell_iterator end() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const Triangulation< dim, spacedim > & get_triangulation() const
active_cell_iterator begin_active(const unsigned int level=0) const
active_cell_iterator end_active(const unsigned int level) const
types::global_dof_index n_dofs() const
cell_iterator begin(const unsigned int level=0) const
unsigned int n_dofs_per_cell() const
unsigned int n_components() const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
std::vector< unsigned int > mg_target_component
ComponentMask mg_component_mask
std::vector< unsigned int > target_component
void build_matrices(const DoFHandler< dim, spacedim > &dof, const DoFHandler< dim, spacedim > &mg_dof)
std::vector< types::global_dof_index > component_start
std::vector< std::shared_ptr< BlockSparseMatrix< double > > > prolongation_matrices
std::vector< std::vector< std::pair< types::global_dof_index, unsigned int > > > copy_to_and_from_indices
std::vector< std::vector< types::global_dof_index > > mg_component_start
std::vector< std::shared_ptr< BlockSparsityPattern > > prolongation_sparsities
void build(const DoFHandler< dim, spacedim > &dof_handler)
std::vector< std::vector< types::global_dof_index > > sizes
std::vector< std::set< types::global_dof_index > > boundary_indices
void do_copy_to_mg(const DoFHandler< dim, spacedim > &mg_dof, MGLevelObject< Vector< number > > &dst, const InVector &src) const
void build_matrices(const DoFHandler< dim, spacedim > &dof, const DoFHandler< dim, spacedim > &mg_dof, unsigned int selected, unsigned int mg_selected, const std::vector< unsigned int > &target_component=std::vector< unsigned int >(), const std::vector< unsigned int > &mg_target_component=std::vector< unsigned int >(), const std::vector< std::set< types::global_dof_index > > &boundary_indices=std::vector< std::set< types::global_dof_index > >())
void build(const DoFHandler< dim, spacedim > &dof, unsigned int selected, unsigned int mg_selected, const std::vector< unsigned int > &target_component=std::vector< unsigned int >(), const std::vector< unsigned int > &mg_target_component=std::vector< unsigned int >(), const std::vector< std::set< types::global_dof_index > > &boundary_indices=std::vector< std::set< types::global_dof_index > >())
unsigned int n_levels() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
std::enable_if< IsBlockVector< VectorType >::value, unsignedint >::type n_blocks(const VectorType &vector)
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
VectorType::value_type * end(VectorType &V)
VectorType::value_type * begin(VectorType &V)
const types::global_dof_index invalid_dof_index