boundary_constraints[level].reinit(
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler, level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_refinement_edge_indices(level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_refinement_edge_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
boundary_constraints[level].close();
}
AffineConstraints<double> level_constraints(
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler, level));
- level_constraints.add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ level_constraints.add_constraint(dof_index, {}, 0.);
level_constraints.close();
typename MatrixFree<dim, float>::AdditionalData additional_data;
0) &&
!constraints_hanging_nodes.is_constrained(index_z))
{
- all_constraints.add_line(index_z);
- all_constraints.set_inhomogeneity(index_z,
- undeformed_gap);
+ all_constraints.add_constraint(index_z,
+ {},
+ undeformed_gap);
distributed_solution(index_z) = undeformed_gap;
active_set.add_index(index_z);
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler,
level));
- level_constraints.add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ level_constraints.add_constraint(dof_index, {}, 0.);
level_constraints.close();
typename MatrixFree<dim, float>::AdditionalData additional_data;
boundary_constraints[level].reinit(
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler, level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_refinement_edge_indices(level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_refinement_edge_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
boundary_constraints[level].close();
}
// this here by marking the first pressure dof, which has index n_u as a
// constrained dof.
if (solver_type == SolverType::UMFPACK)
- constraints.add_line(n_u);
+ constraints.add_constraint(n_u, {}, 0.);
constraints.close();
}
triangulation.n_levels());
for (unsigned int level = 0; level < triangulation.n_levels(); ++level)
{
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_refinement_edge_indices(level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_refinement_edge_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
boundary_constraints[level].close();
- IndexSet idx = mg_constrained_dofs.get_refinement_edge_indices(level) &
- mg_constrained_dofs.get_boundary_indices(level);
+ const IndexSet idx =
+ mg_constrained_dofs.get_refinement_edge_indices(level) &
+ mg_constrained_dofs.get_boundary_indices(level);
- boundary_interface_constraints[level].add_lines(idx);
+ for (const types::global_dof_index dof_index : idx)
+ boundary_interface_constraints[level].add_constraint(dof_index,
+ {},
+ 0.);
boundary_interface_constraints[level].close();
}
boundary_constraints[level].reinit(
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler, level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_refinement_edge_indices(level));
- boundary_constraints[level].add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_refinement_edge_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ boundary_constraints[level].add_constraint(dof_index, {}, 0.);
boundary_constraints[level].close();
}
AffineConstraints<double> level_constraints(
dof_handler.locally_owned_mg_dofs(level),
DoFTools::extract_locally_relevant_level_dofs(dof_handler, level));
- level_constraints.add_lines(
- mg_constrained_dofs.get_boundary_indices(level));
+
+ for (const types::global_dof_index dof_index :
+ mg_constrained_dofs.get_boundary_indices(level))
+ level_constraints.add_constraint(dof_index, {}, 0.);
level_constraints.close();
typename MatrixFree<dim, float>::AdditionalData additional_data;
types::global_dof_index last_density_dof =
density_dofs.nth_index_in_set(density_dofs.n_elements() - 1);
constraints.clear();
- constraints.add_line(last_density_dof);
- for (unsigned int i = 0; i < density_dofs.n_elements() - 1; ++i)
- constraints.add_entry(last_density_dof,
- density_dofs.nth_index_in_set(i),
- -1);
- constraints.set_inhomogeneity(last_density_dof, 0);
-
+ {
+ std::vector<std::pair<types::global_dof_index, double>>
+ constraint_entries;
+ constraint_entries.reserve(density_dofs.n_elements() - 1);
+ for (const types::global_dof_index dof_index : density_dofs)
+ if (dof_index != last_density_dof)
+ constraint_entries.emplace_back(dof_index, -1.);
+
+ constraints.add_constraint(last_density_dof, constraint_entries, 0.);
+ }
constraints.close();
// We can now finally create the sparsity pattern for the
// 3) refill this constraint matrix
for (const auto &line : temporal_constraint_matrix)
- {
- // ... line
- this->add_line(line.index);
-
- // ... inhomogeneity
- if (line.inhomogeneity != number())
- this->set_inhomogeneity(line.index, line.inhomogeneity);
-
- // ... entries
- if (!line.entries.empty())
- this->add_entries(line.index, line.entries);
- }
+ this->add_constraint(line.index, line.entries, line.inhomogeneity);
#ifdef DEBUG
Assert(this->is_consistent_in_parallel(
for (const ConstraintLine &line : lines)
if (mask.is_element(line.index))
{
- const size_type row = mask.index_within_set(line.index);
- view.add_line(row);
+#ifdef DEBUG
for (const std::pair<size_type, number> &entry : line.entries)
{
Assert(
"In creating a view of an AffineConstraints "
"object, the constraint on degree of freedom " +
std::to_string(line.index) + " (which corresponds to the " +
- std::to_string(row) +
+ std::to_string(mask.index_within_set(line.index)) +
"th degree of freedom selected in the mask) "
"is constrained against degree of freedom " +
std::to_string(entry.first) +
", but this degree of freedom is not listed in the mask and "
"consequently cannot be transcribed into the index space "
"of the output object."));
- view.add_entry(row,
- mask.index_within_set(entry.first),
- entry.second);
}
- view.set_inhomogeneity(row, line.inhomogeneity);
+#endif
+
+ std::vector<std::pair<size_type, number>> translated_entries =
+ line.entries;
+ for (auto &entry : translated_entries)
+ entry.first = mask.index_within_set(entry.first);
+
+ view.add_constraint(mask.index_within_set(line.index),
+ translated_entries,
+ line.inhomogeneity);
}
view.close();
!level_constraints[l].is_constrained(dofs_2[i]) &&
!level_constraints[l].is_constrained(dofs_1[i]))
{
- level_constraints[l].add_line(dofs_2[i]);
- level_constraints[l].add_entry(dofs_2[i], dofs_1[i], 1.);
+ level_constraints[l].add_constraint(dofs_2[i],
+ {{dofs_1[i], 1.}},
+ 0.);
}
}
level_constraints[l].close();
// merge constraints
if (add_boundary_indices && this->have_boundary_indices())
- constraints.add_lines(this->get_boundary_indices(level));
+ for (const auto i : this->get_boundary_indices(level))
+ constraints.add_constraint(i, {}, 0.);
if (add_refinement_edge_indices)
- constraints.add_lines(this->get_refinement_edge_indices(level));
+ for (const auto i : this->get_refinement_edge_indices(level))
+ constraints.add_constraint(i, {}, 0.);
if (add_level_constraints)
constraints.merge(this->get_level_constraints(level),
if (!constraints.is_constrained(dof_indices.dof_indices[0]) &&
constraints.can_store_line(dof_indices.dof_indices[0]))
{
- constraints.add_line(dof_indices.dof_indices[0]);
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity / constraining_vector[0]) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity / constraining_vector[0] :
+ 0);
if (std::fabs(constraining_vector[1] /
constraining_vector[0]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[0],
- dof_indices.dof_indices[1],
- -constraining_vector[1] /
- constraining_vector[0]);
-
- if (std::fabs(inhomogeneity / constraining_vector[0]) >
- std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(
- dof_indices.dof_indices[0],
- inhomogeneity / constraining_vector[0]);
+ constraints.add_constraint(dof_indices.dof_indices[0],
+ {{dof_indices.dof_indices[1],
+ -constraining_vector[1] /
+ constraining_vector[0]}},
+ normalized_inhomogeneity);
+ else
+ constraints.add_constraint(dof_indices.dof_indices[0],
+ {},
+ normalized_inhomogeneity);
}
}
else
if (!constraints.is_constrained(dof_indices.dof_indices[1]) &&
constraints.can_store_line(dof_indices.dof_indices[1]))
{
- constraints.add_line(dof_indices.dof_indices[1]);
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity / constraining_vector[1]) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity / constraining_vector[1] :
+ 0);
if (std::fabs(constraining_vector[0] /
constraining_vector[1]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[1],
- dof_indices.dof_indices[0],
- -constraining_vector[0] /
- constraining_vector[1]);
-
- if (std::fabs(inhomogeneity / constraining_vector[1]) >
- std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(
- dof_indices.dof_indices[1],
- inhomogeneity / constraining_vector[1]);
+ constraints.add_constraint(dof_indices.dof_indices[1],
+ {{dof_indices.dof_indices[0],
+ -constraining_vector[0] /
+ constraining_vector[1]}},
+ normalized_inhomogeneity);
+ else
+ constraints.add_constraint(dof_indices.dof_indices[1],
+ {},
+ normalized_inhomogeneity);
}
}
break;
case 3:
{
+ // Store constraints. There are at most 2 entries per
+ // constraint, so a boost::container::small_vector with
+ // static size 2 will do just fine.
+ //
+ // TODO: This could use std::inplace_vector once available (in
+ // C++26?)
+ boost::container::
+ small_vector<std::pair<types::global_dof_index, double>, 2>
+ constraint_entries;
+
if ((std::fabs(constraining_vector[0]) >=
std::fabs(constraining_vector[1]) + 1e-10) &&
(std::fabs(constraining_vector[0]) >=
if (!constraints.is_constrained(dof_indices.dof_indices[0]) &&
constraints.can_store_line(dof_indices.dof_indices[0]))
{
- constraints.add_line(dof_indices.dof_indices[0]);
-
if (std::fabs(constraining_vector[1] /
constraining_vector[0]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[0],
- dof_indices.dof_indices[1],
- -constraining_vector[1] /
- constraining_vector[0]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[1],
+ -constraining_vector[1] / constraining_vector[0]);
if (std::fabs(constraining_vector[2] /
constraining_vector[0]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[0],
- dof_indices.dof_indices[2],
- -constraining_vector[2] /
- constraining_vector[0]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[2],
+ -constraining_vector[2] / constraining_vector[0]);
- if (std::fabs(inhomogeneity / constraining_vector[0]) >
- std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(
- dof_indices.dof_indices[0],
- inhomogeneity / constraining_vector[0]);
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity / constraining_vector[0]) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity / constraining_vector[0] :
+ 0);
+
+ constraints.add_constraint(dof_indices.dof_indices[0],
+ constraint_entries,
+ normalized_inhomogeneity);
}
}
else if ((std::fabs(constraining_vector[1]) + 1e-10 >=
if (!constraints.is_constrained(dof_indices.dof_indices[1]) &&
constraints.can_store_line(dof_indices.dof_indices[1]))
{
- constraints.add_line(dof_indices.dof_indices[1]);
-
if (std::fabs(constraining_vector[0] /
constraining_vector[1]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[1],
- dof_indices.dof_indices[0],
- -constraining_vector[0] /
- constraining_vector[1]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[0],
+ -constraining_vector[0] / constraining_vector[1]);
if (std::fabs(constraining_vector[2] /
constraining_vector[1]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[1],
- dof_indices.dof_indices[2],
- -constraining_vector[2] /
- constraining_vector[1]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[2],
+ -constraining_vector[2] / constraining_vector[1]);
- if (std::fabs(inhomogeneity / constraining_vector[1]) >
- std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(
- dof_indices.dof_indices[1],
- inhomogeneity / constraining_vector[1]);
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity / constraining_vector[1]) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity / constraining_vector[1] :
+ 0);
+
+ constraints.add_constraint(dof_indices.dof_indices[1],
+ constraint_entries,
+ normalized_inhomogeneity);
}
}
else
if (!constraints.is_constrained(dof_indices.dof_indices[2]) &&
constraints.can_store_line(dof_indices.dof_indices[2]))
{
- constraints.add_line(dof_indices.dof_indices[2]);
-
if (std::fabs(constraining_vector[0] /
constraining_vector[2]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[2],
- dof_indices.dof_indices[0],
- -constraining_vector[0] /
- constraining_vector[2]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[0],
+ -constraining_vector[0] / constraining_vector[2]);
if (std::fabs(constraining_vector[1] /
constraining_vector[2]) >
std::numeric_limits<double>::epsilon())
- constraints.add_entry(dof_indices.dof_indices[2],
- dof_indices.dof_indices[1],
- -constraining_vector[1] /
- constraining_vector[2]);
+ constraint_entries.emplace_back(
+ dof_indices.dof_indices[1],
+ -constraining_vector[1] / constraining_vector[2]);
- if (std::fabs(inhomogeneity / constraining_vector[2]) >
- std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(
- dof_indices.dof_indices[2],
- inhomogeneity / constraining_vector[2]);
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity / constraining_vector[2]) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity / constraining_vector[2] :
+ 0);
+
+ constraints.add_constraint(dof_indices.dof_indices[2],
+ constraint_entries,
+ normalized_inhomogeneity);
}
}
if (!constraints.is_constrained(dof_indices.dof_indices[d]) &&
constraints.can_store_line(dof_indices.dof_indices[d]))
{
- constraints.add_line(dof_indices.dof_indices[d]);
-
- if (std::fabs(tangent_vector[d] /
- tangent_vector[largest_component]) >
- std::numeric_limits<double>::epsilon())
- constraints.add_entry(
- dof_indices.dof_indices[d],
- dof_indices.dof_indices[largest_component],
- tangent_vector[d] / tangent_vector[largest_component]);
-
const double inhomogeneity =
(b_values(d) * tangent_vector[largest_component] -
b_values(largest_component) * tangent_vector[d]) /
tangent_vector[largest_component];
- if (std::fabs(inhomogeneity) >
+ const double normalized_inhomogeneity =
+ (std::fabs(inhomogeneity) >
+ std::numeric_limits<double>::epsilon() ?
+ inhomogeneity :
+ 0);
+
+ if (std::fabs(tangent_vector[d] /
+ tangent_vector[largest_component]) >
std::numeric_limits<double>::epsilon())
- constraints.set_inhomogeneity(dof_indices.dof_indices[d],
- inhomogeneity);
+ constraints.add_constraint(
+ dof_indices.dof_indices[d],
+ {{dof_indices.dof_indices[largest_component],
+ tangent_vector[d] / tangent_vector[largest_component]}},
+ normalized_inhomogeneity);
+ else
+ constraints.add_constraint(dof_indices.dof_indices[d],
+ {},
+ normalized_inhomogeneity);
}
}
Assert(primary_dofs[col] != numbers::invalid_dof_index,
ExcInternalError());
- std::vector<
- std::pair<typename AffineConstraints<number2>::size_type, number2>>
+ // Build constraints in a vector of pairs that can be
+ // arbitrarily large, but that holds up to 25 elements without
+ // external memory allocation. This is good enough for hanging
+ // node constraints of Q4 elements in 3d, so covers most
+ // common cases.
+ boost::container::small_vector<
+ std::pair<typename AffineConstraints<number2>::size_type, number2>,
+ 25>
entries;
entries.reserve(n_primary_dofs);
for (unsigned int row = 0; row != n_dependent_dofs; ++row)
std::vector<types::global_dof_index> dofs_on_mother;
std::vector<types::global_dof_index> dofs_on_children;
+ // Build constraints in a vector of pairs that can be
+ // arbitrarily large, but that holds up to 25 elements without
+ // external memory allocation. This is good enough for hanging
+ // node constraints of Q4 elements in 3d, so covers most
+ // common cases.
+ boost::container::small_vector<
+ std::pair<typename AffineConstraints<number>::size_type, number>,
+ 25>
+ constraint_entries;
+
// loop over all lines; only on lines there can be constraints. We do so
// by looping over all active cells and checking whether any of the faces
// are refined which can only be from the neighboring cell because this
// note that even though we may visit a face twice if the neighboring
// cells are equally refined, we can only visit each face with hanging
// nodes once
- typename DoFHandler<dim_, spacedim>::active_cell_iterator
- cell = dof_handler.begin_active(),
- endc = dof_handler.end();
- for (; cell != endc; ++cell)
+ for (const auto &cell : dof_handler.active_cell_iterators())
{
// artificial cells can at best neighbor ghost cells, but we're not
// interested in these interfaces
for (unsigned int row = 0; row != dofs_on_children.size();
++row)
{
- constraints.add_line(dofs_on_children[row]);
+ constraint_entries.clear();
+ constraint_entries.reserve(dofs_on_mother.size());
for (unsigned int i = 0; i != dofs_on_mother.size(); ++i)
- constraints.add_entry(dofs_on_children[row],
- dofs_on_mother[i],
- fe.constraints()(row, i));
+ constraint_entries.emplace_back(dofs_on_mother[i],
+ fe.constraints()(row, i));
- constraints.set_inhomogeneity(dofs_on_children[row], 0.);
+ constraints.add_constraint(dofs_on_children[row],
+ constraint_entries,
+ 0.);
}
}
else
std::vector<types::global_dof_index> dofs_on_mother;
std::vector<types::global_dof_index> dofs_on_children;
+ // Build constraints in a vector of pairs that can be
+ // arbitrarily large, but that holds up to 25 elements without
+ // external memory allocation. This is good enough for hanging
+ // node constraints of Q4 elements in 3d, so covers most
+ // common cases.
+ boost::container::small_vector<
+ std::pair<typename AffineConstraints<number>::size_type, number>,
+ 25>
+ constraint_entries;
+
// loop over all quads; only on quads there can be constraints. We do so
// by looping over all active cells and checking whether any of the faces
// are refined which can only be from the neighboring cell because this
// note that even though we may visit a face twice if the neighboring
// cells are equally refined, we can only visit each face with hanging
// nodes once
- typename DoFHandler<dim_, spacedim>::active_cell_iterator
- cell = dof_handler.begin_active(),
- endc = dof_handler.end();
- for (; cell != endc; ++cell)
+ for (const auto &cell : dof_handler.active_cell_iterators())
{
// artificial cells can at best neighbor ghost cells, but we're not
// interested in these interfaces
for (unsigned int row = 0; row != dofs_on_children.size();
++row)
{
- constraints.add_line(dofs_on_children[row]);
+ constraint_entries.clear();
+ constraint_entries.reserve(dofs_on_mother.size());
for (unsigned int i = 0; i != dofs_on_mother.size(); ++i)
- constraints.add_entry(dofs_on_children[row],
- dofs_on_mother[i],
- fe.constraints()(row, i));
+ constraint_entries.emplace_back(dofs_on_mother[i],
+ fe.constraints()(row, i));
- constraints.set_inhomogeneity(dofs_on_children[row], 0.);
+ constraints.add_constraint(dofs_on_children[row],
+ constraint_entries,
+ 0.);
}
}
else
cell_to_rotated_face_index[cell_index] = i;
}
+ // Build constraints in a vector of pairs that can be
+ // arbitrarily large, but that holds up to 25 elements without
+ // external memory allocation. This is good enough for hanging
+ // node constraints of Q4 elements in 3d, so covers most
+ // common cases.
+ boost::container::small_vector<
+ std::pair<typename AffineConstraints<number>::size_type, number>,
+ 25>
+ constraint_entries;
+
//
// Loop over all dofs on face 2 and constrain them against all
// matching dofs on face 1:
//
-
for (unsigned int i = 0; i < dofs_per_face; ++i)
{
// Obey the component mask
if (affine_constraints.is_constrained(dofs_2[i]))
continue;
- // Enter the constraint piece by piece:
- affine_constraints.add_line(dofs_2[i]);
+ constraint_entries.clear();
+ constraint_entries.reserve(dofs_per_face);
for (unsigned int jj = 0; jj < dofs_per_face; ++jj)
{
jj, 0, face_orientation, face_flip, face_rotation)];
if (std::abs(transformation(i, jj)) > eps)
- affine_constraints.add_entry(dofs_2[i],
- dofs_1[j],
- transformation(i, jj));
+ constraint_entries.emplace_back(dofs_1[j],
+ transformation(i, jj));
}
+ // Enter the constraint::
+ affine_constraints.add_constraint(dofs_2[i],
+ constraint_entries,
+ 0.);
+
+
// Continue with next dof.
continue;
}
if (constraints_are_cyclic)
{
if (std::abs(cycle_constraint_factor - number(1.)) > eps)
- affine_constraints.add_line(dof_left);
+ affine_constraints.add_constraint(dof_left, {}, 0.);
}
else
{
- affine_constraints.add_line(dof_left);
- affine_constraints.add_entry(dof_left,
- dof_right,
- constraint_factor);
+ affine_constraints.add_constraint(
+ dof_left, {{dof_right, constraint_factor}}, 0.);
// The number 1e10 in the assert below is arbitrary. If the
// absolute value of constraint_factor is too large, then probably
// the absolute value of periodicity_factor is too large or too
// otherwise enter all constraints
- constraints.add_line(global_dof);
-
constraint_line.clear();
for (types::global_dof_index row = first_used_row;
row < n_coarse_dofs;
constraint_line.emplace_back(representants[row], j->second);
}
- constraints.add_entries(global_dof, constraint_line);
+ constraints.add_constraint(global_dof, constraint_line, 0.);
}
}
std::vector<types::global_dof_index> cell_dofs;
cell_dofs.reserve(dof.get_fe_collection().max_dofs_per_cell());
- typename DoFHandler<dim, spacedim>::active_cell_iterator
- cell = dof.begin_active(),
- endc = dof.end();
- for (; cell != endc; ++cell)
+ // In looping over faces, we will encounter some DoFs multiple
+ // times (namely, the ones on vertices and (in 3d) edges shared
+ // between multiple boundary faces. Keep track of which DoFs we
+ // have already encountered, so that we do not have to consider
+ // them a second time.
+ std::set<types::global_dof_index> dofs_already_treated;
+
+ for (const auto &cell : dof.active_cell_iterators())
if (!cell->is_artificial() && cell->at_boundary())
{
const FiniteElement<dim, spacedim> &fe = cell->get_fe();
// enter those dofs into the list that match the component
// signature.
for (const types::global_dof_index face_dof : face_dofs)
- {
- // Find out if a dof has a contribution in this component,
- // and if so, add it to the list
- const std::vector<types::global_dof_index>::iterator
- it_index_on_cell = std::find(cell_dofs.begin(),
- cell_dofs.end(),
- face_dof);
- Assert(it_index_on_cell != cell_dofs.end(),
- ExcInvalidIterator());
- const unsigned int index_on_cell =
- std::distance(cell_dofs.begin(), it_index_on_cell);
- const ComponentMask &nonzero_component_array =
- cell->get_fe().get_nonzero_components(index_on_cell);
- bool nonzero = false;
- for (unsigned int c = 0; c < n_components; ++c)
- if (nonzero_component_array[c] && component_mask[c])
- {
- nonzero = true;
- break;
- }
-
- if (nonzero)
- zero_boundary_constraints.add_line(face_dof);
- }
+ if (dofs_already_treated.find(face_dof) ==
+ dofs_already_treated.end())
+ {
+ // Find out if a dof has a contribution in this
+ // component, and if so, add it to the list
+ const std::vector<types::global_dof_index>::iterator
+ it_index_on_cell = std::find(cell_dofs.begin(),
+ cell_dofs.end(),
+ face_dof);
+ Assert(it_index_on_cell != cell_dofs.end(),
+ ExcInvalidIterator());
+ const unsigned int index_on_cell =
+ std::distance(cell_dofs.begin(), it_index_on_cell);
+ const ComponentMask &nonzero_component_array =
+ cell->get_fe().get_nonzero_components(index_on_cell);
+
+ bool nonzero = false;
+ for (unsigned int c = 0; c < n_components; ++c)
+ if (nonzero_component_array[c] && component_mask[c])
+ {
+ nonzero = true;
+ break;
+ }
+
+ if (nonzero)
+ zero_boundary_constraints.add_constraint(face_dof,
+ {},
+ 0.);
+ // We already dealt with this DoF. Make sure we
+ // don't touch it again.
+ dofs_already_treated.insert(face_dof);
+ }
}
}
}
if (map_iter != map_end)
for (unsigned int i = 0; i < dim; ++i)
{
- constraints[i].add_line(cell->vertex_dof_index(vertex_no, 0));
- constraints[i].set_inhomogeneity(
+ constraints[i].add_constraint(
cell->vertex_dof_index(vertex_no, 0),
+ {},
(solve_for_absolute_positions ?
map_iter->second(i) :
map_iter->second(i) - vertex_point[i]));