Reference documentation for deal.II version Git 2449c74b32 2021-01-19 22:49:31 +0100
\(\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_dgp_nonparametric.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2002 - 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_dgp_nonparametric_h
17 #define dealii_fe_dgp_nonparametric_h
18 
19 #include <deal.II/base/config.h>
20 
23 
24 #include <deal.II/fe/fe.h>
25 #include <deal.II/fe/mapping.h>
26 
28 
31 
268 template <int dim, int spacedim = dim>
269 class FE_DGPNonparametric : public FiniteElement<dim, spacedim>
270 {
271 public:
275  FE_DGPNonparametric(const unsigned int k);
276 
282  virtual std::string
283  get_name() const override;
284 
285  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
286  clone() const override;
287 
288  // for documentation, see the FiniteElement base class
289  virtual UpdateFlags
290  requires_update_flags(const UpdateFlags update_flags) const override;
291 
302  virtual double
303  shape_value(const unsigned int i, const Point<dim> &p) const override;
304 
315  virtual double
316  shape_value_component(const unsigned int i,
317  const Point<dim> & p,
318  const unsigned int component) const override;
319 
330  virtual Tensor<1, dim>
331  shape_grad(const unsigned int i, const Point<dim> &p) const override;
332 
343  virtual Tensor<1, dim>
344  shape_grad_component(const unsigned int i,
345  const Point<dim> & p,
346  const unsigned int component) const override;
347 
358  virtual Tensor<2, dim>
359  shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
360 
371  virtual Tensor<2, dim>
372  shape_grad_grad_component(const unsigned int i,
373  const Point<dim> & p,
374  const unsigned int component) const override;
375 
380  unsigned int
381  get_degree() const;
382 
394  virtual void
397  const unsigned int face_no = 0) const override;
398 
410  virtual void
412  const FiniteElement<dim, spacedim> &source,
413  const unsigned int subface,
414  FullMatrix<double> & matrix,
415  const unsigned int face_no = 0) const override;
416 
440  virtual std::vector<std::pair<unsigned int, unsigned int>>
442  const FiniteElement<dim, spacedim> &fe_other) const override;
443 
451  virtual std::vector<std::pair<unsigned int, unsigned int>>
453  const FiniteElement<dim, spacedim> &fe_other) const override;
454 
462  virtual std::vector<std::pair<unsigned int, unsigned int>>
464  const unsigned int face_no = 0) const override;
465 
474  virtual bool
475  hp_constraints_are_implemented() const override;
476 
482  const unsigned int codim = 0) const override final;
483 
492  virtual bool
493  has_support_on_face(const unsigned int shape_index,
494  const unsigned int face_index) const override;
495 
504  virtual std::size_t
505  memory_consumption() const override;
506 
507 protected:
512  virtual std::unique_ptr<
514  get_data(
515  const UpdateFlags update_flags,
516  const Mapping<dim, spacedim> &mapping,
517  const Quadrature<dim> & quadrature,
519  spacedim>
520  &output_data) const override;
521 
522  virtual void
524  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
525  const CellSimilarity::Similarity cell_similarity,
526  const Quadrature<dim> & quadrature,
527  const Mapping<dim, spacedim> & mapping,
528  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
529  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
530  spacedim>
531  & mapping_data,
532  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
534  spacedim>
535  &output_data) const override;
536 
538 
539  virtual void
541  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
542  const unsigned int face_no,
543  const hp::QCollection<dim - 1> & quadrature,
544  const Mapping<dim, spacedim> & mapping,
545  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
546  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
547  spacedim>
548  & mapping_data,
549  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
551  spacedim>
552  &output_data) const override;
553 
554  virtual void
556  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
557  const unsigned int face_no,
558  const unsigned int sub_no,
559  const Quadrature<dim - 1> & quadrature,
560  const Mapping<dim, spacedim> & mapping,
561  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
562  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
563  spacedim>
564  & mapping_data,
565  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
567  spacedim>
568  &output_data) const override;
569 
570 private:
577  static std::vector<unsigned int>
578  get_dpo_vector(const unsigned int degree);
579 
584 
585  // Allow access from other dimensions.
586  template <int, int>
587  friend class FE_DGPNonparametric;
588 };
589 
593 
594 #endif
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
friend class FE_DGPNonparametric
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Contents is actually a matrix.
virtual std::size_t memory_consumption() const override
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
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
const unsigned int degree
Definition: fe_base.h:435
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
const PolynomialSpace< dim > polynomial_space
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
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
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual bool hp_constraints_are_implemented() const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const override
unsigned int get_degree() const
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:303
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
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:380
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:379
virtual std::string get_name() const override
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
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
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override