Reference documentation for deal.II version Git 92c1895bf8 2020-06-04 17:48:28 -0400
\(\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\}}\)
scratch_data.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2019 - 2020 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 
17 
18 #include <memory>
19 
21 
22 namespace MeshWorker
23 {
24  template <int dim, int spacedim>
26  const Mapping<dim, spacedim> & mapping,
28  const Quadrature<dim> & quadrature,
29  const UpdateFlags & update_flags,
30  const Quadrature<dim - 1> & face_quadrature,
31  const UpdateFlags & face_update_flags)
32  : mapping(&mapping)
33  , fe(&fe)
34  , cell_quadrature(quadrature)
35  , face_quadrature(face_quadrature)
36  , cell_update_flags(update_flags)
37  , neighbor_cell_update_flags(update_flags)
38  , face_update_flags(face_update_flags)
39  , neighbor_face_update_flags(face_update_flags)
40  , local_dof_indices(fe.dofs_per_cell)
41  , neighbor_dof_indices(fe.dofs_per_cell)
42  {}
43 
44 
45 
46  template <int dim, int spacedim>
50  const Quadrature<dim> & quadrature,
51  const UpdateFlags & update_flags,
52  const UpdateFlags & neighbor_update_flags,
56  : mapping(&mapping)
57  , fe(&fe)
58  , cell_quadrature(quadrature)
59  , face_quadrature(face_quadrature)
60  , cell_update_flags(update_flags)
61  , neighbor_cell_update_flags(neighbor_update_flags)
62  , face_update_flags(face_update_flags)
63  , neighbor_face_update_flags(neighbor_face_update_flags)
64  , local_dof_indices(fe.dofs_per_cell)
65  , neighbor_dof_indices(fe.dofs_per_cell)
66  {}
67 
68 
69 
70  template <int dim, int spacedim>
73  const Quadrature<dim> & quadrature,
74  const UpdateFlags & update_flags,
77  : ScratchData(StaticMappingQ1<dim, spacedim>::mapping,
78  fe,
79  quadrature,
80  update_flags,
81  face_quadrature,
82  face_update_flags)
83  {}
84 
85 
86 
87  template <int dim, int spacedim>
90  const Quadrature<dim> & quadrature,
91  const UpdateFlags & update_flags,
92  const UpdateFlags & neighbor_update_flags,
96  : ScratchData(StaticMappingQ1<dim, spacedim>::mapping,
97  fe,
98  quadrature,
99  update_flags,
100  neighbor_update_flags,
101  face_quadrature,
102  face_update_flags,
103  neighbor_face_update_flags)
104  {}
105 
106 
107 
108  template <int dim, int spacedim>
110  const ScratchData<dim, spacedim> &scratch)
111  : mapping(scratch.mapping)
112  , fe(scratch.fe)
123  {}
124 
125 
126 
127  template <int dim, int spacedim>
131  {
132  if (!fe_values)
133  fe_values = std::make_unique<FEValues<dim, spacedim>>(*mapping,
134  *fe,
137 
138  fe_values->reinit(cell);
139  cell->get_dof_indices(local_dof_indices);
141  return *fe_values;
142  }
143 
144 
145 
146  template <int dim, int spacedim>
150  const unsigned int face_no)
151  {
152  if (!fe_face_values)
153  fe_face_values = std::make_unique<FEFaceValues<dim, spacedim>>(
155 
156  fe_face_values->reinit(cell, face_no);
157  cell->get_dof_indices(local_dof_indices);
159  return *fe_face_values;
160  }
161 
162 
163 
164  template <int dim, int spacedim>
168  const unsigned int face_no,
169  const unsigned int subface_no)
170  {
171  if (subface_no != numbers::invalid_unsigned_int)
172  {
173  if (!fe_subface_values)
174  fe_subface_values = std::make_unique<FESubfaceValues<dim, spacedim>>(
176  fe_subface_values->reinit(cell, face_no, subface_no);
177  cell->get_dof_indices(local_dof_indices);
178 
180  return *fe_subface_values;
181  }
182  else
183  return reinit(cell, face_no);
184  }
185 
186 
187 
188  template <int dim, int spacedim>
192  {
193  if (!neighbor_fe_values)
194  neighbor_fe_values = std::make_unique<FEValues<dim, spacedim>>(
196 
197  neighbor_fe_values->reinit(cell);
198  cell->get_dof_indices(neighbor_dof_indices);
200  return *neighbor_fe_values;
201  }
202 
203 
204 
205  template <int dim, int spacedim>
209  const unsigned int face_no)
210  {
212  neighbor_fe_face_values = std::make_unique<FEFaceValues<dim, spacedim>>(
214  neighbor_fe_face_values->reinit(cell, face_no);
215  cell->get_dof_indices(neighbor_dof_indices);
217  return *neighbor_fe_face_values;
218  }
219 
220 
221 
222  template <int dim, int spacedim>
226  const unsigned int face_no,
227  const unsigned int subface_no)
228  {
229  if (subface_no != numbers::invalid_unsigned_int)
230  {
233  std::make_unique<FESubfaceValues<dim, spacedim>>(
235  neighbor_fe_subface_values->reinit(cell, face_no, subface_no);
236  cell->get_dof_indices(neighbor_dof_indices);
239  }
240  else
241  return reinit_neighbor(cell, face_no);
242  }
243 
244 
245 
246  template <int dim, int spacedim>
249  {
250  Assert(current_fe_values != nullptr,
251  ExcMessage("You have to initialize the cache using one of the "
252  "reinit functions first!"));
253  return *current_fe_values;
254  }
255 
256 
257 
258  template <int dim, int spacedim>
261  {
263  ExcMessage("You have to initialize the cache using one of the "
264  "reinit functions first!"));
266  }
267 
268 
269 
270  template <int dim, int spacedim>
271  const std::vector<Point<spacedim>> &
273  {
274  return get_current_fe_values().get_quadrature_points();
275  }
276 
277 
278 
279  template <int dim, int spacedim>
280  const std::vector<double> &
282  {
283  return get_current_fe_values().get_JxW_values();
284  }
285 
286 
287 
288  template <int dim, int spacedim>
289  const std::vector<double> &
291  {
292  return get_current_neighbor_fe_values().get_JxW_values();
293  }
294 
295 
296 
297  template <int dim, int spacedim>
298  const std::vector<Tensor<1, spacedim>> &
300  {
301  return get_current_fe_values().get_normal_vectors();
302  }
303 
304 
305 
306  template <int dim, int spacedim>
307  const std::vector<Tensor<1, spacedim>> &
309  {
310  return get_current_neighbor_fe_values().get_normal_vectors();
311  }
312 
313 
314 
315  template <int dim, int spacedim>
316  const std::vector<types::global_dof_index> &
318  {
319  return local_dof_indices;
320  }
321 
322 
323 
324  template <int dim, int spacedim>
325  const std::vector<types::global_dof_index> &
327  {
328  return neighbor_dof_indices;
329  }
330 
331 
332 
333  template <int dim, int spacedim>
336  {
337  return user_data_storage;
338  }
339 
340 
341 
342  template <int dim, int spacedim>
343  const Mapping<dim, spacedim> &
345  {
346  return *mapping;
347  }
348 
349 } // namespace MeshWorker
351 
352 // Explicit instantiations
354 namespace MeshWorker
355 {
356 #include "scratch_data.inst"
357 }
static const unsigned int invalid_unsigned_int
Definition: types.h:191
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:25
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:1408
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
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:361
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
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:360
const FEValuesBase< dim, spacedim > & get_current_fe_values() const
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:329
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