00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__fe_raviart_thomas_h
00013 #define __deal2__fe_raviart_thomas_h
00014
00015 #include <deal.II/base/config.h>
00016 #include <deal.II/base/table.h>
00017 #include <deal.II/base/polynomials_raviart_thomas.h>
00018 #include <deal.II/base/polynomial.h>
00019 #include <deal.II/base/tensor_product_polynomials.h>
00020 #include <deal.II/base/geometry_info.h>
00021 #include <deal.II/fe/fe.h>
00022 #include <deal.II/fe/fe_poly_tensor.h>
00023
00024 #include <vector>
00025
00026 DEAL_II_NAMESPACE_OPEN
00027
00028 template <int dim, int spacedim> class MappingQ;
00029
00030
00033
00100 template <int dim>
00101 class FE_RaviartThomas
00102 :
00103 public FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim>
00104 {
00105 public:
00110 FE_RaviartThomas (const unsigned int p);
00111
00121 virtual std::string get_name () const;
00122
00123
00133 virtual bool has_support_on_face (const unsigned int shape_index,
00134 const unsigned int face_index) const;
00135
00136 virtual void interpolate(std::vector<double>& local_dofs,
00137 const std::vector<double>& values) const;
00138 virtual void interpolate(std::vector<double>& local_dofs,
00139 const std::vector<Vector<double> >& values,
00140 unsigned int offset = 0) const;
00141 virtual void interpolate(
00142 std::vector<double>& local_dofs,
00143 const VectorSlice<const std::vector<std::vector<double> > >& values) const;
00144 virtual std::size_t memory_consumption () const;
00145 virtual FiniteElement<dim> * clone() const;
00146
00147 private:
00158 static std::vector<unsigned int>
00159 get_dpo_vector (const unsigned int degree);
00160
00171 void initialize_support_points (const unsigned int rt_degree);
00172
00184 void initialize_restriction ();
00185
00194 class InternalData : public FiniteElement<dim>::InternalDataBase
00195 {
00196 public:
00220 std::vector<std::vector<Tensor<1,dim> > > shape_values;
00221
00240 std::vector<std::vector<Tensor<2,dim> > > shape_gradients;
00241 };
00242
00256 Table<2, double> boundary_weights;
00268 Table<3, double> interior_weights;
00269
00274 template <int dim1> friend class FE_RaviartThomas;
00275 };
00276
00277
00278
00321 template <int dim>
00322 class FE_RaviartThomasNodal
00323 :
00324 public FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim>
00325 {
00326 public:
00331 FE_RaviartThomasNodal (const unsigned int p);
00332
00342 virtual std::string get_name () const;
00343
00344 virtual FiniteElement<dim>* clone () const;
00345
00346 virtual void interpolate(std::vector<double>& local_dofs,
00347 const std::vector<double>& values) const;
00348 virtual void interpolate(std::vector<double>& local_dofs,
00349 const std::vector<Vector<double> >& values,
00350 unsigned int offset = 0) const;
00351 virtual void interpolate(
00352 std::vector<double>& local_dofs,
00353 const VectorSlice<const std::vector<std::vector<double> > >& values) const;
00354
00355
00356 virtual void get_face_interpolation_matrix (const FiniteElement<dim> &source,
00357 FullMatrix<double> &matrix) const;
00358
00359 virtual void get_subface_interpolation_matrix (const FiniteElement<dim> &source,
00360 const unsigned int subface,
00361 FullMatrix<double> &matrix) const;
00362 virtual bool hp_constraints_are_implemented () const;
00363
00364 virtual std::vector<std::pair<unsigned int, unsigned int> >
00365 hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
00366
00367 virtual std::vector<std::pair<unsigned int, unsigned int> >
00368 hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
00369
00370 virtual std::vector<std::pair<unsigned int, unsigned int> >
00371 hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
00372
00373 virtual FiniteElementDomination::Domination
00374 compare_for_face_domination (const FiniteElement<dim> &fe_other) const;
00375
00376 private:
00387 static std::vector<unsigned int>
00388 get_dpo_vector (const unsigned int degree);
00389
00397 static std::vector<bool>
00398 get_ria_vector (const unsigned int degree);
00408 virtual bool has_support_on_face (const unsigned int shape_index,
00409 const unsigned int face_index) const;
00417 void initialize_support_points (const unsigned int rt_degree);
00418 };
00419
00420
00423
00424
00425 #ifndef DOXYGEN
00426
00427 template <>
00428 void
00429 FE_RaviartThomas<1>::initialize_restriction();
00430
00431 #endif // DOXYGEN
00432
00433 DEAL_II_NAMESPACE_CLOSE
00434
00435 #endif
00436