From: Ralf Hartmann Date: Thu, 15 Mar 2001 18:59:06 +0000 (+0000) Subject: Change private, protected, public of some functions. Remove InternalDataBase argument... X-Git-Tag: v8.0.0~19551 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fc87b67a05f6c3420aa20e35776536afccd26325;p=dealii.git Change private, protected, public of some functions. Remove InternalDataBase argument of transform_unit_to_real_cell function; new implementation of transform_unit_to_real_cell_internal function. New internal helper function: get_cell_data. git-svn-id: https://svn.dealii.org/trunk@4221 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/mapping_q1.h b/deal.II/deal.II/include/fe/mapping_q1.h index bd52d90ab5..5050857235 100644 --- a/deal.II/deal.II/include/fe/mapping_q1.h +++ b/deal.II/deal.II/include/fe/mapping_q1.h @@ -17,6 +17,7 @@ #include #include + /** * Mapping of general quadrilateral/hexahedra by d-linear shape * functions. @@ -36,72 +37,30 @@ template class MappingQ1 : public Mapping { public: -//TODO: why make the following functions public? they are only helpful for fevalues and maybe the finite elements? - /** - * Implementation of the interface in - * @ref{Mapping}. - */ - virtual Mapping::InternalDataBase* - get_data (const UpdateFlags, - const Quadrature& quadrature) const; - /** - * Implementation of the interface in - * @ref{Mapping}. - */ - virtual Mapping::InternalDataBase* - get_face_data (const UpdateFlags flags, - const Quadrature& quadrature) const; - - /** - * Implementation of the interface in - * @ref{Mapping}. - */ - virtual Mapping::InternalDataBase* - get_subface_data (const UpdateFlags flags, - const Quadrature& quadrature) const; - - /** - * Implementation of the interface in - * @ref{Mapping}. - */ - virtual void - fill_fe_values (const DoFHandler::cell_iterator &cell, - const Quadrature& quadrature, - Mapping::InternalDataBase &mapping_data, - std::vector > &quadrature_points, - std::vector &JxW_values) const ; - - /** - * Implementation of the interface in - * @ref{Mapping}. + * Transforms the point @p{p} on + * the unit cell to the point + * @p{p_real} on the real cell + * @p{cell} and returns @p{p_real}. */ - virtual void - fill_fe_face_values (const typename DoFHandler::cell_iterator &cell, - const unsigned int face_no, - const Quadrature& quadrature, - typename Mapping::InternalDataBase &mapping_data, - std::vector > &quadrature_points, - std::vector &JxW_values, - std::vector > &boundary_form, - std::vector > &normal_vectors) const ; - + virtual Point transform_unit_to_real_cell ( + const typename Triangulation::cell_iterator cell, + const Point &p) const; + /** - * Implementation of the interface in - * @ref{Mapping}. + * Transforms the point @p{p} on + * the real cell to the point + * @p{p_unit} on the unit cell + * @p{cell} and returns @p{p_unit}. + * + * Uses Newton iteration and the + * @p{transform_unit_to_real_cell} + * function. */ - virtual void - fill_fe_subface_values (const typename DoFHandler::cell_iterator &cell, - const unsigned int face_no, - const unsigned int sub_no, - const Quadrature& quadrature, - typename Mapping::InternalDataBase &mapping_data, - std::vector > &quadrature_points, - std::vector &JxW_values, - std::vector > &boundary_form, - std::vector > &normal_vectors) const ; - - + virtual Point transform_real_to_unit_cell ( + const typename Triangulation::cell_iterator cell, + const Point &p) const; + /** * Implementation of the interface in * @ref{Mapping}. @@ -142,31 +101,6 @@ class MappingQ1 : public Mapping const Mapping::InternalDataBase &mapping_data, const unsigned int src_offset) const; - /** - * Transforms the point @p{p} on - * the unit cell to the point - * @p{p_real} on the real cell - * @p{cell} and returns @p{p_real}. - */ -//TODO: document meaning of mdata argument - virtual Point transform_unit_to_real_cell ( - const typename Triangulation::cell_iterator cell, - const Point &p, - const typename Mapping::InternalDataBase *const mdata=0) const; - - /** - * Transforms the point @p{p} on - * the real cell to the point - * @p{p_unit} on the unit cell - * @p{cell} and returns @p{p_unit}. - * - * Uses Newton iteration and the - * @p{transform_unit_to_real_cell} - * function. - */ - virtual Point transform_real_to_unit_cell ( - const typename Triangulation::cell_iterator cell, - const Point &p) const; /** * Implementation of the interface in @@ -179,13 +113,49 @@ class MappingQ1 : public Mapping * @ref{Mapping}. */ virtual UpdateFlags update_each (const UpdateFlags) const; + + protected: /** - * Exception + * Implementation of the interface in + * @ref{Mapping}. */ - DeclException0 (ExcInvalidData); + virtual void + fill_fe_values (const DoFHandler::cell_iterator &cell, + const Quadrature& quadrature, + Mapping::InternalDataBase &mapping_data, + std::vector > &quadrature_points, + std::vector &JxW_values) const ; - protected: + /** + * Implementation of the interface in + * @ref{Mapping}. + */ + virtual void + fill_fe_face_values (const typename DoFHandler::cell_iterator &cell, + const unsigned int face_no, + const Quadrature& quadrature, + typename Mapping::InternalDataBase &mapping_data, + std::vector > &quadrature_points, + std::vector &JxW_values, + std::vector > &boundary_form, + std::vector > &normal_vectors) const ; + + /** + * Implementation of the interface in + * @ref{Mapping}. + */ + virtual void + fill_fe_subface_values (const typename DoFHandler::cell_iterator &cell, + const unsigned int face_no, + const unsigned int sub_no, + const Quadrature& quadrature, + typename Mapping::InternalDataBase &mapping_data, + std::vector > &quadrature_points, + std::vector &JxW_values, + std::vector > &boundary_form, + std::vector > &normal_vectors) const ; + /** * Implementation of the * covariant transformation. @@ -340,30 +310,6 @@ class MappingQ1 : public Mapping unsigned int n_shape_functions; }; - /** - * Do the computations for the - * @p{get_face_data} - * functions. Here, the data - * vectors of @p{InternalData} - * are reinitialized to proper - * size and shape values and - * derivatives are - * computed. Furthermore - * @p{unit_tangential} vectors of - * the face are computed. - */ - void compute_face_data (const UpdateFlags flags, - const Quadrature &quadrature, - const unsigned int n_orig_q_points, - InternalData &data) const; - - /** - * Mapping between tensor product - * ordering and real ordering of - * the vertices. - */ - static const unsigned int vertex_mapping[1< */ void compute_shapes (const std::vector > &unit_points, InternalData &data) const; - + /** - * Do the computations for the @p{get_data} - * functions. Here, the data - * vectors of @p{InternalData} are - * reinitialized to proper size and - * shape values are computed. + * Do the computations for the + * @p{get_data} functions. Here, + * the data vectors of + * @p{InternalData} are + * reinitialized to proper size + * and shape values are computed. */ void compute_data (const UpdateFlags flags, const Quadrature& quadrature, const unsigned int n_orig_q_points, InternalData& data) const; + + /** + * Do the computations for the + * @p{get_face_data} + * functions. Here, the data + * vectors of @p{InternalData} + * are reinitialized to proper + * size and shape values and + * derivatives are + * computed. Furthermore + * @p{unit_tangential} vectors of + * the face are computed. + */ + void compute_face_data (const UpdateFlags flags, + const Quadrature &quadrature, + const unsigned int n_orig_q_points, + InternalData &data) const; /** * Do the computation for the @@ -416,14 +380,101 @@ class MappingQ1 : public Mapping std::vector &JxW_values, std::vector > &boundary_form, std::vector > &normal_vectors) const; - + /** * Compute shape values and/or * derivatives. */ virtual void compute_shapes_virtual (const std::vector > &unit_points, InternalData &data) const; + + /** + * Mapping between tensor product + * ordering and real ordering of + * the vertices. + */ + static const unsigned int vertex_mapping[1< transform_unit_to_real_cell_internal ( + const typename Triangulation::cell_iterator cell, + const Point &p, + const InternalData &m_data) const; + + /** + * Returns an @p{InternalData} + * whose data vectors are resized + * corresponding to the + * @p{update_flags} and a + * one-point + * quadrature. Furthermore the + * @p{InternalData} stores the + * mapping support points of the + * given @p{cell}. + * + * This function is called by + * @p{transform_unit_to_real_cell} + * and by + * @p{transform_real_to_unit_cell}. + * The resulting @p{InternalData} + * is given to the + * @p{transform_unit_to_real_internal} + * function. + */ + InternalData* get_cell_data(const typename Triangulation::cell_iterator cell, + const UpdateFlags update_flags) const; + + /** + * Implementation of the interface in + * @ref{Mapping}. + */ + virtual Mapping::InternalDataBase* + get_data (const UpdateFlags, + const Quadrature& quadrature) const; + + /** + * Implementation of the interface in + * @ref{Mapping}. + */ + virtual Mapping::InternalDataBase* + get_face_data (const UpdateFlags flags, + const Quadrature& quadrature) const; + + /** + * Implementation of the interface in + * @ref{Mapping}. + */ + virtual Mapping::InternalDataBase* + get_subface_data (const UpdateFlags flags, + const Quadrature& quadrature) const; /** * Computes the support points of