Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fe_dgp_nonparametric.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_dgp_nonparametric.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2002 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__fe_dgp_nonparametric_h
18 #define __deal2__fe_dgp_nonparametric_h
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/polynomial.h>
22 #include <deal.II/base/polynomial_space.h>
23 #include <deal.II/fe/fe.h>
24 
26 
27 template <int dim> class PolynomialSpace;
28 template <int dim, int spacedim> class MappingQ;
29 
30 
33 
57 template <int dim, int spacedim=dim>
58 class FE_DGPNonparametric : public FiniteElement<dim,spacedim>
59 {
60 public:
65  FE_DGPNonparametric (const unsigned int k);
66 
76  virtual std::string get_name () const;
77 
87  virtual double shape_value (const unsigned int i,
88  const Point<dim> &p) const;
89 
108  virtual double shape_value_component (const unsigned int i,
109  const Point<dim> &p,
110  const unsigned int component) const;
111 
121  virtual Tensor<1,dim> shape_grad (const unsigned int i,
122  const Point<dim> &p) const;
123 
142  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
143  const Point<dim> &p,
144  const unsigned int component) const;
145 
156  virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
157  const Point<dim> &p) const;
158 
177  virtual Tensor<2,dim> shape_grad_grad_component (const unsigned int i,
178  const Point<dim> &p,
179  const unsigned int component) const;
180 
187  unsigned int get_degree () const;
188 
210  virtual void
212  FullMatrix<double> &matrix) const;
213 
235  virtual void
237  const unsigned int subface,
238  FullMatrix<double> &matrix) const;
239 
265  virtual
266  std::vector<std::pair<unsigned int, unsigned int> >
268 
276  virtual
277  std::vector<std::pair<unsigned int, unsigned int> >
278  hp_line_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
279 
287  virtual
288  std::vector<std::pair<unsigned int, unsigned int> >
289  hp_quad_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;
290 
304  virtual bool hp_constraints_are_implemented () const;
305 
318  virtual
321 
339  virtual bool has_support_on_face (const unsigned int shape_index,
340  const unsigned int face_index) const;
341 
353  virtual std::size_t memory_consumption () const;
354 
355 
356 private:
367  struct Matrices
368  {
376 
386  static const unsigned int n_embedding_matrices;
387 
393 
400  static const unsigned int n_projection_matrices;
401  };
402 
403 protected:
404 
412  virtual FiniteElement<dim,spacedim> *clone() const;
413 
419  virtual
421  get_data (const UpdateFlags,
422  const Mapping<dim,spacedim> &mapping,
423  const Quadrature<dim> &quadrature) const ;
424 
430  virtual void
431  fill_fe_values (const Mapping<dim,spacedim> &mapping,
432  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
433  const Quadrature<dim> &quadrature,
434  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
435  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
437  CellSimilarity::Similarity &cell_similarity) const;
438 
444  virtual void
446  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
447  const unsigned int face_no,
448  const Quadrature<dim-1> &quadrature,
449  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
450  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
451  FEValuesData<dim,spacedim> &data) const ;
452 
458  virtual void
460  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
461  const unsigned int face_no,
462  const unsigned int sub_no,
463  const Quadrature<dim-1> &quadrature,
464  typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
465  typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
466  FEValuesData<dim,spacedim> &data) const ;
467 
468 private:
469 
480  static std::vector<unsigned int> get_dpo_vector (const unsigned int degree);
481 
501  virtual UpdateFlags update_once (const UpdateFlags flags) const;
502 
519  virtual UpdateFlags update_each (const UpdateFlags flags) const;
520 
524  const unsigned int degree;
525 
532 
541  class InternalData : public FiniteElement<dim,spacedim>::InternalDataBase
542  {
543  public:
544  // have some scratch arrays
545  std::vector<double> values;
546  std::vector<Tensor<1,dim> > grads;
547  std::vector<Tensor<2,dim> > grad_grads;
548  };
549 
553  template <int, int> friend class FE_DGPNonparametric;
554 
560  template <int, int> friend class MappingQ;
561 // friend class MappingQ<dim>;
562 };
563 
566 #ifndef DOXYGEN
567 
568 // declaration of explicit specializations of member variables, if the
569 // compiler allows us to do that (the standard says we must)
570 #ifndef DEAL_II_MEMBER_VAR_SPECIALIZATION_BUG
571 template <>
573 
574 template <>
576 
577 template <>
579 
580 template <>
582 
583 template <>
585 
586 template <>
588 
589 template <>
591 
592 template <>
594 
595 template <>
597 
598 template <>
600 
601 template <>
603 
604 template <>
606 #endif
607 
608 #endif // DOXYGEN
609 
610 DEAL_II_NAMESPACE_CLOSE
611 
612 #endif
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) const
friend class FE_DGPNonparametric
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
const PolynomialSpace< dim > polynomial_space
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual void fill_fe_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, CellSimilarity::Similarity &cell_similarity) const
virtual UpdateFlags update_once(const UpdateFlags flags) const
unsigned int get_degree() const
static const unsigned int n_embedding_matrices
static const double *const embedding[][GeometryInfo< dim >::max_children_per_cell]
static const double *const projection_matrices[][GeometryInfo< dim >::max_children_per_cell]
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
UpdateFlags
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual FiniteElement< dim, spacedim > * clone() const
virtual bool hp_constraints_are_implemented() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual void fill_fe_face_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual UpdateFlags update_each(const UpdateFlags flags) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
virtual std::string get_name() const
static const unsigned int n_projection_matrices
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
virtual void fill_fe_subface_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
virtual std::size_t memory_consumption() const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const unsigned int degree
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const