Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
mapping_cartesian.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping_cartesian.h 30036 2013-07-18 16:55:32Z maier @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_cartesian_h
18 #define __deal2__mapping_cartesian_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/table.h>
23 #include <cmath>
24 #include <deal.II/fe/mapping.h>
25 
27 
30 
46 template <int dim, int spacedim=dim>
47 class MappingCartesian : public Mapping<dim,spacedim>
48 {
49 public:
50  virtual
52  get_data (const UpdateFlags,
53  const Quadrature<dim> &quadrature) const;
54 
55  virtual
57  get_face_data (const UpdateFlags flags,
58  const Quadrature<dim-1>& quadrature) const;
59 
60  virtual
62  get_subface_data (const UpdateFlags flags,
63  const Quadrature<dim-1>& quadrature) const;
64 
65  virtual void
66  fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
67  const Quadrature<dim> &quadrature,
68  typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
69  std::vector<Point<spacedim> > &quadrature_points,
70  std::vector<double> &JxW_values,
71  std::vector<DerivativeForm<1,dim,spacedim> > &jacobians,
72  std::vector<DerivativeForm<2,dim,spacedim> > &jacobian_grads,
73  std::vector<DerivativeForm<1,spacedim,dim> > &inverse_jacobians,
74  std::vector<Point<spacedim> > &,
75  CellSimilarity::Similarity &cell_similarity) const ;
76 
77 
78  virtual void
79  fill_fe_face_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
80  const unsigned int face_no,
81  const Quadrature<dim-1>& quadrature,
82  typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
83  std::vector<Point<dim> > &quadrature_points,
84  std::vector<double> &JxW_values,
85  std::vector<Tensor<1,dim> > &boundary_form,
86  std::vector<Point<spacedim> > &normal_vectors) const ;
87  virtual void
88  fill_fe_subface_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
89  const unsigned int face_no,
90  const unsigned int sub_no,
91  const Quadrature<dim-1>& quadrature,
92  typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
93  std::vector<Point<dim> > &quadrature_points,
94  std::vector<double> &JxW_values,
95  std::vector<Tensor<1,dim> > &boundary_form,
96  std::vector<Point<spacedim> > &normal_vectors) const ;
97 
98  virtual void
99  transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
100  VectorSlice<std::vector<Tensor<1,spacedim> > > output,
101  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
102  const MappingType type) const;
103 
104  virtual void
105  transform (const VectorSlice<const std::vector<DerivativeForm<1, dim,spacedim> > > input,
106  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
107  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
108  const MappingType type) const;
109 
110 
111  virtual
112  void
113  transform (const VectorSlice<const std::vector<Tensor<2, dim> > > input,
114  VectorSlice<std::vector<Tensor<2,spacedim> > > output,
115  const typename Mapping<dim,spacedim>::InternalDataBase &internal,
116  const MappingType type) const;
117 
118  virtual Point<spacedim>
120  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
121  const Point<dim> &p) const;
122 
133  virtual Point<dim>
135  const typename Triangulation<dim,spacedim>::cell_iterator &cell,
136  const Point<spacedim> &p) const;
137 
138 
144  virtual
145  Mapping<dim, spacedim> *clone () const;
146 
152  bool preserves_vertex_locations () const;
153 
154 protected:
159  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
160  {
161  public:
165  InternalData (const Quadrature<dim> &quadrature);
166 
173  virtual std::size_t memory_consumption () const;
174 
182 
187 
193  std::vector<Point<dim> > quadrature_points;
194  };
195 
200  void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
201  const unsigned int face_no,
202  const unsigned int sub_no,
203  const CellSimilarity::Similarity cell_similarity,
204  InternalData &data,
205  std::vector<Point<dim> > &quadrature_points,
206  std::vector<Point<dim> > &normal_vectors) const;
207 
208 private:
209  virtual UpdateFlags update_once (const UpdateFlags) const;
210  virtual UpdateFlags update_each (const UpdateFlags) const;
211 
218 };
219 
222 /* -------------- declaration of explicit specializations ------------- */
223 
224 #ifndef DOXYGEN
225 
226 template <int dim, int spacedim>
227 inline
228 bool
230 {
231  return true;
232 }
233 
234 #endif // DOXYGEN
235 
236 DEAL_II_NAMESPACE_CLOSE
237 
238 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
void compute_fill(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const CellSimilarity::Similarity cell_similarity, InternalData &data, std::vector< Point< dim > > &quadrature_points, std::vector< Point< dim > > &normal_vectors) const
MappingType
Definition: mapping.h:58
virtual UpdateFlags update_once(const UpdateFlags) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
std::vector< Point< dim > > quadrature_points
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const
UpdateFlags
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const
bool preserves_vertex_locations() const
static const unsigned int invalid_face_number
virtual Mapping< dim, spacedim > * clone() const
InternalData(const Quadrature< dim > &quadrature)
virtual UpdateFlags update_each(const UpdateFlags) const
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const
virtual std::size_t memory_consumption() const