34 #include <deal.II/multigrid/mg_transfer_block.templates.h>
52 template <
int dim,
typename number,
int spacedim>
54 reinit_vector_by_blocks(
55 const ::DoFHandler<dim, spacedim> & dof_handler,
57 const std::vector<bool> & sel,
58 std::vector<std::vector<types::global_dof_index>> &ndofs)
60 std::vector<bool> selected = sel;
62 const unsigned int n_selected =
63 std::accumulate(selected.begin(), selected.end(), 0u);
65 if (ndofs.size() == 0)
67 std::vector<std::vector<types::global_dof_index>> new_dofs(
68 dof_handler.get_triangulation().n_levels(),
69 std::vector<types::global_dof_index>(selected.size()));
74 for (
unsigned int level = v.min_level();
level <= v.max_level(); ++
level)
76 v[
level].reinit(n_selected, 0);
78 for (
unsigned int i = 0;
86 v[
level].collect_sizes();
98 template <
int dim,
typename number,
int spacedim>
100 reinit_vector_by_blocks(
101 const ::DoFHandler<dim, spacedim> & dof_handler,
103 const unsigned int selected_block,
104 std::vector<std::vector<types::global_dof_index>> &ndofs)
106 const unsigned int n_blocks = dof_handler.get_fe().n_blocks();
109 std::vector<bool> selected(
n_blocks,
false);
110 selected[selected_block] =
true;
112 if (ndofs.size() == 0)
114 std::vector<std::vector<types::global_dof_index>> new_dofs(
115 dof_handler.get_triangulation().n_levels(),
116 std::vector<types::global_dof_index>(selected.size()));
121 for (
unsigned int level = v.min_level();
level <= v.max_level(); ++
level)
123 v[
level].reinit(ndofs[
level][selected_block]);
129 template <
typename number>
130 template <
int dim,
typename number2,
int spacedim>
137 reinit_vector_by_blocks(dof_handler, dst, selected_block, sizes);
146 for (IT i = copy_indices[selected_block][
level].
begin();
147 i != copy_indices[selected_block][
level].end();
149 dst[
level](i->second) = src.
block(selected_block)(i->first);
155 template <
typename number>
156 template <
int dim,
typename number2,
int spacedim>
161 const Vector<number2> & src)
const
163 reinit_vector_by_blocks(dof_handler, dst, selected_block, sizes);
171 for (IT i = copy_indices[selected_block][
level].
begin();
172 i != copy_indices[selected_block][
level].end();
174 dst[
level](i->second) = src(i->first);
180 template <
typename number>
181 template <
int dim,
typename number2,
int spacedim>
188 reinit_vector_by_blocks(dof_handler, dst, selected, sizes);
193 for (
unsigned int block = 0; block < selected.size(); ++block)
195 for (IT i = copy_indices[block][
level].
begin();
196 i != copy_indices[block][
level].end();
198 dst[
level].block(mg_block[block])(i->second) =
199 src.
block(block)(i->first);
205 template <
int dim,
int spacedim>
211 build(mg_dof_handler);
216 template <
int dim,
int spacedim>
233 for (
unsigned int i = 0; i <
n_blocks; ++i)
241 sizes.resize(n_levels, std::vector<types::global_dof_index>(fe.
n_blocks()));
254 level_block_start = k;
276 for (
unsigned int block = 0; block <
n_blocks; ++block)
297 for (
unsigned int i = 0; i < n_levels - 1; ++i)
306 std::vector<types::global_dof_index> dof_indices_parent(dofs_per_cell);
307 std::vector<types::global_dof_index> dof_indices_child(dofs_per_cell);
328 for (
unsigned int i = 0; i <
n_blocks; ++i)
329 for (
unsigned int j = 0; j <
n_blocks; ++j)
343 if (cell->has_children())
345 cell->get_mg_dof_indices(dof_indices_parent);
347 Assert(cell->n_children() ==
350 for (
unsigned int child = 0; child < cell->n_children(); ++child)
356 dof_handler.
get_fe().get_prolongation_matrix(
357 child, cell->refinement_case());
359 cell->child(child)->get_mg_dof_indices(dof_indices_child);
364 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
365 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
366 if (prolongation(i, j) != 0)
368 const unsigned int icomp =
370 const unsigned int jcomp =
372 if ((icomp == jcomp) &&
selected[icomp])
374 dof_indices_child[i], dof_indices_parent[j]);
386 if (cell->has_children())
388 cell->get_mg_dof_indices(dof_indices_parent);
390 Assert(cell->n_children() ==
393 for (
unsigned int child = 0; child < cell->n_children(); ++child)
399 dof_handler.
get_fe().get_prolongation_matrix(
400 child, cell->refinement_case());
402 cell->child(child)->get_mg_dof_indices(dof_indices_child);
406 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
407 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
408 if (prolongation(i, j) != 0)
410 const unsigned int icomp =
412 const unsigned int jcomp =
414 if ((icomp == jcomp) &&
selected[icomp])
416 dof_indices_child[i],
417 dof_indices_parent[j],
429 std::vector<types::global_dof_index> constrain_indices;
430 std::vector<std::vector<bool>> constraints_per_block(
n_blocks);
442 constrain_indices.resize(0);
447 for (; dof != endd; ++dof)
448 constrain_indices[*dof] = 1;
450 unsigned int index = 0;
451 for (
unsigned int block = 0; block <
n_blocks; ++block)
455 constraints_per_block[block].resize(0);
456 constraints_per_block[block].resize(n_dofs,
false);
458 constraints_per_block[block][i] =
459 (constrain_indices[index] == 1);
465 ->block(block, block)
469 for (; start_row != end_row; ++start_row)
471 if (constraints_per_block[block][start_row->column()])
472 start_row->value() = 0.;
480 template <
typename number>
481 template <
int dim,
int spacedim>
488 build(mg_dof, select);
491 template <
typename number>
492 template <
int dim,
int spacedim>
501 selected_block = select;
503 selected[select] =
true;
507 std::vector<types::global_dof_index> temp_copy_indices;
508 std::vector<types::global_dof_index> global_dof_indices(fe.
dofs_per_cell);
509 std::vector<types::global_dof_index> level_dof_indices(fe.
dofs_per_cell);
519 temp_copy_indices.resize(0);
520 temp_copy_indices.resize(sizes[
level][selected_block],
525 for (; level_cell != level_end; ++level_cell)
531 level_cell->get_dof_indices(global_dof_indices);
532 level_cell->get_mg_dof_indices(level_dof_indices);
539 if (mg_constrained_dofs !=
nullptr)
541 if (!mg_constrained_dofs->at_refinement_edge(
542 level, level_dof_indices[i]))
543 temp_copy_indices[level_dof_indices[i] -
544 mg_block_start[
level][block]] =
545 global_dof_indices[i] - block_start[block];
548 temp_copy_indices[level_dof_indices[i] -
549 mg_block_start[
level][block]] =
550 global_dof_indices[i] - block_start[block];
562 std::count_if(temp_copy_indices.begin(),
563 temp_copy_indices.end(),
565 return index != numbers::invalid_dof_index;
567 copy_indices[selected_block][
level].resize(n_active_dofs);
571 copy_indices[selected_block][
level][counter++] =
572 std::pair<types::global_dof_index, unsigned int>(
573 temp_copy_indices[i], i);
579 template <
typename number>
580 template <
int dim,
int spacedim>
585 const std::vector<bool> & sel)
590 template <
typename number>
591 template <
int dim,
int spacedim>
594 const std::vector<bool> & sel)
605 if (selected.size() == 0)
606 selected = std::vector<bool>(
n_blocks,
true);
610 std::vector<std::vector<types::global_dof_index>> temp_copy_indices(
n_blocks);
611 std::vector<types::global_dof_index> global_dof_indices(fe.
dofs_per_cell);
612 std::vector<types::global_dof_index> level_dof_indices(fe.
dofs_per_cell);
621 for (
unsigned int block = 0; block <
n_blocks; ++block)
624 temp_copy_indices[block].resize(0);
625 temp_copy_indices[block].resize(sizes[
level][block],
631 for (; level_cell != level_end; ++level_cell)
637 level_cell->get_dof_indices(global_dof_indices);
638 level_cell->get_mg_dof_indices(level_dof_indices);
644 temp_copy_indices[block][level_dof_indices[i] -
645 mg_block_start[
level][block]] =
646 global_dof_indices[i] - block_start[block];
650 for (
unsigned int block = 0; block <
n_blocks; ++block)
654 std::count_if(temp_copy_indices[block].
begin(),
655 temp_copy_indices[block].
end(),
659 copy_indices[block][
level].resize(n_active_dofs);
662 i < temp_copy_indices[block].size();
665 copy_indices[block][
level][counter++] =
666 std::pair<types::global_dof_index, unsigned int>(
667 temp_copy_indices[block][i], i);
676 #include "mg_transfer_block.inst"