00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__fe_base_h
00013 #define __deal2__fe_base_h
00014
00015 #include <deal.II/base/config.h>
00016 #include <deal.II/base/exceptions.h>
00017 #include <deal.II/base/subscriptor.h>
00018 #include <deal.II/base/point.h>
00019 #include <deal.II/base/tensor.h>
00020 #include <deal.II/base/table.h>
00021 #include <deal.II/base/vector_slice.h>
00022 #include <deal.II/base/geometry_info.h>
00023 #include <deal.II/lac/full_matrix.h>
00024 #include <deal.II/fe/fe_update_flags.h>
00025 #include <deal.II/fe/mapping.h>
00026
00027 #include <string>
00028 #include <vector>
00029
00030 DEAL_II_NAMESPACE_OPEN
00031
00032 template<int dim, int spacedim> class FESystem;
00033
00034
00039 namespace FiniteElementDomination
00040 {
00128 enum Domination
00129 {
00130 this_element_dominates,
00131 other_element_dominates,
00132 neither_element_dominates,
00133 either_element_can_dominate,
00134 no_requirements
00135 };
00136
00137
00166 inline Domination operator & (const Domination d1,
00167 const Domination d2);
00168 }
00169
00170
00180 template <int dim>
00181 class FiniteElementData
00182 {
00183 public:
00238 enum Conformity
00239 {
00245 unknown = 0x00,
00246
00251 L2 = 0x01,
00252
00261 Hcurl = 0x02,
00262
00271 Hdiv = 0x04,
00272
00279 H1 = Hcurl | Hdiv,
00280
00288 H2 = 0x0e
00289 };
00290
00296 static const unsigned int dimension = dim;
00297
00302 const unsigned int dofs_per_vertex;
00303
00309 const unsigned int dofs_per_line;
00310
00318 const unsigned int dofs_per_quad;
00319
00327 const unsigned int dofs_per_hex;
00328
00332 const unsigned int first_line_index;
00333
00337 const unsigned int first_quad_index;
00338
00342 const unsigned int first_hex_index;
00343
00347 const unsigned int first_face_line_index;
00348
00352 const unsigned int first_face_quad_index;
00353
00363 const unsigned int dofs_per_face;
00364
00374 const unsigned int dofs_per_cell;
00375
00394 const unsigned int components;
00395
00401 const unsigned int degree;
00402
00406 const Conformity conforming_space;
00407
00417 BlockIndices block_indices_data;
00418
00427 FiniteElementData ();
00428
00453 FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
00454 const unsigned int n_components,
00455 const unsigned int degree,
00456 const Conformity conformity = unknown,
00457 const unsigned int n_blocks = numbers::invalid_unsigned_int);
00458
00462 unsigned int n_dofs_per_vertex () const;
00463
00469 unsigned int n_dofs_per_line () const;
00470
00476 unsigned int n_dofs_per_quad () const;
00477
00483 unsigned int n_dofs_per_hex () const;
00484
00491 unsigned int n_dofs_per_face () const;
00492
00499 unsigned int n_dofs_per_cell () const;
00500
00517 template <int structdim>
00518 unsigned int n_dofs_per_object () const;
00519
00523 unsigned int n_components () const;
00524
00528 unsigned int n_blocks () const;
00529
00533 const BlockIndices& block_indices() const;
00534
00551 bool is_primitive () const;
00552
00562 unsigned int tensor_degree () const;
00563
00575 bool conforms (const Conformity) const;
00576
00584 unsigned int face_to_cell_index (const unsigned int face_index,
00585 const unsigned int face,
00586 const bool face_orientation = true,
00587 const bool face_flip = false,
00588 const bool face_rotation = false) const;
00589
00628 unsigned int face_to_equivalent_cell_index (const unsigned int index) const;
00629
00633 bool operator == (const FiniteElementData &) const;
00634
00635 protected:
00636
00645 void set_primitivity(const bool value);
00646
00647 private:
00657 bool cached_primitivity;
00658 };
00659
00660
00661
00662
00663
00664
00665 #ifndef DOXYGEN
00666
00667 namespace FiniteElementDomination
00668 {
00669 inline
00670 Domination operator & (const Domination d1,
00671 const Domination d2)
00672 {
00673
00674
00675
00676
00677
00678
00679
00680
00681
00682
00683
00684 switch (d1)
00685 {
00686 case this_element_dominates:
00687 if ((d2 == this_element_dominates) ||
00688 (d2 == either_element_can_dominate) ||
00689 (d2 == no_requirements))
00690 return this_element_dominates;
00691 else
00692 return neither_element_dominates;
00693
00694 case other_element_dominates:
00695 if ((d2 == other_element_dominates) ||
00696 (d2 == either_element_can_dominate) ||
00697 (d2 == no_requirements))
00698 return other_element_dominates;
00699 else
00700 return neither_element_dominates;
00701
00702 case neither_element_dominates:
00703 return neither_element_dominates;
00704
00705 case either_element_can_dominate:
00706 if (d2 == no_requirements)
00707 return either_element_can_dominate;
00708 else
00709 return d2;
00710
00711 case no_requirements:
00712 return d2;
00713
00714 default:
00715
00716 Assert (false, ExcInternalError());
00717 }
00718
00719 return neither_element_dominates;
00720 }
00721 }
00722
00723
00724 template <int dim>
00725 inline
00726 unsigned int
00727 FiniteElementData<dim>::n_dofs_per_vertex () const
00728 {
00729 return dofs_per_vertex;
00730 }
00731
00732
00733
00734 template <int dim>
00735 inline
00736 unsigned int
00737 FiniteElementData<dim>::n_dofs_per_line () const
00738 {
00739 return dofs_per_line;
00740 }
00741
00742
00743
00744 template <int dim>
00745 inline
00746 unsigned int
00747 FiniteElementData<dim>::n_dofs_per_quad () const
00748 {
00749 return dofs_per_quad;
00750 }
00751
00752
00753
00754 template <int dim>
00755 inline
00756 unsigned int
00757 FiniteElementData<dim>::n_dofs_per_hex () const
00758 {
00759 return dofs_per_hex;
00760 }
00761
00762
00763
00764 template <int dim>
00765 inline
00766 unsigned int
00767 FiniteElementData<dim>::n_dofs_per_face () const
00768 {
00769 return dofs_per_face;
00770 }
00771
00772
00773
00774 template <int dim>
00775 inline
00776 unsigned int
00777 FiniteElementData<dim>::n_dofs_per_cell () const
00778 {
00779 return dofs_per_cell;
00780 }
00781
00782
00783
00784 template <int dim>
00785 template <int structdim>
00786 inline
00787 unsigned int
00788 FiniteElementData<dim>::n_dofs_per_object () const
00789 {
00790 switch (structdim)
00791 {
00792 case 0:
00793 return dofs_per_vertex;
00794 case 1:
00795 return dofs_per_line;
00796 case 2:
00797 return dofs_per_quad;
00798 case 3:
00799 return dofs_per_hex;
00800 default:
00801 Assert (false, ExcInternalError());
00802 }
00803 return numbers::invalid_unsigned_int;
00804 }
00805
00806
00807
00808 template <int dim>
00809 inline
00810 unsigned int
00811 FiniteElementData<dim>::n_components () const
00812 {
00813 return components;
00814 }
00815
00816
00817
00818 template <int dim>
00819 inline
00820 bool
00821 FiniteElementData<dim>::is_primitive () const
00822 {
00823 return cached_primitivity;
00824 }
00825
00826
00827 template <int dim>
00828 inline
00829 void
00830 FiniteElementData<dim>::set_primitivity (const bool value)
00831 {
00832 cached_primitivity = value;
00833 }
00834
00835
00836 template <int dim>
00837 inline
00838 const BlockIndices&
00839 FiniteElementData<dim>::block_indices () const
00840 {
00841 return block_indices_data;
00842 }
00843
00844
00845
00846 template <int dim>
00847 inline
00848 unsigned int
00849 FiniteElementData<dim>::n_blocks () const
00850 {
00851 return block_indices_data.size();
00852 }
00853
00854
00855
00856 template <int dim>
00857 inline
00858 unsigned int
00859 FiniteElementData<dim>::tensor_degree () const
00860 {
00861 return degree;
00862 }
00863
00864
00865 template <int dim>
00866 inline
00867 bool
00868 FiniteElementData<dim>::conforms (const Conformity space) const
00869 {
00870 return ((space & conforming_space) == space);
00871 }
00872
00873
00874
00875 template <>
00876 inline
00877 unsigned int
00878 FiniteElementData<1>::
00879 face_to_equivalent_cell_index (const unsigned int index) const
00880 {
00881 Assert (index < dofs_per_face,
00882 ExcIndexRange (index, 0, dofs_per_face));
00883 return index;
00884 }
00885
00886
00887
00888 template <>
00889 inline
00890 unsigned int
00891 FiniteElementData<2>::
00892 face_to_equivalent_cell_index (const unsigned int index) const
00893 {
00894 Assert (index < dofs_per_face,
00895 ExcIndexRange (index, 0, dofs_per_face));
00896
00897
00898 if (index < this->dofs_per_vertex)
00899 return index;
00900 else if (index < 2*this->dofs_per_vertex)
00901 return index + this->dofs_per_vertex;
00902 else
00903
00904
00905
00906 return index + 2*this->dofs_per_vertex;
00907 }
00908
00909
00910
00911
00912 template <>
00913 inline
00914 unsigned int
00915 FiniteElementData<3>::
00916 face_to_equivalent_cell_index (const unsigned int index) const
00917 {
00918
00919
00920
00921 return face_to_cell_index(index, 0, true);
00922 }
00923
00924 #endif // DOXYGEN
00925
00926
00927 DEAL_II_NAMESPACE_CLOSE
00928
00929 #endif
00930