57 std::vector<unsigned int> &,
69 std::vector<unsigned int> &,
82 std::vector<unsigned int> &,
93 std::vector<unsigned int> &,
103 template <
int dim,
int spacedim>
106 const unsigned int level,
107 std::vector<unsigned int> & row_lengths,
115 std::fill(row_lengths.begin(), row_lengths.end(), 0);
118 std::vector<bool> old_flags;
128 std::vector<types::global_dof_index> cell_indices;
129 std::vector<types::global_dof_index> neighbor_indices;
142 cell->get_mg_dof_indices(cell_indices);
166 const unsigned int face_no = 0;
171 unsigned int increment =
174 row_lengths[cell_indices[i++]] += increment;
191 row_lengths[cell_indices[i++]] += increment;
200 row_lengths[cell_indices[i++]] += increment;
205 row_lengths[cell_indices[i++]] += increment;
220 bool level_boundary = cell->at_boundary(iface);
224 neighbor = cell->neighbor(iface);
225 if (
static_cast<unsigned int>(neighbor->level()) !=
level)
226 level_boundary =
true;
231 for (
unsigned int local_dof = 0;
234 row_lengths[cell_indices[local_dof]] +=
253 const unsigned int dof_increment =
255 for (
unsigned int local_dof = 0;
258 row_lengths[cell_indices[local_dof]] += dof_increment;
263 if (face->user_flag_set())
265 face->set_user_flag();
279 neighbor->get_mg_dof_indices(neighbor_indices);
282 row_lengths[cell_indices[local_dof]] +=
286 row_lengths[neighbor_indices[local_dof]] +=
295 template <
int dim,
int spacedim>
298 const unsigned int level,
299 std::vector<unsigned int> & row_lengths,
308 std::fill(row_lengths.begin(), row_lengths.end(), 0);
311 std::vector<bool> old_flags;
321 std::vector<types::global_dof_index> cell_indices;
322 std::vector<types::global_dof_index> neighbor_indices;
328 std::vector<Table<2, DoFTools::Coupling>> couple_cell;
329 std::vector<Table<2, DoFTools::Coupling>> couple_face;
343 const unsigned int fe_index = cell->active_fe_index();
347 const unsigned int face_no = 0;
363 cell->get_mg_dof_indices(cell_indices);
386 unsigned int increment;
399 row_lengths[cell_indices[i]] += increment;
424 ((dim > 1) ? (dim - 1) :
427 row_lengths[cell_indices[i]] += increment;
444 ((dim > 2) ? (dim - 2) :
447 row_lengths[cell_indices[i]] += increment;
466 row_lengths[cell_indices[i]] += increment;
484 bool level_boundary = cell->at_boundary(iface);
488 neighbor = cell->neighbor(iface);
489 if (
static_cast<unsigned int>(neighbor->level()) !=
level)
490 level_boundary =
true;
495 for (
unsigned int local_dof = 0;
498 row_lengths[cell_indices[local_dof]] +=
518 for (
unsigned int local_dof = 0;
521 if (couple_face[fe_index](
525 const unsigned int dof_increment =
528 row_lengths[cell_indices[local_dof]] += dof_increment;
534 if (face->user_flag_set())
536 face->set_user_flag();
558 neighbor->get_mg_dof_indices(neighbor_indices);
562 for (
unsigned int local_dof = 0;
565 if (couple_cell[fe_index](
568 row_lengths[cell_indices[local_dof]] +=
573 for (
unsigned int local_dof = 0;
576 if (couple_cell[fe_index](
579 row_lengths[neighbor_indices[local_dof]] +=
590 typename SparsityPatternType,
594 SparsityPatternType & sparsity,
595 const unsigned int level,
597 const bool keep_constrained_dofs)
602 Assert(sparsity.n_rows() == n_dofs,
604 Assert(sparsity.n_cols() == n_dofs,
608 std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
611 for (; cell != endc; ++cell)
614 cell->level_subdomain_id() ==
617 cell->get_mg_dof_indices(dofs_on_this_cell);
620 keep_constrained_dofs);
626 template <
int dim,
typename SparsityPatternType,
int spacedim>
629 SparsityPatternType & sparsity,
630 const unsigned int level)
635 Assert(sparsity.n_rows() == n_dofs,
637 Assert(sparsity.n_cols() == n_dofs,
641 std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
642 std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
645 for (; cell != endc; ++cell)
647 if (!cell->is_locally_owned_on_level())
650 cell->get_mg_dof_indices(dofs_on_this_cell);
652 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
653 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
654 sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
659 bool use_face =
false;
660 if ((!cell->at_boundary(face)) &&
661 (
static_cast<unsigned int>(cell->neighbor_level(face)) ==
664 else if (cell->has_periodic_neighbor(face) &&
665 (
static_cast<unsigned int>(
666 cell->periodic_neighbor_level(face)) ==
level))
672 cell->neighbor_or_periodic_neighbor(face);
673 neighbor->get_mg_dof_indices(dofs_on_other_cell);
677 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
679 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
681 sparsity.add(dofs_on_this_cell[i],
682 dofs_on_other_cell[j]);
685 if (neighbor->is_locally_owned_on_level() ==
false)
686 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
687 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
689 sparsity.add(dofs_on_other_cell[i],
690 dofs_on_other_cell[j]);
691 sparsity.add(dofs_on_other_cell[i],
692 dofs_on_this_cell[j]);
701 template <
int dim,
typename SparsityPatternType,
int spacedim>
704 SparsityPatternType & sparsity,
705 const unsigned int level)
717 Assert(sparsity.n_rows() == coarse_dofs,
719 Assert(sparsity.n_cols() == fine_dofs,
723 std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
724 std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
727 for (; cell != endc; ++cell)
729 if (!cell->is_locally_owned_on_level())
732 cell->get_mg_dof_indices(dofs_on_this_cell);
737 bool use_face =
false;
738 if ((!cell->at_boundary(face)) &&
739 (
static_cast<unsigned int>(cell->neighbor_level(face)) !=
742 else if (cell->has_periodic_neighbor(face) &&
743 (
static_cast<unsigned int>(
744 cell->periodic_neighbor_level(face)) !=
level))
750 cell->neighbor_or_periodic_neighbor(face);
751 neighbor->get_mg_dof_indices(dofs_on_other_cell);
753 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
755 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
757 sparsity.add(dofs_on_other_cell[i],
758 dofs_on_this_cell[j]);
759 sparsity.add(dofs_on_other_cell[j],
760 dofs_on_this_cell[i]);
770 template <
int dim,
typename SparsityPatternType,
int spacedim>
773 SparsityPatternType & sparsity,
774 const unsigned int level,
784 Assert(sparsity.n_rows() == n_dofs,
786 Assert(sparsity.n_cols() == n_dofs,
788 Assert(int_mask.n_rows() == n_comp,
790 Assert(int_mask.n_cols() == n_comp,
792 Assert(flux_mask.n_rows() == n_comp,
794 Assert(flux_mask.n_cols() == n_comp,
798 std::vector<types::global_dof_index> dofs_on_this_cell(total_dofs);
799 std::vector<types::global_dof_index> dofs_on_other_cell(total_dofs);
812 for (
unsigned int i = 0; i < total_dofs; ++i)
824 std::vector<bool> user_flags;
829 for (; cell != endc; ++cell)
831 if (!cell->is_locally_owned_on_level())
834 cell->get_mg_dof_indices(dofs_on_this_cell);
836 for (
unsigned int i = 0; i < total_dofs; ++i)
837 for (
unsigned int j = 0; j < total_dofs; ++j)
839 sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
846 if (cell_face->user_flag_set())
849 if (cell->at_boundary(face) && !cell->has_periodic_neighbor(face))
851 for (
unsigned int i = 0; i < total_dofs; ++i)
853 const bool i_non_zero_i = support_on_face(i, face);
854 for (
unsigned int j = 0; j < total_dofs; ++j)
856 const bool j_non_zero_i = support_on_face(j, face);
859 sparsity.add(dofs_on_this_cell[i],
860 dofs_on_this_cell[j]);
862 i_non_zero_i && j_non_zero_i)
863 sparsity.add(dofs_on_this_cell[i],
864 dofs_on_this_cell[j]);
871 cell->neighbor_or_periodic_neighbor(face);
873 if (neighbor->level() < cell->level())
876 unsigned int neighbor_face =
877 cell->has_periodic_neighbor(face) ?
878 cell->periodic_neighbor_of_periodic_neighbor(face) :
879 cell->neighbor_of_neighbor(face);
881 neighbor->get_mg_dof_indices(dofs_on_other_cell);
882 for (
unsigned int i = 0; i < total_dofs; ++i)
884 const bool i_non_zero_i = support_on_face(i, face);
885 const bool i_non_zero_e = support_on_face(i, neighbor_face);
886 for (
unsigned int j = 0; j < total_dofs; ++j)
888 const bool j_non_zero_i = support_on_face(j, face);
889 const bool j_non_zero_e =
890 support_on_face(j, neighbor_face);
893 sparsity.add(dofs_on_this_cell[i],
894 dofs_on_other_cell[j]);
895 sparsity.add(dofs_on_other_cell[i],
896 dofs_on_this_cell[j]);
897 sparsity.add(dofs_on_this_cell[i],
898 dofs_on_this_cell[j]);
899 sparsity.add(dofs_on_other_cell[i],
900 dofs_on_other_cell[j]);
904 if (i_non_zero_i && j_non_zero_e)
905 sparsity.add(dofs_on_this_cell[i],
906 dofs_on_other_cell[j]);
907 if (i_non_zero_e && j_non_zero_i)
908 sparsity.add(dofs_on_other_cell[i],
909 dofs_on_this_cell[j]);
910 if (i_non_zero_i && j_non_zero_i)
911 sparsity.add(dofs_on_this_cell[i],
912 dofs_on_this_cell[j]);
913 if (i_non_zero_e && j_non_zero_e)
914 sparsity.add(dofs_on_other_cell[i],
915 dofs_on_other_cell[j]);
920 sparsity.add(dofs_on_this_cell[j],
921 dofs_on_other_cell[i]);
922 sparsity.add(dofs_on_other_cell[j],
923 dofs_on_this_cell[i]);
924 sparsity.add(dofs_on_this_cell[j],
925 dofs_on_this_cell[i]);
926 sparsity.add(dofs_on_other_cell[j],
927 dofs_on_other_cell[i]);
931 if (j_non_zero_i && i_non_zero_e)
932 sparsity.add(dofs_on_this_cell[j],
933 dofs_on_other_cell[i]);
934 if (j_non_zero_e && i_non_zero_i)
935 sparsity.add(dofs_on_other_cell[j],
936 dofs_on_this_cell[i]);
937 if (j_non_zero_i && i_non_zero_i)
938 sparsity.add(dofs_on_this_cell[j],
939 dofs_on_this_cell[i]);
940 if (j_non_zero_e && i_non_zero_e)
941 sparsity.add(dofs_on_other_cell[j],
942 dofs_on_other_cell[i]);
946 neighbor->face(neighbor_face)->set_user_flag();
953 .load_user_flags(user_flags);
958 template <
int dim,
typename SparsityPatternType,
int spacedim>
961 SparsityPatternType & sparsity,
962 const unsigned int level,
979 Assert(sparsity.n_rows() == coarse_dofs,
981 Assert(sparsity.n_cols() == fine_dofs,
983 Assert(flux_mask.n_rows() == n_comp,
985 Assert(flux_mask.n_cols() == n_comp,
989 std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
990 std::vector<types::global_dof_index> dofs_on_other_cell(dofs_per_cell);
1000 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1004 for (; cell != endc; ++cell)
1006 if (!cell->is_locally_owned_on_level())
1009 cell->get_mg_dof_indices(dofs_on_this_cell);
1014 bool use_face =
false;
1015 if ((!cell->at_boundary(face)) &&
1016 (
static_cast<unsigned int>(cell->neighbor_level(face)) !=
1019 else if (cell->has_periodic_neighbor(face) &&
1020 (
static_cast<unsigned int>(
1021 cell->periodic_neighbor_level(face)) !=
level))
1027 cell->neighbor_or_periodic_neighbor(face);
1028 neighbor->get_mg_dof_indices(dofs_on_other_cell);
1030 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1032 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
1036 sparsity.add(dofs_on_other_cell[i],
1037 dofs_on_this_cell[j]);
1038 sparsity.add(dofs_on_other_cell[j],
1039 dofs_on_this_cell[i]);
1050 template <
int dim,
int spacedim,
typename SparsityPatternType>
1054 SparsityPatternType & sparsity,
1055 const unsigned int level)
1060 Assert(sparsity.n_rows() == n_dofs,
1062 Assert(sparsity.n_cols() == n_dofs,
1066 std::vector<types::global_dof_index> dofs_on_this_cell(dofs_per_cell);
1069 for (; cell != endc; ++cell)
1070 if (cell->is_locally_owned_on_level())
1072 cell->get_mg_dof_indices(dofs_on_this_cell);
1073 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1074 for (
unsigned int j = 0; j < dofs_per_cell; ++j)
1076 level, dofs_on_this_cell[i], dofs_on_this_cell[j]))
1077 sparsity.add(dofs_on_this_cell[i], dofs_on_this_cell[j]);
1083 template <
int dim,
int spacedim>
1087 std::vector<std::vector<types::global_dof_index>> &result,
1089 std::vector<unsigned int> target_component)
1093 const unsigned int nlevels =
1096 Assert(result.size() == nlevels,
1099 if (target_component.size() == 0)
1101 target_component.resize(n_components);
1102 for (
unsigned int i = 0; i < n_components; ++i)
1103 target_component[i] = i;
1106 Assert(target_component.size() == n_components,
1109 for (
unsigned int l = 0;
l < nlevels; ++
l)
1111 result[
l].resize(n_components);
1112 std::fill(result[
l].
begin(), result[
l].
end(), 0
U);
1118 if (n_components == 1)
1120 result[
l][0] = dof_handler.
n_dofs(
l);
1127 std::vector<std::vector<bool>> dofs_in_component(
1128 n_components, std::vector<bool>(dof_handler.
n_dofs(
l),
false));
1129 std::vector<ComponentMask> component_select(n_components);
1131 for (
unsigned int i = 0; i < n_components; ++i)
1133 void (*fun_ptr)(
const unsigned int level,
1136 std::vector<bool> &) =
1137 &DoFTools::extract_level_dofs<dim, spacedim>;
1139 std::vector<bool> tmp(n_components,
false);
1146 component_select[i],
1147 dofs_in_component[i]);
1152 unsigned int component = 0;
1161 for (
unsigned int dd = 0; dd <
d; ++dd)
1164 result[
l][target_component[component]] +=
1165 std::count(dofs_in_component[component].
begin(),
1166 dofs_in_component[component].
end(),
1174 std::accumulate(result[
l].begin(),
1185 template <
int dim,
int spacedim>
1189 std::vector<std::vector<types::global_dof_index>> &dofs_per_block,
1190 std::vector<unsigned int> target_block)
1194 const unsigned int n_levels =
1199 for (
unsigned int l = 0;
l < n_levels; ++
l)
1200 std::fill(dofs_per_block[
l].
begin(), dofs_per_block[
l].
end(), 0
U);
1204 if (target_block.size() == 0)
1207 for (
unsigned int i = 0; i <
n_blocks; ++i)
1208 target_block[i] = i;
1213 const unsigned int max_block =
1214 *std::max_element(target_block.begin(), target_block.end());
1215 const unsigned int n_target_blocks = max_block + 1;
1216 (void)n_target_blocks;
1218 for (
unsigned int l = 0;
l < n_levels; ++
l)
1227 for (
unsigned int l = 0;
l < n_levels; ++
l)
1228 dofs_per_block[
l][0] = dof_handler.
n_dofs(
l);
1234 for (
unsigned int l = 0;
l < n_levels; ++
l)
1236 std::vector<std::vector<bool>> dofs_in_block(
1238 std::vector<BlockMask> block_select(
n_blocks);
1240 for (
unsigned int i = 0; i <
n_blocks; ++i)
1242 void (*fun_ptr)(
const unsigned int level,
1245 std::vector<bool> &) =
1246 &DoFTools::extract_level_dofs<dim, spacedim>;
1248 std::vector<bool> tmp(
n_blocks,
false);
1250 block_select[i] = tmp;
1253 fun_ptr,
l, dof_handler, block_select[i], dofs_in_block[i]);
1258 for (
unsigned int block = 0; block < fe.
n_blocks(); ++block)
1259 dofs_per_block[
l][target_block[block]] +=
1260 std::count(dofs_in_block[block].
begin(),
1261 dofs_in_block[block].
end(),
1268 template <
int dim,
int spacedim>
1274 std::vector<std::set<types::global_dof_index>> &boundary_indices,
1281 std::set<types::boundary_id> boundary_ids;
1282 for (
const auto &boundary_function : function_map)
1283 boundary_ids.insert(boundary_function.first);
1285 std::vector<IndexSet> boundary_indexset;
1288 boundary_indices[i].
insert(boundary_indexset[i].
begin(),
1289 boundary_indexset[i].
end());
1293 template <
int dim,
int spacedim>
1298 std::vector<IndexSet> &boundary_indices,
1305 std::set<types::boundary_id> boundary_ids;
1306 for (
const auto &boundary_function : function_map)
1307 boundary_ids.insert(boundary_function.first);
1314 template <
int dim,
int spacedim>
1317 const std::set<types::boundary_id> &boundary_ids,
1318 std::vector<IndexSet> & boundary_indices,
1324 if (boundary_ids.size() == 0)
1328 if (boundary_indices[i].size() == 0)
1332 const bool fe_is_system = (n_components != 1);
1334 std::vector<types::global_dof_index> local_dofs;
1338 std::vector<std::vector<types::global_dof_index>> dofs_by_level(
1352 const unsigned int level = cell->level();
1355 if (cell->at_boundary(face_no) ==
true)
1358 cell->face(face_no);
1361 if (boundary_ids.find(bi) != boundary_ids.end())
1364 face->get_mg_dof_indices(
level, local_dofs);
1376 "It's probably worthwhile to select at least one component."));
1384 if (cell->at_boundary(face_no) ==
false)
1388 const unsigned int level = cell->level();
1391 cell->face(face_no);
1393 face->boundary_id();
1394 if (boundary_ids.find(boundary_component) != boundary_ids.end())
1397 for (
unsigned int i = 0;
1398 i < cell->get_fe().n_dofs_per_cell();
1402 cell->get_fe().get_nonzero_components(i);
1406 bool selected =
false;
1407 for (
unsigned int c = 0; c < n_components; ++c)
1408 if (nonzero_component_array[c] ==
true &&
1409 component_mask[c] ==
true)
1415 for (
unsigned int c = 0; c < n_components; ++c)
1417 nonzero_component_array[c] ==
false ||
1418 component_mask[c] ==
true,
1420 "You are using a non-primitive FiniteElement "
1421 "and try to constrain just some of its components!"));
1427 face->get_mg_dof_indices(
level, local_dofs);
1430 for (
unsigned int i = 0; i < local_dofs.size(); ++i)
1432 unsigned int component =
1444 cell->get_fe().get_nonzero_components(i);
1445 for (
unsigned int c = 0; c < n_components; ++c)
1446 if (nonzero_component_array[c] ==
true)
1454 if (component_mask[component] ==
true)
1455 dofs_by_level[
level].push_back(local_dofs[i]);
1469 boundary_indices[
level].add_indices(
1478 template <
int dim,
int spacedim>
1481 std::vector<IndexSet> & interface_dofs)
1483 Assert(interface_dofs.size() ==
1486 interface_dofs.size(),
1489 std::vector<std::vector<types::global_dof_index>> tmp_interface_dofs(
1490 interface_dofs.size());
1496 std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
1498 std::vector<bool> cell_dofs(dofs_per_cell,
false);
1501 endc = mg_dof_handler.
end();
1503 for (; cell != endc; ++cell)
1512 bool has_coarser_neighbor =
false;
1514 std::fill(cell_dofs.begin(), cell_dofs.end(),
false);
1519 cell->face(face_nr);
1520 if (!face->at_boundary() || cell->has_periodic_neighbor(face_nr))
1524 cell->neighbor_or_periodic_neighbor(face_nr);
1531 neighbor->level_subdomain_id() ==
1536 if (neighbor->level() < cell->level())
1542 has_coarser_neighbor =
true;
1547 if (has_coarser_neighbor ==
false)
1550 const unsigned int level = cell->level();
1551 cell->get_mg_dof_indices(local_dof_indices);
1553 for (
unsigned int i = 0; i < dofs_per_cell; ++i)
1556 tmp_interface_dofs[
level].push_back(local_dof_indices[i]);
1560 for (
unsigned int l = 0;
1564 interface_dofs[
l].
clear();
1565 std::sort(tmp_interface_dofs[
l].
begin(), tmp_interface_dofs[
l].
end());
1566 interface_dofs[
l].add_indices(tmp_interface_dofs[
l].
begin(),
1567 std::unique(tmp_interface_dofs[
l].
begin(),
1568 tmp_interface_dofs[
l].
end()));
1569 interface_dofs[
l].compress();
1575 template <
int dim,
int spacedim>
1588 for (; cell != endc; ++cell)
1589 if (cell->is_locally_owned())
1591 min_level = cell->level();
1595 unsigned int global_min = min_level;
1610 template <
int dim,
int spacedim>
1624 tr->is_multilevel_hierarchy_constructed(),
1626 "We can only compute the workload imbalance if the multilevel hierarchy has been constructed!"));
1628 const unsigned int n_proc =
1630 const unsigned int n_global_levels = tr->n_global_levels();
1640 for (
int lvl = n_global_levels - 1; lvl >= 0; --lvl)
1645 for (
const auto &cell : tr->cell_iterators_on_level(lvl))
1646 if (cell->is_locally_owned_on_level())
1647 ++n_owned_cells_on_lvl;
1651 tr->get_communicator());
1653 total_cells_in_hierarchy +=
1655 tr->get_communicator());
1658 const double ideal_work =
1659 total_cells_in_hierarchy /
static_cast<double>(n_proc);
1669#include "mg_tools.inst"
void add_entries_local_to_global(const std::vector< size_type > &local_dof_indices, SparsityPatternType &sparsity_pattern, const bool keep_constrained_entries=true, const Table< 2, bool > &dof_mask=Table< 2, bool >()) const
unsigned int n_selected_components(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
cell_iterator end() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
cell_iterator begin(const unsigned int level=0) const
unsigned int get_first_line_index() const
unsigned int n_dofs_per_cell() const
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_blocks() const
unsigned int n_components() const
ReferenceCell reference_cell() const
unsigned int get_first_hex_index() const
std::pair< unsigned int, types::global_dof_index > system_to_block_index(const unsigned int component) const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
bool is_primitive() const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
unsigned int element_multiplicity(const unsigned int index) const
unsigned int n_base_elements() const
std::pair< unsigned int, unsigned int > face_system_to_component_index(const unsigned int index, const unsigned int face_no=0) const
types::global_dof_index first_block_of_base(const unsigned int b) const
bool is_interface_matrix_entry(const unsigned int level, const types::global_dof_index i, const types::global_dof_index j) const
virtual types::subdomain_id locally_owned_subdomain() const
void save_user_flags(std::ostream &out) const
unsigned int n_levels() const
cell_iterator end() const
virtual unsigned int n_global_levels() const
void load_user_flags(std::istream &in)
active_cell_iterator begin_active(const unsigned int level=0) const
unsigned int max_dofs_per_face() const
unsigned int n_components() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators() const
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ActiveSelector::cell_iterator cell_iterator
Task< RT > new_task(const std::function< RT()> &function)
std::enable_if< IsBlockVector< VectorType >::value, unsignedint >::type n_blocks(const VectorType &vector)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(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)
T min(const T &t, const MPI_Comm &mpi_communicator)
T sum(const T &t, const MPI_Comm &mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
const types::subdomain_id artificial_subdomain_id
const types::subdomain_id invalid_subdomain_id
static const unsigned int invalid_unsigned_int
const types::global_dof_index invalid_dof_index