From: Wolfgang Bangerth Date: Tue, 30 Jun 1998 14:15:55 +0000 (+0000) Subject: Implement cubic (1d and 2d) and quartic (2d) elements. Add doc snippet. X-Git-Tag: v8.0.0~22833 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c1603fad46844dc6069ace2aada653bb800e4085;p=dealii.git Implement cubic (1d and 2d) and quartic (2d) elements. Add doc snippet. git-svn-id: https://svn.dealii.org/trunk@421 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index fa843db075..f648eb7f0e 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -454,6 +454,14 @@ struct FiniteElementBase : public FiniteElementData { * since wrong results will most likely lead to internal errors through * the #Assert# mechanism, but the first places will lead to undiscovered * errors if not thought of properly. + * + * This assumption also comes into play when computing the constraints of + * hanging nodes. If functions not located on a certain face vanish on + * that face (they do for Lagrange elements), then the distribution of + * constrained nodes happens with the face nodes on the large call. If + * the assumption does not hold, then the distribution has to happen + * with all nodes on the small and the large cells. This is not + * implemented in the #DoFHandler# class as of now. * \end{itemize} * * @author Wolfgang Bangerth, 1998 diff --git a/deal.II/deal.II/include/fe/fe_lib.lagrange.h b/deal.II/deal.II/include/fe/fe_lib.lagrange.h index 8d54ce4c61..a0b3fda342 100644 --- a/deal.II/deal.II/include/fe/fe_lib.lagrange.h +++ b/deal.II/deal.II/include/fe/fe_lib.lagrange.h @@ -297,6 +297,313 @@ class FEQuadraticSub : public FiniteElement { +/** + * Define a (bi-, tri-, etc)cubic finite element in #dim# space dimensions. + * A linear (subparametric) mapping from the unit cell + * to the real cell is implemented. + * + * The numbering of degrees of freedom in one spatial dimension is as follows: + * \begin{verbatim} + * 0--2--3--1 + * \end{verbatim} + * + * The numbering of degrees of freedom in two spatial dimension is as follows: + * \begin{verbatim} + * 3--8--9--2 + * | | + * 11 15 14 7 + * | | + * 10 12 13 6 + * | | + * 0--4--5--1 + * \end{verbatim} + * Note the reverse ordering of degrees of freedom on the left and upper + * line and the counterclockwise numbering of the interior degrees of + * freedom. + */ +template +class FECubicSub : public FiniteElement { + public: + /** + * Constructor + */ + FECubicSub (); + + /** + * Return the value of the #i#th shape + * function at point #p# on the unit cell. + */ + virtual double shape_value(const unsigned int i, + const Point& p) const; + + /** + * Return the gradient of the #i#th shape + * function at point #p# on the unit cell. + */ + virtual Point shape_grad(const unsigned int i, + const Point& p) const; + + /** + * Refer to the base class for detailed + * information on this function. + * + * For one dimensional elements, this + * function simply passes through to + * the one implemented in the base class. + */ + virtual void fill_fe_values (const DoFHandler::cell_iterator &cell, + const vector > &unit_points, + vector &jacobians, + const bool compute_jacobians, + vector > &ansatz_points, + const bool compute_ansatz_points, + vector > &q_points, + const bool compute_q_points, + const Boundary &boundary) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_ansatz_points (const DoFHandler::cell_iterator &cell, + const Boundary &boundary, + vector > &ansatz_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_face_ansatz_points (const DoFHandler::face_iterator &face, + const Boundary &boundary, + vector > &ansatz_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_face_jacobians (const DoFHandler::face_iterator &face, + const Boundary &boundary, + const vector > &unit_points, + vector &face_jacobi_determinants) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_subface_jacobians (const DoFHandler::face_iterator &face, + const unsigned int subface_no, + const vector > &unit_points, + vector &face_jacobi_determinants) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_normal_vectors (const DoFHandler::cell_iterator &cell, + const unsigned int face_no, + const Boundary &boundary, + const vector > &unit_points, + vector > &normal_vectors) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_normal_vectors (const DoFHandler::cell_iterator &cell, + const unsigned int subface_no, + const unsigned int face_no, + const vector > &unit_points, + vector > &normal_vectors) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_local_mass_matrix (const DoFHandler::cell_iterator &cell, + const Boundary &boundary, + dFMatrix &local_mass_matrix) const; + + private: + /** + * Return the value of the #i#th shape + * function at point #p# on the unit cell. + * Here, the (bi-)linear basis functions + * are meant, which are used for the + * computation of the transformation from + * unit cell to real space cell. + */ + double linear_shape_value(const unsigned int i, + const Point& p) const; + + /** + * Return the gradient of the #i#th shape + * function at point #p# on the unit cell. + * Here, the (bi-)linear basis functions + * are meant, which are used for the + * computation of the transformation from + * unit cell to real space cell. + */ + Point linear_shape_grad(const unsigned int i, + const Point& p) const; +}; + + + +/** + * Define a (bi-, tri-, etc)cubic finite element in #dim# space dimensions. + * A linear (subparametric) mapping from the unit cell + * to the real cell is implemented. + * + * The numbering of degrees of freedom in one spatial dimension is as follows: + * \begin{verbatim} + * 0--2--3--4--1 + * \end{verbatim} + * + * The numbering of degrees of freedom in two spatial dimension is as follows: + * \begin{verbatim} + * 3--10-11-12-2 + * | | + * 15 19 22 18 9 + * | | + * 14 23 24 21 8 + * | | + * 13 16 20 17 7 + * | | + * 0--4--5--6--1 + * \end{verbatim} + * Note the reverse ordering of degrees of freedom on the left and upper + * line and the numbering of the interior degrees of + * freedom. + */ +template +class FEQuarticSub : public FiniteElement { + public: + /** + * Constructor + */ + FEQuarticSub (); + + /** + * Return the value of the #i#th shape + * function at point #p# on the unit cell. + */ + virtual double shape_value(const unsigned int i, + const Point& p) const; + + /** + * Return the gradient of the #i#th shape + * function at point #p# on the unit cell. + */ + virtual Point shape_grad(const unsigned int i, + const Point& p) const; + + /** + * Refer to the base class for detailed + * information on this function. + * + * For one dimensional elements, this + * function simply passes through to + * the one implemented in the base class. + */ + virtual void fill_fe_values (const DoFHandler::cell_iterator &cell, + const vector > &unit_points, + vector &jacobians, + const bool compute_jacobians, + vector > &ansatz_points, + const bool compute_ansatz_points, + vector > &q_points, + const bool compute_q_points, + const Boundary &boundary) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_ansatz_points (const DoFHandler::cell_iterator &cell, + const Boundary &boundary, + vector > &ansatz_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_face_ansatz_points (const DoFHandler::face_iterator &face, + const Boundary &boundary, + vector > &ansatz_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_face_jacobians (const DoFHandler::face_iterator &face, + const Boundary &boundary, + const vector > &unit_points, + vector &face_jacobi_determinants) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_subface_jacobians (const DoFHandler::face_iterator &face, + const unsigned int subface_no, + const vector > &unit_points, + vector &face_jacobi_determinants) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_normal_vectors (const DoFHandler::cell_iterator &cell, + const unsigned int face_no, + const Boundary &boundary, + const vector > &unit_points, + vector > &normal_vectors) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_normal_vectors (const DoFHandler::cell_iterator &cell, + const unsigned int subface_no, + const unsigned int face_no, + const vector > &unit_points, + vector > &normal_vectors) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_local_mass_matrix (const DoFHandler::cell_iterator &cell, + const Boundary &boundary, + dFMatrix &local_mass_matrix) const; + + private: + /** + * Return the value of the #i#th shape + * function at point #p# on the unit cell. + * Here, the (bi-)linear basis functions + * are meant, which are used for the + * computation of the transformation from + * unit cell to real space cell. + */ + double linear_shape_value(const unsigned int i, + const Point& p) const; + + /** + * Return the gradient of the #i#th shape + * function at point #p# on the unit cell. + * Here, the (bi-)linear basis functions + * are meant, which are used for the + * computation of the transformation from + * unit cell to real space cell. + */ + Point linear_shape_grad(const unsigned int i, + const Point& p) const; +}; + + + + /*---------------------------- fe_lib.h ---------------------------*/ /* end of #ifndef __fe_lib_H */