17 #ifndef __deal2__mapping_h
18 #define __deal2__mapping_h
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>
35 template <
int dim,
int spacedim>
class FEValues;
225 template <
int dim,
int spacedim=dim>
439 std::vector<Tensor<2,spacedim> > support_point_gradients;
448 std::vector<Tensor<2,spacedim> > support_point_inverse_gradients;
608 const unsigned int offset,
617 const
unsigned int offset,
626 const
unsigned int offset,
636 const
unsigned int offset,
645 const
unsigned int index,
655 const
unsigned int index,
665 const
unsigned int index,
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 "
775 const
Quadrature<dim-1>& quadrature) const = 0;
786 const
Quadrature<dim-1>& quadrature) const = 0;
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,
863 const
unsigned int face_no,
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;
876 const
unsigned int face_no,
877 const
unsigned int sub_no,
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;
892 friend class
FEValues<dim,spacedim>;
902 template <
int dim,
int spacedim>
919 template <
int dim,
int spacedim>
929 template <
int dim,
int spacedim>
939 template <
int dim,
int spacedim>
943 const unsigned int index,
947 return internal.support_point_values[index];
951 template <
int dim,
int spacedim>
955 const unsigned int index,
959 return internal.support_point_gradients[index];
963 template <
int dim,
int spacedim>
967 const unsigned int index,
971 return internal.support_point_inverse_gradients[index];
977 DEAL_II_NAMESPACE_CLOSE
virtual void clear_first_cell()
The mapping used for BDM elements.
virtual std::size_t memory_consumption() const
bool is_first_cell() const
Mapping of the gradient of a covariant vector field (see Mapping::transform() for details) ...
virtual InternalDataBase * get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const =0
#define AssertIndexRange(index, range)
The mapping used for Nedelec elements.
std::vector< double > volume_elements
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const =0
virtual UpdateFlags update_once(const UpdateFlags) const =0
Iterator is invalid, probably due to an error.
virtual InternalDataBase * get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0
Mapping of the gradient of a contravariant vector field (see Mapping::transform() for details) ...
#define Assert(cond, exc)
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
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<< ".")
DeclException0(ExcInvalidData)
virtual ~InternalDataBase()
The mapping used for Raviart-Thomas elements.
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
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
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
virtual InternalDataBase * get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const =0
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
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
const Point< spacedim > & support_point_value(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
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
const Tensor< 2, spacedim > & support_point_gradient(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
const Tensor< 2, spacedim > & support_point_inverse_gradient(const unsigned int index, const typename Mapping< dim, spacedim >::InternalDataBase &internal) const
Covariant mapping (see Mapping::transform() for details)
virtual bool preserves_vertex_locations() const =0
virtual Mapping< dim, spacedim > * clone() const =0
::ExceptionBase & ExcInternalError()
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const =0
Contravariant mapping (see Mapping::transform() for details)
std::vector< Point< spacedim > > support_point_values
virtual UpdateFlags update_each(const UpdateFlags) const =0