Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
mapping_info.h
1 // ---------------------------------------------------------------------
2 // @f$Id: mapping_info.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2011 - 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 
18 #ifndef __deal2__matrix_free_mapping_info_h
19 #define __deal2__matrix_free_mapping_info_h
20 
21 
23 #include <deal.II/base/vectorization.h>
24 #include <deal.II/hp/q_collection.h>
25 #include <deal.II/fe/fe.h>
26 #include <deal.II/fe/mapping.h>
27 #include <deal.II/matrix_free/helper_functions.h>
28 
29 #include <memory>
30 
31 
33 
34 
35 namespace internal
36 {
37  namespace MatrixFreeFunctions
38  {
45  template <int dim, typename Number>
46  struct MappingInfo
47  {
52  static const std::size_t n_cell_type_bits = 2;
53 
58  static const unsigned int n_cell_types = 1U<<n_cell_type_bits;
59 
63  MappingInfo();
64 
73  void initialize (const ::Triangulation<dim> &tria,
74  const std::vector<std::pair<unsigned int,unsigned int> > &cells,
75  const std::vector<unsigned int> &active_fe_index,
76  const Mapping<dim> &mapping,
77  const std::vector<::hp::QCollection<1> > &quad,
78  const UpdateFlags update_flags);
79 
85  compute_update_flags (const UpdateFlags update_flags,
86  const std::vector<::hp::QCollection<1> > &quad) const;
87 
91  CellType get_cell_type (const unsigned int cell_chunk_no) const;
92 
96  unsigned int get_cell_data_index (const unsigned int cell_chunk_no) const;
97 
101  void clear ();
102 
106  std::size_t memory_consumption() const;
107 
112  template <typename STREAM>
113  void print_memory_consumption(STREAM &out,
114  const SizeInfo &size_info) const;
115 
124  std::vector<unsigned int> cell_type;
125 
140 
155 
162  {
167  std::vector<unsigned int> rowstart_jacobians;
168 
175 
181 
192 
201  AlignedVector<Tensor<1,(dim>1?dim*(dim-1)/2:1),
203 
210  std::vector<unsigned int> rowstart_q_points;
211 
217 
223 
230 
234  std::vector<unsigned int> n_q_points;
235 
242  std::vector<unsigned int> n_q_points_face;
243 
247  std::vector<AlignedVector<VectorizedArray<Number> > > quadrature_weights;
248 
254  std::vector<unsigned int> quad_index_conversion;
255 
262  unsigned int
263  quad_index_from_n_q_points (const unsigned int n_q_points) const;
264 
265 
270  template <typename STREAM>
271  void print_memory_consumption(STREAM &out,
272  const SizeInfo &size_info) const;
273 
277  std::size_t memory_consumption () const;
278  };
279 
284  std::vector<MappingInfoDependent> mapping_data_gen;
285 
290 
295 
300 
304  struct CellData
305  {
306  CellData (const double jac_size);
307  void resize (const unsigned int size);
308 
313  const double jac_size;
314  };
315 
319  void evaluate_on_cell (const ::Triangulation<dim> &tria,
320  const std::pair<unsigned int,unsigned int> *cells,
321  const unsigned int cell,
322  const unsigned int my_q,
323  CellType (&cell_t_prev)[VectorizedArray<Number>::n_array_elements],
324  CellType (&cell_t)[VectorizedArray<Number>::n_array_elements],
325  FEValues<dim,dim> &fe_values,
326  CellData &cell_data) const;
327  };
328 
329 
330 
331  /* ------------------- inline functions ----------------------------- */
332 
333  template <int dim, typename Number>
334  inline
335  unsigned int
337  quad_index_from_n_q_points (const unsigned int n_q_points) const
338  {
339  for (unsigned int i=0; i<quad_index_conversion.size(); ++i)
340  if (n_q_points == quad_index_conversion[i])
341  return i;
342  return 0;
343  }
344 
345 
346 
347  template <int dim, typename Number>
348  inline
349  CellType
350  MappingInfo<dim,Number>::get_cell_type (const unsigned int cell_no) const
351  {
352  AssertIndexRange (cell_no, cell_type.size());
353  CellType enum_cell_type = (CellType)(cell_type[cell_no] % n_cell_types);
354  Assert(enum_cell_type != undefined, ExcInternalError());
355  return enum_cell_type;
356  }
357 
358 
359 
360  template <int dim, typename Number>
361  inline
362  unsigned int
363  MappingInfo<dim,Number>::get_cell_data_index (const unsigned int cell_no) const
364  {
365  AssertIndexRange (cell_no, cell_type.size());
366  return cell_type[cell_no] >> n_cell_type_bits;
367  }
368 
369  } // end of namespace MatrixFreeFunctions
370 } // end of namespace internal
371 
372 DEAL_II_NAMESPACE_CLOSE
373 
374 #endif
AlignedVector< std::pair< Tensor< 1, dim, VectorizedArray< Number > >, VectorizedArray< Number > > > cartesian_data
Definition: mapping_info.h:139
AlignedVector< VectorizedArray< Number > > JxW_values
Definition: mapping_info.h:180
AlignedVector< std::pair< Tensor< 2, dim, VectorizedArray< Number > >, VectorizedArray< Number > > > affine_data
Definition: mapping_info.h:154
#define AssertIndexRange(index, range)
Definition: exceptions.h:888
std::vector< MappingInfoDependent > mapping_data_gen
Definition: mapping_info.h:284
CellType get_cell_type(const unsigned int cell_chunk_no) const
Definition: mapping_info.h:350
void print_memory_consumption(STREAM &out, const SizeInfo &size_info) const
UpdateFlags compute_update_flags(const UpdateFlags update_flags, const std::vector<::hp::QCollection< 1 > > &quad) const
std::vector< unsigned int > cell_type
Definition: mapping_info.h:124
#define Assert(cond, exc)
Definition: exceptions.h:299
UpdateFlags
AlignedVector< Tensor< 1,(dim >1?dim *(dim-1)/2:1), Tensor< 1, dim, VectorizedArray< Number > > > > jacobians_grad_upper
Definition: mapping_info.h:202
Definition: tensor.h:26
void print_memory_consumption(STREAM &out, const SizeInfo &size_info) const
unsigned int quad_index_from_n_q_points(const unsigned int n_q_points) const
Definition: mapping_info.h:337
unsigned int get_cell_data_index(const unsigned int cell_chunk_no) const
Definition: mapping_info.h:363
static const unsigned int n_cell_types
Definition: mapping_info.h:58
AlignedVector< Point< dim, VectorizedArray< Number > > > quadrature_points
Definition: mapping_info.h:216
void initialize(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< unsigned int > &active_fe_index, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 > > &quad, const UpdateFlags update_flags)
static const std::size_t n_cell_type_bits
Definition: mapping_info.h:52
::ExceptionBase & ExcInternalError()
std::vector< AlignedVector< VectorizedArray< Number > > > quadrature_weights
Definition: mapping_info.h:247
void evaluate_on_cell(const ::Triangulation< dim > &tria, const std::pair< unsigned int, unsigned int > *cells, const unsigned int cell, const unsigned int my_q, CellType(&cell_t_prev)[VectorizedArray< Number >::n_array_elements], CellType(&cell_t)[VectorizedArray< Number >::n_array_elements], FEValues< dim, dim > &fe_values, CellData &cell_data) const
AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > jacobians
Definition: mapping_info.h:174
AlignedVector< Tensor< 2, dim, VectorizedArray< Number > > > jacobians_grad_diag
Definition: mapping_info.h:191