From: Guido Kanschat Date: Fri, 9 Jun 2000 22:11:44 +0000 (+0000) Subject: continuing on FEDG_Pk X-Git-Tag: v8.0.0~20409 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1962417917375f106cbd5e1da27490a95cd024ac;p=dealii.git continuing on FEDG_Pk git-svn-id: https://svn.dealii.org/trunk@3003 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_lib.dgp.h b/deal.II/deal.II/include/fe/fe_lib.dgp.h index 2e407a24fa..cc3f4a637f 100644 --- a/deal.II/deal.II/include/fe/fe_lib.dgp.h +++ b/deal.II/deal.II/include/fe/fe_lib.dgp.h @@ -19,21 +19,16 @@ #include + /** - * Isoparametric Q1 finite element in #dim# space dimensions. - * - * The linear, isoparametric mapping from a point $\vec \xi$ on the unit cell - * to a point $\vec x$ on the real cell is defined as - * $$ \vec x(\vec \xi) = \sum_j {\vec p_j} N_j(\xi) $$ - * where $\vec p_j$ is the vector to the $j$th corner point of the cell in - * real space and $N_j(\vec \xi)$ is the value of the basis function associated - * with the $j$th corner point, on the unit cell at point $\vec \xi$. The sum - * over $j$ runs over all corner points. + * Discontinuous P1-element on hypercubes. * - * The number of degrees of freedom equal the number of the respective vertex - * of the cell + * This is the implementation of a linear (sic) polynomial space on a + * d-dimensional hypercube. The shape functions are the first @p{d+1} + * of @p{1,x,y,z}. Later on, these should be exchanged for mutually + * orthogonal, preferably by changing the unit cell to $[-1,1]^d$. * - * @author Wolfgang Bangerth, 1998, 1999 + * @author Guido Kanschat, 2000 */ template class FEDG_P1 : public FEQ1Mapping @@ -119,86 +114,15 @@ class FEDG_P1 : public FEQ1Mapping /** - * Subparametric Q2 finite element in #dim# space dimensions. - * A Q1 mapping from the unit cell - * to the real cell is implemented. - * - * The numbering of the degrees of freedom is as follows: - * \begin{itemize} - * \item 1D case: - * \begin{verbatim} - * 0---2---1 - * \end{verbatim} - * - * \item 2D case: - * \begin{verbatim} - * 3---6---2 - * | | - * 7 8 5 - * | | - * 0---4---1 - * \end{verbatim} + * Discontinuous P2-element on hypercubes. * - * \item 3D case: - * \begin{verbatim} - * 7--14---6 7--14---6 - * /| | / /| - * 19 | 13 19 1813 - * / 15 | / / | - * 3 | | 3---10--2 | - * | 4--12---5 | | 5 - * | / / | 9 / - * 11 16 17 11 | 17 - * |/ / | |/ - * 0---8---1 0---8---1 + * This is the implementation of a linear (sic) polynomial space on a + * d-dimensional hypercube. The shape functions are those of + * @p{1,x,y,z, x*x, x*y, x*z, y*y, y*z, z*z} applying to the space + * dimension. Later on, these should be exchanged for mutually + * orthogonal, preferably by changing the unit cell to $[-1,1]^d$. * - * *-------* *-------* - * /| | / /| - * / | 21 | / 24 / | - * / | | / / | - * * | | *-------* | - * |25 *-------* | |23 * - * | / / | 20 | / - * | / 22 / | | / - * |/ / | |/ - * *-------* *-------* - * \end{verbatim} - * The center vertex has number 26. - * - * The respective coordinate values of the support points of the degrees - * of freedom are as follows: - * \begin{itemize} - * \item Index 0: #[0, 0, 0]#; - * \item Index 1: #[1, 0, 0]#; - * \item Index 2: #[1, 0, 1]#; - * \item Index 3: #[0, 0, 1]#; - * \item Index 4: #[0, 1, 0]#; - * \item Index 5: #[1, 1, 0]#; - * \item Index 6: #[1, 1, 1]#; - * \item Index 7: #[0, 1, 1]#; - * \item Index 8: #[1/2, 0, 0]#; - * \item Index 9: #[1, 0, 1/2]#; - * \item Index 10: # [1/2, 0, 1]#; - * \item Index 11: # [0, 0, 1/2]#; - * \item Index 12: # [1/2, 1, 0]#; - * \item Index 13: # [1, 1, 1/2]#; - * \item Index 14: # [1/2, 1, 1]#; - * \item Index 15: # [0, 1, 1/2]#; - * \item Index 16: # [0, 1/2, 0]#; - * \item Index 17: # [1, 1/2, 0]#; - * \item Index 18: # [1, 1/2, 1]#; - * \item Index 19: # [0, 1/2, 1]#; - * \item Index 20: # [1/2, 0, 1/2]#; - * \item Index 21: # [1/2, 1, 1/2]#; - * \item Index 22: # [1/2, 1/2, 0]#; - * \item Index 23: # [1, 1/2, 1/2]#; - * \item Index 24: # [1/2, 1/2, 1]#; - * \item Index 25: # [0, 1/2, 1/2]#; - * \item Index 26: # [1/2, 1/2, 1/2]#; - * \end{itemize} - * \end{itemize} - * - * @author Wolfgang Bangerth, 1998, 1999 + * @author Guido Kanschat, 2000 */ template class FEDG_P2 : public FEQ1Mapping @@ -282,30 +206,15 @@ class FEDG_P2 : public FEQ1Mapping /** - * Subparametric Q3 finite element in #dim# space dimensions. - * A Q1 mapping from the unit cell - * to the real cell is implemented. + * Discontinuous P3-element on hypercubes. * - * The numbering of degrees of freedom in one spatial dimension is as follows: - * \begin{verbatim} - * 0--2--3--1 - * \end{verbatim} + * This is the implementation of a linear (sic) polynomial space on a + * d-dimensional hypercube. The shape functions are the basis + * polynomials spanning the space of cubic polynomials. Later on, + * they should be exchanged for mutually orthogonal, preferably by + * changing the unit cell to $[-1,1]^d$. * - * 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. - * - * @author Wolfgang Bangerth, 1998 + * @author Guido Kanschat, 2000 */ template class FEDG_P3 : public FEQ1Mapping @@ -390,32 +299,15 @@ class FEDG_P3 : public FEQ1Mapping /** - * Subparametric Q4 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} + * Discontinuous P3-element on hypercubes. * - * 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. + * This is the implementation of a linear (sic) polynomial space on a + * d-dimensional hypercube. The shape functions are the basis + * polynomials spanning the space of cubic polynomials. Later on, + * they should be exchanged for mutually orthogonal, preferably by + * changing the unit cell to $[-1,1]^d$. * - * @author Wolfgang Bangerth, 1998 + * @author Guido Kanschat, 2000 */ template class FEDG_P4 : public FEQ1Mapping diff --git a/deal.II/deal.II/source/fe/fe_lib.dgp2.cc b/deal.II/deal.II/source/fe/fe_lib.dgp2.cc index 3168fc4d57..555f5fd5c7 100644 --- a/deal.II/deal.II/source/fe/fe_lib.dgp2.cc +++ b/deal.II/deal.II/source/fe/fe_lib.dgp2.cc @@ -237,6 +237,7 @@ FEDG_P2<3>::FEDG_P2 () : FEQ1Mapping<3> (0, 0, 0, 4, 1, vector (1, true)) { + Assert(false, ExcNotImplemented ()); // initialize_matrices (); }; diff --git a/deal.II/deal.II/source/fe/fe_lib.dgp3.cc b/deal.II/deal.II/source/fe/fe_lib.dgp3.cc new file mode 100644 index 0000000000..d4ca02548e --- /dev/null +++ b/deal.II/deal.II/source/fe/fe_lib.dgp3.cc @@ -0,0 +1,375 @@ +//---------------------------- $RCSFile$ --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 1998, 1999, 2000 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- $RCSFile$ --------------------------- + + +#include +#include +#include +#include +#include + + +// declare explicit specializations before use: +template <> void FEDG_P3::initialize_matrices (); + + +#if deal_II_dimension == 1 + +template <> +FEDG_P3<1>::FEDG_P3 () : + FEQ1Mapping<1> (0, 4, 0, 0, 1, + vector (1, true)) +{ +// initialize_matrices (); +}; + + +template <> +void FEDG_P3<1>::initialize_matrices () +{ + Assert(false, ExcNotImplemented()); +}; + + +template <> +double +FEDG_P3<1>::shape_value(const unsigned int i, + const Point<1> &p) const +{ + Assert((i +inline +Tensor<1,1> +FEDG_P3<1>::shape_grad(const unsigned int i, + const Point<1>&p) const +{ + Assert((i, so we + // still construct it as that. it should + // make no difference in practice, + // however + switch (i) + { + case 0: return Point<1>(-1.); + case 1: return Point<1>(1.); + case 2: return Point<1>(2.*p(0)); + case 2: return Point<1>(3.*p(0)*p(0)); + + } + return Point<1>(); +}; + + +template <> +inline +Tensor<2,1> +FEDG_P3<1>::shape_grad_grad (const unsigned int i, + const Point<1> &) const +{ + Assert(false, ExcNotImplemented()); + Assert((i(); +}; + + +template <> +void FEDG_P3<1>::get_unit_support_points (vector > &support_points) const +{ + FiniteElement<1>::get_unit_support_points (support_points); +}; + + +template <> +void FEDG_P3<1>::get_support_points (const DoFHandler<1>::cell_iterator &cell, + vector > &support_points) const +{ + FiniteElement<1>::get_support_points (cell, support_points); +}; + + +template <> +void FEDG_P3<1>::get_face_support_points (const DoFHandler<1>::face_iterator &, + vector > &) const +{ + Assert (false, ExcInternalError()); +}; + + +template <> +void FEDG_P3<1>::get_local_mass_matrix (const DoFHandler<1>::cell_iterator &cell, + FullMatrix &local_mass_matrix) const +{ + Assert(false, ExcNotImplemented()); + Assert (local_mass_matrix.n() == dofs_per_cell, + ExcWrongFieldDimension(local_mass_matrix.n(),dofs_per_cell)); + Assert (local_mass_matrix.m() == dofs_per_cell, + ExcWrongFieldDimension(local_mass_matrix.m(),dofs_per_cell)); +}; + +#endif + + +#if deal_II_dimension == 2 + +template <> +FEDG_P3<2>::FEDG_P3 () : + FEQ1Mapping<2> (0, 0, 10, 0, 1, + vector (1, true)) +{ +// initialize_matrices (); +}; + + +template <> +void FEDG_P3<2>::initialize_matrices () +{ + Assert(false, ExcNotImplemented()); +}; + + +template <> +inline +double +FEDG_P3<2>::shape_value (const unsigned int i, + const Point<2>& p) const +{ + Assert((i +inline +Tensor<1,2> +FEDG_P3<2>::shape_grad (const unsigned int i, + const Point<2>& p) const +{ + Assert((i, so we + // still construct it as that. it should + // make no difference in practice, + // however + switch (i) + { + case 0: return Point<2> (0,0); + case 1: return Point<2> (1,0); + case 2: return Point<2> (0,1); + case 3: return Point<2> (2*p(0),0); + case 4: return Point<2> (p(1),p(0)); + case 5: return Point<2> (0,2*p(1)); + case 6: return Point<2> (3*p(0)*p(0), 0); + case 7: return Point<2> (2*p(0)*p(1), p(0)*p(0)); + case 8: return Point<2> (p(1)*p(1), 2*p(0)*p(1)); + case 9: return Point<2> (0, 3*p(1)*p(1)); + } + return Point<2> (); +}; + + +template <> +inline +Tensor<2,2> +FEDG_P3<2>::shape_grad_grad (const unsigned int i, + const Point<2> &) const +{ + Assert(false, ExcNotImplemented()); + Assert((i(); +}; + + +template <> +void FEDG_P3<2>::get_local_mass_matrix (const DoFHandler<2>::cell_iterator &, + FullMatrix &) const +{ + Assert(false, ExcNotImplemented ()); +}; + + +template <> +void FEDG_P3<2>::get_unit_support_points (vector > &unit_points) const +{ + Assert(false, ExcNotImplemented ()); + Assert (unit_points.size() == dofs_per_cell, + ExcWrongFieldDimension (unit_points.size(), dofs_per_cell)); + + unit_points[0] = Point<2> (.5,.5); + unit_points[1] = Point<2> (1,0); + unit_points[2] = Point<2> (0,1); + unit_points[3] = Point<2> (1,0); + unit_points[4] = Point<2> (0,1); + unit_points[5] = Point<2> (1,1); +}; + + +#endif + + +#if deal_II_dimension == 3 + +template <> +FEDG_P3<3>::FEDG_P3 () : + FEQ1Mapping<3> (0, 0, 0, 4, 1, + vector (1, true)) +{ + Assert(false, ExcNotImplemented ()); +// initialize_matrices (); +}; + + +template <> +void FEDG_P3<3>::initialize_matrices () +{ + Assert(false, ExcNotImplemented()); +}; + + +template <> +inline +double +FEDG_P3<3>::shape_value (const unsigned int i, + const Point<3>& p) const +{ + Assert((i +inline +Tensor<1,3> +FEDG_P3<3>::shape_grad (const unsigned int i, + const Point<3>& p) const +{ + Assert((i, so we + // still construct it as that. it should + // make no difference in practice, + // however + switch (i) + { + case 0: return Point<3>(0,0,0); + case 1: return Point<3>(1,0,0); + case 2: return Point<3>(0,1,0); + case 3: return Point<3>(0,0,1); + } + return Point<3> (); +}; + + +template <> +inline +Tensor<2,3> +FEDG_P3<3>::shape_grad_grad (const unsigned int i, + const Point<3> &p) const +{ + Assert((i return_value; + return return_value; +}; + + +template <> +void FEDG_P3<3>::get_local_mass_matrix (const DoFHandler<3>::cell_iterator &, + FullMatrix &local_mass_matrix) const +{ + Assert (local_mass_matrix.n() == dofs_per_cell, + ExcWrongFieldDimension(local_mass_matrix.n(),dofs_per_cell)); + Assert (local_mass_matrix.m() == dofs_per_cell, + ExcWrongFieldDimension(local_mass_matrix.m(),dofs_per_cell)); + + AssertThrow (false, ExcComputationNotUseful(3)); +}; + + +template <> +void FEDG_P3<3>::get_unit_support_points (vector > &unit_points) const { + Assert (unit_points.size() == dofs_per_cell, + ExcWrongFieldDimension (unit_points.size(), dofs_per_cell)); + + unit_points[0] = Point<3> (.5,.5,.5); + unit_points[1] = Point<3> (1,0,0); + unit_points[2] = Point<3> (0,1,0); + unit_points[3] = Point<3> (0,0,1); +}; + + +#endif + + +template +void +FEDG_P3::get_support_points (const typename DoFHandler::cell_iterator &cell, + vector > &support_points) const +{ + Assert (support_points.size() == dofs_per_cell, + ExcWrongFieldDimension (support_points.size(), dofs_per_cell)); + + for (unsigned int vertex=0; vertex::vertices_per_cell; ++vertex) + support_points[vertex] = cell->vertex(vertex); +}; + + +template +void +FEDG_P3::get_face_support_points (const typename DoFHandler::face_iterator &face, + vector > &support_points) const +{ + Assert ((support_points.size() == dofs_per_face) && + (support_points.size() == GeometryInfo::vertices_per_face), + ExcWrongFieldDimension (support_points.size(), + GeometryInfo::vertices_per_face)); + + for (unsigned int vertex=0; vertexvertex(vertex); +}; + + +// explicit instantiations + +template class FEDG_P3;