Reference documentation for deal.II version Git 170b4c9308 2021-10-26 16:43:28 -0600
\(\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_poly_face.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 2021 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 #ifndef dealii_fe_poly_face_h
17 #define dealii_fe_poly_face_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
24 #include <deal.II/fe/fe.h>
25 
26 #include <memory>
27 
29 
32 
56 template <class PolynomialType,
57  int dim = PolynomialType::dimension + 1,
58  int spacedim = dim>
59 class FE_PolyFace : public FiniteElement<dim, spacedim>
60 {
61 public:
65  FE_PolyFace(const PolynomialType & poly_space,
66  const FiniteElementData<dim> &fe_data,
67  const std::vector<bool> & restriction_is_additive_flags);
68 
73  unsigned int
74  get_degree() const;
75 
76  // for documentation, see the FiniteElement base class
77  virtual UpdateFlags
78  requires_update_flags(const UpdateFlags update_flags) const override;
79 
80 protected:
81  /*
82  * NOTE: The following functions have their definitions inlined into the class
83  * declaration because we otherwise run into a compiler error with MS Visual
84  * Studio.
85  */
86 
87 
88  virtual std::unique_ptr<
91  const UpdateFlags update_flags,
92  const Mapping<dim, spacedim> &mapping,
93  const Quadrature<dim> & quadrature,
95  spacedim>
96  &output_data) const override
97  {
98  (void)update_flags;
99  (void)mapping;
100  (void)quadrature;
101  (void)output_data;
102  return std::make_unique<InternalData>();
103  }
104 
106 
107  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
109  const UpdateFlags update_flags,
110  const Mapping<dim, spacedim> & mapping,
111  const hp::QCollection<dim - 1> &quadrature,
113  spacedim>
114  &output_data) const override
115  {
116  (void)mapping;
117  (void)output_data;
118  AssertDimension(quadrature.size(), 1);
119 
120  // generate a new data object and
121  // initialize some fields
122  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
123  data_ptr = std::make_unique<InternalData>();
124  auto &data = dynamic_cast<InternalData &>(*data_ptr);
125  data.update_each = requires_update_flags(update_flags);
126 
127  const unsigned int n_q_points = quadrature[0].size();
128 
129  // some scratch arrays
130  std::vector<double> values(0);
131  std::vector<Tensor<1, dim - 1>> grads(0);
132  std::vector<Tensor<2, dim - 1>> grad_grads(0);
133  std::vector<Tensor<3, dim - 1>> empty_vector_of_3rd_order_tensors;
134  std::vector<Tensor<4, dim - 1>> empty_vector_of_4th_order_tensors;
135 
136  // initialize fields only if really
137  // necessary. otherwise, don't
138  // allocate memory
139  if (data.update_each & update_values)
140  {
141  values.resize(poly_space.n());
142  data.shape_values.resize(poly_space.n(),
143  std::vector<double>(n_q_points));
144  for (unsigned int i = 0; i < n_q_points; ++i)
145  {
146  poly_space.evaluate(quadrature[0].point(i),
147  values,
148  grads,
149  grad_grads,
150  empty_vector_of_3rd_order_tensors,
151  empty_vector_of_4th_order_tensors);
152 
153  for (unsigned int k = 0; k < poly_space.n(); ++k)
154  data.shape_values[k][i] = values[k];
155  }
156  }
157  // No derivatives of this element
158  // are implemented.
159  if (data.update_each & update_gradients ||
160  data.update_each & update_hessians)
161  {
162  Assert(false, ExcNotImplemented());
163  }
164 
165  return data_ptr;
166  }
167 
168  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
170  const UpdateFlags update_flags,
171  const Mapping<dim, spacedim> &mapping,
172  const Quadrature<dim - 1> & quadrature,
174  spacedim>
175  &output_data) const override
176  {
177  return get_face_data(
178  update_flags,
179  mapping,
181  ReferenceCells::get_hypercube<dim - 1>(), quadrature)),
182  output_data);
183  }
184 
185  virtual void
187  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
188  const CellSimilarity::Similarity cell_similarity,
189  const Quadrature<dim> & quadrature,
190  const Mapping<dim, spacedim> & mapping,
191  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
192  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
193  spacedim>
194  & mapping_data,
195  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
197  spacedim>
198  &output_data) const override;
199 
201 
202  virtual void
204  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
205  const unsigned int face_no,
206  const hp::QCollection<dim - 1> & quadrature,
207  const Mapping<dim, spacedim> & mapping,
208  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
209  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
210  spacedim>
211  & mapping_data,
212  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
214  spacedim>
215  &output_data) const override;
216 
217  virtual void
219  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
220  const unsigned int face_no,
221  const unsigned int sub_no,
222  const Quadrature<dim - 1> & quadrature,
223  const Mapping<dim, spacedim> & mapping,
224  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
225  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
226  spacedim>
227  & mapping_data,
228  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
230  spacedim>
231  &output_data) const override;
232 
239  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
240  {
241  public:
255  std::vector<std::vector<double>> shape_values;
256  };
257 
262  PolynomialType poly_space;
263 };
264 
268 
269 #endif
Shape function values.
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_poly_face.h:90
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
unsigned int size() const
Definition: collection.h:110
std::vector< std::vector< double > > shape_values
Definition: fe_poly_face.h:255
#define Assert(cond, exc)
Definition: exceptions.h:1461
UpdateFlags
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_poly_face.h:169
Abstract base class for mapping classes.
Definition: mapping.h:303
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Second derivatives of shape functions.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:185
PolynomialType poly_space
Definition: fe_poly_face.h:262
Definition: tensor.h:506
unsigned int get_degree() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
Shape function gradients.
static ::ExceptionBase & ExcNotImplemented()
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
FE_PolyFace(const PolynomialType &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags)
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_poly_face.h:108
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2572