From: Wolfgang Bangerth Date: Fri, 22 Jan 2021 00:27:59 +0000 (-0700) Subject: Move the get_*mapping() functions into ReferenceCell::Type. X-Git-Tag: v9.3.0-rc1~556^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0476ea8aaaba14b91c55cd736277f4ddacb2e93a;p=dealii.git Move the get_*mapping() functions into ReferenceCell::Type. --- diff --git a/include/deal.II/fe/fe_interface_values.h b/include/deal.II/fe/fe_interface_values.h index eb2bb05eb8..d565602eee 100644 --- a/include/deal.II/fe/fe_interface_values.h +++ b/include/deal.II/fe/fe_interface_values.h @@ -526,26 +526,26 @@ FEInterfaceValues::FEInterfaceValues( const UpdateFlags update_flags) : n_quadrature_points(quadrature.size()) , internal_fe_face_values( - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags) , internal_fe_subface_values( - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags) , internal_fe_face_values_neighbor( - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags) , internal_fe_subface_values_neighbor( - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags) diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index 88c8c3049a..4e3edf8a92 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -1608,8 +1608,7 @@ namespace FETools ReferenceCell::make_triangulation(reference_cell_type, tr); const auto &mapping = - ReferenceCell::get_default_linear_mapping( - reference_cell_type); + reference_cell_type.template get_default_linear_mapping(); // Choose a Gauss quadrature rule that is exact up to degree 2n-1 const unsigned int degree = @@ -1822,8 +1821,8 @@ namespace FETools const unsigned int degree = fe.degree; const auto &mapping = - ReferenceCell::get_default_linear_mapping( - reference_cell_type); + reference_cell_type + .template get_default_linear_mapping(); const auto &q_fine = ReferenceCell::get_gauss_type_quadrature(reference_cell_type, degree + 1); @@ -2291,8 +2290,8 @@ namespace FETools q_points_coarse[q](j) = q_points_fine[q](j); Quadrature q_coarse(q_points_coarse, fine.get_JxW_values()); FEValues coarse( - ReferenceCell::get_default_linear_mapping( - coarse_cell->reference_cell_type()), + coarse_cell->reference_cell_type() + .template get_default_linear_mapping(), fe, q_coarse, update_values); diff --git a/include/deal.II/fe/mapping_q1.h b/include/deal.II/fe/mapping_q1.h index 266cc7739a..d5f65cb854 100644 --- a/include/deal.II/fe/mapping_q1.h +++ b/include/deal.II/fe/mapping_q1.h @@ -83,7 +83,9 @@ public: * * @note The use of this object should be avoided since it is only applicable * in cases where a mesh consists exclusively of quadrilaterals or hexahedra. - * Use ReferenceCell::get_default_linear_mapping() instead. + * Use + * `ReferenceCell::Type::get_hypercube().get_default_linear_mapping()` + * instead. */ template struct StaticMappingQ1 diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index 17661eeb38..8985886aac 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -176,8 +176,8 @@ namespace GridTools double volume(const Triangulation &tria, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Return an approximation of the diameter of the smallest active cell of a @@ -194,8 +194,8 @@ namespace GridTools minimal_cell_diameter( const Triangulation &triangulation, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Return an approximation of the diameter of the largest active cell of a @@ -212,8 +212,8 @@ namespace GridTools maximal_cell_diameter( const Triangulation &triangulation, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Given a list of vertices (typically obtained using @@ -1044,8 +1044,8 @@ namespace GridTools extract_used_vertices( const Triangulation &container, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Find and return the index of the closest vertex to a given point in the @@ -1872,8 +1872,8 @@ namespace GridTools const typename Triangulation::active_cell_iterator &cell, const Point & position, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Compute a globally unique index for each vertex and hanging node diff --git a/include/deal.II/grid/grid_tools_cache.h b/include/deal.II/grid/grid_tools_cache.h index 99112c91f1..a5b48a8371 100644 --- a/include/deal.II/grid/grid_tools_cache.h +++ b/include/deal.II/grid/grid_tools_cache.h @@ -76,8 +76,8 @@ namespace GridTools */ Cache(const Triangulation &tria, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Destructor. diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index 7bf0fadc60..9c969b6cd4 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -173,6 +173,34 @@ namespace ReferenceCell permute_according_orientation(const std::array &vertices, const unsigned int orientation) const; + /** + * Return a default mapping of degree @p degree matching the current + * reference cell. If this reference cell is a hypercube, then the returned + * mapping is a MappingQGeneric; otherwise, it is an object of type + * MappingFE initialized with Simplex::FE_P (if the reference cell is a + * triangle or tetrahedron), with Simplex::FE_PyramidP (if the reference + * cell is a pyramid), or with Simplex::FE_WedgeP (if the reference cell is + * a wedge). + */ + template + std::unique_ptr> + get_default_mapping(const unsigned int degree) const; + + /** + * Return a default linear mapping matching the current reference cell. + * If this reference cell is a hypercube, then the returned mapping + * is a MappingQ1; otherwise, it is an object of type MappingFE + * initialized with Simplex::FE_P (if the reference cell is a triangle or + * tetrahedron), with Simplex::FE_PyramidP (if the reference cell is a + * pyramid), or with Simplex::FE_WedgeP (if the reference cell is a wedge). + * In other words, the term "linear" in the name of the function has to be + * understood as $d$-linear (i.e., bilinear or trilinear) for some of the + * coordinate directions. + */ + template + const Mapping & + get_default_linear_mapping() const; + /** * Return a text representation of the reference cell represented by the * current object. @@ -693,33 +721,6 @@ namespace ReferenceCell make_triangulation(const Type & reference_cell, Triangulation &tria); - /** - * Return a default mapping of degree @ degree matching the given reference - * cell. If this reference cell is a hypercube, then the returned mapping is a - * MappingQGeneric; otherwise, it is an object of type MappingFE initialized - * with Simplex::FE_P (if the reference cell is a triangle and tetrahedron), - * with Simplex::FE_PyramidP (if the reference cell is a pyramid), or with - * Simplex::FE_WedgeP (if the reference cell is a wedge). - */ - template - std::unique_ptr> - get_default_mapping(const Type &reference_cell, const unsigned int degree); - - /** - * Return a default linear mapping matching the given reference cell. - * If this reference cell is a hypercube, then the returned mapping - * is a MappingQ1; otherwise, it is an object of type MappingFE - * initialized with Simplex::FE_P (if the reference cell is a triangle and - * tetrahedron), with Simplex::FE_PyramidP (if the reference cell is a - * pyramid), or with Simplex::FE_WedgeP (if the reference cell is a wedge). In - * other words, the term "linear" in the name of the function has to be - * understood as $d$-linear (i.e., bilinear or trilinear) for some of the - * coordinate directions. - */ - template - const Mapping & - get_default_linear_mapping(const Type &reference_cell); - /** * Return a default linear mapping that works for the given triangulation. * Internally, this function calls the function above for the reference diff --git a/include/deal.II/numerics/error_estimator.templates.h b/include/deal.II/numerics/error_estimator.templates.h index 9bf3177bba..05aa1373ef 100644 --- a/include/deal.II/numerics/error_estimator.templates.h +++ b/include/deal.II/numerics/error_estimator.templates.h @@ -1055,7 +1055,7 @@ KellyErrorEstimator::estimate( const types::material_id material_id, const Strategy strategy) { - estimate(ReferenceCell::get_default_linear_mapping( + estimate(ReferenceCell::get_default_linear_mapping( dof_handler.get_triangulation()), dof_handler, quadrature, @@ -1126,7 +1126,7 @@ KellyErrorEstimator::estimate( const types::material_id material_id, const Strategy strategy) { - estimate(ReferenceCell::get_default_linear_mapping( + estimate(ReferenceCell::get_default_linear_mapping( dof_handler.get_triangulation()), dof_handler, quadrature, @@ -1372,7 +1372,7 @@ KellyErrorEstimator::estimate( const types::material_id material_id, const Strategy strategy) { - estimate(ReferenceCell::get_default_linear_mapping( + estimate(ReferenceCell::get_default_linear_mapping( dof_handler.get_triangulation()), dof_handler, quadrature, @@ -1407,7 +1407,7 @@ KellyErrorEstimator::estimate( const types::material_id material_id, const Strategy strategy) { - estimate(ReferenceCell::get_default_linear_mapping( + estimate(ReferenceCell::get_default_linear_mapping( dof_handler.get_triangulation()), dof_handler, quadrature, diff --git a/include/deal.II/numerics/vector_tools_constraints.h b/include/deal.II/numerics/vector_tools_constraints.h index 3d53a53a4b..b4e5c5a61f 100644 --- a/include/deal.II/numerics/vector_tools_constraints.h +++ b/include/deal.II/numerics/vector_tools_constraints.h @@ -280,8 +280,8 @@ namespace VectorTools & function_map, AffineConstraints & constraints, const Mapping &mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * This function does the same as the @@ -303,8 +303,8 @@ namespace VectorTools const std::set &boundary_ids, AffineConstraints & constraints, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Compute the constraints that correspond to boundary conditions of the @@ -332,8 +332,8 @@ namespace VectorTools & function_map, AffineConstraints & constraints, const Mapping &mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * Same as above for homogeneous tangential-flux constraints. @@ -351,8 +351,8 @@ namespace VectorTools const std::set &boundary_ids, AffineConstraints & constraints, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); //@} } // namespace VectorTools diff --git a/include/deal.II/numerics/vector_tools_interpolate.templates.h b/include/deal.II/numerics/vector_tools_interpolate.templates.h index 80138b7c19..429ae5d81b 100644 --- a/include/deal.II/numerics/vector_tools_interpolate.templates.h +++ b/include/deal.II/numerics/vector_tools_interpolate.templates.h @@ -653,8 +653,9 @@ namespace VectorTools { const Quadrature quad(fe.get_unit_support_points()); - const auto map_q = ReferenceCell::get_default_mapping( - fe.reference_cell_type(), fe.degree); + const auto map_q = + fe.reference_cell_type().template get_default_mapping( + fe.degree); FEValues fe_v(*map_q, fe, quad, diff --git a/include/deal.II/particles/generators.h b/include/deal.II/particles/generators.h index a5dd84ba75..46ef1118b6 100644 --- a/include/deal.II/particles/generators.h +++ b/include/deal.II/particles/generators.h @@ -69,8 +69,8 @@ namespace Particles const std::vector> & particle_reference_locations, ParticleHandler & particle_handler, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * A function that generates one particle at a random location in cell @p cell and with @@ -108,8 +108,8 @@ namespace Particles const types::particle_index id, std::mt19937 & random_number_generator, const Mapping &mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube())); + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping()); /** * A function that generates particles randomly in the domain with a @@ -164,8 +164,8 @@ namespace Particles const types::particle_index n_particles_to_create, ParticleHandler & particle_handler, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube()), + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping(), const unsigned int random_number_seed = 5432); @@ -211,8 +211,8 @@ namespace Particles & global_bounding_boxes, ParticleHandler &particle_handler, const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube()), + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping(), const ComponentMask & components = ComponentMask(), const std::vector> &properties = {}); @@ -249,16 +249,15 @@ namespace Particles */ template void - quadrature_points( - const Triangulation &triangulation, - const Quadrature & quadrature, - const std::vector>> - & global_bounding_boxes, - ParticleHandler &particle_handler, - const Mapping & mapping = - ReferenceCell::get_default_linear_mapping( - ReferenceCell::Type::get_hypercube()), - const std::vector> &properties = {}); + quadrature_points(const Triangulation &triangulation, + const Quadrature & quadrature, + const std::vector>> + & global_bounding_boxes, + ParticleHandler &particle_handler, + const Mapping & mapping = + ReferenceCell::Type::get_hypercube() + .template get_default_linear_mapping(), + const std::vector> &properties = {}); } // namespace Generators } // namespace Particles diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 4f5348d568..d6f1b76beb 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -4388,8 +4388,8 @@ FEValues::FEValues(const FiniteElement &fe, q.size(), fe.n_dofs_per_cell(), update_default, - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe) , quadrature(q) { @@ -4714,8 +4714,8 @@ FEFaceValues::FEFaceValues( : FEFaceValuesBase( fe.n_dofs_per_cell(), update_flags, - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature) { @@ -4946,8 +4946,8 @@ FESubfaceValues::FESubfaceValues( : FEFaceValuesBase( fe.n_dofs_per_cell(), update_flags, - ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature) { diff --git a/source/grid/reference_cell.cc b/source/grid/reference_cell.cc index 7eeebe9f07..2b6c2a634e 100644 --- a/source/grid/reference_cell.cc +++ b/source/grid/reference_cell.cc @@ -179,19 +179,19 @@ namespace ReferenceCell template std::unique_ptr> - get_default_mapping(const Type &reference_cell, const unsigned int degree) + Type::get_default_mapping(const unsigned int degree) const { - AssertDimension(dim, reference_cell.get_dimension()); + AssertDimension(dim, get_dimension()); - if (reference_cell == Type::get_hypercube()) + if (*this == Type::get_hypercube()) return std::make_unique>(degree); - else if (reference_cell == Type::Tri || reference_cell == Type::Tet) + else if (*this == Type::Tri || *this == Type::Tet) return std::make_unique>( Simplex::FE_P(degree)); - else if (reference_cell == Type::Pyramid) + else if (*this == Type::Pyramid) return std::make_unique>( Simplex::FE_PyramidP(degree)); - else if (reference_cell == Type::Wedge) + else if (*this == Type::Wedge) return std::make_unique>( Simplex::FE_WedgeP(degree)); else @@ -203,29 +203,30 @@ namespace ReferenceCell } + template const Mapping & - get_default_linear_mapping(const Type &reference_cell) + Type::get_default_linear_mapping() const { - AssertDimension(dim, reference_cell.get_dimension()); + AssertDimension(dim, get_dimension()); - if (reference_cell == Type::get_hypercube()) + if (*this == Type::get_hypercube()) { return StaticMappingQ1::mapping; } - else if (reference_cell == Type::Tri || reference_cell == Type::Tet) + else if (*this == Type::Tri || *this == Type::Tet) { static const MappingFE mapping( Simplex::FE_P(1)); return mapping; } - else if (reference_cell == Type::Pyramid) + else if (*this == Type::Pyramid) { static const MappingFE mapping( Simplex::FE_PyramidP(1)); return mapping; } - else if (reference_cell == Type::Wedge) + else if (*this == Type::Wedge) { static const MappingFE mapping( Simplex::FE_WedgeP(1)); @@ -255,8 +256,8 @@ namespace ReferenceCell "cells of the triangulation. The triangulation you are " "passing to this function uses multiple cell types.")); - return get_default_linear_mapping( - reference_cell_types.front()); + return reference_cell_types.front() + .template get_default_linear_mapping(); } diff --git a/source/grid/reference_cell.inst.in b/source/grid/reference_cell.inst.in index e5b8045cdf..b232ba517a 100644 --- a/source/grid/reference_cell.inst.in +++ b/source/grid/reference_cell.inst.in @@ -22,10 +22,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) template std::unique_ptr< Mapping> - get_default_mapping(const Type &reference_cell, const unsigned int degree); + Type::get_default_mapping(const unsigned int degree) const; template const Mapping - &get_default_linear_mapping(const Type &reference_cell); + &Type::get_default_linear_mapping() const; template const Mapping &get_default_linear_mapping( diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 89a6c94cfb..955b9131ab 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -9611,9 +9611,9 @@ namespace internal // to check it, transform to the unit cell // with a linear mapping const Point new_unit = - ReferenceCell::get_default_linear_mapping( - cell->reference_cell_type()) + cell->reference_cell_type() + .template get_default_linear_mapping() .transform_real_to_unit_cell(cell, new_bound); // Now, we have to calculate the distance from diff --git a/source/grid/tria_accessor.cc b/source/grid/tria_accessor.cc index bf6f98c58f..4a79123df3 100644 --- a/source/grid/tria_accessor.cc +++ b/source/grid/tria_accessor.cc @@ -1857,8 +1857,8 @@ CellAccessor<3>::point_inside(const Point<3> &p) const { const TriaRawIterator> cell_iterator(*this); return (GeometryInfo::is_inside_unit_cell( - ReferenceCell::get_default_linear_mapping( - reference_cell_type()) + reference_cell_type() + .template get_default_linear_mapping() .transform_real_to_unit_cell(cell_iterator, p))); } catch (const Mapping::ExcTransformationFailed &) diff --git a/source/meshworker/scratch_data.cc b/source/meshworker/scratch_data.cc index 5fbd2d2077..99e4138660 100644 --- a/source/meshworker/scratch_data.cc +++ b/source/meshworker/scratch_data.cc @@ -74,8 +74,8 @@ namespace MeshWorker const UpdateFlags & update_flags, const Quadrature & face_quadrature, const UpdateFlags & face_update_flags) - : ScratchData(ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + : ScratchData(fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags, @@ -94,8 +94,8 @@ namespace MeshWorker const Quadrature & face_quadrature, const UpdateFlags & face_update_flags, const UpdateFlags & neighbor_face_update_flags) - : ScratchData(ReferenceCell::get_default_linear_mapping( - fe.reference_cell_type()), + : ScratchData(fe.reference_cell_type() + .template get_default_linear_mapping(), fe, quadrature, update_flags, diff --git a/source/numerics/derivative_approximation.cc b/source/numerics/derivative_approximation.cc index 0af8face47..fb2e29998b 100644 --- a/source/numerics/derivative_approximation.cc +++ b/source/numerics/derivative_approximation.cc @@ -1101,8 +1101,8 @@ namespace DerivativeApproximation { // just call the respective function with Q1 mapping approximate_derivative_tensor( - ReferenceCell::get_default_linear_mapping( - cell->reference_cell_type()), + cell->reference_cell_type() + .template get_default_linear_mapping(), dof, solution, cell,