1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2012 - 2018 by the deal.II authors 4 // 5 // This file is part of the deal.II library. 6 // 7 // The deal.II library is free software; you can use it, redistribute 8 // it, and/or modify it under the terms of the GNU Lesser General 9 // Public License as published by the Free Software Foundation; either 10 // version 2.1 of the License, or (at your option) any later version. 11 // The full text of the license can be found in the file LICENSE.md at 12 // the top level directory of deal.II. 13 // 14 // --------------------------------------------------------------------- 15 16 17 #ifndef dealii_fe_q_dg0_h 18 #define dealii_fe_q_dg0_h 19 20 #include <deal.II/base/config.h> 21 22 #include <deal.II/base/tensor_product_polynomials_const.h> 23 24 #include <deal.II/fe/fe_q_base.h> 25 26 DEAL_II_NAMESPACE_OPEN 27 28 29 /*!@addtogroup fe */ 30 /*@{*/ 31 32 /** 33 * Implementation of a scalar Lagrange finite element @p Qp+DG0 that yields 34 * the finite element space of continuous, piecewise polynomials of degree @p 35 * p in each coordinate direction plus the space of locally constant 36 * functions. This class is realized using tensor product polynomials based on 37 * equidistant or given support points. 38 * 39 * The standard constructor of this class takes the degree @p p of this finite 40 * element. Alternatively, it can take a quadrature formula @p points defining 41 * the support points of the Lagrange interpolation in one coordinate 42 * direction. 43 * 44 * For more information about the <tt>spacedim</tt> template parameter check 45 * the documentation of FiniteElement or the one of Triangulation. 46 * 47 * For more information regarding this element see: Boffi, D., et al. "Local 48 * Mass Conservation of Stokes Finite Elements." Journal of Scientific 49 * Computing (2012): 1-18. 50 * 51 * <h3>Implementation</h3> 52 * 53 * The constructor creates a TensorProductPolynomials object that includes the 54 * tensor product of @p LagrangeEquidistant polynomials of degree @p p plus 55 * the locally constant function. This @p TensorProductPolynomialsConst object 56 * provides all values and derivatives of the shape functions. In case a 57 * quadrature rule is given, the constructor creates a 58 * TensorProductPolynomialsConst object that includes the tensor product of @p 59 * Lagrange polynomials with the support points from @p points and a locally 60 * constant function. 61 * 62 * Furthermore the constructor fills the @p interface_constrains, the @p 63 * prolongation (embedding) and the @p restriction matrices. 64 * 65 * <h3>Numbering of the degrees of freedom (DoFs)</h3> 66 * 67 * The original ordering of the shape functions represented by the 68 * TensorProductPolynomialsConst is a tensor product numbering. However, the 69 * shape functions on a cell are renumbered beginning with the shape functions 70 * whose support points are at the vertices, then on the line, on the quads, 71 * and finally (for 3d) on the hexes. Finally there is a support point for the 72 * discontinuous shape function in the middle of the cell. To be explicit, 73 * these numberings are listed in the following: 74 * 75 * <h4>Q1 elements</h4> 76 * <ul> 77 * <li> 1D case: 78 * @verbatim 79 * 0---2---1 80 * @endverbatim 81 * 82 * <li> 2D case: 83 * @verbatim 84 * 2-------3 85 * | | 86 * | 5 | 87 * | | 88 * 0-------1 89 * @endverbatim 90 * 91 * <li> 3D case: 92 * @verbatim 93 * 6-------7 6-------7 94 * /| | / /| 95 * / | | / / | 96 * / | | / / | 97 * 4 | 8 | 4-------5 | 98 * | 2-------3 | | 3 99 * | / / | | / 100 * | / / | | / 101 * |/ / | |/ 102 * 0-------1 0-------1 103 * @endverbatim 104 * 105 * The respective coordinate values of the support points of the degrees of 106 * freedom are as follows: 107 * <ul> 108 * <li> Index 0: <tt>[ 0, 0, 0]</tt>; 109 * <li> Index 1: <tt>[ 1, 0, 0]</tt>; 110 * <li> Index 2: <tt>[ 0, 1, 0]</tt>; 111 * <li> Index 3: <tt>[ 1, 1, 0]</tt>; 112 * <li> Index 4: <tt>[ 0, 0, 1]</tt>; 113 * <li> Index 5: <tt>[ 1, 0, 1]</tt>; 114 * <li> Index 6: <tt>[ 0, 1, 1]</tt>; 115 * <li> Index 7: <tt>[ 1, 1, 1]</tt>; 116 * <li> Index 8: <tt>[1/2, 1/2, 1/2]</tt>; 117 * </ul> 118 * </ul> 119 * <h4>Q2 elements</h4> 120 * <ul> 121 * <li> 1D case: 122 * @verbatim 123 * 0---2---1 124 * @endverbatim 125 * Index 3 has the same coordinates as index 2 126 * 127 * <li> 2D case: 128 * @verbatim 129 * 2---7---3 130 * | | 131 * 4 8 5 132 * | | 133 * 0---6---1 134 * @endverbatim 135 * Index 9 has the same coordinates as index 2 136 * 137 * <li> 3D case: 138 * @verbatim 139 * 6--15---7 6--15---7 140 * /| | / /| 141 * 12 | 19 12 1319 142 * / 18 | / / | 143 * 4 | | 4---14--5 | 144 * | 2---11--3 | | 3 145 * | / / | 17 / 146 * 16 8 9 16 | 9 147 * |/ / | |/ 148 * 0---10--1 0---8---1 149 * 150 * *-------* *-------* 151 * /| | / /| 152 * / | 23 | / 25 / | 153 * / | | / / | 154 * * | | *-------* | 155 * |20 *-------* | |21 * 156 * | / / | 22 | / 157 * | / 24 / | | / 158 * |/ / | |/ 159 * *-------* *-------* 160 * @endverbatim 161 * The center vertices have number 26 and 27. 162 * 163 * The respective coordinate values of the support points of the degrees of 164 * freedom are as follows: 165 * <ul> 166 * <li> Index 0: <tt>[0, 0, 0]</tt>; 167 * <li> Index 1: <tt>[1, 0, 0]</tt>; 168 * <li> Index 2: <tt>[0, 1, 0]</tt>; 169 * <li> Index 3: <tt>[1, 1, 0]</tt>; 170 * <li> Index 4: <tt>[0, 0, 1]</tt>; 171 * <li> Index 5: <tt>[1, 0, 1]</tt>; 172 * <li> Index 6: <tt>[0, 1, 1]</tt>; 173 * <li> Index 7: <tt>[1, 1, 1]</tt>; 174 * <li> Index 8: <tt>[0, 1/2, 0]</tt>; 175 * <li> Index 9: <tt>[1, 1/2, 0]</tt>; 176 * <li> Index 10: <tt>[1/2, 0, 0]</tt>; 177 * <li> Index 11: <tt>[1/2, 1, 0]</tt>; 178 * <li> Index 12: <tt>[0, 1/2, 1]</tt>; 179 * <li> Index 13: <tt>[1, 1/2, 1]</tt>; 180 * <li> Index 14: <tt>[1/2, 0, 1]</tt>; 181 * <li> Index 15: <tt>[1/2, 1, 1]</tt>; 182 * <li> Index 16: <tt>[0, 0, 1/2]</tt>; 183 * <li> Index 17: <tt>[1, 0, 1/2]</tt>; 184 * <li> Index 18: <tt>[0, 1, 1/2]</tt>; 185 * <li> Index 19: <tt>[1, 1, 1/2]</tt>; 186 * <li> Index 20: <tt>[0, 1/2, 1/2]</tt>; 187 * <li> Index 21: <tt>[1, 1/2, 1/2]</tt>; 188 * <li> Index 22: <tt>[1/2, 0, 1/2]</tt>; 189 * <li> Index 23: <tt>[1/2, 1, 1/2]</tt>; 190 * <li> Index 24: <tt>[1/2, 1/2, 0]</tt>; 191 * <li> Index 25: <tt>[1/2, 1/2, 1]</tt>; 192 * <li> Index 26: <tt>[1/2, 1/2, 1/2]</tt>; 193 * <li> Index 27: <tt>[1/2, 1/2, 1/2]</tt>; 194 * </ul> 195 * </ul> 196 * <h4>Q3 elements</h4> 197 * <ul> 198 * <li> 1D case: 199 * @verbatim 200 * 0--2-4-3--1 201 * @endverbatim 202 * 203 * <li> 2D case: 204 * @verbatim 205 * 2--10-11-3 206 * | | 207 * 5 14 15 7 208 * | 16 | 209 * 4 12 13 6 210 * | | 211 * 0--8--9--1 212 * @endverbatim 213 * </ul> 214 * <h4>Q4 elements</h4> 215 * <ul> 216 * <li> 1D case: 217 * @verbatim 218 * 0--2--3--4--1 219 * @endverbatim 220 * Index 5 has the same coordinates as index 3 221 * 222 * <li> 2D case: 223 * @verbatim 224 * 2--13-14-15-3 225 * | | 226 * 6 22 23 24 9 227 * | | 228 * 5 19 20 21 8 229 * | | 230 * 4 16 17 18 7 231 * | | 232 * 0--10-11-12-1 233 * @endverbatim 234 * Index 21 has the same coordinates as index 20 235 * </ul> 236 */ 237 template <int dim, int spacedim = dim> 238 class FE_Q_DG0 239 : public FE_Q_Base<TensorProductPolynomialsConst<dim>, dim, spacedim> 240 { 241 public: 242 /** 243 * Constructor for tensor product polynomials of degree @p p plus locally 244 * constant functions. 245 */ 246 FE_Q_DG0(const unsigned int p); 247 248 /** 249 * Constructor for tensor product polynomials with support points @p points 250 * plus locally constant functions based on a one-dimensional quadrature 251 * formula. The degree of the finite element is <tt>points.size()-1</tt>. 252 * Note that the first point has to be 0 and the last one 1. 253 */ 254 FE_Q_DG0(const Quadrature<1> &points); 255 256 /** 257 * Return a string that uniquely identifies a finite element. This class 258 * returns <tt>FE_Q_DG0<dim>(degree)</tt>, with @p dim and @p degree 259 * replaced by appropriate values. 260 */ 261 virtual std::string 262 get_name() const override; 263 264 // documentation inherited from the base class 265 virtual void 266 convert_generalized_support_point_values_to_dof_values( 267 const std::vector<Vector<double>> &support_point_values, 268 std::vector<double> & nodal_values) const override; 269 270 /** 271 * Return the matrix interpolating from the given finite element to the 272 * present one. The size of the matrix is then @p dofs_per_cell times 273 * <tt>source.n_dofs_per_cell()</tt>. 274 * 275 * These matrices are only available if the source element is also a @p 276 * FE_Q_DG0 element. Otherwise, an exception of type 277 * FiniteElement<dim,spacedim>::ExcInterpolationNotImplemented is thrown. 278 */ 279 virtual void 280 get_interpolation_matrix(const FiniteElement<dim, spacedim> &source, 281 FullMatrix<double> &matrix) const override; 282 283 284 /** 285 * This function returns @p true, if the shape function @p shape_index has 286 * non-zero function values somewhere on the face @p face_index. 287 */ 288 virtual bool 289 has_support_on_face(const unsigned int shape_index, 290 const unsigned int face_index) const override; 291 292 /** 293 * Return a list of constant modes of the element. For this element, there 294 * are two constant modes despite the element is scalar: The first constant 295 * mode is all ones for the usual FE_Q basis and the second one only using 296 * the discontinuous part. 297 */ 298 virtual std::pair<Table<2, bool>, std::vector<unsigned int>> 299 get_constant_modes() const override; 300 301 virtual std::unique_ptr<FiniteElement<dim, spacedim>> 302 clone() const override; 303 304 /** 305 * @copydoc FiniteElement::compare_for_domination() 306 */ 307 virtual FiniteElementDomination::Domination 308 compare_for_domination(const FiniteElement<dim, spacedim> &fe_other, 309 const unsigned int codim = 0) const override final; 310 311 private: 312 /** 313 * Return the restriction_is_additive flags. Only the last component is 314 * true. 315 */ 316 static std::vector<bool> 317 get_riaf_vector(const unsigned int degree); 318 319 /** 320 * Only for internal use. Its full name is @p get_dofs_per_object_vector 321 * function and it creates the @p dofs_per_object vector that is needed 322 * within the constructor to be passed to the constructor of @p 323 * FiniteElementData. 324 */ 325 static std::vector<unsigned int> 326 get_dpo_vector(const unsigned int degree); 327 }; 328 329 330 331 /*@}*/ 332 333 334 DEAL_II_NAMESPACE_CLOSE 335 336 #endif 337