Reference documentation for deal.II version 9.2.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_trace.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2020 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 #include <deal.II/base/config.h>
17 
22 
23 #include <deal.II/fe/fe_nothing.h>
25 #include <deal.II/fe/fe_q.h>
26 #include <deal.II/fe/fe_tools.h>
27 #include <deal.II/fe/fe_trace.h>
28 
29 #include <fstream>
30 #include <iostream>
31 #include <sstream>
32 
34 
35 
36 
37 template <int dim, int spacedim>
38 FE_TraceQ<dim, spacedim>::FE_TraceQ(const unsigned int degree)
39  : FE_PolyFace<TensorProductPolynomials<dim - 1>, dim, spacedim>(
40  TensorProductPolynomials<dim - 1>(
42  QGaussLobatto<1>(degree + 1).get_points())),
43  FiniteElementData<dim>(get_dpo_vector(degree),
44  1,
45  degree,
46  FiniteElementData<dim>::L2),
47  std::vector<bool>(1, true))
48  , fe_q(degree)
49 {
50  Assert(degree > 0,
51  ExcMessage("FE_Trace can only be used for polynomial degrees "
52  "greater than zero"));
54  FETools::hierarchic_to_lexicographic_numbering<dim - 1>(degree));
55 
56  // Initialize face support points
57  this->unit_face_support_points = fe_q.get_unit_face_support_points();
58 
59  // initialize unit support points (this makes it possible to assign initial
60  // values to FE_TraceQ). Note that we simply take the points of fe_q but
61  // skip the last ones which are associated with the interior of FE_Q.
62  this->unit_support_points.resize(this->dofs_per_cell);
63  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
64  this->unit_support_points[i] = fe_q.get_unit_support_points()[i];
65 
66  // Initialize constraint matrices
67  this->interface_constraints = fe_q.constraints();
68 }
69 
70 
71 
72 template <int dim, int spacedim>
73 std::unique_ptr<FiniteElement<dim, spacedim>>
75 {
76  return std_cxx14::make_unique<FE_TraceQ<dim, spacedim>>(this->degree);
77 }
78 
79 
80 
81 template <int dim, int spacedim>
82 std::string
84 {
85  // note that the FETools::get_fe_by_name function depends on the
86  // particular format of the string this function returns, so they have to be
87  // kept in synch
88 
89  std::ostringstream namebuf;
90  namebuf << "FE_TraceQ<" << Utilities::dim_string(dim, spacedim) << ">("
91  << this->degree << ")";
92 
93  return namebuf.str();
94 }
95 
96 
97 
98 template <int dim, int spacedim>
99 bool
101  const unsigned int shape_index,
102  const unsigned int face_index) const
103 {
104  AssertIndexRange(shape_index, this->dofs_per_cell);
106 
107  // FE_TraceQ shares the numbering of elemental degrees of freedom with FE_Q
108  // except for the missing interior ones (quad dofs in 2D and hex dofs in
109  // 3D). Therefore, it is safe to ask fe_q for the corresponding
110  // information. The assertion 'shape_index < this->dofs_per_cell' will make
111  // sure that we only access the trace dofs.
112  return fe_q.has_support_on_face(shape_index, face_index);
113 }
114 
115 
116 
117 template <int dim, int spacedim>
118 std::pair<Table<2, bool>, std::vector<unsigned int>>
120 {
121  Table<2, bool> constant_modes(1, this->dofs_per_cell);
122  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
123  constant_modes(0, i) = true;
124  return std::pair<Table<2, bool>, std::vector<unsigned int>>(
125  constant_modes, std::vector<unsigned int>(1, 0));
126 }
127 
128 template <int dim, int spacedim>
129 void
132  const std::vector<Vector<double>> &support_point_values,
133  std::vector<double> & nodal_values) const
134 {
135  AssertDimension(support_point_values.size(),
136  this->get_unit_support_points().size());
137  AssertDimension(support_point_values.size(), nodal_values.size());
138  AssertDimension(this->dofs_per_cell, nodal_values.size());
139 
140  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
141  {
142  AssertDimension(support_point_values[i].size(), 1);
143 
144  nodal_values[i] = support_point_values[i](0);
145  }
146 }
147 
148 
149 template <int dim, int spacedim>
150 std::vector<unsigned int>
152 {
153  // This constructs FE_TraceQ in exactly the same way as FE_Q except for the
154  // interior degrees of freedom that are not present here (line in 1D, quad
155  // in 2D, hex in 3D).
156  AssertThrow(deg > 0, ExcMessage("FE_TraceQ needs to be of degree > 0."));
157  std::vector<unsigned int> dpo(dim + 1, 1U);
158  dpo[dim] = 0;
159  dpo[0] = 1;
160  for (unsigned int i = 1; i < dim; ++i)
161  dpo[i] = dpo[i - 1] * (deg - 1);
162  return dpo;
163 }
164 
165 
166 
167 template <int dim, int spacedim>
168 bool
170 {
171  return fe_q.hp_constraints_are_implemented();
172 }
173 
174 
175 template <int dim, int spacedim>
178  const FiniteElement<dim, spacedim> &fe_other,
179  const unsigned int codim) const
180 {
181  Assert(codim <= dim, ExcImpossibleInDim(dim));
182  (void)codim;
183 
184  // vertex/line/face/cell domination
185  // --------------------------------
186  if (const FE_TraceQ<dim, spacedim> *fe_traceq_other =
187  dynamic_cast<const FE_TraceQ<dim, spacedim> *>(&fe_other))
188  {
189  if (this->degree < fe_traceq_other->degree)
191  else if (this->degree == fe_traceq_other->degree)
193  else
195  }
196  else if (const FE_Nothing<dim> *fe_nothing =
197  dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
198  {
199  if (fe_nothing->is_dominating())
201  else
202  // the FE_Nothing has no degrees of freedom and it is typically used
203  // in a context where we don't require any continuity along the
204  // interface
206  }
207 
208  Assert(false, ExcNotImplemented());
210 }
211 
212 
213 
214 template <int dim, int spacedim>
215 void
217  const FiniteElement<dim, spacedim> &source_fe,
218  FullMatrix<double> & interpolation_matrix) const
219 {
220  get_subface_interpolation_matrix(source_fe,
222  interpolation_matrix);
223 }
224 
225 
226 
227 template <int dim, int spacedim>
228 void
230  const FiniteElement<dim, spacedim> &x_source_fe,
231  const unsigned int subface,
232  FullMatrix<double> & interpolation_matrix) const
233 {
234  // this is the code from FE_FaceQ
235  Assert(interpolation_matrix.n() == this->dofs_per_face,
236  ExcDimensionMismatch(interpolation_matrix.n(), this->dofs_per_face));
237  Assert(interpolation_matrix.m() == x_source_fe.dofs_per_face,
238  ExcDimensionMismatch(interpolation_matrix.m(),
239  x_source_fe.dofs_per_face));
240 
241  // see if source is a FaceQ element
242  if (const FE_TraceQ<dim, spacedim> *source_fe =
243  dynamic_cast<const FE_TraceQ<dim, spacedim> *>(&x_source_fe))
244  {
245  fe_q.get_subface_interpolation_matrix(source_fe->fe_q,
246  subface,
247  interpolation_matrix);
248  }
249  else if (dynamic_cast<const FE_Nothing<dim> *>(&x_source_fe) != nullptr)
250  {
251  // nothing to do here, the FE_Nothing has no degrees of freedom anyway
252  }
253  else
254  AssertThrow(
255  false,
256  (typename FiniteElement<dim,
257  spacedim>::ExcInterpolationNotImplemented()));
258 }
259 
260 
261 
262 template <int spacedim>
263 FE_TraceQ<1, spacedim>::FE_TraceQ(const unsigned int degree)
264  : FE_FaceQ<1, spacedim>(degree)
265 {}
266 
267 
268 
269 template <int spacedim>
270 std::string
272 {
273  // note that the FETools::get_fe_by_name function depends on the
274  // particular format of the string this function returns, so they have to be
275  // kept in synch
276  std::ostringstream namebuf;
277  namebuf << "FE_TraceQ<" << Utilities::dim_string(1, spacedim) << ">("
278  << this->degree << ")";
279 
280  return namebuf.str();
281 }
282 
283 
284 
285 // explicit instantiations
286 #include "fe_trace.inst"
287 
288 
FE_TraceQ::fe_q
FE_Q< dim, spacedim > fe_q
Definition: fe_trace.h:138
tensor_product_polynomials.h
FE_Nothing< dim >
FE_PolyFace< TensorProductPolynomials< dim - 1 >, dim, dim >::poly_space
TensorProductPolynomials< dim - 1 > poly_space
Definition: fe_poly_face.h:251
FE_TraceQ::convert_generalized_support_point_values_to_dof_values
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_trace.cc:131
FE_PolyFace
Definition: fe_poly_face.h:61
FE_TraceQ::hp_constraints_are_implemented
virtual bool hp_constraints_are_implemented() const override
Definition: fe_trace.cc:169
FE_TraceQ::compare_for_domination
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_trace.cc:177
Polynomials
Definition: polynomial.h:41
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
FiniteElementData
Definition: fe_base.h:147
FiniteElement::interface_constraints
FullMatrix< double > interface_constraints
Definition: fe.h:2440
FullMatrix::m
size_type m() const
FE_TraceQ::get_dpo_vector
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
Definition: fe_trace.cc:151
LAPACKSupport::U
static const char U
Definition: lapack_support.h:167
GeometryInfo
Definition: geometry_info.h:1224
AssertIndexRange
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
FullMatrix::n
size_type n() const
FE_FaceQ
Definition: fe_face.h:57
quadrature_lib.h
fe_nothing.h
Polynomials::generate_complete_Lagrange_basis
std::vector< Polynomial< double > > generate_complete_Lagrange_basis(const std::vector< Point< 1 >> &points)
Definition: polynomial.cc:834
FE_TraceQ::get_name
virtual std::string get_name() const override
Definition: fe_trace.cc:83
FiniteElementDomination::either_element_can_dominate
@ either_element_can_dominate
Definition: fe_base.h:100
Table
Definition: table.h:699
fe_poly_face.h
FiniteElementData::degree
const unsigned int degree
Definition: fe_base.h:298
bool
FiniteElementDomination::no_requirements
@ no_requirements
Definition: fe_base.h:104
FiniteElementDomination::neither_element_dominates
@ neither_element_dominates
Definition: fe_base.h:96
FiniteElementDomination::other_element_dominates
@ other_element_dominates
Definition: fe_base.h:92
FE_TraceQ::clone
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_trace.cc:74
FE_TraceQ::get_subface_interpolation_matrix
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const override
Definition: fe_trace.cc:229
FiniteElementDomination::Domination
Domination
Definition: fe_base.h:83
FE_TraceQ
Definition: fe_trace.h:46
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
FiniteElementData::dofs_per_cell
const unsigned int dofs_per_cell
Definition: fe_base.h:282
Utilities::dim_string
std::string dim_string(const int dim, const int spacedim)
Definition: utilities.cc:559
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
FE_TraceQ::FE_TraceQ
FE_TraceQ(unsigned int p)
Definition: fe_trace.cc:38
FE_TraceQ::get_constant_modes
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition: fe_trace.cc:119
TensorProductPolynomials
Definition: tensor_product_polynomials.h:74
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1579
fe_q.h
QGaussLobatto
Definition: quadrature_lib.h:76
FiniteElement
Definition: fe.h:648
fe_tools.h
FiniteElementData::dofs_per_face
const unsigned int dofs_per_face
Definition: fe_base.h:275
LocalIntegrators::L2::L2
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
Definition: l2.h:171
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
FiniteElement::unit_face_support_points
std::vector< Point< dim - 1 > > unit_face_support_points
Definition: fe.h:2459
StandardExceptions::ExcImpossibleInDim
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:191
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
FE_TraceQ::has_support_on_face
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition: fe_trace.cc:100
FiniteElementDomination::this_element_dominates
@ this_element_dominates
Definition: fe_base.h:88
config.h
FullMatrix< double >
memory.h
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
fe_trace.h
Vector< double >
TensorProductPolynomials::set_numbering
void set_numbering(const std::vector< unsigned int > &renumber)
Definition: tensor_product_polynomials.cc:117
AssertThrow
#define AssertThrow(cond, exc)
Definition: exceptions.h:1531
qprojector.h
FE_TraceQ::get_face_interpolation_matrix
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override
Definition: fe_trace.cc:216
FiniteElement::unit_support_points
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2452