From 9756336788d08d560e4a3ef0447746c34f768e53 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sat, 4 May 2024 14:38:21 +0200 Subject: [PATCH] Add template argument to VectorTools::compute_no_normal_flux_constraints() --- .../numerics/vector_tools_constraints.h | 4 +- .../vector_tools_constraints.templates.h | 93 ++++++++++--------- .../numerics/vector_tools_constraints.inst.in | 8 ++ 3 files changed, 57 insertions(+), 48 deletions(-) diff --git a/include/deal.II/numerics/vector_tools_constraints.h b/include/deal.II/numerics/vector_tools_constraints.h index 28d7ffa778..33e683e135 100644 --- a/include/deal.II/numerics/vector_tools_constraints.h +++ b/include/deal.II/numerics/vector_tools_constraints.h @@ -343,13 +343,13 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template void compute_no_normal_flux_constraints( const DoFHandler &dof_handler, const unsigned int first_vector_component, const std::set &boundary_ids, - AffineConstraints &constraints, + AffineConstraints &constraints, const Mapping &mapping = (ReferenceCells::get_hypercube() #ifndef _MSC_VER diff --git a/include/deal.II/numerics/vector_tools_constraints.templates.h b/include/deal.II/numerics/vector_tools_constraints.templates.h index 2e55ab4a35..c0f19b66db 100644 --- a/include/deal.II/numerics/vector_tools_constraints.templates.h +++ b/include/deal.II/numerics/vector_tools_constraints.templates.h @@ -118,12 +118,12 @@ namespace VectorTools * The function does not add constraints if a degree of freedom is already * constrained in the constraints object. */ - template + template void add_constraint(const VectorDoFTuple &dof_indices, const Tensor<1, dim> &constraining_vector, - AffineConstraints &constraints, - const double inhomogeneity = 0) + AffineConstraints &constraints, + const number inhomogeneity = 0) { // choose the DoF that has the largest component in the // constraining_vector as the one to be constrained as this makes the @@ -169,15 +169,15 @@ namespace VectorTools if (!constraints.is_constrained(dof_indices.dof_indices[0]) && constraints.can_store_line(dof_indices.dof_indices[0])) { - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity / constraining_vector[0]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity / constraining_vector[0] : 0); if (std::fabs(constraining_vector[1] / constraining_vector[0]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraints.add_constraint(dof_indices.dof_indices[0], {{dof_indices.dof_indices[1], -constraining_vector[1] / @@ -194,15 +194,15 @@ namespace VectorTools if (!constraints.is_constrained(dof_indices.dof_indices[1]) && constraints.can_store_line(dof_indices.dof_indices[1])) { - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity / constraining_vector[1]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity / constraining_vector[1] : 0); if (std::fabs(constraining_vector[0] / constraining_vector[1]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraints.add_constraint(dof_indices.dof_indices[1], {{dof_indices.dof_indices[0], -constraining_vector[0] / @@ -226,7 +226,7 @@ namespace VectorTools // TODO: This could use std::inplace_vector once available (in // C++26?) boost::container:: - small_vector, 2> + small_vector, 2> constraint_entries; if ((std::fabs(constraining_vector[0]) >= @@ -239,7 +239,7 @@ namespace VectorTools { if (std::fabs(constraining_vector[1] / constraining_vector[0]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[1], -constraining_vector[1] / @@ -247,15 +247,15 @@ namespace VectorTools if (std::fabs(constraining_vector[2] / constraining_vector[0]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[2], -constraining_vector[2] / constraining_vector[0])); - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity / constraining_vector[0]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity / constraining_vector[0] : 0); @@ -274,7 +274,7 @@ namespace VectorTools { if (std::fabs(constraining_vector[0] / constraining_vector[1]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[0], -constraining_vector[0] / @@ -282,15 +282,15 @@ namespace VectorTools if (std::fabs(constraining_vector[2] / constraining_vector[1]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[2], -constraining_vector[2] / constraining_vector[1])); - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity / constraining_vector[1]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity / constraining_vector[1] : 0); @@ -306,7 +306,7 @@ namespace VectorTools { if (std::fabs(constraining_vector[0] / constraining_vector[2]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[0], -constraining_vector[0] / @@ -314,15 +314,15 @@ namespace VectorTools if (std::fabs(constraining_vector[1] / constraining_vector[2]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraint_entries.emplace_back( std::make_pair(dof_indices.dof_indices[1], -constraining_vector[1] / constraining_vector[2])); - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity / constraining_vector[2]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity / constraining_vector[2] : 0); @@ -352,13 +352,13 @@ namespace VectorTools * The function does not add constraints if a degree of freedom is already * constrained in the constraints object. */ - template + template void add_tangentiality_constraints( const VectorDoFTuple &dof_indices, const Tensor<1, dim> &tangent_vector, - AffineConstraints &constraints, - const Vector &b_values = Vector(dim)) + AffineConstraints &constraints, + const Vector &b_values = Vector(dim)) { // choose the DoF that has the // largest component in the @@ -384,20 +384,20 @@ namespace VectorTools if (!constraints.is_constrained(dof_indices.dof_indices[d]) && constraints.can_store_line(dof_indices.dof_indices[d])) { - const double inhomogeneity = + const number inhomogeneity = (b_values(d) * tangent_vector[largest_component] - b_values(largest_component) * tangent_vector[d]) / tangent_vector[largest_component]; - const double normalized_inhomogeneity = + const number normalized_inhomogeneity = (std::fabs(inhomogeneity) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? inhomogeneity : 0); if (std::fabs(tangent_vector[d] / tangent_vector[largest_component]) > - std::numeric_limits::epsilon()) + std::numeric_limits::epsilon()) constraints.add_constraint( dof_indices.dof_indices[d], {{dof_indices.dof_indices[largest_component], @@ -596,13 +596,13 @@ namespace VectorTools * Compute the mappings from vector degrees of freedom to normal vectors @p dof_to_normals_map * and vector degrees of freedom to prescribed normal fluxes @p dof_vector_to_b_values. */ - template + template void map_dofs_to_normal_vectors_and_normal_fluxes( const typename DoFHandler::cell_iterator &cell, const unsigned int first_vector_component, const std::set &boundary_ids, - const std::map *> + const std::map *> &function_map, hp::FEFaceValues &x_fe_face_values, const unsigned int n_dofs, @@ -613,7 +613,7 @@ namespace VectorTools std::pair, typename DoFHandler::cell_iterator>> &dof_to_normals_map, - std::map, Vector> &dof_vector_to_b_values) + std::map, Vector> &dof_vector_to_b_values) { // mapping from (active_fe_index, face_no and local // dof index) to dim vector indices of the same @@ -749,7 +749,7 @@ namespace VectorTools normal_vector /= normal_vector.norm(); const Point &point = fe_values.quadrature_point(i); - Vector b_values(dim); + Vector b_values(dim); function_map.at(*b_id)->vector_value(point, b_values); // now enter the (dofs,(normal_vector,cell)) entry into @@ -781,15 +781,15 @@ namespace VectorTools * compute_nonzero_normal_flux_constraints_on_level() so as to have * separate interfaces for the active and level cells. */ - template + template void compute_nonzero_normal_flux_constraints_active_or_level( const DoFHandler &dof_handler, const unsigned int first_vector_component, const std::set &boundary_ids, - const std::map *> + const std::map *> &function_map, - AffineConstraints &constraints, + AffineConstraints &constraints, const Mapping &mapping, const IndexSet &refinement_edge_indices = IndexSet(), const unsigned int level = numbers::invalid_unsigned_int) @@ -851,7 +851,7 @@ namespace VectorTools VectorDoFTuple, std::pair, typename DoFHandler::cell_iterator>>; - std::map, Vector> dof_vector_to_b_values; + std::map, Vector> dof_vector_to_b_values; DoFToNormalsMap dof_to_normals_map; @@ -1025,8 +1025,8 @@ namespace VectorTools // then construct constraints from this: const VectorDoFTuple &dof_indices = same_dof_range[0]->first; - double normal_value = 0.; - const Vector b_values = + number normal_value = 0.; + const Vector b_values = dof_vector_to_b_values[dof_indices]; for (unsigned int i = 0; i < dim; ++i) normal_value += b_values[i] * normal[i]; @@ -1080,7 +1080,7 @@ namespace VectorTools // ignore dofs already constrained const VectorDoFTuple &dof_indices = same_dof_range[0]->first; - const Vector b_values = + const Vector b_values = dof_vector_to_b_values[dof_indices]; for (unsigned int i = 0; i < dim; ++i) if (!constraints.is_constrained( @@ -1092,7 +1092,7 @@ namespace VectorTools dof_indices.dof_indices[i], {}, (std::fabs(b_values[i]) > - std::numeric_limits::epsilon() ? + std::numeric_limits::epsilon() ? b_values[i] : 0)); } @@ -1249,7 +1249,7 @@ namespace VectorTools // the vector is parallel to the tangent const VectorDoFTuple &dof_indices = same_dof_range[0]->first; - const Vector b_values = + const Vector b_values = dof_vector_to_b_values[dof_indices]; add_tangentiality_constraints(dof_indices, average_tangent, @@ -1543,17 +1543,18 @@ namespace VectorTools } - template + template void compute_no_normal_flux_constraints( const DoFHandler &dof_handler, const unsigned int first_vector_component, const std::set &boundary_ids, - AffineConstraints &constraints, + AffineConstraints &constraints, const Mapping &mapping) { - Functions::ZeroFunction zero_function(dim); - std::map *> function_map; + Functions::ZeroFunction zero_function(dim); + std::map *> + function_map; for (const types::boundary_id boundary_id : boundary_ids) function_map[boundary_id] = &zero_function; internal::compute_nonzero_normal_flux_constraints_active_or_level( diff --git a/source/numerics/vector_tools_constraints.inst.in b/source/numerics/vector_tools_constraints.inst.in index 76eeb8e975..babc0e82af 100644 --- a/source/numerics/vector_tools_constraints.inst.in +++ b/source/numerics/vector_tools_constraints.inst.in @@ -64,6 +64,14 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) AffineConstraints &constraints, const Mapping &mapping); + template void + compute_no_normal_flux_constraints( + const DoFHandler &dof_handler, + const unsigned int first_vector_component, + const std::set &boundary_ids, + AffineConstraints &constraints, + const Mapping &mapping); + template void compute_no_normal_flux_constraints_on_level( const DoFHandler &dof_handler, -- 2.39.5