Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fe_nedelec.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_nedelec.h 31791 2013-11-25 10:36:38Z felix.gruber @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_nedelec_h
18 #define __deal2__fe_nedelec_h
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/table.h>
22 #include <deal.II/base/tensor.h>
23 #include <deal.II/base/tensor_base.h>
24 #include <deal.II/base/polynomials_nedelec.h>
25 #include <deal.II/base/polynomial.h>
26 #include <deal.II/base/tensor_product_polynomials.h>
27 #include <deal.II/base/geometry_info.h>
28 #include <deal.II/fe/fe.h>
29 #include <deal.II/fe/fe_poly_tensor.h>
30 #include <vector>
31 
33 
34 template <int dim, int spacedim> class MappingQ;
35 
36 
39 
158 template <int dim>
159 class FE_Nedelec : public FE_PolyTensor<PolynomialsNedelec<dim>, dim>
160 {
161 public:
166  FE_Nedelec (const unsigned int p);
167 
177  virtual std::string get_name () const;
178 
179 
184  virtual bool has_support_on_face (const unsigned int shape_index,
185  const unsigned int face_index) const;
186 
197  virtual bool hp_constraints_are_implemented () const;
198 
204  compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
205 
223  virtual std::vector<std::pair<unsigned int, unsigned int> >
224  hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
225 
230  virtual std::vector<std::pair<unsigned int, unsigned int> >
231  hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
232 
237  virtual std::vector<std::pair<unsigned int, unsigned int> >
238  hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
239 
253  virtual void
255  FullMatrix<double> &matrix) const;
256 
270  virtual void
272  const unsigned int subface,
273  FullMatrix<double> &matrix) const;
274 
275  virtual void interpolate (std::vector<double> &local_dofs,
276  const std::vector<double> &values) const;
277 
278  virtual void interpolate (std::vector<double> &local_dofs,
279  const std::vector<Vector<double> > &values,
280  unsigned int offset = 0) const;
281  virtual void interpolate (std::vector<double> &local_dofs,
282  const VectorSlice<const std::vector<std::vector<double> > > &values)
283  const;
284  virtual std::size_t memory_consumption () const;
285  virtual FiniteElement<dim> *clone() const;
286 
287 private:
305  static std::vector<unsigned int>
306  get_dpo_vector (const unsigned int degree, bool dg=false);
307 
318  void initialize_support_points (const unsigned int degree);
319 
331  void initialize_restriction ();
332 
342  {
343  public:
367  std::vector<std::vector<Tensor<1, dim> > > shape_values;
368 
387  std::vector<std::vector<Tensor<2, dim> > > shape_gradients;
388  };
389 
402 
407  template <int dim1> friend class FE_Nedelec;
408 };
409 
410 /* -------------- declaration of explicit specializations ------------- */
411 
412 #ifndef DOXYGEN
413 
414 template <>
415 void
417 
418 #endif // DOXYGEN
419 
422 DEAL_II_NAMESPACE_CLOSE
423 
424 #endif
virtual FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim > &fe_other) const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim > &source, const unsigned int subface, FullMatrix< double > &matrix) const
friend class FE_Nedelec
Definition: fe_nedelec.h:407
virtual bool hp_constraints_are_implemented() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim > &fe_other) const
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree, bool dg=false)
virtual void get_face_interpolation_matrix(const FiniteElement< dim > &source, FullMatrix< double > &matrix) const
virtual std::string get_name() const
std::vector< std::vector< Tensor< 1, dim > > > shape_values
Definition: fe_nedelec.h:367
Table< 2, double > boundary_weights
Definition: fe_nedelec.h:401
const unsigned int degree
Definition: fe_base.h:287
virtual FiniteElement< dim > * clone() const
virtual std::size_t memory_consumption() const
std::vector< std::vector< Tensor< 2, dim > > > shape_gradients
Definition: fe_nedelec.h:387
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const
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 > &fe_other) const
void initialize_restriction()
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim > &fe_other) const
void initialize_support_points(const unsigned int degree)