Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
mapping_q.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping_q.h 30450 2013-08-23 15:48:29Z kronbichler @f$
3 //
4 // Copyright (C) 2001 - 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_q_h
18 #define __deal2__mapping_q_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/table.h>
23 #include <deal.II/fe/mapping_q1.h>
24 #include <deal.II/fe/fe_q.h>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 template <int dim, typename POLY> class TensorProductPolynomials;
29 
30 
33 
55 template <int dim, int spacedim=dim>
56 class MappingQ : public MappingQ1<dim,spacedim>
57 {
58 public:
72  MappingQ (const unsigned int p,
73  const bool use_mapping_q_on_all_cells = false);
74 
80  MappingQ (const MappingQ<dim,spacedim> &mapping);
81 
85  virtual ~MappingQ ();
86 
91  virtual Point<spacedim>
93  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
94  const Point<dim> &p) const;
95 
117  virtual Point<dim>
119  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
120  const Point<spacedim> &p) const;
121 
122  virtual void
123  transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
124  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
125  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
126  const MappingType type) const;
127 
128  virtual void
129  transform (const VectorSlice<const std::vector<DerivativeForm<1, dim, spacedim> > > input,
130  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
131  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
132  const MappingType type) const;
133 
134  virtual
135  void
136  transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
137  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
138  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
139  const MappingType type) const;
140 
145  unsigned int get_degree () const;
146 
151  virtual
152  Mapping<dim,spacedim> *clone () const;
153 
157  class InternalData : public MappingQ1<dim,spacedim>::InternalData
158  {
159  public:
163  InternalData (const unsigned int n_shape_functions);
164 
165 
169  virtual std::size_t memory_consumption () const;
170 
177  std::vector<std::vector<Point<dim> > > unit_normals;
178 
186 
191  };
192 
193 protected:
197  virtual void
199  const Quadrature<dim> &quadrature,
200  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
201  typename std::vector<Point<spacedim> > &quadrature_points,
202  std::vector<double> &JxW_values,
203  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
204  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
205  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
206  std::vector<Point<spacedim> > &cell_normal_vectors,
207  CellSimilarity::Similarity &cell_similarity) const ;
208 
212  virtual void
214  const unsigned int face_no,
215  const Quadrature<dim-1>& quadrature,
216  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
217  typename std::vector<Point<spacedim> > &quadrature_points,
218  std::vector<double> &JxW_values,
219  typename std::vector<Tensor<1,spacedim> > &exterior_form,
220  typename std::vector<Point<spacedim> > &normal_vectors) const ;
221 
225  virtual void
227  const unsigned int face_no,
228  const unsigned int sub_no,
229  const Quadrature<dim-1>& quadrature,
230  typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
231  typename std::vector<Point<spacedim> > &quadrature_points,
232  std::vector<double> &JxW_values,
233  typename std::vector<Tensor<1,spacedim> > &exterior_form,
234  typename std::vector<Point<spacedim> > &normal_vectors) const ;
235 
248  virtual void
250  std::vector<Point<spacedim> > &a) const;
251 
264  virtual void
266  std::vector<Point<spacedim> > &a) const;
267 
268 private:
269 
270  virtual
272  get_data (const UpdateFlags,
273  const Quadrature<dim> &quadrature) const;
274 
275  virtual
277  get_face_data (const UpdateFlags flags,
278  const Quadrature<dim-1>& quadrature) const;
279 
280  virtual
282  get_subface_data (const UpdateFlags flags,
283  const Quadrature<dim-1>& quadrature) const;
284 
288  virtual void
289  compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
290  typename MappingQ1<dim,spacedim>::InternalData &data) const;
291 
303  void
305 
315  void set_laplace_on_hex_vector(Table<2,double> &lohvs) const;
316 
326  void compute_laplace_vector(Table<2,double> &lvs) const;
327 
338  void apply_laplace_vector(const Table<2,double> &lvs,
339  std::vector<Point<spacedim> > &a) const;
340 
344  virtual void compute_mapping_support_points(
345  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
346  std::vector<Point<spacedim> > &a) const;
347 
356  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
357  std::vector<Point<spacedim> > &a) const;
358 
372 
381 
385  DeclException1 (ExcLaplaceVectorNotSet,
386  int,
387  << "laplace_vector not set for degree=" << arg1 << ".");
388 
393  const unsigned int degree;
394 
398  const unsigned int n_inner;
399 
403  const unsigned int n_outer;
404 
410 
414  const unsigned int n_shape_functions;
415 
420  const std::vector<unsigned int> renumber;
421 
427 
435  const FE_Q<dim> feq;
436 
440  template <int,int> friend class MappingQ;
441 };
442 
445 /* -------------- declaration of explicit specializations ------------- */
446 
447 #ifndef DOXYGEN
448 
449 template<> MappingQ<1>::MappingQ (const unsigned int,
450  const bool);
451 template<> MappingQ<1>::~MappingQ ();
452 
453 template<>
454 void MappingQ<1>::compute_shapes_virtual (const std::vector<Point<1> > &unit_points,
455  MappingQ1<1>::InternalData &data) const;
456 template <>
458 
459 template <>
461 
462 template <>
464 template <>
466  std::vector<Point<1> > &) const;
467 template <>
469  std::vector<Point<2> > &) const;
470 template <>
472  std::vector<Point<3> > &) const;
473 
474 template<>
476  std::vector<Point<3> > &a) const;
477 
478 #endif // DOXYGEN
479 
480 DEAL_II_NAMESPACE_CLOSE
481 
482 #endif
void compute_laplace_vector(Table< 2, double > &lvs) 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, 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 > > &exterior_form, typename std::vector< Point< spacedim > > &normal_vectors) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
unsigned int get_degree() 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
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const
InternalData(const unsigned int n_shape_functions)
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
const unsigned int n_outer
Definition: mapping_q.h:403
MappingType
Definition: mapping.h:58
bool use_mapping_q1_on_current_cell
Definition: mapping_q.h:185
virtual void compute_shapes_virtual(const std::vector< Point< dim > > &unit_points, typename MappingQ1< dim, spacedim >::InternalData &data) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const
friend class MappingQ
Definition: mapping_q.h:440
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 > > &exterior_form, typename std::vector< Point< spacedim > > &normal_vectors) const
const bool use_mapping_q_on_all_cells
Definition: mapping_q.h:426
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const
void compute_support_points_laplace(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
const FE_Q< dim > feq
Definition: mapping_q.h:435
const unsigned int degree
Definition: mapping_q.h:393
void set_laplace_on_quad_vector(Table< 2, double > &loqvs) const
const unsigned int n_inner
Definition: mapping_q.h:398
Table< 2, double > laplace_on_hex_vector
Definition: mapping_q.h:380
DeclException1(ExcLaplaceVectorNotSet, int,<< "laplace_vector not set for degree="<< arg1<< ".")
UpdateFlags
const std::vector< unsigned int > renumber
Definition: mapping_q.h:420
void set_laplace_on_hex_vector(Table< 2, double > &lohvs) const
const unsigned int n_shape_functions
Definition: mapping_q.h:414
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
virtual std::size_t memory_consumption() const
Table< 2, double > laplace_on_quad_vector
Definition: mapping_q.h:371
const TensorProductPolynomials< dim > * tensor_pols
Definition: mapping_q.h:409
virtual void compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim > > &a) const
std::vector< std::vector< Point< dim > > > unit_normals
Definition: mapping_q.h:177
virtual ~MappingQ()
unsigned int n_shape_functions
Definition: mapping_q1.h:379
virtual Mapping< dim, spacedim > * clone() const
MappingQ1< dim, spacedim >::InternalData mapping_q1_data
Definition: mapping_q.h:190
void apply_laplace_vector(const Table< 2, double > &lvs, std::vector< Point< spacedim > > &a) const