From: Jaeryun Yim Date: Tue, 27 Sep 2016 13:44:04 +0000 (+0900) Subject: Edit comments for the description. X-Git-Tag: v8.5.0-rc1~624^2~11 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b703fc25143c03552b406ef568b83538562855c3;p=dealii.git Edit comments for the description. --- diff --git a/include/deal.II/fe/fe_p1nc.h b/include/deal.II/fe/fe_p1nc.h index 23f948fc44..b4799077ed 100644 --- a/include/deal.II/fe/fe_p1nc.h +++ b/include/deal.II/fe/fe_p1nc.h @@ -49,32 +49,43 @@ DEAL_II_NAMESPACE_OPEN * 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

+ *

Dice Rule

* Since any function in the P1 nonconforming space is piecewise linear on each element, * 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 (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': + * Thus for the P1 nonconforming element, the function values at midpoints on edges of a cell are important. + * The first attempt to define (local) degrees of freedom (DOFs) on a quadrilateral + * is by using midpoint values of a function as usual nonconforming finite elements. + * + * However, these 4 functionals are not linearly independent + * because a linear function on 2D is uniquely determined by only 3 independent values. + * A simple observation reads that any linear function on a quadrilateral should satisfies the '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. + * position is equal to the sum of those of the other edge pair. + * This is called the 'dice rule' because the number of points on opposite sides of a dice always + * adds up to the same number as well (in the case of dice, to seven). * - * $\phi(m_0) + \phi(m_1) = \phi(m_2) + \phi(m_3)$ + * In formulas, the dice rule is written as $\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$. + * Here, we assume the standard numbering convention for edges used in deal.II + * and described for class GeometryInfo. * * 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. - * It means that the genuine number of (independent) DOFs on a quad is 3, - * and it is same as the dimension of the linear polynomial space in 2D. - * + * It means that the number of independent local functionals on a cell is 3, + * and it is same as the dimension of the linear polynomial space on a cell in 2D. * *

Shape functions

+ * Before introducing the DOFs, we present 4 local shape functions on a cell. + * Due to the dice rule, we need a special construction for shape functions. + * Although the following 4 shape functions are not linearly independent within a cell, + * they are helpful to define the global basis functions which are linearly independent on whole domain. + * Again, we assume the standard numbering for vertices used in deal.II. * * @verbatim * 2---------|---------3 @@ -90,15 +101,14 @@ 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 cell, 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 associated with vertex $v_j$ by $\phi_j$. + * Then the set of 4 shape functions is a partition of unity on a cell: $\sum_{j=0}^{3} \phi_j = 1$. * * 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 local DOFs are defined by the coefficients of the shape functions associated vertices, respectively. + * Although these 4 local DOFs are not linearly independent within a single cell as well, + * this definition is a good start point for the definition of the global DOFs. + * + * We want to emphasize that the shape functions are constructed on each cell, not on the reference cell only. + * Usual finite elements are defined based on a 'parametric' concept. + * That means, a function space for a finite element is defined on one reference cell, and it is transfomed + * into each cell via a mapping from the reference cell. + * However the P1 nonconforming element does not follow such concept. It defines a function space with + * linear shape functions on each cell without any help of a function space on the reference cell. + * In other words, the element is defined in real space, not via a mapping from a reference cell. + * + * Thus this implementation does not have to compute shape values on the reference cell. + * Rather than, the shape values are computed by construction of the shape functions + * on each cell independently. + * + *

DOFs

+ * We have to consider the basis function for the element space in global domain + * because the system of equations which we have to solve at last is for a global system, not local. + * The global basis function associated with a node is defined by a cell-wise composition of + * local shape functions associated with the node on each element. + * And we define a global DOF associated with a node by a coefficient of the basis function associated with that node. + * + * There is a theoretical result about the linear independency of the global basis functions + * depending on the type of the boundary condition we consider. + * + * When the homogeneous Dirichlet boundary condition is given, + * the global basis functions associated with interior nodes are linearly independent. + * And the number of DOFs is equal to the number of interior nodes, + * same as the number of DOFs for the standard bilinear finite element @p Q_1. * - * 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. - * 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. + * When the Neumann boundary condition is given, + * the global basis functions associated with all nodes (including boundary nodes) + * are actually not linearly independent. There exists 1 redundancy. + * Thus in this case, the number of DOFs is equal to the number of all nodes minus 1. * *

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 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. + * For a smooth function, we construct a piecewise linear function which belongs to the element space by + * using its nodal values as DOF values. + * This interpolation is implemented by using appropriate @p unit_support_points. + * + * Note that two nodal values of a smooth function and its interpolant do not coincide in general, + * contrast with ordinary Lagrange finite elements. + * Of course, it is meaningless to refer 'nodal value' because the element space has nonconformity. + * But it is also true even though the single global basis function associated a node is considered + * with a valid expression 'nodal value'. + * For instance, consider the basis function associated with a node. + * Consider two lines representing the level sets for value 0.5 and 0, respectively, by connecting two midpoints. + * Then we cut the quad into two sub-triangles by the diagonal which is placed along those two lines. + * It gives another level set for value 0.25 which coincides with the cutting diagonal. + * Therefore these three level sets are all parallel and it gives the value 0.75 at the base node, not value 1. + * Even though a general quad is given, this is also true. * *

References

- * You can find the paper about the P1NC element at - * http://epubs.siam.org/doi/abs/10.1137/S0036142902404923. + * You can find the paper about the P1NC element + * Park & Sheen (2003). P1-nonconforming quadrilateral finite element methods for second-order elliptic problems. + * SIAM Journal on Numerical Analysis, 41(2), 624-640, + * available at http://epubs.siam.org/doi/abs/10.1137/S0036142902404923. * - **/ + */ class FE_P1NC : public FiniteElement<2,2> { public: /** - * Constructor for nonparametric version of P1 nonconforming element. + * Constructor for the P1 nonconforming element. + * It is only for 2D and codimension = 0. */ FE_P1NC() ; - /** - * Return the name of the class for the element. - */ virtual std::string get_name () const ; - /** - * Return the update flags which are needed. - */ virtual UpdateFlags requires_update_flags (const UpdateFlags flags) const ; - /** - * Copy constructor. - */ virtual FiniteElement<2,2> *clone () const ; /** @@ -232,9 +267,10 @@ private: 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. + * Compute the values of the variables a,b and c which are for the coefficients of + * the standard linear shape function $\phi_j(x,y) = a_j x + b_j y + c_j$ on given cell. + * Since there are 4 local shape functions on each cell, + * each variable is an array consisting of 4 values which are coefficients corresponding to 4 shape functions, respectively. */ static void get_linear_shape (const Triangulation<2,2>::cell_iterator &cell, @@ -244,8 +280,9 @@ private: /** * Do the work which is needed before cellwise data computation. - * Since the basis functions are constructed independently on each cell, + * Since the shape functions are constructed independently on each cell, * the data on the reference cell is not necessary. + * It returns an empty variable type of @ InternalDataBase and updates @ update_flags. */ virtual FiniteElement<2,2>::InternalDataBase * get_data (const UpdateFlags update_flags,