Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_nedelec_sz.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2015 - 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_nedelec_sz_h
17 #define dealii_fe_nedelec_sz_h
18 
19 #include <deal.II/base/derivative_form.h>
20 #include <deal.II/base/polynomials_integrated_legendre_sz.h>
21 #include <deal.II/base/qprojector.h>
22 
23 #include <deal.II/fe/fe.h>
24 #include <deal.II/fe/fe_values.h>
25 #include <deal.II/fe/mapping.h>
26 
27 DEAL_II_NAMESPACE_OPEN
28 
31 
72 template <int dim, int spacedim = dim>
73 class FE_NedelecSZ : public FiniteElement<dim, dim>
74 {
75 public:
76  static_assert(dim == spacedim,
77  "FE_NedelecSZ is only implemented for dim==spacedim!");
78 
82  FE_NedelecSZ(const unsigned int degree);
83 
84  virtual UpdateFlags
85  requires_update_flags(const UpdateFlags update_flags) const override;
86 
87  virtual std::string
88  get_name() const override;
89 
90  virtual std::unique_ptr<FiniteElement<dim, dim>>
91  clone() const override;
92 
97  virtual double
98  shape_value(const unsigned int i, const Point<dim> &p) const override;
99 
103  virtual double
104  shape_value_component(const unsigned int i,
105  const Point<dim> & p,
106  const unsigned int component) const override;
107 
112  virtual Tensor<1, dim>
113  shape_grad(const unsigned int i, const Point<dim> &p) const override;
114 
118  virtual Tensor<1, dim>
119  shape_grad_component(const unsigned int i,
120  const Point<dim> & p,
121  const unsigned int component) const override;
122 
127  virtual Tensor<2, dim>
128  shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
129 
133  virtual Tensor<2, dim>
134  shape_grad_grad_component(const unsigned int i,
135  const Point<dim> & p,
136  const unsigned int component) const override;
137 
144  update_once(const UpdateFlags flags) const;
145 
152  update_each(const UpdateFlags flags) const;
153 
154 protected:
160 
161  virtual std::unique_ptr<
162  typename ::FiniteElement<dim, spacedim>::InternalDataBase>
163  get_data(
164  const UpdateFlags update_flags,
165  const Mapping<dim, spacedim> &mapping,
166  const Quadrature<dim> & quadrature,
168  spacedim>
169  &output_data) const override;
170 
176  virtual void
178  const typename Triangulation<dim, dim>::cell_iterator &cell,
179  const CellSimilarity::Similarity cell_similarity,
180  const Quadrature<dim> & quadrature,
181  const Mapping<dim, dim> & mapping,
182  const typename Mapping<dim, dim>::InternalDataBase & mapping_internal,
183  const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
184  & mapping_data,
185  const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
187  &data) const override;
188 
194  virtual void
196  const typename Triangulation<dim, dim>::cell_iterator &cell,
197  const unsigned int face_no,
198  const Quadrature<dim - 1> & quadrature,
199  const Mapping<dim, dim> & mapping,
200  const typename Mapping<dim, dim>::InternalDataBase & mapping_internal,
201  const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
202  & mapping_data,
203  const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
205  &data) const override;
206 
210  virtual void
212  const typename Triangulation<dim, dim>::cell_iterator &cell,
213  const unsigned int face_no,
214  const unsigned int sub_no,
215  const Quadrature<dim - 1> & quadrature,
216  const Mapping<dim, dim> & mapping,
217  const typename Mapping<dim, dim>::InternalDataBase & mapping_internal,
218  const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
219  & mapping_data,
220  const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
222  &data) const override;
223 
243  class InternalData : public FiniteElement<dim, dim>::InternalDataBase
244  {
245  public:
252  mutable std::vector<std::vector<Tensor<1, dim>>> shape_values;
253 
261  mutable std::vector<std::vector<DerivativeForm<1, dim, dim>>> shape_grads;
262 
278  std::vector<std::vector<std::vector<double>>> sigma_imj_values;
279 
295  std::vector<std::vector<std::vector<double>>> sigma_imj_grads;
296 
308  std::vector<std::vector<double>> edge_sigma_values;
309 
322  std::vector<std::vector<double>> edge_sigma_grads;
323 
339  std::vector<std::vector<double>> edge_lambda_values;
340 
350  std::vector<std::vector<double>> edge_lambda_grads_2d;
351 
361  std::vector<std::vector<std::vector<double>>> edge_lambda_grads_3d;
362 
373  std::vector<std::vector<std::vector<double>>> edge_lambda_gradgrads_3d;
374 
391  std::vector<std::vector<double>> face_lambda_values;
392 
401  std::vector<std::vector<double>> face_lambda_grads;
402  };
403 
404 private:
413  static std::vector<unsigned int>
414  get_dpo_vector(const unsigned int degree);
415 
419  std::vector<Polynomials::Polynomial<double>> IntegratedLegendrePolynomials;
420 
425  void
426  create_polynomials(const unsigned int degree);
427 
431  unsigned int
432  compute_num_dofs(const unsigned int degree) const;
433 
438  void
440  const Quadrature<dim> &quadrature,
441  const InternalData & fedata) const;
442 
447  void
449  const Quadrature<dim> &quadrature,
450  const InternalData & fedata) const;
451 };
452 
453 
454 
457 DEAL_II_NAMESPACE_CLOSE
458 
459 #endif
virtual void fill_fe_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
std::vector< std::vector< double > > edge_lambda_grads_2d
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
UpdateFlags update_once(const UpdateFlags flags) const
std::vector< std::vector< double > > face_lambda_values
void fill_edge_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
MappingType
Definition: mapping.h:61
std::vector< std::vector< double > > edge_sigma_values
const unsigned int degree
Definition: fe_base.h:298
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
void fill_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
MappingType mapping_type
virtual std::string get_name() const override
virtual void fill_fe_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
std::vector< std::vector< DerivativeForm< 1, dim, dim > > > shape_grads
std::vector< std::vector< std::vector< double > > > edge_lambda_gradgrads_3d
UpdateFlags
virtual void fill_fe_subface_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
Abstract base class for mapping classes.
Definition: dof_tools.h:57
std::vector< Polynomials::Polynomial< double > > IntegratedLegendrePolynomials
std::vector< std::vector< Tensor< 1, dim > > > shape_values
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
void create_polynomials(const unsigned int degree)
std::vector< std::vector< double > > edge_lambda_values
UpdateFlags update_each(const UpdateFlags flags) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< std::vector< std::vector< double > > > edge_lambda_grads_3d
FE_NedelecSZ(const unsigned int degree)
std::vector< std::vector< double > > edge_sigma_grads
std::vector< std::vector< std::vector< double > > > sigma_imj_values
std::vector< std::vector< double > > face_lambda_grads
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
std::vector< std::vector< std::vector< double > > > sigma_imj_grads
unsigned int compute_num_dofs(const unsigned int degree) const