Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
fe_tools.h
1 // ---------------------------------------------------------------------
2 // @f$Id: fe_tools.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2000 - 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_tools_H
18 #define __deal2__fe_tools_H
19 
20 
21 
22 #include <deal.II/base/config.h>
24 #include <deal.II/base/geometry_info.h>
25 #include <deal.II/base/tensor.h>
26 #include <deal.II/base/symmetric_tensor.h>
27 
28 #include <vector>
29 #include <string>
30 
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 template <typename number> class FullMatrix;
35 template <typename number> class Vector;
36 template <int dim> class Quadrature;
37 template <int dim, int spacedim> class FiniteElement;
38 template <int dim, int spacedim> class DoFHandler;
39 namespace hp
40 {
41  template <int dim, int spacedim> class DoFHandler;
42 }
43 template <int dim> class FiniteElementData;
44 class ConstraintMatrix;
45 
46 
47 
50 
51 
69 namespace FETools
70 {
81  template <int dim, int spacedim=dim>
83  {
84  public:
89  get (const unsigned int degree) const = 0;
90 
97  get (const Quadrature<1> &quad) const = 0;
101  virtual ~FEFactoryBase();
102  };
103 
115  template <class FE>
116  class FEFactory : public FEFactoryBase<FE::dimension,FE::dimension>
117  {
118  public:
123  get (const unsigned int degree) const;
124 
130  get (const Quadrature<1> &quad) const;
131  };
132 
148  template<int dim, int spacedim>
150  const FiniteElement<dim,spacedim> &fe,
151  std::vector<unsigned int> &renumbering,
152  std::vector<std::vector<unsigned int> > &start_indices);
153 
170  template<int dim, int spacedim>
172  const FiniteElement<dim,spacedim> &fe,
173  std::vector<types::global_dof_index> &renumbering,
174  std::vector<types::global_dof_index> &block_data,
175  bool return_start_indices = true);
176 
190  template <int dim, typename number, int spacedim>
191  void
193  const FiniteElement<dim,spacedim> &fe2,
194  FullMatrix<number> &interpolation_matrix);
195 
207  template <int dim, typename number, int spacedim>
208  void
210  const FiniteElement<dim,spacedim> &fe2,
211  FullMatrix<number> &interpolation_matrix);
212 
223  template <int dim, typename number, int spacedim>
224  void
226  const FiniteElement<dim,spacedim> &fe2,
227  FullMatrix<number> &difference_matrix);
228 
232  template <int dim, typename number, int spacedim>
234  const FiniteElement<dim,spacedim> &fe2,
235  FullMatrix<number> &matrix);
236 
264  template <int dim, int spacedim>
266  const FiniteElement<dim,spacedim> &fe);
267 
303  template <int dim, typename number, int spacedim>
305  std::vector<std::vector<FullMatrix<number> > > &matrices,
306  const bool isotropic_only = false);
307 
321  template <int dim, typename number, int spacedim>
322  void
325  const unsigned int face_coarse,
326  const unsigned int face_fine);
327 
355  template <int dim, typename number, int spacedim>
357  const FiniteElement<dim,spacedim> &fe,
358  std::vector<std::vector<FullMatrix<number> > > &matrices,
359  const bool isotropic_only = false);
360 
446  template <int dim, int spacedim>
447  void
449  const Quadrature<dim> &lhs_quadrature,
450  const Quadrature<dim> &rhs_quadrature,
451  FullMatrix<double> &X);
452 
460  template <int dim, int spacedim>
461  void
463  const Quadrature<dim> &quadrature,
464  FullMatrix<double> &I_q);
465 
479  template <int dim>
480  void
482  const FullMatrix<double> &projection_matrix,
483  const std::vector< Tensor<1, dim > > &vector_of_tensors_at_qp,
484  std::vector< Tensor<1, dim > > &vector_of_tensors_at_nodes);
485 
486 
487 
491  template <int dim>
492  void
494  const FullMatrix<double> &projection_matrix,
495  const std::vector< SymmetricTensor<2, dim > > &vector_of_tensors_at_qp,
496  std::vector< SymmetricTensor<2, dim > > &vector_of_tensors_at_nodes);
497 
498 
499 
500 
510  template <int dim, int spacedim>
511  void
513  const Quadrature<dim-1> &lhs_quadrature,
514  const Quadrature<dim-1> &rhs_quadrature,
515  const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
516  const unsigned int face,
517  FullMatrix<double> &X);
518 
519 
520 
522 
554  template <int dim, int spacedim,
555  template <int,int> class DH1,
556  template <int,int> class DH2,
557  class InVector, class OutVector>
558  void
559  interpolate (const DH1<dim,spacedim> &dof1,
560  const InVector &u1,
561  const DH2<dim,spacedim> &dof2,
562  OutVector &u2);
563 
580  template <int dim, int spacedim,
581  template <int, int> class DH1,
582  template <int, int> class DH2,
583  class InVector, class OutVector>
584  void interpolate (const DH1<dim,spacedim> &dof1,
585  const InVector &u1,
586  const DH2<dim,spacedim> &dof2,
587  const ConstraintMatrix &constraints,
588  OutVector &u2);
589 
603  template <int dim, class InVector, class OutVector, int spacedim>
604  void back_interpolate (const DoFHandler<dim,spacedim> &dof1,
605  const InVector &u1,
606  const FiniteElement<dim,spacedim> &fe2,
607  OutVector &u1_interpolated);
608 
613  template <int dim,
614  template <int> class DH,
615  class InVector, class OutVector, int spacedim>
616  void back_interpolate (const DH<dim> &dof1,
617  const InVector &u1,
618  const FiniteElement<dim,spacedim> &fe2,
619  OutVector &u1_interpolated);
620 
634  template <int dim, class InVector, class OutVector, int spacedim>
635  void back_interpolate (const DoFHandler<dim,spacedim> &dof1,
636  const ConstraintMatrix &constraints1,
637  const InVector &u1,
638  const DoFHandler<dim,spacedim> &dof2,
639  const ConstraintMatrix &constraints2,
640  OutVector &u1_interpolated);
641 
651  template <int dim, class InVector, class OutVector, int spacedim>
653  const InVector &z1,
654  const FiniteElement<dim,spacedim> &fe2,
655  OutVector &z1_difference);
656 
669  template <int dim, class InVector, class OutVector, int spacedim>
671  const ConstraintMatrix &constraints1,
672  const InVector &z1,
673  const DoFHandler<dim,spacedim> &dof2,
674  const ConstraintMatrix &constraints2,
675  OutVector &z1_difference);
676 
677 
678 
687  template <int dim, class InVector, class OutVector, int spacedim>
688  void project_dg (const DoFHandler<dim,spacedim> &dof1,
689  const InVector &u1,
690  const DoFHandler<dim,spacedim> &dof2,
691  OutVector &u2);
692 
712  template <int dim, class InVector, class OutVector, int spacedim>
713  void extrapolate (const DoFHandler<dim,spacedim> &dof1,
714  const InVector &z1,
715  const DoFHandler<dim,spacedim> &dof2,
716  OutVector &z2);
717 
728  template <int dim, class InVector, class OutVector, int spacedim>
729  void extrapolate (const DoFHandler<dim,spacedim> &dof1,
730  const InVector &z1,
731  const DoFHandler<dim,spacedim> &dof2,
732  const ConstraintMatrix &constraints,
733  OutVector &z2);
735 
759  template <int dim>
760  void
762  std::vector<unsigned int> &h2l);
763 
768  template <int dim>
769  std::vector<unsigned int>
771 
777  template <int dim>
778  void
780  std::vector<unsigned int> &l2h);
781 
786  template <int dim>
787  std::vector<unsigned int>
789 
822  template <int dim>
824  get_fe_from_name (const std::string &name);
825 
826 
869  template <int dim, int spacedim>
870  void add_fe_name (const std::string &name,
871  const FEFactoryBase<dim,spacedim> *factory);
872 
882  DeclException1 (ExcInvalidFEName,
883  std::string,
884  << "Can't re-generate a finite element from the string '"
885  << arg1 << "'.");
886 
898  DeclException2 (ExcInvalidFEDimension,
899  char, int,
900  << "The dimension " << arg1
901  << " in the finite element string must match "
902  << "the space dimension "
903  << arg2 << ".");
904 
910  DeclException0 (ExcInvalidFE);
911 
917  DeclException0 (ExcFENotPrimitive);
923  DeclException0 (ExcTriangulationMismatch);
924 
931  DeclException1 (ExcHangingNodesNotAllowed,
932  int,
933  << "You are using continuous elements on a grid with "
934  << "hanging nodes but without providing hanging node "
935  << "constraints. Use the respective function with "
936  << "additional ConstraintMatrix argument(s), instead."
937  << (arg1?"":""));
943  DeclException0 (ExcGridNotRefinedAtLeastOnce);
949  DeclException4 (ExcMatrixDimensionMismatch,
950  int, int, int, int,
951  << "This is a " << arg1 << "x" << arg2 << " matrix, "
952  << "but should be a " << arg3 << "x" << arg4 << " matrix.");
953 
959  DeclException1(ExcLeastSquaresError, double,
960  << "Least squares fit leaves a gap of " << arg1);
961 
967  DeclException2 (ExcNotGreaterThan,
968  int, int,
969  << arg1 << " must be greater than " << arg2);
970 }
971 
972 
973 #ifndef DOXYGEN
974 
975 namespace FETools
976 {
977  template <class FE>
979  FEFactory<FE>::get (const unsigned int degree) const
980  {
981  return new FE(degree);
982  }
983 }
984 
985 #endif
986 
989 DEAL_II_NAMESPACE_CLOSE
990 
991 /*---------------------------- fe_tools.h ---------------------------*/
992 /* end of #ifndef __deal2__fe_tools_H */
993 #endif
994 /*---------------------------- fe_tools.h ---------------------------*/
void project_dg(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const DoFHandler< dim, spacedim > &dof2, OutVector &u2)
void compute_interpolation_to_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, FullMatrix< double > &I_q)
void compute_component_wise(const FiniteElement< dim, spacedim > &fe, std::vector< unsigned int > &renumbering, std::vector< std::vector< unsigned int > > &start_indices)
void compute_block_renumbering(const FiniteElement< dim, spacedim > &fe, std::vector< types::global_dof_index > &renumbering, std::vector< types::global_dof_index > &block_data, bool return_start_indices=true)
void interpolation_difference(const DoFHandler< dim, spacedim > &dof1, const InVector &z1, const FiniteElement< dim, spacedim > &fe2, OutVector &z1_difference)
DeclException2(ExcInvalidFEDimension, char, int,<< "The dimension "<< arg1<< " in the finite element string must match "<< "the space dimension "<< arg2<< ".")
void extrapolate(const DoFHandler< dim, spacedim > &dof1, const InVector &z1, const DoFHandler< dim, spacedim > &dof2, OutVector &z2)
void get_projection_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &matrix)
void get_interpolation_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &interpolation_matrix)
void compute_projection_from_face_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim-1 > &lhs_quadrature, const Quadrature< dim-1 > &rhs_quadrature, const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell, const unsigned int face, FullMatrix< double > &X)
void back_interpolate(const DoFHandler< dim, spacedim > &dof1, const InVector &u1, const FiniteElement< dim, spacedim > &fe2, OutVector &u1_interpolated)
void compute_embedding_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number > > > &matrices, const bool isotropic_only=false)
void compute_projection_from_quadrature_points_matrix(const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &lhs_quadrature, const Quadrature< dim > &rhs_quadrature, FullMatrix< double > &X)
void get_back_interpolation_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &interpolation_matrix)
void compute_face_embedding_matrices(const FiniteElement< dim, spacedim > &fe, FullMatrix< number >(&matrices)[GeometryInfo< dim >::max_children_per_face], const unsigned int face_coarse, const unsigned int face_fine)
void compute_projection_from_quadrature_points(const FullMatrix< double > &projection_matrix, const std::vector< Tensor< 1, dim > > &vector_of_tensors_at_qp, std::vector< Tensor< 1, dim > > &vector_of_tensors_at_nodes)
FiniteElement< dim, dim > * get_fe_from_name(const std::string &name)
void lexicographic_to_hierarchic_numbering(const FiniteElementData< dim > &fe_data, std::vector< unsigned int > &l2h)
DeclException1(ExcInvalidFEName, std::string,<< "Can't re-generate a finite element from the string '"<< arg1<< "'.")
void add_fe_name(const std::string &name, const FEFactoryBase< dim, spacedim > *factory)
DeclException0(ExcInvalidFE)
void get_interpolation_difference_matrix(const FiniteElement< dim, spacedim > &fe1, const FiniteElement< dim, spacedim > &fe2, FullMatrix< number > &difference_matrix)
void interpolate(const DH1< dim, spacedim > &dof1, const InVector &u1, const DH2< dim, spacedim > &dof2, OutVector &u2)
void compute_projection_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number > > > &matrices, const bool isotropic_only=false)
void compute_node_matrix(FullMatrix< double > &M, const FiniteElement< dim, spacedim > &fe)
void hierarchic_to_lexicographic_numbering(const FiniteElementData< dim > &fe_data, std::vector< unsigned int > &h2l)
virtual FiniteElement< FE::dimension, FE::dimension > * get(const unsigned int degree) const
DeclException4(ExcMatrixDimensionMismatch, int, int, int, int,<< "This is a "<< arg1<< "x"<< arg2<< " matrix, "<< "but should be a "<< arg3<< "x"<< arg4<< " matrix.")