Reference documentation for deal.II version 9.2.0
\(\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 
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>
48  const Mapping<dim, spacedim> & mapping,
50  const Quadrature<dim> & quadrature,
51  const UpdateFlags & update_flags,
52  const UpdateFlags & neighbor_update_flags,
53  const Quadrature<dim - 1> & face_quadrature,
54  const UpdateFlags & face_update_flags,
55  const UpdateFlags & neighbor_face_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,
75  const Quadrature<dim - 1> & face_quadrature,
76  const UpdateFlags & face_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,
93  const Quadrature<dim - 1> & face_quadrature,
94  const UpdateFlags & face_update_flags,
95  const UpdateFlags & neighbor_face_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)
113  , cell_quadrature(scratch.cell_quadrature)
114  , face_quadrature(scratch.face_quadrature)
115  , cell_update_flags(scratch.cell_update_flags)
116  , neighbor_cell_update_flags(scratch.neighbor_cell_update_flags)
117  , face_update_flags(scratch.face_update_flags)
118  , neighbor_face_update_flags(scratch.neighbor_face_update_flags)
119  , local_dof_indices(scratch.local_dof_indices)
120  , neighbor_dof_indices(scratch.neighbor_dof_indices)
121  , user_data_storage(scratch.user_data_storage)
122  , internal_data_storage(scratch.internal_data_storage)
123  {}
124 
125 
126 
127  template <int dim, int spacedim>
131  {
132  if (!fe_values)
133  fe_values = std_cxx14::make_unique<FEValues<dim, spacedim>>(
134  *mapping, *fe, cell_quadrature, cell_update_flags);
135 
136  fe_values->reinit(cell);
137  cell->get_dof_indices(local_dof_indices);
138  current_fe_values = fe_values.get();
139  return *fe_values;
140  }
141 
142 
143 
144  template <int dim, int spacedim>
148  const unsigned int face_no)
149  {
150  if (!fe_face_values)
151  fe_face_values = std_cxx14::make_unique<FEFaceValues<dim, spacedim>>(
152  *mapping, *fe, face_quadrature, face_update_flags);
153 
154  fe_face_values->reinit(cell, face_no);
155  cell->get_dof_indices(local_dof_indices);
156  current_fe_values = fe_face_values.get();
157  return *fe_face_values;
158  }
159 
160 
161 
162  template <int dim, int spacedim>
166  const unsigned int face_no,
167  const unsigned int subface_no)
168  {
169  if (subface_no != numbers::invalid_unsigned_int)
170  {
171  if (!fe_subface_values)
172  fe_subface_values =
173  std_cxx14::make_unique<FESubfaceValues<dim, spacedim>>(
174  *mapping, *fe, face_quadrature, face_update_flags);
175  fe_subface_values->reinit(cell, face_no, subface_no);
176  cell->get_dof_indices(local_dof_indices);
177 
178  current_fe_values = fe_subface_values.get();
179  return *fe_subface_values;
180  }
181  else
182  return reinit(cell, face_no);
183  }
184 
185 
186 
187  template <int dim, int spacedim>
191  {
192  if (!neighbor_fe_values)
193  neighbor_fe_values = std_cxx14::make_unique<FEValues<dim, spacedim>>(
194  *mapping, *fe, cell_quadrature, neighbor_cell_update_flags);
195 
196  neighbor_fe_values->reinit(cell);
197  cell->get_dof_indices(neighbor_dof_indices);
198  current_neighbor_fe_values = neighbor_fe_values.get();
199  return *neighbor_fe_values;
200  }
201 
202 
203 
204  template <int dim, int spacedim>
208  const unsigned int face_no)
209  {
210  if (!neighbor_fe_face_values)
211  neighbor_fe_face_values =
212  std_cxx14::make_unique<FEFaceValues<dim, spacedim>>(
213  *mapping, *fe, face_quadrature, neighbor_face_update_flags);
214  neighbor_fe_face_values->reinit(cell, face_no);
215  cell->get_dof_indices(neighbor_dof_indices);
216  current_neighbor_fe_values = neighbor_fe_face_values.get();
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  {
231  if (!neighbor_fe_subface_values)
232  neighbor_fe_subface_values =
233  std_cxx14::make_unique<FESubfaceValues<dim, spacedim>>(
234  *mapping, *fe, face_quadrature, neighbor_face_update_flags);
235  neighbor_fe_subface_values->reinit(cell, face_no, subface_no);
236  cell->get_dof_indices(neighbor_dof_indices);
237  current_neighbor_fe_values = neighbor_fe_subface_values.get();
238  return *neighbor_fe_subface_values;
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  {
262  Assert(current_neighbor_fe_values != nullptr,
263  ExcMessage("You have to initialize the cache using one of the "
264  "reinit functions first!"));
265  return *current_neighbor_fe_values;
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 }
DoFHandler::active_cell_iterator
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:329
StaticMappingQ1
Definition: mapping_q1.h:88
GeneralDataStorage
Definition: general_data_storage.h:56
MeshWorker
Definition: assemble_flags.h:30
MeshWorker::ScratchData::get_neighbor_dof_indices
const std::vector< types::global_dof_index > & get_neighbor_dof_indices() const
Definition: scratch_data.cc:326
MeshWorker::ScratchData::get_normal_vectors
const std::vector< Tensor< 1, spacedim > > & get_normal_vectors() const
Definition: scratch_data.cc:299
MeshWorker::ScratchData::get_neighbor_JxW_values
const std::vector< double > & get_neighbor_JxW_values() const
Definition: scratch_data.cc:290
MeshWorker::ScratchData::get_local_dof_indices
const std::vector< types::global_dof_index > & get_local_dof_indices() const
Definition: scratch_data.cc:317
MeshWorker::ScratchData::get_current_neighbor_fe_values
const FEValuesBase< dim, spacedim > & get_current_neighbor_fe_values() const
Definition: scratch_data.cc:260
FEValues
Definition: fe.h:38
MeshWorker::ScratchData::reinit
const FEValues< dim, spacedim > & reinit(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
Definition: scratch_data.cc:129
MeshWorker::ScratchData
Definition: scratch_data.h:232
Mapping
Abstract base class for mapping classes.
Definition: mapping.h:302
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
internal::reinit
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:621
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
MeshWorker::ScratchData::get_current_fe_values
const FEValuesBase< dim, spacedim > & get_current_fe_values() const
Definition: scratch_data.cc:248
MeshWorker::ScratchData::get_neighbor_normal_vectors
const std::vector< Tensor< 1, spacedim > > & get_neighbor_normal_vectors()
Definition: scratch_data.cc:308
FiniteElement
Definition: fe.h:648
FEValuesBase
Definition: fe.h:36
UpdateFlags
UpdateFlags
Definition: fe_update_flags.h:66
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:191
scratch_data.h
FEFaceValuesBase
Definition: fe_values.h:3735
MeshWorker::ScratchData::get_JxW_values
const std::vector< double > & get_JxW_values() const
Definition: scratch_data.cc:281
MeshWorker::ScratchData::get_quadrature_points
const std::vector< Point< spacedim > > & get_quadrature_points() const
Definition: scratch_data.cc:272
memory.h
FEFaceValues
Definition: fe.h:40
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
MeshWorker::ScratchData::get_general_data_storage
GeneralDataStorage & get_general_data_storage()
Definition: scratch_data.cc:335
MeshWorker::ScratchData::get_mapping
const Mapping< dim, spacedim > & get_mapping() const
Definition: scratch_data.cc:344
Quadrature
Definition: quadrature.h:85
MeshWorker::ScratchData::ScratchData
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
MeshWorker::ScratchData::reinit_neighbor
const FEValues< dim, spacedim > & reinit_neighbor(const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell)
Definition: scratch_data.cc:189