From c8a10a6e4441cdc9721ff273b84b34f72c867945 Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 3 Sep 2016 13:08:50 -0400 Subject: [PATCH] Use consistent dof handler template parameters. --- doc/news/changes.h | 8 + include/deal.II/numerics/vector_tools.h | 163 ++++++------ .../deal.II/numerics/vector_tools.templates.h | 244 +++++++++--------- 3 files changed, 222 insertions(+), 193 deletions(-) diff --git a/doc/news/changes.h b/doc/news/changes.h index 568ec07d8a..a45c9e2c5a 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -38,6 +38,14 @@ inconvenience this causes.

    +
  1. Changed: The template parameter order in many VectorTools functions is now + different; this was done so that the order is the same across similar functions. + This will only effect code that explicitly specifies template parameters for + overloaded VectorTools functions (no known deal.II-based projects do this). +
    + (David Wells, 2016/09/06) +
  2. +
  3. New: There is now the possibility to store information about the time of an output time step within the .visit file created by the DataOutInterface::write_visit_record function. diff --git a/include/deal.II/numerics/vector_tools.h b/include/deal.II/numerics/vector_tools.h index 2dc321075d..4b24acebca 100644 --- a/include/deal.II/numerics/vector_tools.h +++ b/include/deal.II/numerics/vector_tools.h @@ -583,10 +583,11 @@ namespace VectorTools * @todo The @p mapping argument should be replaced by a * hp::MappingCollection in case of a hp::DoFHandler. */ - template class DoFHandlerType> + template class DoFHandlerType> void interpolate (const Mapping &mapping, const DoFHandlerType &dof, - const Function &function, + const Function &function, VectorType &vec, const ComponentMask &component_mask = ComponentMask()); @@ -594,11 +595,13 @@ namespace VectorTools * Call the @p interpolate() function above with * mapping=MappingQGeneric1@@(). */ - template - void interpolate (const DoFHandlerType &dof, - const Function &function, - VectorType &vec, - const ComponentMask &component_mask = ComponentMask()); + template class DoFHandlerType> + void interpolate + (const DoFHandlerType &dof, + const Function &function, + VectorType &vec, + const ComponentMask &component_mask = ComponentMask()); /** * Interpolate different finite element spaces. The interpolation of vector @@ -669,14 +672,16 @@ namespace VectorTools * * @author Valentin Zingan, 2013 */ - template + template class DoFHandlerType> void interpolate_based_on_material_id - (const Mapping &mapping, - const DoFHandlerType &dof_handler, - const std::map *> &function_map, - VectorType &dst, - const ComponentMask &component_mask = ComponentMask()); + (const Mapping &mapping, + const DoFHandlerType &dof_handler, + const std::map *> &function_map, + VectorType &dst, + const ComponentMask &component_mask = ComponentMask()); /** * Gives the interpolation of a @p dof1-function @p u1 to a @p dof2-function @@ -702,9 +707,8 @@ namespace VectorTools * parallel::distributed::Triangulation::no_automatic_repartitioning * flag). */ - template class DoFHandlerType, - typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const DoFHandlerType &dof1, const VectorType &u1, @@ -724,9 +728,8 @@ namespace VectorTools * Without it - due to cellwise interpolation - the resulting output vector * does not necessarily respect continuity requirements at hanging nodes. */ - template class DoFHandlerType, - typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const DoFHandlerType &dof1, const VectorType &u1, @@ -742,9 +745,8 @@ namespace VectorTools * @p intergridmap has to be initialized via InterGridMap::make_mapping * pointing from a source DoFHandler to a destination DoFHandler. */ - template class DoFHandlerType, - typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const InterGridMap > &intergridmap, @@ -903,14 +905,15 @@ namespace VectorTools * * See the general documentation of this namespace for more information. */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const std::map*> &function_map, - std::map &boundary_values, - const ComponentMask &component_mask = ComponentMask()); + (const Mapping &mapping, + const DoFHandlerType &dof, + const std::map*> &function_map, + std::map &boundary_values, + const ComponentMask &component_mask = ComponentMask()); /** * Like the previous function, but take a mapping collection to go with the @@ -934,15 +937,16 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - std::map &boundary_values, - const ComponentMask &component_mask = ComponentMask()); + (const Mapping &mapping, + const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + std::map &boundary_values, + const ComponentMask &component_mask = ComponentMask()); /** * Call the other interpolate_boundary_values() function, see above, with @@ -953,14 +957,15 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - std::map &boundary_values, - const ComponentMask &component_mask = ComponentMask()); + (const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + std::map &boundary_values, + const ComponentMask &component_mask = ComponentMask()); /** @@ -969,13 +974,14 @@ namespace VectorTools * apply as for the previous function, in particular about the use of the * component mask and the requires size of the function object. */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const std::map*> &function_map, - std::map &boundary_values, - const ComponentMask &component_mask = ComponentMask()); + (const DoFHandlerType &dof, + const std::map*> &function_map, + std::map &boundary_values, + const ComponentMask &component_mask = ComponentMask()); /** @@ -1039,14 +1045,15 @@ namespace VectorTools * * @ingroup constraints */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const std::map*> &function_map, - ConstraintMatrix &constraints, - const ComponentMask &component_mask = ComponentMask()); + (const Mapping &mapping, + const DoFHandlerType &dof, + const std::map*> &function_map, + ConstraintMatrix &constraints, + const ComponentMask &component_mask = ComponentMask()); /** * Same function as above, but taking only one pair of boundary indicator @@ -1059,15 +1066,16 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - ConstraintMatrix &constraints, - const ComponentMask &component_mask = ComponentMask()); + (const Mapping &mapping, + const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + ConstraintMatrix &constraints, + const ComponentMask &component_mask = ComponentMask()); /** * Call the other interpolate_boundary_values() function, see above, with @@ -1080,14 +1088,15 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - ConstraintMatrix &constraints, - const ComponentMask &component_mask = ComponentMask()); + (const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + ConstraintMatrix &constraints, + const ComponentMask &component_mask = ComponentMask()); /** @@ -1098,13 +1107,14 @@ namespace VectorTools * * @ingroup constraints */ - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const std::map*> &function_map, - ConstraintMatrix &constraints, - const ComponentMask &component_mask = ComponentMask()); + (const DoFHandlerType &dof, + const std::map*> &function_map, + ConstraintMatrix &constraints, + const ComponentMask &component_mask = ComponentMask()); /** @@ -1762,7 +1772,7 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_nonzero_normal_flux_constraints (const DoFHandlerType &dof_handler, @@ -1783,7 +1793,7 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_no_normal_flux_constraints (const DoFHandlerType &dof_handler, @@ -1808,7 +1818,7 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_nonzero_tangential_flux_constraints (const DoFHandlerType &dof_handler, @@ -1826,7 +1836,7 @@ namespace VectorTools * @see * @ref GlossBoundaryIndicator "Glossary entry on boundary indicators" */ - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_normal_flux_constraints (const DoFHandlerType &dof_handler, @@ -2816,8 +2826,9 @@ namespace VectorTools * * @author Luca Heltai, 2015 */ - template - void get_position_vector(const DoFHandlerType &dh, + template class DoFHandlerType, + typename VectorType> + void get_position_vector(const DoFHandlerType &dh, VectorType &vector, const ComponentMask &mask = ComponentMask()); diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 5e7be2ce0d..3fa704b982 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -75,7 +75,7 @@ DEAL_II_NAMESPACE_OPEN namespace VectorTools { - template class DoFHandlerType> void interpolate (const Mapping &mapping, const DoFHandlerType &dof, @@ -279,14 +279,14 @@ namespace VectorTools } - template - void interpolate (const DoFHandlerType &dof, - const Function &function, - VectorType &vec, - const ComponentMask &component_mask) + template class DoFHandlerType> + void interpolate (const DoFHandlerType &dof, + const Function &function, + VectorType &vec, + const ComponentMask &component_mask) { - interpolate(StaticMappingQ1::mapping, + interpolate(StaticMappingQ1::mapping, dof, function, vec, @@ -348,17 +348,17 @@ namespace VectorTools } - template + template class DoFHandlerType> void interpolate_based_on_material_id - (const Mapping &mapping, - const DoFHandlerType &dof, - const std::map *> &function_map, - VectorType &dst, - const ComponentMask &component_mask) + (const Mapping &mapping, + const DoFHandlerType &dof, + const std::map *> &function_map, + VectorType &dst, + const ComponentMask &component_mask) { typedef typename VectorType::value_type number; - const unsigned int dim = DoFHandlerType::dimension; Assert( component_mask.represents_n_components(dof.get_fe().n_components()), ExcMessage("The number of components in the mask has to be either " @@ -372,7 +372,7 @@ namespace VectorTools ExcMessage("You cannot specify the invalid material indicator " "in your function map.")); - for (typename std::map* > + for (typename std::map* > ::const_iterator iter = function_map.begin(); iter != function_map.end(); @@ -382,13 +382,14 @@ namespace VectorTools ExcDimensionMismatch(dof.get_fe().n_components(), iter->second->n_components) ); } - const hp::FECollection + const hp::FECollection fe(dof.get_fe()); const unsigned int n_components = fe.n_components(); const bool fe_is_system = (n_components != 1); - typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(), - endc = dof.end(); + typename DoFHandlerType::active_cell_iterator + cell = dof.begin_active(), + endc = dof.end(); std::vector< std::vector< Point > > unit_support_points(fe.size()); for (unsigned int fe_index = 0; fe_index < fe.size(); ++fe_index) @@ -433,7 +434,7 @@ namespace VectorTools const unsigned int max_rep_points = *std::max_element(n_rep_points.begin(), n_rep_points.end()); std::vector< types::global_dof_index> dofs_on_cell(fe.max_dofs_per_cell()); - std::vector< Point > rep_points(max_rep_points); + std::vector< Point > rep_points(max_rep_points); std::vector< std::vector > function_values_scalar(fe.size()); std::vector< std::vector< Vector > > function_values_system(fe.size()); @@ -442,11 +443,11 @@ namespace VectorTools for (unsigned int fe_index = 0; fe_index < fe.size(); ++fe_index) support_quadrature.push_back( Quadrature(unit_support_points[fe_index]) ); - hp::MappingCollection mapping_collection(mapping); - hp::FEValues fe_values(mapping_collection, - fe, - support_quadrature, - update_quadrature_points); + hp::MappingCollection mapping_collection(mapping); + hp::FEValues fe_values(mapping_collection, + fe, + support_quadrature, + update_quadrature_points); for ( ; cell != endc; ++cell) if ( cell->is_locally_owned() ) @@ -456,7 +457,7 @@ namespace VectorTools fe_values.reinit(cell); - const std::vector< Point > &support_points + const std::vector< Point > &support_points = fe_values.get_present_fe_values().get_quadrature_points(); rep_points.resize( dofs_of_rep_points[fe_index].size() ); @@ -509,13 +510,12 @@ namespace VectorTools * mapping here because the function we evaluate for the DoFs is zero in * the mapped locations as well as in the original, unmapped locations */ - template + template class DoFHandlerType, + typename number> void - interpolate_zero_boundary_values (const DoFHandlerType &dof_handler, + interpolate_zero_boundary_values (const DoFHandlerType &dof_handler, std::map &boundary_values) { - const unsigned int dim = DoFHandlerType::dimension; - // loop over all boundary faces // to get all dof indices of // dofs on the boundary. note @@ -537,7 +537,7 @@ namespace VectorTools // that is actually wholly on // the boundary, not only by // one line or one vertex - typename DoFHandlerType::active_cell_iterator + typename DoFHandlerType::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end(); std::vector face_dof_indices; @@ -562,7 +562,8 @@ namespace VectorTools - template class DoFHandlerType, typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const DoFHandlerType &dof1, const VectorType &u1, @@ -584,7 +585,8 @@ namespace VectorTools - template class DoFHandlerType, typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const DoFHandlerType &dof1, const VectorType &u1, @@ -622,7 +624,8 @@ namespace VectorTools } - template class DoFHandlerType, typename VectorType> + template class DoFHandlerType> void interpolate_to_different_mesh (const InterGridMap > &intergridmap, @@ -1685,19 +1688,18 @@ namespace VectorTools namespace { - template class M_or_MC> + template class DoFHandlerType, + template class M_or_MC> static inline void do_interpolate_boundary_values - (const M_or_MC &mapping, - const DoFHandlerType &dof, - const std::map*> &function_map, - std::map &boundary_values, - const ComponentMask &component_mask) + (const M_or_MC &mapping, + const DoFHandlerType &dof, + const std::map*> &function_map, + std::map &boundary_values, + const ComponentMask &component_mask) { - const unsigned int dim = DoFHandlerType::dimension; - const unsigned int spacedim=DoFHandlerType::space_dimension; - Assert (component_mask.represents_n_components(dof.get_fe().n_components()), ExcMessage ("The number of components in the mask has to be either " "zero or equal to the number of components in the finite " @@ -1726,7 +1728,8 @@ namespace VectorTools // individual vertices if (dim == 1) { - for (typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(); + for (typename DoFHandlerType::active_cell_iterator + cell = dof.begin_active(); cell != dof.end(); ++cell) for (unsigned int direction=0; direction::faces_per_cell; ++direction) @@ -1734,7 +1737,7 @@ namespace VectorTools && (function_map.find(cell->face(direction)->boundary_id()) != function_map.end())) { - const Function &boundary_function + const Function &boundary_function = *function_map.find(cell->face(direction)->boundary_id())->second; // get the FE corresponding to this cell @@ -1840,8 +1843,9 @@ namespace VectorTools dealii::hp::FEFaceValues x_fe_values (mapping_collection, finite_elements, q_collection, update_quadrature_points); - typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(), - endc = dof.end(); + typename DoFHandlerType::active_cell_iterator + cell = dof.begin_active(), + endc = dof.end(); for (; cell!=endc; ++cell) if (!cell->is_artificial()) for (unsigned int face_no = 0; face_no < GeometryInfo::faces_per_cell; @@ -1869,7 +1873,8 @@ namespace VectorTools "elements")); } - const typename DoFHandlerType::face_iterator face = cell->face(face_no); + const typename DoFHandlerType::face_iterator + face = cell->face(face_no); const types::boundary_id boundary_component = face->boundary_id(); // see if this face is part of the boundaries for which we are @@ -1980,14 +1985,15 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const std::map*> &function_map, - std::map &boundary_values, - const ComponentMask &component_mask_) + (const Mapping &mapping, + const DoFHandlerType &dof, + const std::map*> &function_map, + std::map &boundary_values, + const ComponentMask &component_mask_) { do_interpolate_boundary_values (mapping, dof, function_map, boundary_values, component_mask_); @@ -1995,17 +2001,18 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - std::map &boundary_values, - const ComponentMask &component_mask) + (const Mapping &mapping, + const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + std::map &boundary_values, + const ComponentMask &component_mask) { - std::map*> function_map; + std::map*> function_map; function_map[boundary_component] = &boundary_function; interpolate_boundary_values (mapping, dof, function_map, boundary_values, component_mask); @@ -2027,32 +2034,34 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - std::map &boundary_values, - const ComponentMask &component_mask) + (const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + std::map &boundary_values, + const ComponentMask &component_mask) { - interpolate_boundary_values(StaticMappingQ1::mapping, + interpolate_boundary_values(StaticMappingQ1::mapping, dof, boundary_component, boundary_function, boundary_values, component_mask); } - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const std::map*> &function_map, - std::map &boundary_values, - const ComponentMask &component_mask) + (const DoFHandlerType &dof, + const std::map*> &function_map, + std::map &boundary_values, + const ComponentMask &component_mask) { interpolate_boundary_values - (StaticMappingQ1::mapping, + (StaticMappingQ1::mapping, dof, function_map, boundary_values, component_mask); } @@ -2064,14 +2073,15 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const std::map*> &function_map, - ConstraintMatrix &constraints, - const ComponentMask &component_mask_) + (const Mapping &mapping, + const DoFHandlerType &dof, + const std::map*> &function_map, + ConstraintMatrix &constraints, + const ComponentMask &component_mask_) { std::map boundary_values; interpolate_boundary_values (mapping, dof, function_map, @@ -2093,17 +2103,18 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const Mapping &mapping, - const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - ConstraintMatrix &constraints, - const ComponentMask &component_mask) + (const Mapping &mapping, + const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + ConstraintMatrix &constraints, + const ComponentMask &component_mask) { - std::map*> function_map; + std::map*> function_map; function_map[boundary_component] = &boundary_function; interpolate_boundary_values (mapping, dof, function_map, constraints, component_mask); @@ -2111,33 +2122,35 @@ namespace VectorTools - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const types::boundary_id boundary_component, - const Function &boundary_function, - ConstraintMatrix &constraints, - const ComponentMask &component_mask) + (const DoFHandlerType &dof, + const types::boundary_id boundary_component, + const Function &boundary_function, + ConstraintMatrix &constraints, + const ComponentMask &component_mask) { interpolate_boundary_values - (StaticMappingQ1::mapping, + (StaticMappingQ1::mapping, dof, boundary_component, boundary_function, constraints, component_mask); } - template + template class DoFHandlerType, + typename number> void interpolate_boundary_values - (const DoFHandlerType &dof, - const std::map*> &function_map, - ConstraintMatrix &constraints, - const ComponentMask &component_mask) + (const DoFHandlerType &dof, + const std::map*> &function_map, + ConstraintMatrix &constraints, + const ComponentMask &component_mask) { interpolate_boundary_values - (StaticMappingQ1::mapping, + (StaticMappingQ1::mapping, dof, function_map, constraints, component_mask); } @@ -2225,7 +2238,7 @@ namespace VectorTools void do_project_boundary_values (const M_or_MC &mapping, - const DoFHandlerType &dof, + const DoFHandlerType &dof, const std::map*> &boundary_functions, const Q_or_QC &q, std::map &boundary_values, @@ -5285,7 +5298,7 @@ namespace VectorTools - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_no_normal_flux_constraints (const DoFHandlerType &dof_handler, const unsigned int first_vector_component, @@ -5307,7 +5320,7 @@ namespace VectorTools mapping); } - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_nonzero_normal_flux_constraints (const DoFHandlerType &dof_handler, @@ -5872,7 +5885,7 @@ namespace VectorTools - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_normal_flux_constraints (const DoFHandlerType &dof_handler, const unsigned int first_vector_component, @@ -5894,7 +5907,7 @@ namespace VectorTools mapping); } - template class DoFHandlerType, int spacedim> + template class DoFHandlerType> void compute_nonzero_tangential_flux_constraints (const DoFHandlerType &dof_handler, @@ -7288,18 +7301,15 @@ namespace VectorTools } - template - void get_position_vector(const DoFHandlerType &dh, - VectorType &vector, - const ComponentMask &mask) + template class DoFHandlerType, + typename VectorType> + void get_position_vector(const DoFHandlerType &dh, + VectorType &vector, + const ComponentMask &mask) { AssertDimension(vector.size(), dh.n_dofs()); - - const unsigned int dim=DoFHandlerType::dimension; - const unsigned int spacedim=DoFHandlerType::space_dimension; const FiniteElement &fe = dh.get_fe(); - // Construct default fe_mask; const ComponentMask fe_mask(mask.size() ? mask : ComponentMask(fe.get_nonzero_components(0).size(), true)); @@ -7319,7 +7329,7 @@ namespace VectorTools if ( fe.has_support_points() ) { - typename DoFHandlerType::active_cell_iterator cell; + typename DoFHandlerType::active_cell_iterator cell; const Quadrature quad(fe.get_unit_support_points()); MappingQ map_q(fe.degree); @@ -7372,7 +7382,7 @@ namespace VectorTools // carefully selecting the right components. FESystem feq(FE_Q(degree), spacedim); - DoFHandlerType dhq(dh.get_triangulation()); + DoFHandlerType dhq(dh.get_triangulation()); dhq.distribute_dofs(feq); Vector eulerq(dhq.n_dofs()); const ComponentMask maskq(spacedim, true); -- 2.39.5