Reference documentation for deal.II version 8.1.0
mapping_q1.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping_q1.h 30036 2013-07-18 16:55:32Z maier @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_q1_h
18 #define __deal2__mapping_q1_h
19 
20 
21 #include <deal.II/base/derivative_form.h>
22 #include <deal.II/base/config.h>
23 #include <deal.II/base/table.h>
24 #include <deal.II/base/qprojector.h>
25 #include <deal.II/grid/tria_iterator.h>
26 #include <deal.II/dofs/dof_accessor.h>
27 #include <deal.II/fe/mapping.h>
28 
29 #include <cmath>
30 
32 
35 
36 
56 template <int dim, int spacedim=dim>
57 class MappingQ1 : public Mapping<dim,spacedim>
58 {
59 public:
63  MappingQ1 ();
64 
65  virtual Point<spacedim>
67  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
68  const Point<dim> &p) const;
69 
112  virtual Point<dim>
114  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
115  const Point<spacedim> &p) const;
116 
117  virtual void
118  transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
119  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
120  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
121  const MappingType type) const;
122 
123  virtual void
124  transform (const VectorSlice<const std::vector<DerivativeForm<1, dim,spacedim> > > input,
125  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
126  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
127  const MappingType type) const;
128 
129  virtual
130  void
131  transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
132  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
133  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
134  const MappingType type) const;
135 
136 
137 protected:
143  template < int rank >
144  void
145  transform_fields(const VectorSlice<const std::vector<Tensor<rank,dim> > > input,
146  VectorSlice< std::vector<Tensor<rank,spacedim> > > output,
147  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
148  const MappingType type) const;
152  template < int rank >
153  void
154  transform_gradients(const VectorSlice<const std::vector<Tensor<rank,dim> > > input,
155  VectorSlice< std::vector<Tensor<rank,spacedim> > > output,
156  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
157  const MappingType type) const;
161  template < int rank >
162  void
164  const VectorSlice<const std::vector<DerivativeForm<rank, dim, spacedim> > > input,
166  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
167  const MappingType type) const;
168 
169 public:
170 
171 
172 
173 
179  virtual
180  Mapping<dim,spacedim> *clone () const;
181 
186  class InternalData : public Mapping<dim,spacedim>::InternalDataBase
187  {
188  public:
193  InternalData(const unsigned int n_shape_functions);
194 
202  double shape (const unsigned int qpoint,
203  const unsigned int shape_nr) const;
204 
209  double &shape (const unsigned int qpoint,
210  const unsigned int shape_nr);
211 
217  Tensor<1,dim> derivative (const unsigned int qpoint,
218  const unsigned int shape_nr) const;
219 
225  Tensor<1,dim> &derivative (const unsigned int qpoint,
226  const unsigned int shape_nr);
227 
233  Tensor<2,dim> second_derivative (const unsigned int qpoint,
234  const unsigned int shape_nr) const;
235 
241  Tensor<2,dim> &second_derivative (const unsigned int qpoint,
242  const unsigned int shape_nr);
243 
250  virtual std::size_t memory_consumption () const;
251 
259  std::vector<double> shape_values;
260 
268  std::vector<Tensor<1,dim> > shape_derivatives;
269 
278  std::vector<Tensor<2,dim> > shape_second_derivatives;
279 
299  std::vector<DerivativeForm<1,dim, spacedim > > covariant;
300 
312  std::vector< DerivativeForm<1,dim,spacedim> > contravariant;
313 
336  std::vector<std::vector<Tensor<1,dim> > > unit_tangentials;
337 
341  std::vector<std::vector<Tensor<1,spacedim> > > aux;
342 
348  std::vector<Point<spacedim> > mapping_support_points;
349 
356 
365 
379  unsigned int n_shape_functions;
380  };
381 
389  typedef
392 
397  virtual void
399  const Quadrature<dim> &quadrature,
400  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
401  typename std::vector<Point<spacedim> > &quadrature_points,
402  std::vector<double> &JxW_values,
403  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
404  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
405  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
406  std::vector<Point<spacedim> > &cell_normal_vectors,
407  CellSimilarity::Similarity &cell_similarity) const;
408 
413  virtual void
415  const unsigned int face_no,
416  const Quadrature<dim-1> &quadrature,
417  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
418  typename std::vector<Point<spacedim> > &quadrature_points,
419  std::vector<double> &JxW_values,
420  typename std::vector<Tensor<1,spacedim> > &boundary_form,
421  typename std::vector<Point<spacedim> > &normal_vectors) const ;
422 
427  virtual void
429  const unsigned int face_no,
430  const unsigned int sub_no,
431  const Quadrature<dim-1>& quadrature,
432  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
433  typename std::vector<Point<spacedim> > &quadrature_points,
434  std::vector<double> &JxW_values,
435  typename std::vector<Tensor<1,spacedim> > &boundary_form,
436  typename std::vector<Point<spacedim> > &normal_vectors) const ;
437 
450  void compute_shapes (const std::vector<Point<dim> > &unit_points,
451  InternalData &data) const;
452 
461  void compute_data (const UpdateFlags flags,
462  const Quadrature<dim> &quadrature,
463  const unsigned int n_orig_q_points,
464  InternalData &data) const;
465 
478  void compute_face_data (const UpdateFlags flags,
479  const Quadrature<dim> &quadrature,
480  const unsigned int n_orig_q_points,
481  InternalData &data) const;
482 
487  void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
488  const unsigned int npts,
489  const DataSetDescriptor data_set,
490  const CellSimilarity::Similarity cell_similarity,
491  InternalData &data,
492  std::vector<Point<spacedim> > &quadrature_points) const;
493 
499  const unsigned int face_no,
500  const unsigned int subface_no,
501  const unsigned int npts,
502  const DataSetDescriptor data_set,
503  const std::vector<double> &weights,
504  InternalData &mapping_data,
505  std::vector<Point<spacedim> > &quadrature_points,
506  std::vector<double> &JxW_values,
507  std::vector<Tensor<1,spacedim> > &boundary_form,
508  std::vector<Point<spacedim> > &normal_vectors) const;
509 
514  virtual void compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
515  InternalData &data) const;
516 
545 
576  Point<dim>
578  const Point<spacedim> &p,
579  const Point<dim> &initial_p_unit,
580  InternalData &mdata) const;
581 
586  virtual
587  bool preserves_vertex_locations () const;
588 
589 protected:
590  /* Trick to templatize transform_real_to_unit_cell<dim, dim+1> */
591  template<int dim_>
593  transform_real_to_unit_cell_internal_codim1
594  (const typename Triangulation<dim_,dim_+1>::cell_iterator &cell,
595  const Point<dim_+1> &p,
596  const Point<dim_> &initial_p_unit,
597  InternalData &mdata) const;
598 
638  Point<dim>
639  transform_real_to_unit_cell_initial_guess (const std::vector<Point<spacedim> > &vertex,
640  const Point<spacedim> &p) const;
641 
642 
643 private:
672  virtual UpdateFlags update_once (const UpdateFlags flags) const;
673 
718  virtual UpdateFlags update_each (const UpdateFlags flags) const;
719 
720  virtual
722  get_data (const UpdateFlags,
723  const Quadrature<dim> &quadrature) const;
724 
725  virtual
727  get_face_data (const UpdateFlags flags,
728  const Quadrature<dim-1>& quadrature) const;
729 
730  virtual
732  get_subface_data (const UpdateFlags flags,
733  const Quadrature<dim-1>& quadrature) const;
734 
751  virtual void compute_mapping_support_points(
752  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
753  std::vector<Point<spacedim> > &a) const;
754 
761 };
762 
763 
764 // explicit specializations
765 
766 template<>
767 Point<2>
771  const Point<3> &p,
772  const Point<2> &initial_p_unit,
773  InternalData &mdata) const;
774 
775 template<>
776 Point<1>
780  const Point<2> &p,
781  const Point<1> &initial_p_unit,
782  InternalData &mdata) const;
783 
784 template<>
785 Point<1>
789  const Point<3> &p,
790  const Point<1> &initial_p_unit,
791  InternalData &mdata) const;
792 
793 
800 template <int dim, int spacedim=dim>
802 {
803  static MappingQ1<dim, spacedim> mapping;
804 };
805 
806 
809 /*----------------------------------------------------------------------*/
810 
811 #ifndef DOXYGEN
812 
813 template<int dim, int spacedim>
814 inline
815 double
816 MappingQ1<dim,spacedim>::InternalData::shape (const unsigned int qpoint,
817  const unsigned int shape_nr) const
818 {
819  Assert(qpoint*n_shape_functions + shape_nr < shape_values.size(),
820  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
821  shape_values.size()));
822  return shape_values [qpoint*n_shape_functions + shape_nr];
823 }
824 
825 
826 
827 template<int dim, int spacedim>
828 inline
829 double &
830 MappingQ1<dim,spacedim>::InternalData::shape (const unsigned int qpoint,
831  const unsigned int shape_nr)
832 {
833  Assert(qpoint*n_shape_functions + shape_nr < shape_values.size(),
834  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
835  shape_values.size()));
836  return shape_values [qpoint*n_shape_functions + shape_nr];
837 }
838 
839 
840 template<int dim, int spacedim>
841 inline
843 MappingQ1<dim,spacedim>::InternalData::derivative (const unsigned int qpoint,
844  const unsigned int shape_nr) const
845 {
846  Assert(qpoint*n_shape_functions + shape_nr < shape_derivatives.size(),
847  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
848  shape_derivatives.size()));
849  return shape_derivatives [qpoint*n_shape_functions + shape_nr];
850 }
851 
852 
853 
854 template<int dim, int spacedim>
855 inline
857 MappingQ1<dim,spacedim>::InternalData::derivative (const unsigned int qpoint,
858  const unsigned int shape_nr)
859 {
860  Assert(qpoint*n_shape_functions + shape_nr < shape_derivatives.size(),
861  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
862  shape_derivatives.size()));
863  return shape_derivatives [qpoint*n_shape_functions + shape_nr];
864 }
865 
866 
867 template <int dim, int spacedim>
868 inline
871  const unsigned int shape_nr) const
872 {
873  Assert(qpoint*n_shape_functions + shape_nr < shape_second_derivatives.size(),
874  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
875  shape_second_derivatives.size()));
876  return shape_second_derivatives [qpoint*n_shape_functions + shape_nr];
877 }
878 
879 
880 
881 template <int dim, int spacedim>
882 inline
885  const unsigned int shape_nr)
886 {
887  Assert(qpoint*n_shape_functions + shape_nr < shape_second_derivatives.size(),
888  ExcIndexRange(qpoint*n_shape_functions + shape_nr, 0,
889  shape_second_derivatives.size()));
890  return shape_second_derivatives [qpoint*n_shape_functions + shape_nr];
891 }
892 
893 
894 
895 template <int dim, int spacedim>
896 inline
897 bool
899 {
900  return true;
901 }
902 
903 #endif // DOXYGEN
904 
905 /* -------------- declaration of explicit specializations ------------- */
906 
907 
908 DEAL_II_NAMESPACE_CLOSE
909 
910 #endif
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_q1.h:312
virtual bool preserves_vertex_locations() const
Tensor< 2, dim > second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
MappingType
Definition: mapping.h:58
virtual void compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
std::vector< std::vector< Tensor< 1, dim > > > unit_tangentials
Definition: mapping_q1.h:336
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, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, typename std::vector< Tensor< 1, spacedim > > &boundary_form, typename std::vector< Point< spacedim > > &normal_vectors) const
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const
virtual UpdateFlags update_each(const UpdateFlags flags) const
virtual UpdateFlags update_once(const UpdateFlags flags) const
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_q1.h:348
Point< spacedim > transform_unit_to_real_cell_internal(const InternalData &mdata) const
Point< dim > transform_real_to_unit_cell_initial_guess(const std::vector< Point< spacedim > > &vertex, const Point< spacedim > &p) const
QProjector< dim >::DataSetDescriptor DataSetDescriptor
Definition: mapping_q1.h:391
std::vector< std::vector< Tensor< 1, spacedim > > > aux
Definition: mapping_q1.h:341
InternalData(const unsigned int n_shape_functions)
virtual Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
#define Assert(cond, exc)
Definition: exceptions.h:299
UpdateFlags
void compute_shapes(const std::vector< Point< dim > > &unit_points, InternalData &data) const
void transform_gradients(const VectorSlice< const std::vector< Tensor< rank, dim > > > input, VectorSlice< std::vector< Tensor< rank, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const
std::vector< double > shape_values
Definition: mapping_q1.h:259
virtual std::size_t memory_consumption() const
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_q1.h:355
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_q1.h:299
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
void compute_fill_face(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int npts, const DataSetDescriptor data_set, const std::vector< double > &weights, InternalData &mapping_data, 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
virtual void compute_shapes_virtual(const std::vector< Point< dim > > &unit_points, InternalData &data) const
void compute_face_data(const UpdateFlags flags, const Quadrature< dim > &quadrature, const unsigned int n_orig_q_points, InternalData &data) const
Tensor< 1, dim > derivative(const unsigned int qpoint, const unsigned int shape_nr) const
void transform_differential_forms(const VectorSlice< const std::vector< DerivativeForm< rank, dim, spacedim > > > input, VectorSlice< std::vector< DerivativeForm< rank, spacedim, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const
std::vector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_q1.h:268
std::vector< Tensor< 2, dim > > shape_second_derivatives
Definition: mapping_q1.h:278
void transform_fields(const VectorSlice< const std::vector< Tensor< rank, dim > > > input, VectorSlice< std::vector< Tensor< rank, spacedim > > > output, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const MappingType type) const
Point< dim > transform_real_to_unit_cell_internal(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit, InternalData &mdata) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const
virtual Mapping< dim, spacedim > * clone() const
double shape(const unsigned int qpoint, const unsigned int shape_nr) const
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename 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
unsigned int n_shape_functions
Definition: mapping_q1.h:379
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, typename std::vector< Point< spacedim > > &quadrature_points, std::vector< double > &JxW_values, typename std::vector< Tensor< 1, spacedim > > &boundary_form, typename std::vector< Point< spacedim > > &normal_vectors) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const
void compute_fill(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int npts, const DataSetDescriptor data_set, const CellSimilarity::Similarity cell_similarity, InternalData &data, std::vector< Point< spacedim > > &quadrature_points) const
void compute_data(const UpdateFlags flags, const Quadrature< dim > &quadrature, const unsigned int n_orig_q_points, InternalData &data) const
static const unsigned int n_shape_functions
Definition: mapping_q1.h:760