Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
scratch_data.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2019 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 #include <deal.II/meshworker/scratch_data.h>
17 
18 DEAL_II_NAMESPACE_OPEN
19 
20 namespace MeshWorker
21 {
22  template <int dim, int spacedim>
24  const Mapping<dim, spacedim> & mapping,
26  const Quadrature<dim> & quadrature,
27  const UpdateFlags & update_flags,
28  const Quadrature<dim - 1> & face_quadrature,
29  const UpdateFlags & face_update_flags)
30  : mapping(&mapping)
31  , fe(&fe)
32  , cell_quadrature(quadrature)
33  , face_quadrature(face_quadrature)
34  , cell_update_flags(update_flags)
35  , neighbor_cell_update_flags(update_flags)
36  , face_update_flags(face_update_flags)
37  , neighbor_face_update_flags(face_update_flags)
38  , local_dof_indices(fe.dofs_per_cell)
39  , neighbor_dof_indices(fe.dofs_per_cell)
40  {}
41 
42 
43 
44  template <int dim, int spacedim>
48  const Quadrature<dim> & quadrature,
49  const UpdateFlags & update_flags,
50  const UpdateFlags & neighbor_update_flags,
54  : mapping(&mapping)
55  , fe(&fe)
56  , cell_quadrature(quadrature)
57  , face_quadrature(face_quadrature)
58  , cell_update_flags(update_flags)
59  , neighbor_cell_update_flags(neighbor_update_flags)
60  , face_update_flags(face_update_flags)
61  , neighbor_face_update_flags(neighbor_face_update_flags)
62  , local_dof_indices(fe.dofs_per_cell)
63  , neighbor_dof_indices(fe.dofs_per_cell)
64  {}
65 
66 
67 
68  template <int dim, int spacedim>
71  const Quadrature<dim> & quadrature,
72  const UpdateFlags & update_flags,
75  : ScratchData(StaticMappingQ1<dim, spacedim>::mapping,
76  fe,
77  quadrature,
78  update_flags,
79  face_quadrature,
80  face_update_flags)
81  {}
82 
83 
84 
85  template <int dim, int spacedim>
88  const Quadrature<dim> & quadrature,
89  const UpdateFlags & update_flags,
90  const UpdateFlags & neighbor_update_flags,
94  : ScratchData(StaticMappingQ1<dim, spacedim>::mapping,
95  fe,
96  quadrature,
97  update_flags,
98  neighbor_update_flags,
99  face_quadrature,
100  face_update_flags,
101  neighbor_face_update_flags)
102  {}
103 
104 
105 
106  template <int dim, int spacedim>
108  const ScratchData<dim, spacedim> &scratch)
109  : mapping(scratch.mapping)
110  , fe(scratch.fe)
121  {}
122 
123 
124 
125  template <int dim, int spacedim>
129  {
130  if (!fe_values)
131  fe_values = std_cxx14::make_unique<FEValues<dim, spacedim>>(
133 
134  fe_values->reinit(cell);
135  cell->get_dof_indices(local_dof_indices);
137  return *fe_values;
138  }
139 
140 
141 
142  template <int dim, int spacedim>
146  const unsigned int face_no)
147  {
148  if (!fe_face_values)
149  fe_face_values = std_cxx14::make_unique<FEFaceValues<dim, spacedim>>(
151 
152  fe_face_values->reinit(cell, face_no);
153  cell->get_dof_indices(local_dof_indices);
155  return *fe_face_values;
156  }
157 
158 
159 
160  template <int dim, int spacedim>
164  const unsigned int face_no,
165  const unsigned int subface_no)
166  {
167  if (subface_no != numbers::invalid_unsigned_int)
168  {
169  if (!fe_subface_values)
171  std_cxx14::make_unique<FESubfaceValues<dim, spacedim>>(
173  fe_subface_values->reinit(cell, face_no, subface_no);
174  cell->get_dof_indices(local_dof_indices);
175 
177  return *fe_subface_values;
178  }
179  else
180  return reinit(cell, face_no);
181  }
182 
183 
184 
185  template <int dim, int spacedim>
189  {
190  if (!neighbor_fe_values)
191  neighbor_fe_values = std_cxx14::make_unique<FEValues<dim, spacedim>>(
193 
194  neighbor_fe_values->reinit(cell);
195  cell->get_dof_indices(neighbor_dof_indices);
197  return *neighbor_fe_values;
198  }
199 
200 
201 
202  template <int dim, int spacedim>
206  const unsigned int face_no)
207  {
210  std_cxx14::make_unique<FEFaceValues<dim, spacedim>>(
212  neighbor_fe_face_values->reinit(cell, face_no);
213  cell->get_dof_indices(neighbor_dof_indices);
215  return *neighbor_fe_face_values;
216  }
217 
218 
219 
220  template <int dim, int spacedim>
224  const unsigned int face_no,
225  const unsigned int subface_no)
226  {
227  if (subface_no != numbers::invalid_unsigned_int)
228  {
231  std_cxx14::make_unique<FESubfaceValues<dim, spacedim>>(
233  neighbor_fe_subface_values->reinit(cell, face_no, subface_no);
234  cell->get_dof_indices(neighbor_dof_indices);
237  }
238  else
239  return reinit_neighbor(cell, face_no);
240  }
241 
242 
243 
244  template <int dim, int spacedim>
247  {
248  Assert(current_fe_values != nullptr,
249  ExcMessage("You have to initialize the cache using one of the "
250  "reinit functions first!"));
251  return *current_fe_values;
252  }
253 
254 
255 
256  template <int dim, int spacedim>
259  {
261  ExcMessage("You have to initialize the cache using one of the "
262  "reinit functions first!"));
264  }
265 
266 
267 
268  template <int dim, int spacedim>
269  const std::vector<Point<spacedim>> &
271  {
272  return get_current_fe_values().get_quadrature_points();
273  }
274 
275 
276 
277  template <int dim, int spacedim>
278  const std::vector<double> &
280  {
281  return get_current_fe_values().get_JxW_values();
282  }
283 
284 
285 
286  template <int dim, int spacedim>
287  const std::vector<double> &
289  {
290  return get_current_neighbor_fe_values().get_JxW_values();
291  }
292 
293 
294 
295  template <int dim, int spacedim>
296  const std::vector<Tensor<1, spacedim>> &
298  {
299  return get_current_fe_values().get_normal_vectors();
300  }
301 
302 
303 
304  template <int dim, int spacedim>
305  const std::vector<Tensor<1, spacedim>> &
307  {
308  return get_current_neighbor_fe_values().get_normal_vectors();
309  }
310 
311 
312 
313  template <int dim, int spacedim>
314  const std::vector<types::global_dof_index> &
316  {
317  return local_dof_indices;
318  }
319 
320 
321 
322  template <int dim, int spacedim>
323  const std::vector<types::global_dof_index> &
325  {
326  return neighbor_dof_indices;
327  }
328 
329 
330 
331  template <int dim, int spacedim>
334  {
335  return user_data_storage;
336  }
337 
338 
339 
340  template <int dim, int spacedim>
341  const Mapping<dim, spacedim> &
343  {
344  return *mapping;
345  }
346 
347 } // namespace MeshWorker
348 DEAL_II_NAMESPACE_CLOSE
349 
350 // Explicit instantiations
351 DEAL_II_NAMESPACE_OPEN
352 namespace MeshWorker
353 {
354 #include "scratch_data.inst"
355 }
356 DEAL_II_NAMESPACE_CLOSE
static const unsigned int invalid_unsigned_int
Definition: types.h:187
GeneralDataStorage internal_data_storage
Definition: scratch_data.h:885
SmartPointer< FEValuesBase< dim, spacedim > > current_fe_values
Definition: scratch_data.h:891
std::unique_ptr< FEFaceValues< dim, spacedim > > fe_face_values
Definition: scratch_data.h:845
Quadrature< dim > cell_quadrature
Definition: scratch_data.h:807
SmartPointer< const Mapping< dim, spacedim > > mapping
Definition: scratch_data.h:795
const Mapping< dim, spacedim > & get_mapping() const
std::unique_ptr< FEValues< dim, spacedim > > neighbor_fe_values
Definition: scratch_data.h:855
UpdateFlags neighbor_face_update_flags
Definition: scratch_data.h:835
const std::vector< types::global_dof_index > & get_local_dof_indices() const
const FEValuesBase< dim, spacedim > & get_current_neighbor_fe_values() const
const FEValues< dim, spacedim > & reinit(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
UpdateFlags neighbor_cell_update_flags
Definition: scratch_data.h:823
const std::vector< Tensor< 1, spacedim > > & get_neighbor_normal_vectors()
std::vector< types::global_dof_index > neighbor_dof_indices
Definition: scratch_data.h:875
ScratchData(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags &update_flags, const Quadrature< dim - 1 > &face_quadrature=Quadrature< dim - 1 >(), const UpdateFlags &face_update_flags=update_default)
Definition: scratch_data.cc:23
GeneralDataStorage & get_general_data_storage()
std::unique_ptr< FESubfaceValues< dim, spacedim > > fe_subface_values
Definition: scratch_data.h:850
std::unique_ptr< FESubfaceValues< dim, spacedim > > neighbor_fe_subface_values
Definition: scratch_data.h:865
static ::ExceptionBase & ExcMessage(std::string arg1)
const FEValues< dim, spacedim > & reinit_neighbor(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
#define Assert(cond, exc)
Definition: exceptions.h:1407
UpdateFlags
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Abstract base class for mapping classes.
Definition: mapping.h:302
const std::vector< Point< spacedim > > & get_quadrature_points() const
SmartPointer< FEValuesBase< dim, spacedim > > current_neighbor_fe_values
Definition: scratch_data.h:897
std::unique_ptr< FEFaceValues< dim, spacedim > > neighbor_fe_face_values
Definition: scratch_data.h:860
UpdateFlags face_update_flags
Definition: scratch_data.h:829
GeneralDataStorage user_data_storage
Definition: scratch_data.h:880
UpdateFlags cell_update_flags
Definition: scratch_data.h:818
Quadrature< dim - 1 > face_quadrature
Definition: scratch_data.h:813
std::unique_ptr< FEValues< dim, spacedim > > fe_values
Definition: scratch_data.h:840
const std::vector< types::global_dof_index > & get_neighbor_dof_indices() const
const FEValuesBase< dim, spacedim > & get_current_fe_values() const
Definition: fe.h:38
const std::vector< double > & get_neighbor_JxW_values() const
const std::vector< double > & get_JxW_values() const
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:328
std::vector< types::global_dof_index > local_dof_indices
Definition: scratch_data.h:870
SmartPointer< const FiniteElement< dim, spacedim > > fe
Definition: scratch_data.h:801