Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_dgq.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 2018 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_dgq_h
17 #define dealii_fe_dgq_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/base/tensor_product_polynomials.h>
22 #include <deal.II/base/thread_management.h>
23 
24 #include <deal.II/fe/fe_poly.h>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 template <int dim, int spacedim>
29 class MappingQ;
30 template <int dim>
31 class Quadrature;
32 
35 
108 template <int dim, int spacedim = dim>
109 class FE_DGQ : public FE_Poly<TensorProductPolynomials<dim>, dim, spacedim>
110 {
111 public:
118  FE_DGQ(const unsigned int p);
119 
125  virtual std::string
126  get_name() const override;
127 
137  virtual void
139  FullMatrix<double> &matrix) const override;
140 
152  virtual void
154  FullMatrix<double> &matrix) const override;
155 
167  virtual void
169  const unsigned int subface,
170  FullMatrix<double> &matrix) const override;
171 
189  virtual const FullMatrix<double> &
191  const unsigned int child,
192  const RefinementCase<dim> &refinement_case =
194 
216  virtual const FullMatrix<double> &
218  const unsigned int child,
219  const RefinementCase<dim> &refinement_case =
221 
245  virtual std::vector<std::pair<unsigned int, unsigned int>>
247  const FiniteElement<dim, spacedim> &fe_other) const override;
248 
256  virtual std::vector<std::pair<unsigned int, unsigned int>>
258  const FiniteElement<dim, spacedim> &fe_other) const override;
259 
267  virtual std::vector<std::pair<unsigned int, unsigned int>>
269  const FiniteElement<dim, spacedim> &fe_other) const override;
270 
279  virtual bool
280  hp_constraints_are_implemented() const override;
281 
287  const unsigned int codim = 0) const override final;
288 
297  virtual bool
298  has_support_on_face(const unsigned int shape_index,
299  const unsigned int face_index) const override;
300 
305  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
306  get_constant_modes() const override;
307 
315  virtual void
317  const std::vector<Vector<double>> &support_point_values,
318  std::vector<double> & nodal_values) const override;
319 
328  virtual std::size_t
329  memory_consumption() const override;
330 
331  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
332  clone() const override;
333 
334 protected:
343  FE_DGQ(const std::vector<Polynomials::Polynomial<double>> &polynomials);
344 
345 private:
352  static std::vector<unsigned int>
353  get_dpo_vector(const unsigned int degree);
354 
371  void
372  rotate_indices(std::vector<unsigned int> &indices,
373  const char direction) const;
374 
375  /*
376  * Mutex for protecting initialization of restriction and embedding matrix.
377  */
378  mutable Threads::Mutex mutex;
379 
383  template <int dim1, int spacedim1>
384  friend class FE_DGQ;
385 
389  template <int dim1, int spacedim1>
390  friend class MappingQ;
391 };
392 
393 
394 
411 template <int dim, int spacedim = dim>
412 class FE_DGQArbitraryNodes : public FE_DGQ<dim, spacedim>
413 {
414 public:
421  FE_DGQArbitraryNodes(const Quadrature<1> &points);
422 
428  virtual std::string
429  get_name() const override;
430 
438  virtual void
440  const std::vector<Vector<double>> &support_point_values,
441  std::vector<double> & nodal_values) const override;
442  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
443  clone() const override;
444 };
445 
446 
447 
458 template <int dim, int spacedim = dim>
459 class FE_DGQLegendre : public FE_DGQ<dim, spacedim>
460 {
461 public:
466  FE_DGQLegendre(const unsigned int degree);
467 
473  virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
474  get_constant_modes() const override;
475 
482  virtual std::string
483  get_name() const override;
484 
485  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
486  clone() const override;
487 };
488 
489 
490 
509 template <int dim, int spacedim = dim>
510 class FE_DGQHermite : public FE_DGQ<dim, spacedim>
511 {
512 public:
517  FE_DGQHermite(const unsigned int degree);
518 
525  virtual std::string
526  get_name() const override;
527 
528  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
529  clone() const override;
530 };
531 
532 
535 DEAL_II_NAMESPACE_CLOSE
536 
537 #endif
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition: fe_dgq.cc:995
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_dgq.cc:596
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
Definition: fe_dgq.cc:399
virtual bool hp_constraints_are_implemented() const override
Definition: fe_dgq.cc:559
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_dgq.cc:1048
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_dgq.cc:568
FE_DGQLegendre(const unsigned int degree)
Definition: fe_dgq.cc:986
Definition: fe_dgq.h:109
const unsigned int degree
Definition: fe_base.h:298
virtual std::string get_name() const override
Definition: fe_dgq.cc:849
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override
Definition: fe_dgq.cc:346
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_dgq.cc:610
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_dgq.cc:1019
virtual std::size_t memory_consumption() const override
Definition: fe_dgq.cc:824
virtual std::string get_name() const override
Definition: fe_dgq.cc:1009
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition: fe_dgq.cc:812
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
Definition: fe_dgq.cc:177
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_dgq.cc:143
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const override
Definition: fe_dgq.cc:372
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_dgq.cc:164
void rotate_indices(std::vector< unsigned int > &indices, const char direction) const
Definition: fe_dgq.cc:190
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_dgq.cc:967
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition: fe_dgq.cc:719
FE_DGQArbitraryNodes(const Quadrature< 1 > &points)
Definition: fe_dgq.cc:835
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override
Definition: fe_dgq.cc:267
FE_DGQHermite(const unsigned int degree)
Definition: fe_dgq.cc:1029
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_dgq.cc:946
virtual std::string get_name() const override
Definition: fe_dgq.cc:1038
virtual std::string get_name() const override
Definition: fe_dgq.cc:127
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const override
Definition: fe_dgq.cc:479
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_dgq.cc:582
friend class FE_DGQ
Definition: fe_dgq.h:384