403 return !cell.is_artificial();
407 const unsigned int n_raw_lines =
triangulation.n_raw_lines();
408 this->line_to_cells.resize(n_raw_lines);
416 const unsigned int line_to_children[12][2] = {{0, 2},
430 boost::container::small_vector<std::array<unsigned int, 3>, 6>>
431 line_to_inactive_cells(n_raw_lines);
436 const unsigned int cell_level = cell->
level();
438 for (
unsigned int line = 0; line < GeometryInfo<3>::lines_per_cell;
441 const unsigned int line_idx = cell->
line_index(line);
443 line_to_cells[line_idx].push_back(
446 line_to_inactive_cells[line_idx].push_back(
454 for (
unsigned int line_idx = 0; line_idx < n_raw_lines; ++line_idx)
456 if ((line_to_cells[line_idx].
size() > 0) &&
457 line_to_inactive_cells[line_idx].size() > 0)
463 line_to_inactive_cells[line_idx][0][0],
464 line_to_inactive_cells[line_idx][0][1]);
465 const unsigned int neighbor_line =
466 line_to_inactive_cells[line_idx][0][2];
468 for (
unsigned int c = 0; c < 2; ++c)
471 inactive_cell->child(line_to_children[neighbor_line][c]);
472 const unsigned int child_line_idx =
473 child->line_index(neighbor_line);
478 for (
const auto &cl : line_to_cells[line_idx])
479 line_to_cells[child_line_idx].push_back(cl);
525 const CellIterator &cell)
const
528 if ((dim == 3 && line_to_cells.empty()) ||
529 (cell->reference_cell().is_hyper_cube() ==
false))
532 if (cell->level() == 0)
535 const std::uint16_t subcell =
536 cell->parent()->child_iterator_to_index(cell);
537 const std::uint16_t
subcell_x = (subcell >> 0) & 1;
538 const std::uint16_t
subcell_y = (subcell >> 1) & 1;
539 const std::uint16_t
subcell_z = (subcell >> 2) & 1;
541 std::uint16_t face = 0;
542 std::uint16_t edge = 0;
544 for (
unsigned int direction = 0; direction < dim; ++direction)
546 const auto side = (subcell >> direction) & 1U;
547 const auto face_no = direction * 2 + side;
550 if (cell->at_boundary(face_no))
553 const auto &neighbor = cell->neighbor(face_no);
557 if (neighbor->has_children() || neighbor->is_artificial() ||
558 neighbor->level() == cell->level())
562 if (neighbor->get_fe().n_dofs_per_cell() == 0)
565 face |= 1 << direction;
569 for (
unsigned int direction = 0; direction < dim; ++direction)
570 if (face == 0 || face == (1 << direction))
572 const unsigned int line_no =
579 const unsigned int line_index = cell->line_index(line_no);
581 const auto edge_neighbor =
582 std::find_if(line_to_cells[line_index].begin(),
583 line_to_cells[line_index].end(),
584 [&cell](
const auto &edge_neighbor) {
586 &cell->get_triangulation(),
589 &cell->get_dof_handler());
590 return dof_cell.is_artificial() ==
false &&
591 dof_cell.level() < cell->level() &&
592 dof_cell.
get_fe().n_dofs_per_cell() > 0;
595 if (edge_neighbor == line_to_cells[line_index].end())
598 edge |= 1 << direction;
601 if ((face == 0) && (edge == 0))
604 const std::uint16_t inverted_subcell = (subcell ^ (dim == 2 ? 3 : 7));
607 inverted_subcell + (face << 3) + (edge << 6));
609 return refinement_configuration;
618 const CellIterator &cell,
619 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
620 const std::vector<std::vector<unsigned int>> &lexicographic_mapping,
621 const std::vector<std::vector<bool>> &supported_components,
623 std::vector<types::global_dof_index> &dof_indices)
const
625 if (std::find(supported_components[cell->active_fe_index()].begin(),
626 supported_components[cell->active_fe_index()].end(),
628 supported_components[cell->active_fe_index()].end())
631 const auto &fe = cell->get_fe();
634 std::vector<std::vector<unsigned int>>
635 component_to_system_index_face_array(fe.n_components());
637 for (
unsigned int i = 0; i < fe.n_dofs_per_face(0); ++i)
638 component_to_system_index_face_array
639 [fe.face_system_to_component_index(i, 0).first]
642 std::vector<unsigned int> idx_offset = {0};
644 for (
unsigned int base_element_index = 0;
645 base_element_index < fe.n_base_elements();
646 ++base_element_index)
647 for (
unsigned int c = 0;
648 c < fe.element_multiplicity(base_element_index);
650 idx_offset.push_back(
652 fe.base_element(base_element_index).n_dofs_per_cell());
654 std::vector<types::global_dof_index> neighbor_dofs_all(idx_offset.back());
655 std::vector<types::global_dof_index> neighbor_dofs_all_temp(
657 std::vector<types::global_dof_index> neighbor_dofs_face(
658 fe.n_dofs_per_face(0));
661 const auto get_face_idx = [](
const auto n_dofs_1d,
664 const auto j) ->
unsigned int {
665 const auto direction = face_no / 2;
666 const auto side = face_no % 2;
667 const auto offset = (side == 1) ? (n_dofs_1d - 1) : 0;
670 return (direction == 0) ? (n_dofs_1d * i + offset) :
671 (n_dofs_1d * offset + i);
676 return n_dofs_1d * n_dofs_1d * i + n_dofs_1d * j + offset;
678 return n_dofs_1d * n_dofs_1d * j + n_dofs_1d * offset + i;
680 return n_dofs_1d * n_dofs_1d * offset + n_dofs_1d * i + j;
690 const std::uint16_t kind =
691 static_cast<std::uint16_t
>(refinement_configuration);
692 const std::uint16_t subcell = (kind >> 0) & 7;
693 const std::uint16_t
subcell_x = (subcell >> 0) & 1;
694 const std::uint16_t
subcell_y = (subcell >> 1) & 1;
695 const std::uint16_t
subcell_z = (subcell >> 2) & 1;
696 const std::uint16_t face = (kind >> 3) & 7;
697 const std::uint16_t edge = (kind >> 6) & 7;
699 for (
unsigned int direction = 0; direction < dim; ++direction)
700 if ((face >> direction) & 1U)
702 const auto side = ((subcell >> direction) & 1U) == 0;
703 const auto face_no = direction * 2 + side;
706 cell->neighbor(face_no)
707 ->face(cell->neighbor_face_no(face_no))
708 ->get_dof_indices(neighbor_dofs_face,
709 cell->neighbor(face_no)->active_fe_index());
713 for (
auto &
index : neighbor_dofs_face)
716 for (
unsigned int base_element_index = 0, comp = 0;
717 base_element_index < fe.n_base_elements();
718 ++base_element_index)
719 for (
unsigned int c = 0;
720 c < fe.element_multiplicity(base_element_index);
723 if (supported_components[cell->active_fe_index()][comp] ==
727 const unsigned int n_dofs_1d =
729 .base_element(base_element_index)
732 const unsigned int dofs_per_face =
734 std::vector<types::global_dof_index> neighbor_dofs(
736 const auto lex_face_mapping =
741 for (
unsigned int i = 0; i < dofs_per_face; ++i)
742 neighbor_dofs[i] = neighbor_dofs_face
743 [component_to_system_index_face_array[comp][i]];
750 Assert(cell->combined_face_orientation(face_no) ==
756 orient_face(cell->combined_face_orientation(face_no),
766 for (
unsigned int i = 0, k = 0; i < n_dofs_1d; ++i)
767 for (
unsigned int j = 0; j < (dim == 2 ? 1 : n_dofs_1d);
769 dof_indices[get_face_idx(n_dofs_1d, face_no, i, j) +
771 neighbor_dofs[lex_face_mapping[k]];
776 for (
unsigned int direction = 0; direction < dim; ++direction)
777 if ((edge >> direction) & 1U)
779 const unsigned int line_no =
785 const unsigned int line_index = cell->line(line_no)->index();
787 const auto edge_neighbor =
788 std::find_if(line_to_cells[line_index].begin(),
789 line_to_cells[line_index].end(),
790 [&cell](
const auto &edge_array) {
792 edge_neighbor(&cell->get_triangulation(),
795 return edge_neighbor->is_artificial() ==
false &&
796 edge_neighbor->level() < cell->level();
799 if (edge_neighbor == line_to_cells[line_index].end())
803 &cell->get_triangulation(),
806 &cell->get_dof_handler());
807 const auto local_line_neighbor = (*edge_neighbor)[2];
812 for (
auto &
index : neighbor_dofs_all)
815 for (
unsigned int i = 0; i < neighbor_dofs_all_temp.size(); ++i)
816 neighbor_dofs_all_temp[i] = neighbor_dofs_all
817 [lexicographic_mapping[cell->active_fe_index()][i]];
820 cell->line_orientation(line_no) !=
821 neighbor_cell.line_orientation(local_line_neighbor);
823 for (
unsigned int base_element_index = 0, comp = 0;
824 base_element_index < fe.n_base_elements();
825 ++base_element_index)
826 for (
unsigned int c = 0;
827 c < fe.element_multiplicity(base_element_index);
830 if (supported_components[cell->active_fe_index()][comp] ==
834 const unsigned int n_dofs_1d =
836 .base_element(base_element_index)
840 for (
unsigned int i = 0; i < n_dofs_1d; ++i)
841 dof_indices[line_dof_idx(line_no, i, n_dofs_1d) +
842 idx_offset[comp]] = neighbor_dofs_all_temp
843 [line_dof_idx(local_line_neighbor,
844 flipped ? (n_dofs_1d - 1 - i) : i,
857 const CellIterator &cell,
858 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
859 const std::vector<std::vector<unsigned int>> &lexicographic_mapping,
860 std::vector<types::global_dof_index> &dof_indices,
864 const auto supported_components = compute_supported_components(
865 cell->get_dof_handler().get_fe_collection());
867 if (std::none_of(supported_components.begin(),
868 supported_components.end(),
870 return *std::max_element(a.begin(), a.end());
875 const auto refinement_configuration =
876 compute_refinement_configuration(cell);
882 update_dof_indices(cell,
884 lexicographic_mapping,
885 supported_components,
886 refinement_configuration,
890 for (
unsigned int c = 0; c < supported_components[0].size(); ++c)
891 if (supported_components[cell->active_fe_index()][c])
892 masks[c] = refinement_configuration;
918 const unsigned int n_dofs_1d,
919 std::vector<types::global_dof_index> &dofs)
const
921 const auto [orientation, rotation, flip] =
923 const int n_rotations =
924 rotation || flip ? 4 -
int(rotation) - 2 *
int(flip) : 0;
926 const unsigned int rot_mapping[4] = {2, 0, 3, 1};
930 const unsigned int dofs_per_line = n_dofs_1d - 2;
933 std::vector<types::global_dof_index> copy(dofs.size());
934 for (
int t = 0; t < n_rotations; ++t)
936 std::swap(copy, dofs);
939 for (
unsigned int i = 0; i < 4; ++i)
940 dofs[rot_mapping[i]] = copy[i];
943 unsigned int offset = 4;
944 for (
unsigned int i = 0; i < dofs_per_line; ++i)
948 copy[offset + 2 * dofs_per_line + (dofs_per_line - 1 - i)];
950 dofs[offset + dofs_per_line + i] =
951 copy[offset + 3 * dofs_per_line + (dofs_per_line - 1 - i)];
953 dofs[offset + 2 * dofs_per_line + i] =
954 copy[offset + dofs_per_line + i];
956 dofs[offset + 3 * dofs_per_line + i] = copy[offset + i];
960 offset += 4 * dofs_per_line;
962 for (
unsigned int i = 0; i < dofs_per_line; ++i)
963 for (
unsigned int j = 0; j < dofs_per_line; ++j)
964 dofs[offset + i * dofs_per_line + j] =
965 copy[offset + j * dofs_per_line + (dofs_per_line - 1 - i)];
979 unsigned int offset = 4;
980 for (
unsigned int i = 0; i < dofs_per_line; ++i)
983 dofs[offset + i] = copy[offset + 2 * dofs_per_line + i];
985 dofs[offset + dofs_per_line + i] =
986 copy[offset + 3 * dofs_per_line + i];
988 dofs[offset + 2 * dofs_per_line + i] = copy[offset + i];
990 dofs[offset + 3 * dofs_per_line + i] =
991 copy[offset + dofs_per_line + i];
995 offset += 4 * dofs_per_line;
996 for (
unsigned int i = 0; i < dofs_per_line; ++i)
997 for (
unsigned int j = 0; j < dofs_per_line; ++j)
998 dofs[offset + i * dofs_per_line + j] =
999 copy[offset + j * dofs_per_line + i];