Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_poly_face.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 2019 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/qprojector.h>
21 #include <deal.II/base/std_cxx14/memory.h>
22 
23 #include <deal.II/fe/fe.h>
24 
25 
26 DEAL_II_NAMESPACE_OPEN
27 
30 
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  return std_cxx14::make_unique<InternalData>();
99  }
100 
101  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
103  const UpdateFlags update_flags,
104  const Mapping<dim, spacedim> & /*mapping*/,
105  const Quadrature<dim - 1> &quadrature,
107  spacedim>
108  & /*output_data*/) const override
109  {
110  // generate a new data object and
111  // initialize some fields
112  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
113  data_ptr = std_cxx14::make_unique<InternalData>();
114  auto &data = dynamic_cast<InternalData &>(*data_ptr);
115  data.update_each = requires_update_flags(update_flags);
116 
117  const unsigned int n_q_points = quadrature.size();
118 
119  // some scratch arrays
120  std::vector<double> values(0);
121  std::vector<Tensor<1, dim - 1>> grads(0);
122  std::vector<Tensor<2, dim - 1>> grad_grads(0);
123  std::vector<Tensor<3, dim - 1>> empty_vector_of_3rd_order_tensors;
124  std::vector<Tensor<4, dim - 1>> empty_vector_of_4th_order_tensors;
125 
126  // initialize fields only if really
127  // necessary. otherwise, don't
128  // allocate memory
129  if (data.update_each & update_values)
130  {
131  values.resize(poly_space.n());
132  data.shape_values.resize(poly_space.n(),
133  std::vector<double>(n_q_points));
134  for (unsigned int i = 0; i < n_q_points; ++i)
135  {
136  poly_space.compute(quadrature.point(i),
137  values,
138  grads,
139  grad_grads,
140  empty_vector_of_3rd_order_tensors,
141  empty_vector_of_4th_order_tensors);
142 
143  for (unsigned int k = 0; k < poly_space.n(); ++k)
144  data.shape_values[k][i] = values[k];
145  }
146  }
147  // No derivatives of this element
148  // are implemented.
149  if (data.update_each & update_gradients ||
150  data.update_each & update_hessians)
151  {
152  Assert(false, ExcNotImplemented());
153  }
154 
155  return data_ptr;
156  }
157 
158  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
160  const UpdateFlags update_flags,
161  const Mapping<dim, spacedim> &mapping,
162  const Quadrature<dim - 1> & quadrature,
164  spacedim>
165  &output_data) const override
166  {
167  return get_face_data(update_flags,
168  mapping,
170  quadrature),
171  output_data);
172  }
173 
174  virtual void
175  fill_fe_values(
176  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
177  const CellSimilarity::Similarity cell_similarity,
178  const Quadrature<dim> & quadrature,
179  const Mapping<dim, spacedim> & mapping,
180  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
181  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
182  spacedim>
183  & mapping_data,
184  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
186  spacedim>
187  &output_data) const override;
188 
189  virtual void
190  fill_fe_face_values(
191  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
192  const unsigned int face_no,
193  const Quadrature<dim - 1> & quadrature,
194  const Mapping<dim, spacedim> & mapping,
195  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
196  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
197  spacedim>
198  & mapping_data,
199  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
201  spacedim>
202  &output_data) const override;
203 
204  virtual void
205  fill_fe_subface_values(
206  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
207  const unsigned int face_no,
208  const unsigned int sub_no,
209  const Quadrature<dim - 1> & quadrature,
210  const Mapping<dim, spacedim> & mapping,
211  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
212  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
213  spacedim>
214  & mapping_data,
215  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
217  spacedim>
218  &output_data) const override;
219 
226  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
227  {
228  public:
242  std::vector<std::vector<double>> shape_values;
243  };
244 
249  PolynomialType poly_space;
250 };
251 
254 DEAL_II_NAMESPACE_CLOSE
255 
256 #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, const Mapping< dim, spacedim > &, const Quadrature< dim > &, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &) const override
Definition: fe_poly_face.h:90
std::vector< std::vector< double > > shape_values
Definition: fe_poly_face.h:242
const Point< dim > & point(const unsigned int i) const
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &) const override
Definition: fe_poly_face.h:102
#define Assert(cond, exc)
Definition: exceptions.h:1407
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:159
Abstract base class for mapping classes.
Definition: dof_tools.h:57
Second derivatives of shape functions.
PolynomialType poly_space
Definition: fe_poly_face.h:249
unsigned int size() const
Definition: mpi.h:90
unsigned int get_degree() const
Shape function gradients.
static ::ExceptionBase & ExcNotImplemented()
FE_PolyFace(const PolynomialType &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags)
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2609