Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
mapping.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping.h 31258 2013-10-16 14:41:36Z kronbichler @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__mapping_h
18 #define __deal2__mapping_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/derivative_form.h>
23 #include <deal.II/base/vector_slice.h>
24 #include <deal.II/grid/tria.h>
25 #include <deal.II/dofs/dof_handler.h>
26 #include <deal.II/fe/fe_update_flags.h>
27 
28 #include <cmath>
29 
31 
32 template <int dim> class Quadrature;
33 template <int dim, int spacedim> class FEValuesData;
34 template <int dim, int spacedim> class FEValuesBase;
35 template <int dim, int spacedim> class FEValues;
36 template <int dim, int spacedim> class FEFaceValues;
37 template <int dim, int spacedim> class FESubfaceValues;
38 
59 {
61  mapping_none = 0x0000,
81  mapping_piola = 0x0100,
91 
100  mapping_nedelec = 0x0200,
105 };
106 
107 
225 template <int dim, int spacedim=dim>
226 class Mapping : public Subscriptor
227 {
228 public:
229 
233  virtual ~Mapping ();
234 
241  virtual Point<spacedim>
243  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
244  const Point<dim> &p) const = 0;
245 
284  virtual Point<dim>
286  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
287  const Point<spacedim> &p) const = 0;
288 
327  {
328  private:
333 
334  public:
341  InternalDataBase ();
342 
347  virtual ~InternalDataBase ();
348 
354 
360 
366 
377 
398  bool is_first_cell () const;
399 
408  virtual void clear_first_cell ();
409 
416  virtual std::size_t memory_consumption () const;
417 
424  std::vector<double> volume_elements;
425 
431  std::vector<Point<spacedim> > support_point_values;
432 
433  /*
434  * The Jacobian of the
435  * transformation in the
436  * (generalized) support
437  * points.
438  */
439  std::vector<Tensor<2,spacedim> > support_point_gradients;
440 
441  /*
442  * The inverse of the
443  * Jacobian of the
444  * transformation in the
445  * (generalized) support
446  * points.
447  */
448  std::vector<Tensor<2,spacedim> > support_point_inverse_gradients;
449 
450 
451  private:
457  };
458 
506  virtual
507  void
508  transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
509  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
510  const InternalDataBase &internal,
511  const MappingType type) const = 0;
512 
513 
514 
545  virtual
546  void
547  transform (const VectorSlice<const std::vector< DerivativeForm<1, dim, spacedim> > > input,
548  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
549  const InternalDataBase &internal,
550  const MappingType type) const = 0;
551 
552 
553 
596  virtual
597  void
598  transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
599  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
600  const InternalDataBase &internal,
601  const MappingType type) const = 0;
602 
606  void
607  transform_covariant (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
608  const unsigned int offset,
609  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
610  const InternalDataBase &internal) const DEAL_II_DEPRECATED;
611 
615  void
616  transform_covariant (const VectorSlice<const std::vector<DerivativeForm<1, dim ,spacedim> > > input,
617  const unsigned int offset,
618  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
619  const InternalDataBase &internal) const DEAL_II_DEPRECATED;
620 
624  void
625  transform_contravariant (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
626  const unsigned int offset,
627  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
628  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const DEAL_II_DEPRECATED;
629 
634  void
635  transform_contravariant (const VectorSlice<const std::vector<DerivativeForm<1, dim,spacedim> > > input,
636  const unsigned int offset,
637  const VectorSlice<std::vector<Tensor<2,spacedim> > > output,
638  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const DEAL_II_DEPRECATED;
639 
644  const Point<spacedim> &support_point_value(
645  const unsigned int index,
646  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
647 
654  const Tensor<2,spacedim> &support_point_gradient(
655  const unsigned int index,
656  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
657 
664  const Tensor<2,spacedim> &support_point_inverse_gradient(
665  const unsigned int index,
666  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const;
667 
681  virtual
682  Mapping<dim,spacedim> *clone () const = 0;
683 
699  virtual
700  bool preserves_vertex_locations () const = 0;
701 
705  DeclException0 (ExcInvalidData);
706 
707 
717  DeclException0(ExcTransformationFailed);
718 
726  DeclException3 (ExcDistortedMappedCell,
727  Point<spacedim>, double, int,
728  << "The image of the mapping applied to cell with center ["
729  << arg1 << "] is distorted. The cell geometry or the "
730  << "mapping are invalid, giving a non-positive volume "
731  << "fraction of " << arg2 << " in quadrature point "
732  << arg3 << ".");
733 
734 private:
735 
747  virtual UpdateFlags update_once (const UpdateFlags) const = 0;
748 
756  virtual UpdateFlags update_each (const UpdateFlags) const = 0;
757 
763  virtual InternalDataBase *
764  get_data (const UpdateFlags,
765  const Quadrature<dim> &quadrature) const = 0;
766 
773  virtual InternalDataBase *
774  get_face_data (const UpdateFlags flags,
775  const Quadrature<dim-1>& quadrature) const = 0;
776 
784  virtual InternalDataBase *
785  get_subface_data (const UpdateFlags flags,
786  const Quadrature<dim-1>& quadrature) const = 0;
787 
788 
817  /* virtual void */
818  /* fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell, */
819  /* const Quadrature<dim> &quadrature, */
820  /* InternalDataBase &internal, */
821  /* std::vector<Point<spacedim> > &quadrature_points, */
822  /* std::vector<double> &JxW_values) const = 0; */
823 
829  virtual void
830  fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
831  const Quadrature<dim> &quadrature,
833  std::vector<Point<spacedim> > &quadrature_points,
834  std::vector<double> &JxW_values,
835  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
836  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
837  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
838  std::vector<Point<spacedim> > &cell_normal_vectors,
839  CellSimilarity::Similarity &cell_similarity
840  ) const=0;
841 
842 
843 
861  virtual void
862  fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
863  const unsigned int face_no,
864  const Quadrature<dim-1> &quadrature,
865  InternalDataBase &internal,
866  std::vector<Point<spacedim> > &quadrature_points,
867  std::vector<double> &JxW_values,
868  std::vector<Tensor<1,spacedim> > &boundary_form,
869  std::vector<Point<spacedim> > &normal_vectors) const = 0;
870 
874  virtual void
875  fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
876  const unsigned int face_no,
877  const unsigned int sub_no,
878  const Quadrature<dim-1> &quadrature,
879  InternalDataBase &internal,
880  std::vector<Point<spacedim> > &quadrature_points,
881  std::vector<double> &JxW_values,
882  std::vector<Tensor<1,spacedim> > &boundary_form,
883  std::vector<Point<spacedim> > &normal_vectors) const = 0;
884 
891  friend class FEValuesBase<dim,spacedim>;
892  friend class FEValues<dim,spacedim>;
893  friend class FEFaceValues<dim,spacedim>;
894  friend class FESubfaceValues<dim,spacedim>;
895 };
896 
897 
898 /* ------------------------- inline functions ------------------------- */
899 
900 #ifndef DOXYGEN
901 
902 template <int dim, int spacedim>
903 inline
904 UpdateFlags
906 {
907  if (first_cell)
908  {
909  Assert (update_flags==(update_once|update_each),
910  ExcInternalError());
911  return update_flags;
912  }
913  else
914  return update_each;
915 }
916 
917 
918 
919 template <int dim, int spacedim>
920 inline
921 bool
923 {
924  return first_cell;
925 }
926 
927 
928 
929 template <int dim, int spacedim>
930 inline
931 void
933 {
934  first_cell = false;
935 }
936 
937 
938 
939 template <int dim, int spacedim>
940 inline
941 const Point<spacedim> &
943  const unsigned int index,
944  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const
945 {
946  AssertIndexRange(index, internal.support_point_values.size());
947  return internal.support_point_values[index];
948 }
949 
950 
951 template <int dim, int spacedim>
952 inline
953 const Tensor<2,spacedim> &
955  const unsigned int index,
956  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const
957 {
958  AssertIndexRange(index, internal.support_point_gradients.size());
959  return internal.support_point_gradients[index];
960 }
961 
962 
963 template <int dim, int spacedim>
964 inline
965 const Tensor<2,spacedim> &
967  const unsigned int index,
968  const typename Mapping<dim,spacedim>::InternalDataBase &internal) const
969 {
970  AssertIndexRange(index, internal.support_point_inverse_gradients.size());
971  return internal.support_point_inverse_gradients[index];
972 }
973 
974 
975 #endif // DOXYGEN
976 
977 DEAL_II_NAMESPACE_CLOSE
978 
979 #endif
virtual void clear_first_cell()
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const =0
virtual UpdateFlags update_once(const UpdateFlags) const =0
virtual ~Mapping()
The mapping used for BDM elements.
Definition: mapping.h:104
virtual std::size_t memory_consumption() const
bool is_first_cell() const
MappingType
Definition: mapping.h:58
Mapping of the gradient of a covariant vector field (see Mapping::transform() for details) ...
Definition: mapping.h:67
#define AssertIndexRange(index, range)
Definition: exceptions.h:888
STL namespace.
The mapping used for Nedelec elements.
Definition: mapping.h:100
std::vector< double > volume_elements
Definition: mapping.h:424
DeclException0(ExcInvalidData)
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< DerivativeForm< 1, dim, spacedim > > &jacobians, std::vector< DerivativeForm< 2, dim, spacedim > > &jacobian_grads, std::vector< DerivativeForm< 1, spacedim, dim > > &inverse_jacobians, std::vector< Point< spacedim > > &cell_normal_vectors, CellSimilarity::Similarity &cell_similarity) const =0
const Tensor< 2, spacedim > & support_point_gradient(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
Iterator is invalid, probably due to an error.
virtual void transform(const VectorSlice< const std::vector< Tensor< 1, dim > > > input, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const InternalDataBase &internal, const MappingType type) const =0
No mapping.
Definition: mapping.h:61
Mapping of the gradient of a contravariant vector field (see Mapping::transform() for details) ...
Definition: mapping.h:69
#define Assert(cond, exc)
Definition: exceptions.h:299
UpdateFlags
UpdateFlags current_update_flags() const
DeclException3(ExcDistortedMappedCell, Point< spacedim >, double, int,<< "The image of the mapping applied to cell with center ["<< arg1<< "] is distorted. The cell geometry or the "<< "mapping are invalid, giving a non-positive volume "<< "fraction of "<< arg2<< " in quadrature point "<< arg3<< ".")
UpdateFlags update_flags
Definition: mapping.h:353
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< Tensor< 1, spacedim > > &boundary_form, std::vector< Point< spacedim > > &normal_vectors) const =0
The mapping used for Raviart-Thomas elements.
Definition: mapping.h:102
virtual Mapping< dim, spacedim > * clone() const =0
const Tensor< 2, spacedim > & support_point_inverse_gradient(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
virtual InternalDataBase * get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0
void transform_contravariant(const VectorSlice< const std::vector< Tensor< 1, dim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const DEAL_II_DEPRECATED
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const =0
virtual InternalDataBase * get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0
virtual bool preserves_vertex_locations() const =0
const Point< spacedim > & support_point_value(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, InternalDataBase &internal, std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, std::vector< Tensor< 1, spacedim > > &boundary_form, std::vector< Point< spacedim > > &normal_vectors) const =0
Covariant mapping (see Mapping::transform() for details)
Definition: mapping.h:63
virtual InternalDataBase * get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const =0
UpdateFlags update_once
Definition: mapping.h:359
::ExceptionBase & ExcInternalError()
Contravariant mapping (see Mapping::transform() for details)
Definition: mapping.h:65
void transform_covariant(const VectorSlice< const std::vector< Tensor< 1, dim > > > input, const unsigned int offset, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const InternalDataBase &internal) const DEAL_II_DEPRECATED
virtual UpdateFlags update_each(const UpdateFlags) const =0
std::vector< Point< spacedim > > support_point_values
Definition: mapping.h:431
UpdateFlags update_each
Definition: mapping.h:365