Reference documentation for deal.II version GIT 472adc65c1 2022-06-27 15:55:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
mapping_data_on_the_fly.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2022 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 
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 
34 
35 
37 
38 
39 namespace internal
40 {
41  namespace MatrixFreeFunctions
42  {
59  template <int dim, typename Number>
61  {
62  public:
66  MappingDataOnTheFly() = default;
67 
74  MappingDataOnTheFly(const Mapping<dim> & mapping,
75  const Quadrature<1> &quadrature,
76  const UpdateFlags update_flags);
77 
83  MappingDataOnTheFly(const Quadrature<1> &quadrature,
84  const UpdateFlags update_flags);
85 
89  void
90  reinit(typename ::Triangulation<dim>::cell_iterator cell);
91 
96  bool
97  is_initialized() const;
98 
102  typename ::Triangulation<dim>::cell_iterator
103  get_cell() const;
104 
111  const ::FEValues<dim> &
112  get_fe_values() const;
113 
121  get_data_storage() const;
122 
132 
136  const Quadrature<1> &
137  get_quadrature() const;
138 
139  private:
144  typename ::Triangulation<dim>::cell_iterator present_cell;
145 
151 
155  std::unique_ptr<::FEValues<dim>> fe_values;
156 
161 
167  };
168 
169 
170  /*-------------------------- Inline functions ---------------------------*/
171 
172  template <int dim, typename Number>
174  const Mapping<dim> & mapping,
175  const Quadrature<1> &quadrature,
176  const UpdateFlags update_flags)
177  : fe_values(std::make_unique<::FEValues<dim>>(
178  mapping,
179  fe_dummy,
180  Quadrature<dim>(quadrature),
181  MappingInfoStorage<dim, dim, Number>::compute_update_flags(
182  update_flags)))
183  , quadrature_1d(quadrature)
184  {
186  mapping_info_storage.descriptor[0].initialize(quadrature);
188  mapping_info_storage.JxW_values.resize(fe_values->n_quadrature_points);
189  mapping_info_storage.jacobians[0].resize(fe_values->n_quadrature_points);
190  if (update_flags & update_quadrature_points)
191  {
194  fe_values->n_quadrature_points);
195  }
196  if (fe_values->get_update_flags() & update_normal_vectors)
197  {
199  fe_values->n_quadrature_points);
201  fe_values->n_quadrature_points);
202  }
203  Assert(!(fe_values->get_update_flags() & update_jacobian_grads),
205  }
206 
207 
208 
209  template <int dim, typename Number>
211  const Quadrature<1> &quadrature,
212  const UpdateFlags update_flags)
213  : MappingDataOnTheFly(::StaticMappingQ1<dim, dim>::mapping,
214  quadrature,
215  update_flags)
216  {}
217 
218 
219 
220  template <int dim, typename Number>
221  inline void
223  typename ::Triangulation<dim>::cell_iterator cell)
224  {
225  if (present_cell == cell)
226  return;
227  present_cell = cell;
228  fe_values->reinit(present_cell);
229  for (unsigned int q = 0; q < fe_values->get_quadrature().size(); ++q)
230  {
231  if (fe_values->get_update_flags() & update_JxW_values)
232  mapping_info_storage.JxW_values[q] = fe_values->JxW(q);
233  if (fe_values->get_update_flags() & update_jacobians)
234  {
235  Tensor<2, dim> jac = fe_values->jacobian(q);
236  jac = invert(transpose(jac));
237  for (unsigned int d = 0; d < dim; ++d)
238  for (unsigned int e = 0; e < dim; ++e)
239  mapping_info_storage.jacobians[0][q][d][e] = jac[d][e];
240  }
241  if (fe_values->get_update_flags() & update_quadrature_points)
242  for (unsigned int d = 0; d < dim; ++d)
243  mapping_info_storage.quadrature_points[q][d] =
244  fe_values->quadrature_point(q)[d];
245  if (fe_values->get_update_flags() & update_normal_vectors)
246  {
247  for (unsigned int d = 0; d < dim; ++d)
248  mapping_info_storage.normal_vectors[q][d] =
249  fe_values->normal_vector(q)[d];
250  mapping_info_storage.normals_times_jacobians[0][q] =
251  mapping_info_storage.normal_vectors[q] *
252  mapping_info_storage.jacobians[0][q];
253  }
254  }
255  }
256 
257 
258 
259  template <int dim, typename Number>
260  inline bool
262  {
263  return present_cell !=
264  typename ::Triangulation<dim>::cell_iterator();
265  }
266 
267 
268 
269  template <int dim, typename Number>
270  inline typename ::Triangulation<dim>::cell_iterator
272  {
273  return fe_values->get_cell();
274  }
275 
276 
277 
278  template <int dim, typename Number>
279  inline const ::FEValues<dim> &
281  {
282  return *fe_values;
283  }
284 
285 
286 
287  template <int dim, typename Number>
290  {
291  return mapping_info_storage;
292  }
293 
294 
295 
296  template <int dim, typename Number>
299  {
300  return mapping_info_storage;
301  }
302 
303 
304 
305  template <int dim, typename Number>
306  inline const Quadrature<1> &
308  {
309  return quadrature_1d;
310  }
311 
312 
313  } // end of namespace MatrixFreeFunctions
314 } // end of namespace internal
315 
316 
318 
319 #endif
void resize(const size_type new_size)
Abstract base class for mapping classes.
Definition: mapping.h:311
const MappingInfoStorage< dim, dim, Number > & get_data_storage() const
void reinit(typename ::Triangulation< dim >::cell_iterator cell)
typename ::Triangulation< dim >::cell_iterator present_cell
MappingInfoStorage< dim, dim, Number > mapping_info_storage
typename ::Triangulation< dim >::cell_iterator get_cell() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
UpdateFlags
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_jacobians
Volume element.
@ update_quadrature_points
Transformed quadrature points.
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNotImplemented()
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
std::vector< QuadratureDescriptor > descriptor
AlignedVector< Tensor< 1, spacedim, Number > > normal_vectors
std::array< AlignedVector< Tensor< 1, spacedim, Number > >, 2 > normals_times_jacobians
AlignedVector< Point< spacedim, Number > > quadrature_points
std::array< AlignedVector< Tensor< 2, spacedim, Number > >, 2 > jacobians
constexpr SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)