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_bubbles_h 18 #define dealii_fe_q_bubbles_h 19 20 #include <deal.II/base/config.h> 21 22 #include <deal.II/base/tensor_product_polynomials_bubbles.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 Q_p^+ that yields the 34 * finite element space of continuous, piecewise polynomials of degree @p p in 35 * each coordinate direction plus some bubble enrichment space spanned by 36 * $(2x_j-1)^{p-1}\prod_{i=0}^{dim-1}(x_i(1-x_i))$. Therefore the highest 37 * polynomial degree is $p+1$. This class is realized using tensor product 38 * polynomials based on equidistant or given support points. 39 * 40 * The standard constructor of this class takes the degree @p p of this finite 41 * element. Alternatively, it can take a quadrature formula @p points defining 42 * the support points of the Lagrange interpolation in one coordinate 43 * direction. 44 * 45 * For more information about the <tt>spacedim</tt> template parameter check 46 * the documentation of FiniteElement or the one of Triangulation. 47 * 48 * Due to the fact that the enrichments are small almost everywhere for large 49 * p, the condition number for the mass and stiffness matrix fastly 50 * increaseses with increasing p. Below you see a comparison with 51 * FE_Q(QGaussLobatto(p+1)) for dim=1. 52 * 53 * <p ALIGN="center"> 54 * @image html fe_q_bubbles_conditioning.png 55 * </p> 56 * 57 * Therefore, this element should be used with care for $p>3$. 58 * 59 * <h3>Implementation</h3> 60 * 61 * The constructor creates a TensorProductPolynomials object that includes the 62 * tensor product of @p LagrangeEquidistant polynomials of degree @p p plus 63 * the bubble enrichments. This @p TensorProductPolynomialsBubbles object 64 * provides all values and derivatives of the shape functions. In case a 65 * quadrature rule is given, the constructor creates a 66 * TensorProductPolynomialsBubbles object that includes the tensor product of 67 * @p Lagrange polynomials with the support points from @p points and the 68 * bubble enrichments as defined above. 69 * 70 * Furthermore the constructor fills the @p interface_constrains, the @p 71 * prolongation (embedding) and the @p restriction matrices. 72 * 73 * <h3>Numbering of the degrees of freedom (DoFs)</h3> 74 * 75 * The original ordering of the shape functions represented by the 76 * TensorProductPolynomialsBubbles is a tensor product numbering. However, the 77 * shape functions on a cell are renumbered beginning with the shape functions 78 * whose support points are at the vertices, then on the line, on the quads, 79 * and finally (for 3d) on the hexes. Finally, there are support points for 80 * the bubble enrichments in the middle of the cell. 81 */ 82 template <int dim, int spacedim = dim> 83 class FE_Q_Bubbles 84 : public FE_Q_Base<TensorProductPolynomialsBubbles<dim>, dim, spacedim> 85 { 86 public: 87 /** 88 * Constructor for tensor product polynomials of degree @p p plus bubble 89 * enrichments 90 */ 91 FE_Q_Bubbles(const unsigned int p); 92 93 /** 94 * Constructor for tensor product polynomials with support points @p points 95 * plus bubble enrichments based on a one-dimensional quadrature formula. 96 * The degree of the finite element is <tt>points.size()</tt>. Note that the 97 * first point has to be 0 and the last one 1. 98 */ 99 FE_Q_Bubbles(const Quadrature<1> &points); 100 101 /** 102 * Return a string that uniquely identifies a finite element. This class 103 * returns <tt>FE_Q_Bubbles<dim>(degree)</tt>, with @p dim and @p degree 104 * replaced by appropriate values. 105 */ 106 virtual std::string 107 get_name() const override; 108 109 // documentation inherited from the base class 110 virtual void 111 convert_generalized_support_point_values_to_dof_values( 112 const std::vector<Vector<double>> &support_point_values, 113 std::vector<double> & nodal_values) const override; 114 115 /** 116 * Return the matrix interpolating from the given finite element to the 117 * present one. The size of the matrix is then @p dofs_per_cell times 118 * <tt>source.n_dofs_per_cell()</tt>. 119 * 120 * These matrices are only available if the source element is also a @p 121 * FE_Q_Bubbles element. Otherwise, an exception of type 122 * FiniteElement<dim,spacedim>::ExcInterpolationNotImplemented is thrown. 123 */ 124 virtual void 125 get_interpolation_matrix(const FiniteElement<dim, spacedim> &source, 126 FullMatrix<double> &matrix) const override; 127 128 virtual const FullMatrix<double> & 129 get_prolongation_matrix( 130 const unsigned int child, 131 const RefinementCase<dim> &refinement_case) const override; 132 133 virtual const FullMatrix<double> & 134 get_restriction_matrix( 135 const unsigned int child, 136 const RefinementCase<dim> &refinement_case) const override; 137 138 /** 139 * Check for non-zero values on a face. 140 * 141 * This function returns @p true, if the shape function @p shape_index has 142 * non-zero values on the face @p face_index. 143 * 144 * Implementation of the interface in FiniteElement 145 */ 146 virtual bool 147 has_support_on_face(const unsigned int shape_index, 148 const unsigned int face_index) const override; 149 150 virtual std::unique_ptr<FiniteElement<dim, spacedim>> 151 clone() const override; 152 153 /** 154 * @copydoc FiniteElement::compare_for_domination() 155 */ 156 virtual FiniteElementDomination::Domination 157 compare_for_domination(const FiniteElement<dim, spacedim> &fe_other, 158 const unsigned int codim = 0) const override final; 159 160 private: 161 /** 162 * Return the restriction_is_additive flags. Only the last components for 163 * the bubble enrichments are true. 164 */ 165 static std::vector<bool> 166 get_riaf_vector(const unsigned int degree); 167 168 /** 169 * Only for internal use. Its full name is @p get_dofs_per_object_vector 170 * function and it creates the @p dofs_per_object vector that is needed 171 * within the constructor to be passed to the constructor of @p 172 * FiniteElementData. 173 */ 174 static std::vector<unsigned int> 175 get_dpo_vector(const unsigned int degree); 176 177 /** 178 * Number of additional bubble functions 179 */ 180 const unsigned int n_bubbles; 181 }; 182 183 184 185 /*@}*/ 186 187 188 DEAL_II_NAMESPACE_CLOSE 189 190 #endif 191