From 6f6b331353ff253fe3d6b01e47f20531e10b3d39 Mon Sep 17 00:00:00 2001 From: hartmann Date: Wed, 16 Dec 1998 13:35:26 +0000 Subject: [PATCH] Check in initial release of DG elements. git-svn-id: https://svn.dealii.org/trunk@714 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_lib.dg.h | 267 +++++++++++++++++++++++++ 1 file changed, 267 insertions(+) create mode 100644 deal.II/deal.II/include/fe/fe_lib.dg.h diff --git a/deal.II/deal.II/include/fe/fe_lib.dg.h b/deal.II/deal.II/include/fe/fe_lib.dg.h new file mode 100644 index 0000000000..db1ab4c585 --- /dev/null +++ b/deal.II/deal.II/include/fe/fe_lib.dg.h @@ -0,0 +1,267 @@ +/*---------------------------- fe_lib.dg.h ---------------------------*/ +/* $Id$ */ +/* Ralf Hartmann, University of Heidelberg, Dez 98 */ +#ifndef __fe_lib_dg_H +#define __fe_lib_dg_H +/*---------------------------- fe_lib.dg.h ---------------------------*/ + + +#include + + +/** + * Define a constant discontinuous finite element in #dim# + * space dimensions, along with (bi-, tri-)linear (therefore isoparametric) + * transforms from the unit cell to the real cell. + * @author Ralf Hartmann, 1998 + */ + +template +class FEDGConstant : public FELinearMapping { + public: + /** + * Constructor + */ + FEDGConstant (); + + /** + * 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 Tensor<1,dim> shape_grad(const unsigned int i, + const Point& p) const; + + /** + * Return the tensor of second derivatives + * of the #i#th shape function at + * point #p# on the unit cell. + * + * For linear elements, all second + * derivatives on the unit cell are zero. + */ + virtual Tensor<2,dim> shape_grad_grad (const unsigned int i, + const Point &p) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_unit_support_points (vector > &support_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_support_points (const DoFHandler::cell_iterator &cell, + const Boundary &boundary, + vector > &support_points) const; + + /** + * Refer to the base class for detailed + * information on this function. + */ + virtual void get_face_support_points (const DoFHandler::face_iterator &face, + const Boundary &boundary, + vector > &support_points) 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; + + /** + * Return a readonly reference to the + * matrix which describes the transfer of a + * child with the given number to the + * mother cell. See the #restriction# array + * for more information. + * + * This function returns an error since the + * correct use of the restriction + * matrices is not yet finally decided + * about. + */ + const dFMatrix & restrict (const unsigned int) const; + + /** + * Exception + */ + DeclException0 (ExcNotImplemented); +}; + + + +/** + * Define a (bi-, tri-, etc)linear finite element in #dim# space dimensions, + * along with (bi-, tri-)linear (therefore isoparametric) transforms from the + * unit cell to the real cell allowing discontinuous Galerkin methods. + * + * + * This class is derived from and provides substantially the same + * as the #FELinear# class. The only difference is the new constructor that + * calls #FELinear::FELinear(const int)#, the protected constructor of the + * #FELinear# class using a #FiniteElement# with no dofs in the vertices and + * $2^d$ dofs per cell. As now the cells do not share any vertex-dof with + * a neighboring cell the $2^d$ dofs per cell can be choosen independently not + * needing any constraints and allowing the use of discontinuous Galerkin + * methods. Although the basis functions now are not longer associated + * with the vertices but with the cell they retain their shape. As already + * explained no constraint matrices needed to be implemented. + * To use this element you need to think about the jump terms in your + * weak formulation of your discontinuous Galerkin scheme. + * @author Ralf Hartmann, 1998 + */ + + +template +class FEDGLinear : public FELinear{ + public: + /** + * Constructor + */ + FEDGLinear(); + /** + * This function returns an error since the + * correct use of the restriction + * matrices is not yet finally decided + * about. + */ + const dFMatrix & restrict (const unsigned int) const; +}; + + + +/** + * Define a (bi-, tri-, etc)quadratic finite element in #dim# space dimensions, + * along with (bi-, tri-)quadratic (therefore isoparametric) transforms from the + * unit cell to the real cell allowing discontinuous Galerkin methods. + * + * This class is derived from and provides substantially the same + * as the #FEQuadraticSub# class. The only difference is the new constructor that + * calls #FEQuadraticSub::FEQuadraticSub(const int)#, the protected constructor of the + * #FEQuadraticSub# class using a #FiniteElement# with no dofs in the vertices, no dofs on the lines and + * $3^d$ dofs per cell. As now the cells do not share any vertex-dof with + * a neighboring cell the $3^d$ dofs per cell can be choosen independently not + * needing any constraints and allowing the use of discontinuous Galerkin + * methods. Although the basis functions now are not longer associated + * with the vertices but with the cell they retain their shape. As already + * explained no constraint matrices needed to be implemented. + * To use this element you need to think about the jump terms in your + * weak formulation of your discontinuous Galerkin scheme. + * @author Ralf Hartmann, 1998 + */ + + +template +class FEDGQuadraticSub : public FEQuadraticSub{ + public: + /** + * Constructor + */ + FEDGQuadraticSub(); + /** + * This function returns an error since the + * correct use of the restriction + * matrices is not yet finally decided + * about. + */ + const dFMatrix & restrict (const unsigned int) const; +}; + + + + +/** + * Define a (bi-, tri-, etc)cubic finite element in #dim# space dimensions, + * along with (bi-, tri-)cubic (therefore isoparametric) transforms from the + * unit cell to the real cell allowing discontinuous Galerkin methods. + * + * This class is derived from and provides substantially the same + * as the #FECubicSub# class. The only difference is the new constructor that + * calls #FECubicSub::FECubicSub(const int)#, the protected constructor of the + * #FECubicSub# class using a #FiniteElement# with no dofs in the vertices, no dofs on the lines and + * $4^d$ dofs per cell. As now the cells do not share any vertex-dof with + * a neighboring cell the $4^d$ dofs per cell can be choosen independently not + * needing any constraints and allowing the use of discontinuous Galerkin + * methods. Although the basis functions now are not longer associated + * with the vertices but with the cell they retain their shape. As already + * explained no constraint matrices needed to be implemented. + * To use this element you need to think about the jump terms in your + * weak formulation of your discontinuous Galerkin scheme. + * @author Ralf Hartmann, 1998 + */ + + +template +class FEDGCubicSub : public FECubicSub{ + public: + /** + * Constructor + */ + FEDGCubicSub(); + /** + * This function returns an error since the + * correct use of the restriction + * matrices is not yet finally decided + * about. + */ + const dFMatrix & restrict (const unsigned int) const; +}; + + + +/** + * Define a (bi-, tri-, etc)quartic finite element in #dim# space dimensions, + * along with (bi-, tri-)quartic (therefore isoparametric) transforms from the + * unit cell to the real cell allowing discontinuous Galerkin methods. + * + * This class is derived from and provides substantially the same + * as the #FEQuarticSub# class. The only difference is the new constructor that + * calls #FEQuarticSub::FEQuarticSub(const int)#, the protected constructor of the + * #FEQuarticSub# class using a #FiniteElement# with no dofs in the vertices, no dofs on the lines and + * $5^d$ dofs per cell. As now the cells do not share any vertex-dof with + * a neighboring cell the $5^d$ dofs per cell can be choosen independently not + * needing any constraints and allowing the use of discontinuous Galerkin + * methods. Although the basis functions now are not longer associated + * with the vertices but with the cell they retain their shape. As already + * explained no constraint matrices needed to be implemented. + * To use this element you need to think about the jump terms in your + * weak formulation of your discontinuous Galerkin scheme. + * @author Ralf Hartmann, 1998 + */ + + +template +class FEDGQuarticSub : public FEQuarticSub{ + public: + /** + * Constructor + */ + FEDGQuarticSub(); + /** + * This function returns an error since the + * correct use of the restriction + * matrices is not yet finally decided + * about. + */ + const dFMatrix & restrict (const unsigned int) const; +}; + + + + + +/*---------------------------- fe_lib.dg.h ---------------------------*/ +/* end of #ifndef __fe_lib_dg_H */ +#endif +/*---------------------------- fe_lib.dg.h ---------------------------*/ -- 2.39.5