From: Niklas Wik Date: Sat, 4 Jun 2022 15:10:36 +0000 (+0200) Subject: Add typename and instantiation X-Git-Tag: v9.4.0-rc1~69^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2d20b7de4cbebe1bc3b3cc5d180c96bded3d6710;p=dealii.git Add typename and instantiation --- diff --git a/include/deal.II/numerics/vector_tools_boundary.h b/include/deal.II/numerics/vector_tools_boundary.h index 55d07a5ef1..0d139e1b25 100644 --- a/include/deal.II/numerics/vector_tools_boundary.h +++ b/include/deal.II/numerics/vector_tools_boundary.h @@ -718,15 +718,15 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template void project_boundary_values_div_conforming( - const DoFHandler & dof_handler, - const unsigned int first_vector_component, - const Function &boundary_function, - const types::boundary_id boundary_component, - AffineConstraints & constraints, - const Mapping & mapping); + const DoFHandler & dof_handler, + const unsigned int first_vector_component, + const Function &boundary_function, + const types::boundary_id boundary_component, + AffineConstraints & constraints, + const Mapping & mapping); /** * Same as above for the hp-namespace. @@ -736,14 +736,14 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template void project_boundary_values_div_conforming( const DoFHandler & dof_handler, const unsigned int first_vector_component, - const Function & boundary_function, + const Function & boundary_function, const types::boundary_id boundary_component, - AffineConstraints & constraints, + AffineConstraints & constraints, const hp::MappingCollection &mapping_collection = hp::StaticMappingQ1::mapping_collection); diff --git a/include/deal.II/numerics/vector_tools_boundary.templates.h b/include/deal.II/numerics/vector_tools_boundary.templates.h index ec71f01e86..b84f413ac0 100644 --- a/include/deal.II/numerics/vector_tools_boundary.templates.h +++ b/include/deal.II/numerics/vector_tools_boundary.templates.h @@ -2149,16 +2149,16 @@ namespace VectorTools { // This function computes the projection of the boundary function on the // boundary in 2d. - template + template void compute_face_projection_div_conforming( const cell_iterator & cell, const unsigned int face, const FEFaceValues<2> & fe_values, const unsigned int first_vector_component, - const Function<2> & boundary_function, + const Function<2, number2> & boundary_function, const std::vector> &jacobians, - AffineConstraints & constraints) + AffineConstraints & constraints) { // Compute the integral over the product of the normal components of // the boundary function times the normal components of the shape @@ -2171,9 +2171,9 @@ namespace VectorTools 1, 0, 0}; - std::vector> values(fe_values.n_quadrature_points, - Vector(2)); - Vector dof_values(fe.n_dofs_per_face(face)); + std::vector> values(fe_values.n_quadrature_points, + Vector(2)); + Vector dof_values(fe.n_dofs_per_face(face)); // Get the values of the boundary function at the quadrature points. { @@ -2186,7 +2186,7 @@ namespace VectorTools for (unsigned int q_point = 0; q_point < fe_values.n_quadrature_points; ++q_point) { - double tmp = 0.0; + number2 tmp = 0.0; for (unsigned int d = 0; d < 2; ++d) tmp += normals[q_point][d] * values[q_point](d); @@ -2235,32 +2235,35 @@ namespace VectorTools } // dummy implementation of above function for all other dimensions - template + template void compute_face_projection_div_conforming( const cell_iterator &, const unsigned int, const FEFaceValues &, const unsigned int, - const Function &, + const Function &, const std::vector> &, - AffineConstraints &) + AffineConstraints &) { Assert(false, ExcNotImplemented()); } // This function computes the projection of the boundary function on the // boundary in 3d. - template + template void compute_face_projection_div_conforming( const cell_iterator & cell, const unsigned int face, const FEFaceValues<3> & fe_values, const unsigned int first_vector_component, - const Function<3> & boundary_function, + const Function<3, number2> & boundary_function, const std::vector> &jacobians, - std::vector & dof_values, + std::vector & dof_values, std::vector & projected_dofs) { // Compute the intergral over the product of the normal components of @@ -2272,9 +2275,9 @@ namespace VectorTools const unsigned int face_coordinate_directions[GeometryInfo<3>::faces_per_cell][2] = { {1, 2}, {1, 2}, {2, 0}, {2, 0}, {0, 1}, {0, 1}}; - std::vector> values(fe_values.n_quadrature_points, - Vector(3)); - Vector dof_values_local(fe.n_dofs_per_face(face)); + std::vector> values(fe_values.n_quadrature_points, + Vector(3)); + Vector dof_values_local(fe.n_dofs_per_face(face)); { const std::vector> &quadrature_points = @@ -2286,7 +2289,7 @@ namespace VectorTools for (unsigned int q_point = 0; q_point < fe_values.n_quadrature_points; ++q_point) { - double tmp = 0.0; + number2 tmp = 0.0; for (unsigned int d = 0; d < 3; ++d) tmp += normals[q_point][d] * values[q_point](d); @@ -2342,16 +2345,19 @@ namespace VectorTools // dummy implementation of above // function for all other // dimensions - template + template void compute_face_projection_div_conforming( const cell_iterator &, const unsigned int, const FEFaceValues &, const unsigned int, - const Function &, + const Function &, const std::vector> &, - std::vector &, + std::vector &, std::vector &) { Assert(false, ExcNotImplemented()); @@ -2359,15 +2365,15 @@ namespace VectorTools } // namespace internals - template + template void project_boundary_values_div_conforming( - const DoFHandler & dof_handler, - const unsigned int first_vector_component, - const Function & boundary_function, - const types::boundary_id boundary_component, - AffineConstraints &constraints, - const Mapping & mapping) + const DoFHandler & dof_handler, + const unsigned int first_vector_component, + const Function &boundary_function, + const types::boundary_id boundary_component, + AffineConstraints & constraints, + const Mapping & mapping) { const unsigned int spacedim = dim; // Interpolate the normal components @@ -2430,7 +2436,9 @@ namespace VectorTools { AssertThrow( dynamic_cast *>(&fe) != - nullptr, + nullptr || + dynamic_cast *>( + &fe) != nullptr, typename FiniteElement< dim>::ExcInterpolationNotImplemented()); } @@ -2485,7 +2493,9 @@ namespace VectorTools { AssertThrow( dynamic_cast *>(&fe) != - nullptr, + nullptr || + dynamic_cast *>( + &fe) != nullptr, typename FiniteElement< dim>::ExcInterpolationNotImplemented()); } @@ -2529,14 +2539,14 @@ namespace VectorTools } - template + template void project_boundary_values_div_conforming( const DoFHandler & dof_handler, const unsigned int first_vector_component, - const Function & boundary_function, + const Function & boundary_function, const types::boundary_id boundary_component, - AffineConstraints & constraints, + AffineConstraints & constraints, const hp::MappingCollection &mapping_collection) { const unsigned int spacedim = dim; @@ -2618,7 +2628,7 @@ namespace VectorTools case 3: { const unsigned int n_dofs = dof_handler.n_dofs(); - std::vector dof_values(n_dofs); + std::vector dof_values(n_dofs); std::vector projected_dofs(n_dofs); for (unsigned int dof = 0; dof < n_dofs; ++dof) diff --git a/source/numerics/vector_tools_boundary.inst.in b/source/numerics/vector_tools_boundary.inst.in index 1801335127..0604787d39 100644 --- a/source/numerics/vector_tools_boundary.inst.in +++ b/source/numerics/vector_tools_boundary.inst.in @@ -268,21 +268,37 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) AffineConstraints> &, const hp::MappingCollection &); # endif +# endif +#endif + \} + } + +for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; + number : REAL_SCALARS; + number2 : REAL_SCALARS) + { + namespace VectorTools + \{ +#if deal_II_dimension == deal_II_space_dimension + +# if deal_II_dimension != 1 + template void project_boundary_values_div_conforming( const DoFHandler &, const unsigned int, - const Function &, + const Function &, const types::boundary_id, - AffineConstraints &, + AffineConstraints &, const Mapping &); + template void project_boundary_values_div_conforming( const DoFHandler &, const unsigned int, - const Function &, + const Function &, const types::boundary_id, - AffineConstraints &, + AffineConstraints &, const hp::MappingCollection &); # endif #endif