Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_q.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 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 #include <deal.II/base/quadrature_lib.h>
18 #include <deal.II/base/std_cxx14/memory.h>
19 
20 #include <deal.II/fe/fe_dgq.h>
21 #include <deal.II/fe/fe_nothing.h>
22 #include <deal.II/fe/fe_q.h>
23 
24 #include <deal.II/lac/vector.h>
25 
26 #include <sstream>
27 #include <vector>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 
32 namespace internal
33 {
34  namespace FE_Q
35  {
36  namespace
37  {
38  std::vector<Point<1>>
39  get_QGaussLobatto_points(const unsigned int degree)
40  {
41  if (degree > 0)
42  return QGaussLobatto<1>(degree + 1).get_points();
43  else
44  {
45  using FEQ = ::FE_Q_Base<TensorProductPolynomials<1>, 1, 1>;
46  AssertThrow(false, FEQ::ExcFEQCannotHaveDegree0());
47  }
48  return std::vector<Point<1>>();
49  }
50  } // namespace
51  } // namespace FE_Q
52 } // namespace internal
53 
54 
55 
56 template <int dim, int spacedim>
57 FE_Q<dim, spacedim>::FE_Q(const unsigned int degree)
58  : FE_Q_Base<TensorProductPolynomials<dim>, dim, spacedim>(
60  Polynomials::generate_complete_Lagrange_basis(
61  internal::FE_Q::get_QGaussLobatto_points(degree))),
62  FiniteElementData<dim>(this->get_dpo_vector(degree),
63  1,
64  degree,
65  FiniteElementData<dim>::H1),
66  std::vector<bool>(1, false))
67 {
68  this->initialize(internal::FE_Q::get_QGaussLobatto_points(degree));
69 }
70 
71 
72 
73 template <int dim, int spacedim>
75  : FE_Q_Base<TensorProductPolynomials<dim>, dim, spacedim>(
77  Polynomials::generate_complete_Lagrange_basis(points.get_points())),
78  FiniteElementData<dim>(this->get_dpo_vector(points.size() - 1),
79  1,
80  points.size() - 1,
81  FiniteElementData<dim>::H1),
82  std::vector<bool>(1, false))
83 {
84  this->initialize(points.get_points());
85 }
86 
87 
88 
89 template <int dim, int spacedim>
90 std::string
92 {
93  // note that the FETools::get_fe_by_name function depends on the
94  // particular format of the string this function returns, so they have to be
95  // kept in synch
96 
97  std::ostringstream namebuf;
98  bool equidistant = true;
99  std::vector<double> points(this->degree + 1);
100 
101  // Decode the support points in one coordinate direction.
102  std::vector<unsigned int> lexicographic =
103  this->poly_space.get_numbering_inverse();
104  for (unsigned int j = 0; j <= this->degree; j++)
105  points[j] = this->unit_support_points[lexicographic[j]][0];
106 
107  // Check whether the support points are equidistant.
108  for (unsigned int j = 0; j <= this->degree; j++)
109  if (std::fabs(points[j] - static_cast<double>(j) / this->degree) > 1e-15)
110  {
111  equidistant = false;
112  break;
113  }
114 
115  if (equidistant == true)
116  {
117  if (this->degree > 2)
118  namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim)
119  << ">(QIterated(QTrapez()," << this->degree << "))";
120  else
121  namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim) << ">("
122  << this->degree << ")";
123  }
124  else
125  {
126  // Check whether the support points come from QGaussLobatto.
127  const QGaussLobatto<1> points_gl(this->degree + 1);
128  bool gauss_lobatto = true;
129  for (unsigned int j = 0; j <= this->degree; j++)
130  if (points[j] != points_gl.point(j)(0))
131  {
132  gauss_lobatto = false;
133  break;
134  }
135  if (gauss_lobatto == true)
136  namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim) << ">("
137  << this->degree << ")";
138  else
139  namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim)
140  << ">(QUnknownNodes(" << this->degree << "))";
141  }
142  return namebuf.str();
143 }
144 
145 
146 
147 template <int dim, int spacedim>
148 void
150  const std::vector<Vector<double>> &support_point_values,
151  std::vector<double> & nodal_values) const
152 {
153  AssertDimension(support_point_values.size(),
154  this->get_unit_support_points().size());
155  AssertDimension(support_point_values.size(), nodal_values.size());
156  AssertDimension(this->dofs_per_cell, nodal_values.size());
157 
158  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
159  {
160  AssertDimension(support_point_values[i].size(), 1);
161 
162  nodal_values[i] = support_point_values[i](0);
163  }
164 }
165 
166 
167 
168 template <int dim, int spacedim>
169 std::unique_ptr<FiniteElement<dim, spacedim>>
171 {
172  return std_cxx14::make_unique<FE_Q<dim, spacedim>>(*this);
173 }
174 
175 
176 
177 template <int dim, int spacedim>
180  const FiniteElement<dim, spacedim> &fe_other,
181  const unsigned int codim) const
182 {
183  Assert(codim <= dim, ExcImpossibleInDim(dim));
184 
185  // vertex/line/face domination
186  // (if fe_other is derived from FE_DGQ)
187  // ------------------------------------
188  if (codim > 0)
189  if (dynamic_cast<const FE_DGQ<dim, spacedim> *>(&fe_other) != nullptr)
190  // there are no requirements between continuous and discontinuous elements
192 
193  // vertex/line/face domination
194  // (if fe_other is not derived from FE_DGQ)
195  // & cell domination
196  // ----------------------------------------
197  if (const FE_Q<dim, spacedim> *fe_q_other =
198  dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other))
199  {
200  if (this->degree < fe_q_other->degree)
202  else if (this->degree == fe_q_other->degree)
204  else
206  }
207  else if (const FE_Nothing<dim, spacedim> *fe_nothing =
208  dynamic_cast<const FE_Nothing<dim, spacedim> *>(&fe_other))
209  {
210  if (fe_nothing->is_dominating())
212  else
213  // the FE_Nothing has no degrees of freedom and it is typically used
214  // in a context where we don't require any continuity along the
215  // interface
217  }
218 
219  Assert(false, ExcNotImplemented());
221 }
222 
223 
224 // explicit instantiations
225 #include "fe_q.inst"
226 
227 DEAL_II_NAMESPACE_CLOSE
virtual std::string get_name() const override
Definition: fe_q.cc:91
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1567
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_q.cc:170
FE_Q(const unsigned int p)
Definition: fe_q.cc:57
const std::vector< Point< dim > > & get_points() const
Definition: fe_dgq.h:109
const unsigned int degree
Definition: fe_base.h:298
const Point< dim > & point(const unsigned int i) const
STL namespace.
#define AssertThrow(cond, exc)
Definition: exceptions.h:1519
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const override
Definition: fe_q.cc:149
Definition: fe_q.h:554
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
std::string dim_string(const int dim, const int spacedim)
Definition: utilities.cc:457
void initialize(const std::vector< Point< 1 >> &support_points_1d)
Definition: fe_q_base.cc:452
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_q.cc:179
static ::ExceptionBase & ExcNotImplemented()