Reference documentation for deal.II version 8.1.0
fe_raviart_thomas.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_raviart_thomas.h 31893 2013-12-05 03:00:41Z heister @f$
3 //
4 // Copyright (C) 2003 - 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_raviart_thomas_h
18 #define __deal2__fe_raviart_thomas_h
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/table.h>
22 #include <deal.II/base/polynomials_raviart_thomas.h>
23 #include <deal.II/base/polynomial.h>
24 #include <deal.II/base/tensor_product_polynomials.h>
25 #include <deal.II/base/geometry_info.h>
26 #include <deal.II/fe/fe.h>
27 #include <deal.II/fe/fe_poly_tensor.h>
28 
29 #include <vector>
30 
32 
33 template <int dim, int spacedim> class MappingQ;
34 
35 
38 
105 template <int dim>
107  :
108  public FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim>
109 {
110 public:
115  FE_RaviartThomas (const unsigned int p);
116 
126  virtual std::string get_name () const;
127 
128 
138  virtual bool has_support_on_face (const unsigned int shape_index,
139  const unsigned int face_index) const;
140 
141  virtual void interpolate(std::vector<double> &local_dofs,
142  const std::vector<double> &values) const;
143  virtual void interpolate(std::vector<double> &local_dofs,
144  const std::vector<Vector<double> > &values,
145  unsigned int offset = 0) const;
146  virtual void interpolate(
147  std::vector<double> &local_dofs,
148  const VectorSlice<const std::vector<std::vector<double> > > &values) const;
149  virtual std::size_t memory_consumption () const;
150  virtual FiniteElement<dim> *clone() const;
151 
152 private:
163  static std::vector<unsigned int>
164  get_dpo_vector (const unsigned int degree);
165 
176  void initialize_support_points (const unsigned int rt_degree);
177 
189  void initialize_restriction ();
190 
199  class InternalData : public FiniteElement<dim>::InternalDataBase
200  {
201  public:
225  std::vector<std::vector<Tensor<1,dim> > > shape_values;
226 
245  std::vector<std::vector<Tensor<2,dim> > > shape_gradients;
246  };
247 
277 
282  template <int dim1> friend class FE_RaviartThomas;
283 };
284 
285 
286 
329 template <int dim>
331  :
332  public FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim>
333 {
334 public:
339  FE_RaviartThomasNodal (const unsigned int p);
340 
350  virtual std::string get_name () const;
351 
352  virtual FiniteElement<dim> *clone () const;
353 
354  virtual void interpolate(std::vector<double> &local_dofs,
355  const std::vector<double> &values) const;
356  virtual void interpolate(std::vector<double> &local_dofs,
357  const std::vector<Vector<double> > &values,
358  unsigned int offset = 0) const;
359  virtual void interpolate(
360  std::vector<double> &local_dofs,
361  const VectorSlice<const std::vector<std::vector<double> > > &values) const;
362 
363 
364  virtual void get_face_interpolation_matrix (const FiniteElement<dim> &source,
365  FullMatrix<double> &matrix) const;
366 
367  virtual void get_subface_interpolation_matrix (const FiniteElement<dim> &source,
368  const unsigned int subface,
369  FullMatrix<double> &matrix) const;
370  virtual bool hp_constraints_are_implemented () const;
371 
372  virtual std::vector<std::pair<unsigned int, unsigned int> >
373  hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
374 
375  virtual std::vector<std::pair<unsigned int, unsigned int> >
376  hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
377 
378  virtual std::vector<std::pair<unsigned int, unsigned int> >
379  hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
380 
382  compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
383 
384 private:
395  static std::vector<unsigned int>
396  get_dpo_vector (const unsigned int degree);
397 
405  static std::vector<bool>
406  get_ria_vector (const unsigned int degree);
416  virtual bool has_support_on_face (const unsigned int shape_index,
417  const unsigned int face_index) const;
428  void initialize_support_points (const unsigned int rt_degree);
429 };
430 
431 
434 /* -------------- declaration of explicit specializations ------------- */
435 
436 #ifndef DOXYGEN
437 
438 template <>
439 void
441 
442 #endif // DOXYGEN
443 
444 DEAL_II_NAMESPACE_CLOSE
445 
446 #endif
virtual std::size_t memory_consumption() const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
std::vector< std::vector< Tensor< 1, dim > > > shape_values
virtual FiniteElement< dim > * clone() const
FE_RaviartThomasNodal(const unsigned int p)
const unsigned int degree
Definition: fe_base.h:287
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
friend class FE_RaviartThomas
std::vector< std::vector< Tensor< 2, dim > > > shape_gradients
void initialize_support_points(const unsigned int rt_degree)
void initialize_support_points(const unsigned int rt_degree)
static std::vector< bool > get_ria_vector(const unsigned int degree)
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const
void initialize_restriction()
virtual FiniteElement< dim > * clone() const
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const
virtual bool hp_constraints_are_implemented() const
virtual std::string get_name() const
Table< 3, double > interior_weights
virtual std::string get_name() const
Table< 2, double > boundary_weights