Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mapping_data_on_the_fly.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_matrix_free_mapping_data_on_the_fly_h
18 #define dealii_matrix_free_mapping_data_on_the_fly_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #include <deal.II/base/aligned_vector.h>
24 #include <deal.II/base/exceptions.h>
25 #include <deal.II/base/subscriptor.h>
26 #include <deal.II/base/vectorization.h>
27 
28 #include <deal.II/fe/fe_nothing.h>
29 #include <deal.II/fe/fe_values.h>
30 #include <deal.II/fe/mapping_q1.h>
31 
32 #include <deal.II/matrix_free/mapping_info.h>
33 #include <deal.II/matrix_free/shape_info.h>
34 
35 
36 DEAL_II_NAMESPACE_OPEN
37 
38 
39 namespace internal
40 {
41  namespace MatrixFreeFunctions
42  {
61  template <int dim, typename Number = double>
63  {
64  public:
71  MappingDataOnTheFly(const Mapping<dim> & mapping,
72  const Quadrature<1> &quadrature,
73  const UpdateFlags update_flags);
74 
80  MappingDataOnTheFly(const Quadrature<1> &quadrature,
81  const UpdateFlags update_flags);
82 
86  void
87  reinit(typename ::Triangulation<dim>::cell_iterator cell);
88 
93  bool
94  is_initialized() const;
95 
99  typename ::Triangulation<dim>::cell_iterator
100  get_cell() const;
101 
108  const ::FEValues<dim> &
109  get_fe_values() const;
110 
118  get_data_storage() const;
119 
123  const Quadrature<1> &
124  get_quadrature() const;
125 
126  private:
131  typename ::Triangulation<dim>::cell_iterator present_cell;
132 
138 
143 
148 
154  };
155 
156 
157  /*-------------------------- Inline functions ---------------------------*/
158 
159  template <int dim, typename Number>
161  const Mapping<dim> & mapping,
162  const Quadrature<1> &quadrature,
163  const UpdateFlags update_flags)
164  : fe_values(mapping,
165  fe_dummy,
166  Quadrature<dim>(quadrature),
167  MappingInfo<dim, Number>::compute_update_flags(update_flags))
168  , quadrature_1d(quadrature)
169  {
171  mapping_info_storage.descriptor[0].initialize(quadrature);
175  if (update_flags & update_quadrature_points)
176  {
180  }
182  {
187  }
190  }
191 
192 
193 
194  template <int dim, typename Number>
196  const Quadrature<1> &quadrature,
197  const UpdateFlags update_flags)
198  : MappingDataOnTheFly(::StaticMappingQ1<dim, dim>::mapping,
199  quadrature,
200  update_flags)
201  {}
202 
203 
204 
205  template <int dim, typename Number>
206  inline void
208  typename ::Triangulation<dim>::cell_iterator cell)
209  {
210  if (present_cell == cell)
211  return;
212  present_cell = cell;
213  fe_values.reinit(present_cell);
214  for (unsigned int q = 0; q < fe_values.get_quadrature().size(); ++q)
215  {
216  if (fe_values.get_update_flags() & update_JxW_values)
217  mapping_info_storage.JxW_values[q] = fe_values.JxW(q);
218  if (fe_values.get_update_flags() & update_jacobians)
219  {
220  Tensor<2, dim> jac = fe_values.jacobian(q);
221  jac = invert(transpose(jac));
222  for (unsigned int d = 0; d < dim; ++d)
223  for (unsigned int e = 0; e < dim; ++e)
224  mapping_info_storage.jacobians[0][q][d][e] = jac[d][e];
225  }
226  if (fe_values.get_update_flags() & update_quadrature_points)
227  for (unsigned int d = 0; d < dim; ++d)
228  mapping_info_storage.quadrature_points[q][d] =
229  fe_values.quadrature_point(q)[d];
230  if (fe_values.get_update_flags() & update_normal_vectors)
231  {
232  for (unsigned int d = 0; d < dim; ++d)
233  mapping_info_storage.normal_vectors[q][d] =
234  fe_values.normal_vector(q)[d];
235  mapping_info_storage.normals_times_jacobians[0][q] =
236  mapping_info_storage.normal_vectors[q] *
237  mapping_info_storage.jacobians[0][q];
238  }
239  }
240  }
241 
242 
243 
244  template <int dim, typename Number>
245  inline bool
247  {
248  return present_cell !=
249  typename ::Triangulation<dim>::cell_iterator();
250  }
251 
252 
253 
254  template <int dim, typename Number>
255  inline typename ::Triangulation<dim>::cell_iterator
257  {
258  return fe_values.get_cell();
259  }
260 
261 
262 
263  template <int dim, typename Number>
264  inline const ::FEValues<dim> &
266  {
267  return fe_values;
268  }
269 
270 
271 
272  template <int dim, typename Number>
275  {
276  return mapping_info_storage;
277  }
278 
279 
280 
281  template <int dim, typename Number>
282  inline const Quadrature<1> &
284  {
285  return quadrature_1d;
286  }
287 
288 
289  } // end of namespace MatrixFreeFunctions
290 } // end of namespace internal
291 
292 
293 DEAL_II_NAMESPACE_CLOSE
294 
295 #endif
Transformed quadrature weights.
AlignedVector< Tensor< 1, spacedim, VectorizedArray< Number > > > normals_times_jacobians[2]
Definition: mapping_info.h:241
typename ::Triangulation< dim >::cell_iterator present_cell
AlignedVector< unsigned int > data_index_offsets
Definition: mapping_info.h:181
AlignedVector< unsigned int > quadrature_point_offsets
Definition: mapping_info.h:249
Volume element.
const MappingInfoStorage< dim, dim, Number > & get_data_storage() const
UpdateFlags get_update_flags() const
Transformed quadrature points.
void resize(const size_type size_in)
std::vector< QuadratureDescriptor > descriptor
Definition: mapping_info.h:173
SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)
typename ::Triangulation< dim >::cell_iterator get_cell() const
#define Assert(cond, exc)
Definition: exceptions.h:1407
UpdateFlags
SymmetricTensor< rank_, dim, Number > transpose(const SymmetricTensor< rank_, dim, Number > &t)
AlignedVector< Point< spacedim, VectorizedArray< Number > > > quadrature_points
Definition: mapping_info.h:258
Gradient of volume element.
const unsigned int n_quadrature_points
Definition: fe_values.h:2099
MappingInfoStorage< dim, dim, Number > mapping_info_storage
AlignedVector< Tensor< 1, spacedim, VectorizedArray< Number > > > normal_vectors
Definition: mapping_info.h:198
Normal vectors.
static ::ExceptionBase & ExcNotImplemented()
AlignedVector< Tensor< 2, spacedim, VectorizedArray< Number > > > jacobians[2]
Definition: mapping_info.h:211
AlignedVector< VectorizedArray< Number > > JxW_values
Definition: mapping_info.h:190
void reinit(typename ::Triangulation< dim >::cell_iterator cell)
MappingDataOnTheFly(const Mapping< dim > &mapping, const Quadrature< 1 > &quadrature, const UpdateFlags update_flags)