From f921dd7b1ffff881b03cdbf8be7debcc893639a1 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 8 Aug 2001 08:31:19 +0000 Subject: [PATCH] Remove the kludge that the QProjector happened to become a real class derived from Quadrature. Rather, make the respective constructor call two functions in their own that return a Quadrature object. git-svn-id: https://svn.dealii.org/trunk@4869 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/quadrature.h | 111 ++++++++++------- deal.II/base/source/quadrature.cc | 116 ++++++++++-------- deal.II/deal.II/source/fe/fe.cc | 8 +- .../deal.II/source/fe/mapping_cartesian.cc | 8 +- deal.II/deal.II/source/fe/mapping_q.cc | 4 +- deal.II/deal.II/source/fe/mapping_q1.cc | 12 +- deal.II/deal.II/source/grid/grid_out.cc | 6 +- tests/base/quadrature_test.cc | 4 +- 8 files changed, 155 insertions(+), 114 deletions(-) diff --git a/deal.II/base/include/base/quadrature.h b/deal.II/base/include/base/quadrature.h index 53802bce32..56bb445fbc 100644 --- a/deal.II/base/include/base/quadrature.h +++ b/deal.II/base/include/base/quadrature.h @@ -276,55 +276,30 @@ class QIterated : public Quadrature * for a description of the orientation of the different faces. */ template -class QProjector : public Quadrature +class QProjector { public: -//TODO:[GK] remove this constructor again and return to old state with static-only functions. move this function to something in the finite element classes which is the only place where it is used. - /** - * Constructor for a quadrature rule on all (sub)faces. - * The quadrature rule - * @p{quadrature} is applied to - * each face or subface, according - * to the @{sub} flag. - * - * The weights of the new rule are - * replications of the original - * weights. This is not a proper - * handling, but it is - * consistent with later use. If - * the (sub)face rule is applied to - * the unity function, the result - * is the number of (sub)faces. - */ - QProjector (const Quadrature& quadrature, - const bool sub); - + /** - * Compute the quadrature points on the - * cell if the given quadrature formula - * is used on face @p{face_no}. For further - * details, see the general doc for - * this class. + * Compute the quadrature points + * on the cell if the given + * quadrature formula is used on + * face @p{face_no}. For further + * details, see the general doc + * for this class. */ static void project_to_face (const Quadrature &quadrature, const unsigned int face_no, typename std::vector > &q_points); - /** - * Projection to all faces. - * Generate a formula that integrates - * over all faces at the same time. - */ - static void project_to_faces (const Quadrature &quadrature, - typename std::vector > &q_points); - /** - * Compute the quadrature points on the - * cell if the given quadrature formula - * is used on face @p{face_no}, subface - * number @p{subface_no}. For further - * details, see the general doc for - * this class. + * Compute the quadrature points + * on the cell if the given + * quadrature formula is used on + * face @p{face_no}, subface + * number @p{subface_no}. For + * further details, see the + * general doc for this class. */ static void project_to_subface (const Quadrature &quadrature, const unsigned int face_no, @@ -332,14 +307,46 @@ class QProjector : public Quadrature typename std::vector > &q_points); /** - * Projection to all child faces. - * Project to the children of all - * faces at the same time. The - * ordering is first by face, - * then by subface + * Take a face quadrature formula + * and generate a cell quadrature + * formula from it where the + * quadrature points of the given + * argument are projected on all + * faces. + * + * The weights of the new rule + * are replications of the + * original weights. This is not + * a proper handling, in that the + * sum of weights does not equal + * one, but it is consistent with + * the use of this function, + * namely to generate sets of + * face quadrature points on a + * cell, one set of which will + * then be selected at each + * time. This is used in the + * @ref{FEFaceValues} class, + * where we initialize the + * values, derivatives, etc on + * all faces at once, while + * selecting the data of one + * particular face only happens + * later. + */ + static Quadrature + project_to_all_faces (const Quadrature &quadrature); + + /** + * This function is alike the + * previous one, but projects the + * given face quadrature formula + * to the subfaces of a cell, + * i.e. to the children of the + * faces of the unit cell. */ - static void project_to_subfaces (const Quadrature &quadrature, - typename std::vector > &q_points); + static Quadrature + project_to_all_subfaces (const Quadrature &quadrature); }; @@ -374,6 +381,11 @@ template <> void QProjector<3>::project_to_face (const Quadrature<2> &quadrature, const unsigned int face_no, std::vector > &q_points); + +template <> +Quadrature<1> QProjector<1>::project_to_all_faces (const Quadrature<0> &quadrature); + + template <> void QProjector<1>::project_to_subface (const Quadrature<0> &, const unsigned int, @@ -389,6 +401,11 @@ void QProjector<3>::project_to_subface (const Quadrature<2> &quadrature, const unsigned int face_no, const unsigned int subface_no, std::vector > &q_points); + +template <> +Quadrature<1> QProjector<1>::project_to_all_subfaces (const Quadrature<0> &quadrature); + + template <> bool QIterated<1>::uses_both_endpoints (const Quadrature<1> &base_quadrature); diff --git a/deal.II/base/source/quadrature.cc b/deal.II/base/source/quadrature.cc index 18ae378a93..62ab8d2fd9 100644 --- a/deal.II/base/source/quadrature.cc +++ b/deal.II/base/source/quadrature.cc @@ -222,27 +222,6 @@ Quadrature::memory_consumption () const //----------------------------------------------------------------------// -template -QProjector::QProjector (const Quadrature &quadrature, - const bool sub) - : - Quadrature (quadrature.n_quadrature_points - * GeometryInfo::faces_per_cell - * (sub ? GeometryInfo::subfaces_per_face : 1)) -{ - if (sub) - project_to_subfaces (quadrature, quadrature_points); - else - project_to_faces (quadrature, quadrature_points); - - const unsigned int n = GeometryInfo::faces_per_cell - * (sub ? GeometryInfo::subfaces_per_face : 1); - unsigned int k=0; - for (unsigned int i=0; i void QProjector<1>::project_to_face (const Quadrature<0> &, @@ -576,49 +555,88 @@ void QProjector<3>::project_to_subface (const Quadrature<2> &quadrature, +template <> +Quadrature<1> +QProjector<1>::project_to_all_faces (const Quadrature<0> &) +{ + Assert (false, ExcImpossibleInDim(1)); + return Quadrature<1>(0); +}; + + + template -void -QProjector::project_to_faces (const Quadrature &quadrature, - typename std::vector > &q_points) +Quadrature +QProjector::project_to_all_faces (const Quadrature &quadrature) { - unsigned int npt = quadrature.n_quadrature_points; - unsigned int nf = GeometryInfo::faces_per_cell; - - q_points.resize (npt*nf); - std::vector > help(npt); + const unsigned int n_points = quadrature.n_quadrature_points, + n_faces = GeometryInfo::faces_per_cell; + + // first fix quadrature points + typename std::vector > q_points (n_points * n_faces); + std::vector > help(n_points); - unsigned k=0; - for (unsigned int i=0;i weights (n_points * n_faces); + for (unsigned int face=0; face(q_points, weights); } +template <> +Quadrature<1> +QProjector<1>::project_to_all_subfaces (const Quadrature<0> &) +{ + Assert (false, ExcImpossibleInDim(1)); + return Quadrature<1>(0); +}; + + template -void -QProjector::project_to_subfaces (const Quadrature &quadrature, - typename std::vector > &q_points) +Quadrature +QProjector::project_to_all_subfaces (const Quadrature &quadrature) { - unsigned int npt = quadrature.n_quadrature_points; - unsigned int nf = GeometryInfo::faces_per_cell; - unsigned int nc = GeometryInfo::subfaces_per_face; + const unsigned int n_points = quadrature.n_quadrature_points, + n_faces = GeometryInfo::faces_per_cell, + subfaces_per_face = GeometryInfo::subfaces_per_face; - q_points.resize (npt*nf*nc); - std::vector > help(npt); + // first fix quadrature points + typename std::vector > q_points (n_points * n_faces * subfaces_per_face); + std::vector > help(n_points); - unsigned k=0; - for (unsigned int i=0;i weights (n_points * n_faces * subfaces_per_face); + for (unsigned int face=0; face(q_points, weights); } diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index d7b070aa13..8b16e29b7a 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -424,8 +424,8 @@ FiniteElement::get_face_data (const UpdateFlags flags, const Mapping &mapping, const Quadrature &quadrature) const { - QProjector q(quadrature, false); - return get_data (flags, mapping, q); + return get_data (flags, mapping, + QProjector::project_to_all_faces(quadrature)); } @@ -436,8 +436,8 @@ FiniteElement::get_subface_data (const UpdateFlags flags, const Mapping &mapping, const Quadrature &quadrature) const { - QProjector q(quadrature, true); - return get_data (flags, mapping, q); + return get_data (flags, mapping, + QProjector::project_to_all_subfaces(quadrature)); } diff --git a/deal.II/deal.II/source/fe/mapping_cartesian.cc b/deal.II/deal.II/source/fe/mapping_cartesian.cc index 4516ad50e4..45de9489d5 100644 --- a/deal.II/deal.II/source/fe/mapping_cartesian.cc +++ b/deal.II/deal.II/source/fe/mapping_cartesian.cc @@ -99,8 +99,8 @@ typename Mapping::InternalDataBase * MappingCartesian::get_face_data (const UpdateFlags update_flags, const Quadrature& quadrature) const { - QProjector q (quadrature, false); - InternalData* data = new InternalData (q); + InternalData* data + = new InternalData (QProjector::project_to_all_faces(quadrature)); data->update_once = update_once(update_flags); data->update_each = update_each(update_flags); @@ -116,8 +116,8 @@ typename Mapping::InternalDataBase * MappingCartesian::get_subface_data (const UpdateFlags update_flags, const Quadrature &quadrature) const { - QProjector q (quadrature, true); - InternalData* data = new InternalData (q); + InternalData* data + = new InternalData (QProjector::project_to_all_subfaces(quadrature)); data->update_once = update_once(update_flags); data->update_each = update_each(update_flags); diff --git a/deal.II/deal.II/source/fe/mapping_q.cc b/deal.II/deal.II/source/fe/mapping_q.cc index ae159135f5..926440c822 100644 --- a/deal.II/deal.II/source/fe/mapping_q.cc +++ b/deal.II/deal.II/source/fe/mapping_q.cc @@ -221,7 +221,7 @@ MappingQ::get_face_data (const UpdateFlags update_flags, const Quadrature& quadrature) const { InternalData *data = new InternalData(n_shape_functions); - QProjector q (quadrature, false); + Quadrature q (QProjector::project_to_all_faces(quadrature)); compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data); if (!use_mapping_q_on_all_cells) @@ -239,7 +239,7 @@ MappingQ::get_subface_data (const UpdateFlags update_flags, const Quadrature& quadrature) const { InternalData *data = new InternalData(n_shape_functions); - QProjector q (quadrature, true); + Quadrature q (QProjector::project_to_all_subfaces(quadrature)); compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data); if (!use_mapping_q_on_all_cells) diff --git a/deal.II/deal.II/source/fe/mapping_q1.cc b/deal.II/deal.II/source/fe/mapping_q1.cc index 2b5837bea7..1a0e665d11 100644 --- a/deal.II/deal.II/source/fe/mapping_q1.cc +++ b/deal.II/deal.II/source/fe/mapping_q1.cc @@ -462,8 +462,10 @@ MappingQ1::get_face_data (const UpdateFlags update_flags, const Quadrature &quadrature) const { InternalData* data = new InternalData(n_shape_functions); - QProjector q (quadrature, false); - compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data); + compute_face_data (update_flags, + QProjector::project_to_all_faces(quadrature), + quadrature.n_quadrature_points, + *data); return data; } @@ -476,8 +478,10 @@ MappingQ1::get_subface_data (const UpdateFlags update_flags, const Quadrature& quadrature) const { InternalData* data = new InternalData(n_shape_functions); - QProjector q (quadrature, true); - compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data); + compute_face_data (update_flags, + QProjector::project_to_all_subfaces(quadrature), + quadrature.n_quadrature_points, + *data); return data; } diff --git a/deal.II/deal.II/source/grid/grid_out.cc b/deal.II/deal.II/source/grid/grid_out.cc index 9a02433b87..58424c7641 100644 --- a/deal.II/deal.II/source/grid/grid_out.cc +++ b/deal.II/deal.II/source/grid/grid_out.cc @@ -279,7 +279,7 @@ void GridOut::write_gnuplot (const Triangulation &tria, // quadrature formula which will be // used to probe boundary points at // curved faces - QProjector *q_projector=0; + Quadrature *q_projector=0; if (mapping!=0) { typename std::vector > boundary_points(n_points); @@ -287,7 +287,7 @@ void GridOut::write_gnuplot (const Triangulation &tria, boundary_points[i](0)= 1.*(i+1)/(n_points+1); Quadrature quadrature(boundary_points); - q_projector = new QProjector (quadrature, false); + q_projector = new Quadrature (QProjector::project_to_all_faces(quadrature)); } for (; cell!=endc; ++cell) @@ -602,7 +602,7 @@ void GridOut::write_eps (const Triangulation &tria, boundary_points[i](0) = 1.*(i+1)/(n_points+1); Quadrature quadrature (boundary_points); - QProjector q_projector (quadrature, false); + Quadrature q_projector (QProjector::project_to_all_faces(quadrature)); // next loop over all // boundary faces and diff --git a/tests/base/quadrature_test.cc b/tests/base/quadrature_test.cc index c7d3d887ad..2ab4caa994 100644 --- a/tests/base/quadrature_test.cc +++ b/tests/base/quadrature_test.cc @@ -118,7 +118,9 @@ check_faces (const std::vector*>& quadratures, const bool sub) for (unsigned int n=0; n quadrature(*quadratures[n], sub); + Quadrature quadrature (sub == false? + QProjector::project_to_all_faces(*quadratures[n]) : + QProjector::project_to_all_subfaces(*quadratures[n])); const std::vector > &points=quadrature.get_points(); const std::vector &weights=quadrature.get_weights(); -- 2.39.5