00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__mapping_h
00013 #define __deal2__mapping_h
00014
00015
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/derivative_form.h>
00018 #include <deal.II/base/vector_slice.h>
00019 #include <deal.II/grid/tria.h>
00020 #include <deal.II/dofs/dof_handler.h>
00021 #include <deal.II/fe/fe_update_flags.h>
00022
00023 #include <cmath>
00024
00025 DEAL_II_NAMESPACE_OPEN
00026
00027 template <int dim> class Quadrature;
00028 template <int dim, int spacedim> class FEValuesData;
00029 template <int dim, int spacedim> class FEValuesBase;
00030 template <int dim, int spacedim> class FEValues;
00031 template <int dim, int spacedim> class FEFaceValues;
00032 template <int dim, int spacedim> class FESubfaceValues;
00033
00051 enum MappingType
00052 {
00054 mapping_none = 0x0000,
00056 mapping_covariant = 0x0001,
00058 mapping_contravariant = 0x0002,
00060 mapping_covariant_gradient = 0x0003,
00062 mapping_contravariant_gradient = 0x0004,
00074 mapping_piola = 0x0100,
00082 mapping_piola_gradient = 0x0101,
00084
00093 mapping_nedelec = 0x0200,
00095 mapping_raviart_thomas = 0x0300,
00097 mapping_bdm = mapping_raviart_thomas
00098 };
00099
00100
00218 template <int dim, int spacedim=dim>
00219 class Mapping : public Subscriptor
00220 {
00221 public:
00222
00226 virtual ~Mapping ();
00227
00234 virtual Point<spacedim>
00235 transform_unit_to_real_cell (
00236 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00237 const Point<dim> &p) const = 0;
00238
00252 virtual Point<dim>
00253 transform_real_to_unit_cell (
00254 const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00255 const Point<spacedim> &p) const = 0;
00256
00294 class InternalDataBase: public Subscriptor
00295 {
00296 private:
00300 InternalDataBase (const InternalDataBase&);
00301
00302 public:
00309 InternalDataBase ();
00310
00315 virtual ~InternalDataBase ();
00316
00321 UpdateFlags update_flags;
00322
00327 UpdateFlags update_once;
00328
00333 UpdateFlags update_each;
00334
00344 UpdateFlags current_update_flags() const;
00345
00366 bool is_first_cell () const;
00367
00376 virtual void clear_first_cell ();
00377
00384 virtual std::size_t memory_consumption () const;
00385
00392 std::vector<double> volume_elements;
00393
00399 std::vector<Point<spacedim> > support_point_values;
00400
00401
00402
00403
00404
00405
00406
00407 std::vector<Tensor<2,spacedim> > support_point_gradients;
00408
00409
00410
00411
00412
00413
00414
00415
00416 std::vector<Tensor<2,spacedim> > support_point_inverse_gradients;
00417
00418
00419 private:
00424 bool first_cell;
00425 };
00426
00474 virtual
00475 void
00476 transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
00477 VectorSlice<std::vector<Tensor<1,spacedim> > > output,
00478 const InternalDataBase &internal,
00479 const MappingType type) const = 0;
00480
00481
00482
00513 virtual
00514 void
00515 transform (const VectorSlice<const std::vector< DerivativeForm<1, dim, spacedim> > > input,
00516 VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00517 const InternalDataBase &internal,
00518 const MappingType type) const = 0;
00519
00520
00521
00565 virtual
00566 void
00567 transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
00568 VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00569 const InternalDataBase &internal,
00570 const MappingType type) const = 0;
00571
00572
00573
00574
00575
00576
00577
00578
00579
00583 void
00584 transform_covariant (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
00585 const unsigned int offset,
00586 VectorSlice<std::vector<Tensor<1,spacedim> > > output,
00587 const InternalDataBase &internal) const;
00588
00592 void
00593 transform_covariant (const VectorSlice<const std::vector<DerivativeForm<1, dim ,spacedim> > > input,
00594 const unsigned int offset,
00595 VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00596 const InternalDataBase &internal) const;
00597
00601 void
00602 transform_contravariant (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
00603 const unsigned int offset,
00604 VectorSlice<std::vector<Tensor<1,spacedim> > > output,
00605 const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
00606
00611 void
00612 transform_contravariant (const VectorSlice<const std::vector<DerivativeForm<1, dim,spacedim> > > input,
00613 const unsigned int offset,
00614 const VectorSlice<std::vector<Tensor<2,spacedim> > > output,
00615 const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
00616
00621 const Point<spacedim>& support_point_value(
00622 const unsigned int index,
00623 const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
00624
00631 const Tensor<2,spacedim>& support_point_gradient(
00632 const unsigned int index,
00633 const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
00634
00641 const Tensor<2,spacedim>& support_point_inverse_gradient(
00642 const unsigned int index,
00643 const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
00644
00658 virtual
00659 Mapping<dim,spacedim> * clone () const = 0;
00660
00676 virtual
00677 bool preserves_vertex_locations () const = 0;
00678
00682 DeclException0 (ExcInvalidData);
00683
00684 private:
00685
00697 virtual UpdateFlags update_once (const UpdateFlags) const = 0;
00698
00706 virtual UpdateFlags update_each (const UpdateFlags) const = 0;
00707
00713 virtual InternalDataBase*
00714 get_data (const UpdateFlags,
00715 const Quadrature<dim>& quadrature) const = 0;
00716
00723 virtual InternalDataBase*
00724 get_face_data (const UpdateFlags flags,
00725 const Quadrature<dim-1>& quadrature) const = 0;
00726
00734 virtual InternalDataBase*
00735 get_subface_data (const UpdateFlags flags,
00736 const Quadrature<dim-1>& quadrature) const = 0;
00737
00738
00767
00768
00769
00770
00771
00772
00773
00779 virtual void
00780 fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00781 const Quadrature<dim> &quadrature,
00782 InternalDataBase &internal,
00783 std::vector<Point<spacedim> > &quadrature_points,
00784 std::vector<double> &JxW_values,
00785 std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
00786 std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
00787 std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
00788 std::vector<Point<spacedim> > &cell_normal_vectors,
00789 CellSimilarity::Similarity &cell_similarity
00790 ) const=0;
00791
00792
00793
00811 virtual void
00812 fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00813 const unsigned int face_no,
00814 const Quadrature<dim-1> &quadrature,
00815 InternalDataBase &internal,
00816 std::vector<Point<spacedim> > &quadrature_points,
00817 std::vector<double> &JxW_values,
00818 std::vector<Tensor<1,spacedim> > &boundary_form,
00819 std::vector<Point<spacedim> > &normal_vectors) const = 0;
00820
00824 virtual void
00825 fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
00826 const unsigned int face_no,
00827 const unsigned int sub_no,
00828 const Quadrature<dim-1> &quadrature,
00829 InternalDataBase &internal,
00830 std::vector<Point<spacedim> > &quadrature_points,
00831 std::vector<double> &JxW_values,
00832 std::vector<Tensor<1,spacedim> > &boundary_form,
00833 std::vector<Point<spacedim> > &normal_vectors) const = 0;
00834
00841 friend class FEValuesBase<dim,spacedim>;
00842 friend class FEValues<dim,spacedim>;
00843 friend class FEFaceValues<dim,spacedim>;
00844 friend class FESubfaceValues<dim,spacedim>;
00845 };
00846
00847
00848
00849
00850 #ifndef DOXYGEN
00851
00852 template <int dim, int spacedim>
00853 inline
00854 UpdateFlags
00855 Mapping<dim,spacedim>::InternalDataBase::current_update_flags () const
00856 {
00857 if (first_cell)
00858 {
00859 Assert (update_flags==(update_once|update_each),
00860 ExcInternalError());
00861 return update_flags;
00862 }
00863 else
00864 return update_each;
00865 }
00866
00867
00868
00869 template <int dim, int spacedim>
00870 inline
00871 bool
00872 Mapping<dim,spacedim>::InternalDataBase::is_first_cell () const
00873 {
00874 return first_cell;
00875 }
00876
00877
00878
00879 template <int dim, int spacedim>
00880 inline
00881 void
00882 Mapping<dim,spacedim>::InternalDataBase::clear_first_cell ()
00883 {
00884 first_cell = false;
00885 }
00886
00887
00888
00889 template <int dim, int spacedim>
00890 inline
00891 const Point<spacedim>&
00892 Mapping<dim,spacedim>::support_point_value(
00893 const unsigned int index,
00894 const typename Mapping<dim,spacedim>::InternalDataBase& internal) const
00895 {
00896 AssertIndexRange(index, internal.support_point_values.size());
00897 return internal.support_point_values[index];
00898 }
00899
00900
00901 template <int dim, int spacedim>
00902 inline
00903 const Tensor<2,spacedim>&
00904 Mapping<dim,spacedim>::support_point_gradient(
00905 const unsigned int index,
00906 const typename Mapping<dim,spacedim>::InternalDataBase& internal) const
00907 {
00908 AssertIndexRange(index, internal.support_point_gradients.size());
00909 return internal.support_point_gradients[index];
00910 }
00911
00912
00913 template <int dim, int spacedim>
00914 inline
00915 const Tensor<2,spacedim>&
00916 Mapping<dim,spacedim>::support_point_inverse_gradient(
00917 const unsigned int index,
00918 const typename Mapping<dim,spacedim>::InternalDataBase& internal) const
00919 {
00920 AssertIndexRange(index, internal.support_point_inverse_gradients.size());
00921 return internal.support_point_inverse_gradients[index];
00922 }
00923
00924
00925 #endif // DOXYGEN
00926
00927 DEAL_II_NAMESPACE_CLOSE
00928
00929 #endif
00930