From: Jaeryun Yim Date: Tue, 13 Sep 2016 06:18:56 +0000 (+0900) Subject: Edit comments, adjust indents and spacing. X-Git-Tag: v8.5.0-rc1~624^2~12 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=92e277b3d42d8422b6e25c186db2857f6b97dc84;p=dealii.git Edit comments, adjust indents and spacing. --- diff --git a/include/deal.II/fe/fe_p1nc.h b/include/deal.II/fe/fe_p1nc.h index cf210d7ef1..23f948fc44 100644 --- a/include/deal.II/fe/fe_p1nc.h +++ b/include/deal.II/fe/fe_p1nc.h @@ -31,41 +31,51 @@ DEAL_II_NAMESPACE_OPEN /** * Implementation for the scalar version of the P1 nonconforming finite * element, a piecewise linear element on quadrilaterals in 2D. + * This implementation is only for 2D and codimension = 0. * - * Unlike any continuous conforming finite element belongs to H^1_0, + * Unlike any continuous conforming finite element which belongs to $H^1_0$, * the P1 nonconforming element does not enforce the continuity across edges. * But it requires the continuity just in integral sense: - * any function in the space should have the same integral values on two sides of the common edge shared by two adjacent elements. - - * Thus each function in the nonconforming element space can be discontinuous, not included in H^1_0, as functions in Discontinuous Galerkin (DG) finite element spaces. - * Although any function in DG element space also has nonconformity, it is completely discontinuous across edges without any relation. - * It is a reason why usual weak formulations for DG schemes contain additional penalty terms for jump across edges to control discontinuity. - * However nonconforming elements usually do not need additional terms in their weak formulations due to the continuity in integral on edges. - + * any function in the space should have the same integral values + * on two sides of the common edge shared by two adjacent elements. + * + * Thus each function in the nonconforming element space can be discontinuous, + * not included in $H^1_0$, as functions in Discontinuous Galerkin (DG) finite + * element spaces. + * Although any function in DG space also has nonconformity, + * it is completely discontinuous across edges without any relation. + * This is a reason why usual weak formulations for DG schemes contain + * additional penalty terms for jump across edges to control discontinuity. + * However nonconforming elements usually do not need additional terms + * in their weak formulations due to the continuity in integral on edges. + * *

DOFs and Dice Rule

* Since any function in the P1 nonconforming space is piecewise linear on each element, - * the function value at the mipoint of each edge is same to the mean value on the edge. + * the function value at the midpoint of each edge is same as the mean value on the edge. * Thus the continuity of the integral value across each edge is equivalent to * the continuity of the midpoint value of each edge in this case. * - * The degrees of freedom (DOFs) on a quadrilateral are defined by midpoint values on edges. - * But these four DOFs are not independent in fact. - * The simple observation reads that any linear function on a quadrilateral - * satisfies 'dice rule': the sum of two function values at two midpoints of the edge pair on opposite + * The (local) degrees of freedom (DOFs) on a quadrilateral are defined by + * midpoint values on edges. + * But these four (local) DOFs are not independent, in fact. + * A simple observation reads that any linear function on a quadrilateral satisfies 'dice rule': + * the sum of two function values at two midpoints of the edge pair on opposite * position is equal to the sum of those of the another edge pair. * - * \phi(m_0) + \phi(m_1) = \phi(m_2) + \phi(m_3). + * $\phi(m_0) + \phi(m_1) = \phi(m_2) + \phi(m_3)$ + * for all $\phi$ in the function space where $m_j$ is the midpoint of the edge $e_j$. * * Conversely if 4 values at midpoints satisfying the dice rule are just given, * then there always exists the unique linear function which coincides with 4 midpoints values. * - * Due to the dice rule, three values at any three midpoints can determine the last value at the last midpoint. + * Due to the dice rule, three values at any three midpoints can determine + * the last value at the last midpoint. * It means that the genuine number of (independent) DOFs on a quad is 3, - * and it is the same number to the dimension of the linear polynomial space in 2D. - - + * and it is same as the dimension of the linear polynomial space in 2D. + * + * *

Shape functions

- + * * @verbatim * 2---------|---------3 * | | @@ -79,18 +89,19 @@ DEAL_II_NAMESPACE_OPEN * | | * 0---------|---------1 * @endverbatim - - * For each vertex v_j of given quad, there are two edges of which v_j is one of end points. + * + * For each vertex $v_j$ of given quad, there are two edges of which $v_j$ is one of end points. * Consider a linear function such that 0.5 value at two midpoints of such edges, * and 0.0 at two midpoints of other edges. * Note that the set of these values satisfies the dice rule which is described above. - * We denote such a function assoicated with vertex v_j by \phi_j. - - * The canonical (local) basis functions are given as any three shape functions of - * the following four linear functions: - + * We denote such a function associated with vertex $v_j$ by $\phi_j$. + * + * The following figures represent $\phi_j$, $j=0,\cdots,3$ with its values at midpoints. + * Canonical (local) basis functions are given by any three shape functions among + * the following four linear functions. + * * - + * * Note that above shape functions are constructed on each cell, not on the reference cell only. - * @p get_linear_shape computes the coefficients for shape functions when @p fill_fe_values is called on each cell. - - * The (global) basis function associated with a node is defined by the composition of + * @p get_linear_shape computes the coefficients for shape functions when @p fill_fe_values + * is called on each cell. + * + * The basis function associated with a node in global scope is defined by the composition of * (local) basis functions associated with the node on each element. - * When a problem with homogeneous Dirichlet boundary condition is dealt, - * the total number of DOFs is equal to the number of interior nodes, as the standard bilinear finite element @p Q_1. - + * And DOF associated with the node represents the coefficient of the basis function associated with the node. + * In this frame, all DOFs are truly independent. + * When a problem with homogeneous Dirichlet boundary condition is considered, the total number of DOFs + * is equal to the number of interior nodes, as the number of DOFs which the standard bilinear + * finite element @p Q_1 has. + * *

Unit support points

- * Contrast with ordinary Lagrange finite elements, DOF value with respect to the P1 nonconforming element at given node does not coincide with the function value at that node. - * For instance, the (global) basis function associated with a node has 0.75 at that node, not 1.0. - * Thus we need an interpolation operator which maps any smooth function into a function with proper DOF values in the P1 element space. - * One natural interpolant associated with given smooth function is the linear function whose midpoint value at each edge is defined by - * the average of two values at endpoints of the edge. + * Contrast with ordinary Lagrange finite elements, DOF value with respect to the P1 nonconforming + * element at given node does not coincide with the function value at that node. + * For instance, the basis function associated with a node has value 0.75 at that node, not 1.0. + * Thus we need an interpolation operator which maps any smooth function into a function + * with proper DOF values in the P1 element space. + * One natural interpolant associated with given smooth function is the linear function whose midpoint + * value at each edge is defined by the average of two values at endpoints of the edge. * It provides appropriate weights used in @p unit_support_points. - + * *

References

* You can find the paper about the P1NC element at * http://epubs.siam.org/doi/abs/10.1137/S0036142902404923. - + * **/ class FE_P1NC : public FiniteElement<2,2> @@ -191,7 +208,7 @@ public: /** * Return the update flags which are needed. */ - virtual UpdateFlags requires_update_flags (const UpdateFlags flags) const ; + virtual UpdateFlags requires_update_flags (const UpdateFlags flags) const ; /** * Copy constructor. @@ -207,33 +224,24 @@ public: private: - static - std::vector - get_nonzero_component(); - + static std::vector get_nonzero_component(); /** * Return the vector consists of the numbers of degrees of freedom per objects. */ - static - std::vector - get_dpo_vector (); - + static std::vector get_dpo_vector (); /** * Compute the linear shape functions phi(x,y) = ax + by + c * such that each midpoint value on two connecting edges is a half, * and two other midpoint values are all zero. */ - static - void + static void get_linear_shape (const Triangulation<2,2>::cell_iterator &cell, std::vector &a, std::vector &b, std::vector &c); - - /** * Do the work which is needed before cellwise data computation. * Since the basis functions are constructed independently on each cell, @@ -245,14 +253,10 @@ private: const Quadrature<2> &, dealii::internal::FEValues::FiniteElementRelatedData<2,2> &output_data) const ; - - - /** * Compute the data on the current cell. */ - virtual - void + virtual void fill_fe_values (const Triangulation<2,2>::cell_iterator &cell, const CellSimilarity::Similarity , const Quadrature<2> &quadrature, @@ -262,13 +266,10 @@ private: const FiniteElement<2,2>::InternalDataBase &fe_internal, internal::FEValues::FiniteElementRelatedData<2,2> &output_data) const; - - /** * Compute the data on the face of the current cell. */ - virtual - void + virtual void fill_fe_face_values (const Triangulation<2,2>::cell_iterator &cell, const unsigned int face_no, const Quadrature<1> &quadrature, @@ -278,14 +279,10 @@ private: const InternalDataBase &fe_internal, dealii::internal::FEValues::FiniteElementRelatedData<2,2> &output_data) const; - - - /** * Compute the data on the subface of the current cell. */ - virtual - void + virtual void fill_fe_subface_values (const Triangulation<2,2>::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, @@ -296,8 +293,6 @@ private: const InternalDataBase &fe_internal, dealii::internal::FEValues::FiniteElementRelatedData<2,2> &output_data) const; - - /** * Create the constraints matrix for hanging edges. */ diff --git a/source/fe/fe_p1nc.cc b/source/fe/fe_p1nc.cc index 8b85842e83..cacaf1c7ff 100644 --- a/source/fe/fe_p1nc.cc +++ b/source/fe/fe_p1nc.cc @@ -34,7 +34,6 @@ FE_P1NC::FE_P1NC() unit_face_support_points[0][0] = 0.0 ; unit_face_support_points[1][0] = 1.0 ; - // initialize constraints matrix initialize_constraints () ; } @@ -47,7 +46,8 @@ std::string FE_P1NC::get_name () const } -UpdateFlags FE_P1NC::requires_update_flags (const UpdateFlags flags) const + +UpdateFlags FE_P1NC::requires_update_flags (const UpdateFlags flags) const { UpdateFlags out = update_default; @@ -59,21 +59,23 @@ UpdateFlags FE_P1NC::requires_update_flags (const UpdateFlags flags) const out |= update_cell_normal_vectors | update_JxW_values; if (flags & update_hessians) out |= update_hessians; - if (flags & update_hessians) - out |= update_hessians; return out; } + FiniteElement<2,2> *FE_P1NC::clone () const { return new FE_P1NC(*this); } + + FE_P1NC::~FE_P1NC () {} + std::vector FE_P1NC::get_nonzero_component() { @@ -85,6 +87,7 @@ FE_P1NC::get_nonzero_component() } + std::vector FE_P1NC::get_dpo_vector () { @@ -96,6 +99,7 @@ FE_P1NC::get_dpo_vector () } + void FE_P1NC::get_linear_shape (const Triangulation<2,2>::cell_iterator &cell, std::vector &a, @@ -122,8 +126,6 @@ FE_P1NC::get_linear_shape (const Triangulation<2,2>::cell_iterator &cell, cpt(0) = (mpt[0](0) + mpt[1](0) + mpt[2](0) + mpt[3](0))/4.0 ; cpt(1) = (mpt[0](1) + mpt[1](1) + mpt[2](1) + mpt[3](1))/4.0 ; - - double det ; det = (mpt[0](0)-mpt[1](0))*(mpt[2](1)-mpt[3](1)) - (mpt[2](0)-mpt[3](0))*(mpt[0](1)-mpt[1](1)) ; @@ -141,12 +143,10 @@ FE_P1NC::get_linear_shape (const Triangulation<2,2>::cell_iterator &cell, c[1] = 0.25 - cpt(0)*a[1] - cpt(1)*b[1] ; c[2] = 0.25 - cpt(0)*a[2] - cpt(1)*b[2] ; c[3] = 0.25 - cpt(0)*a[3] - cpt(1)*b[3] ; - } - FiniteElement<2,2>::InternalDataBase * FE_P1NC::get_data (const UpdateFlags update_flags, const Mapping<2,2> &, @@ -158,10 +158,6 @@ FE_P1NC::get_data (const UpdateFlags update_flags, data->update_each = requires_update_flags(update_flags); return data; - - - - } @@ -208,7 +204,6 @@ FE_P1NC::fill_fe_values (const Triangulation<2,2>::cell_iterator &cell output_data.shape_gradients[k][i] = grads[k]; } } - } // hessian @@ -242,13 +237,12 @@ FE_P1NC::fill_fe_values (const Triangulation<2,2>::cell_iterator &cell realquadrature = mapping.transform_unit_to_real_cell(cell, quadrature.point(i)) ; values[k] = a[k]*realquadrature(0) + b[k]*realquadrature(1) + c[k] ; output_data.shape_values[k][i] = values[k]; - } } - } + void FE_P1NC::fill_fe_face_values (const Triangulation<2,2>::cell_iterator &cell, const unsigned int face_no, @@ -320,18 +314,12 @@ FE_P1NC::fill_fe_face_values (const Triangulation<2,2>::cell_iterator grads[k][1] = b[k] ; output_data.shape_gradients[k][i] = grads[k]; } - } - } - } - - - void FE_P1NC::fill_fe_subface_values (const Triangulation<2,2>::cell_iterator &cell, const unsigned int face_no, @@ -376,7 +364,6 @@ FE_P1NC::fill_fe_subface_values (const Triangulation<2,2>::cell_iterator output_data.shape_gradients[k][i] = grads[k]; } } - } // When this function is called for computation of facial jump residual, @@ -404,13 +391,12 @@ FE_P1NC::fill_fe_subface_values (const Triangulation<2,2>::cell_iterator grads[k][1] = b[k] ; output_data.shape_gradients[k][i] = grads[k]; } - } - } } + void FE_P1NC::initialize_constraints () { @@ -423,16 +409,10 @@ void FE_P1NC::initialize_constraints () interface_constraints .TableBase<2,double>::reinit (interface_constraints_size()); - interface_constraints(0,0) = 0.5 ; interface_constraints(0,1) = 0.5 ; } - - - - - DEAL_II_NAMESPACE_CLOSE