Reference documentation for deal.II version Git f40be01994 2020-04-09 07:13:12 +0200
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_values.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2018 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/fe/mapping_q1.h>
17 
18 #include <deal.II/hp/fe_values.h>
19 
20 DEAL_II_NAMESPACE_OPEN
21 
22 namespace hp
23 {
24  // -------------------------- FEValuesBase -------------------------
25 
26  template <int dim, int q_dim, class FEValuesType>
29  & mapping_collection,
31  const QCollection<q_dim> & q_collection,
32  const UpdateFlags update_flags)
33  : fe_collection(&fe_collection)
34  , mapping_collection(&mapping_collection)
35  , q_collection(q_collection)
36  , fe_values_table(fe_collection.size(),
37  mapping_collection.size(),
38  q_collection.size())
39  , present_fe_values_index(numbers::invalid_unsigned_int,
40  numbers::invalid_unsigned_int,
41  numbers::invalid_unsigned_int)
42  , update_flags(update_flags)
43  {}
44 
45 
46  template <int dim, int q_dim, class FEValuesType>
49  const QCollection<q_dim> & q_collection,
51  : fe_collection(&fe_collection)
52  , mapping_collection(
53  &::hp::StaticMappingQ1<dim, FEValuesType::space_dimension>::
54  mapping_collection)
55  , q_collection(q_collection)
56  , fe_values_table(fe_collection.size(), 1, q_collection.size())
57  , present_fe_values_index(numbers::invalid_unsigned_int,
58  numbers::invalid_unsigned_int,
59  numbers::invalid_unsigned_int)
60  , update_flags(update_flags)
61  {}
62 
63 
64 
65  template <int dim, int q_dim, class FEValuesType>
66  FEValuesType &
68  const unsigned int fe_index,
69  const unsigned int mapping_index,
70  const unsigned int q_index)
71  {
72  AssertIndexRange(fe_index, fe_collection->size());
73  AssertIndexRange(mapping_index, mapping_collection->size());
74  AssertIndexRange(q_index, q_collection.size());
75 
76 
77  // set the triple of indices
78  // that we want to work with
79  present_fe_values_index = TableIndices<3>(fe_index, mapping_index, q_index);
80 
81  // first check whether we
82  // already have an object for
83  // this particular combination
84  // of indices
85  if (fe_values_table(present_fe_values_index).get() == nullptr)
86  fe_values_table(present_fe_values_index) =
87  std::make_shared<FEValuesType>((*mapping_collection)[mapping_index],
88  (*fe_collection)[fe_index],
89  q_collection[q_index],
90  update_flags);
91 
92  // now there definitely is one!
93  return *fe_values_table(present_fe_values_index);
94  }
95 } // namespace hp
96 
97 
98 namespace hp
99 {
100  // -------------------------- FEValues -------------------------
101 
102 
103  template <int dim, int spacedim>
106  const FECollection<dim, spacedim> & fe_collection,
107  const QCollection<dim> & q_collection,
109  : hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>(mapping,
110  fe_collection,
111  q_collection,
112  update_flags)
113  {}
114 
115 
116  template <int dim, int spacedim>
121  : hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>(fe_collection,
122  q_collection,
123  update_flags)
124  {}
125 
126 
127  template <int dim, int spacedim>
128  template <typename DoFHandlerType, bool lda>
129  void
132  const unsigned int q_index,
133  const unsigned int mapping_index,
134  const unsigned int fe_index)
135  {
136  // determine which indices we
137  // should actually use
138  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
139  real_fe_index = fe_index;
140 
141  if (real_q_index == numbers::invalid_unsigned_int)
142  {
143  if (this->q_collection.size() > 1)
144  real_q_index = cell->active_fe_index();
145  else
146  real_q_index = 0;
147  }
148 
149  if (real_mapping_index == numbers::invalid_unsigned_int)
150  {
151  if (this->mapping_collection->size() > 1)
152  real_mapping_index = cell->active_fe_index();
153  else
154  real_mapping_index = 0;
155  }
156 
157  if (real_fe_index == numbers::invalid_unsigned_int)
158  real_fe_index = cell->active_fe_index();
159 
160  // some checks
161  AssertIndexRange(real_q_index, this->q_collection.size());
162  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
163  AssertIndexRange(real_fe_index, this->fe_collection->size());
164 
165  // now finally actually get the
166  // corresponding object and
167  // initialize it
168  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
169  .reinit(cell);
170  }
171 
172 
173 
174  template <int dim, int spacedim>
175  void
177  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
178  const unsigned int q_index,
179  const unsigned int mapping_index,
180  const unsigned int fe_index)
181  {
182  // determine which indices we
183  // should actually use
184  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
185  real_fe_index = fe_index;
186 
187  if (real_q_index == numbers::invalid_unsigned_int)
188  real_q_index = 0;
189 
190  if (real_mapping_index == numbers::invalid_unsigned_int)
191  real_mapping_index = 0;
192 
193  if (real_fe_index == numbers::invalid_unsigned_int)
194  real_fe_index = 0;
195 
196  // some checks
197  AssertIndexRange(real_q_index, this->q_collection.size());
198  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
199  AssertIndexRange(real_fe_index, this->fe_collection->size());
200 
201  // now finally actually get the
202  // corresponding object and
203  // initialize it
204  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
205  .reinit(cell);
206  }
207 
208 
209  // -------------------------- FEFaceValues -------------------------
210 
211 
212  template <int dim, int spacedim>
218  : hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>(
219  mapping,
220  fe_collection,
221  q_collection,
222  update_flags)
223  {}
224 
225 
226  template <int dim, int spacedim>
231  : hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>(
232  fe_collection,
233  q_collection,
234  update_flags)
235  {}
236 
237 
238  template <int dim, int spacedim>
239  template <typename DoFHandlerType, bool lda>
240  void
243  const unsigned int face_no,
244  const unsigned int q_index,
245  const unsigned int mapping_index,
246  const unsigned int fe_index)
247  {
248  // determine which indices we
249  // should actually use
250  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
251  real_fe_index = fe_index;
252 
253  if (real_q_index == numbers::invalid_unsigned_int)
254  {
255  if (this->q_collection.size() > 1)
256  real_q_index = cell->active_fe_index();
257  else
258  real_q_index = 0;
259  }
260 
261  if (real_mapping_index == numbers::invalid_unsigned_int)
262  {
263  if (this->mapping_collection->size() > 1)
264  real_mapping_index = cell->active_fe_index();
265  else
266  real_mapping_index = 0;
267  }
268 
269  if (real_fe_index == numbers::invalid_unsigned_int)
270  real_fe_index = cell->active_fe_index();
271 
272  // some checks
273  AssertIndexRange(real_q_index, this->q_collection.size());
274  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
275  AssertIndexRange(real_fe_index, this->fe_collection->size());
276 
277  // now finally actually get the
278  // corresponding object and
279  // initialize it
280  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
281  .reinit(cell, face_no);
282  }
283 
284 
285 
286  template <int dim, int spacedim>
287  void
289  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
290  const unsigned int face_no,
291  const unsigned int q_index,
292  const unsigned int mapping_index,
293  const unsigned int fe_index)
294  {
295  // determine which indices we
296  // should actually use
297  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
298  real_fe_index = fe_index;
299 
300  if (real_q_index == numbers::invalid_unsigned_int)
301  real_q_index = 0;
302 
303  if (real_mapping_index == numbers::invalid_unsigned_int)
304  real_mapping_index = 0;
305 
306  if (real_fe_index == numbers::invalid_unsigned_int)
307  real_fe_index = 0;
308 
309  // some checks
310  AssertIndexRange(real_q_index, this->q_collection.size());
311  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
312  AssertIndexRange(real_fe_index, this->fe_collection->size());
313 
314  // now finally actually get the
315  // corresponding object and
316  // initialize it
317  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
318  .reinit(cell, face_no);
319  }
320 
321 
322  // -------------------------- FESubfaceValues -------------------------
323 
324 
325  template <int dim, int spacedim>
331  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
332  mapping,
333  fe_collection,
334  q_collection,
335  update_flags)
336  {}
337 
338 
339  template <int dim, int spacedim>
341  const hp::FECollection<dim, spacedim> &fe_collection,
342  const hp::QCollection<dim - 1> & q_collection,
343  const UpdateFlags update_flags)
344  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
345  fe_collection,
346  q_collection,
347  update_flags)
348  {}
349 
350 
351  template <int dim, int spacedim>
352  template <typename DoFHandlerType, bool lda>
353  void
356  const unsigned int face_no,
357  const unsigned int subface_no,
358  const unsigned int q_index,
359  const unsigned int mapping_index,
360  const unsigned int fe_index)
361  {
362  // determine which indices we
363  // should actually use
364  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
365  real_fe_index = fe_index;
366 
367  if (real_q_index == numbers::invalid_unsigned_int)
368  {
369  if (this->q_collection.size() > 1)
370  real_q_index = cell->active_fe_index();
371  else
372  real_q_index = 0;
373  }
374 
375  if (real_mapping_index == numbers::invalid_unsigned_int)
376  {
377  if (this->mapping_collection->size() > 1)
378  real_mapping_index = cell->active_fe_index();
379  else
380  real_mapping_index = 0;
381  }
382 
383  if (real_fe_index == numbers::invalid_unsigned_int)
384  real_fe_index = cell->active_fe_index();
385 
386  // some checks
387  AssertIndexRange(real_q_index, this->q_collection.size());
388  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
389  AssertIndexRange(real_fe_index, this->fe_collection->size());
390 
391  // now finally actually get the
392  // corresponding object and
393  // initialize it
394  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
395  .reinit(cell, face_no, subface_no);
396  }
397 
398 
399 
400  template <int dim, int spacedim>
401  void
403  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
404  const unsigned int face_no,
405  const unsigned int subface_no,
406  const unsigned int q_index,
407  const unsigned int mapping_index,
408  const unsigned int fe_index)
409  {
410  // determine which indices we
411  // should actually use
412  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
413  real_fe_index = fe_index;
414 
415  if (real_q_index == numbers::invalid_unsigned_int)
416  real_q_index = 0;
417 
418  if (real_mapping_index == numbers::invalid_unsigned_int)
419  real_mapping_index = 0;
420 
421  if (real_fe_index == numbers::invalid_unsigned_int)
422  real_fe_index = 0;
423 
424  // some checks
425  AssertIndexRange(real_q_index, this->q_collection.size());
426  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
427  AssertIndexRange(real_fe_index, this->fe_collection->size());
428 
429  // now finally actually get the
430  // corresponding object and
431  // initialize it
432  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
433  .reinit(cell, face_no, subface_no);
434  }
435 } // namespace hp
436 
437 
438 // explicit instantiations
439 #include "fe_values.inst"
440 
441 
442 DEAL_II_NAMESPACE_CLOSE
static const unsigned int invalid_unsigned_int
Definition: types.h:190
FEValues(const MappingCollection< dim, spacedim > &mapping_collection, const FECollection< dim, spacedim > &fe_collection, const QCollection< dim > &q_collection, const UpdateFlags update_flags)
Definition: fe_values.cc:104
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType, lda >> cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
Definition: fe_values.cc:354
FEValuesBase(const unsigned int n_q_points, const unsigned int dofs_per_cell, const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe)
Definition: fe_values.cc:3076
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
const SmartPointer< const MappingCollection< dim, ::FEValues< dim, spacedim > ::space_dimension >, FEValuesBase< dim, q_dim, ::FEValues< dim, spacedim > > > mapping_collection
Definition: fe_values.h:147
static const unsigned int space_dimension
Definition: fe_values.h:2095
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3387
UpdateFlags
FEFaceValues(const hp::MappingCollection< dim, spacedim > &mapping_collection, const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &q_collection, const UpdateFlags update_flags)
Definition: fe_values.cc:213
unsigned int size() const
Definition: q_collection.h:176
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType, lda >> cell, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
Definition: fe_values.cc:130
Definition: hp.h:117
const SmartPointer< const FECollection< dim, ::FEValues< dim, spacedim > ::space_dimension >, FEValuesBase< dim, q_dim, ::FEValues< dim, spacedim > > > fe_collection
Definition: fe_values.h:139
::FEValues< dim, spacedim > & select_fe_values(const unsigned int fe_index, const unsigned int mapping_index, const unsigned int q_index)
Definition: fe_values.cc:67
FEValuesBase(const MappingCollection< dim, FEValuesType::space_dimension > &mapping_collection, const FECollection< dim, FEValuesType::space_dimension > &fe_collection, const QCollection< q_dim > &q_collection, const UpdateFlags update_flags)
Definition: fe_values.cc:27
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4390
FESubfaceValues(const hp::MappingCollection< dim, spacedim > &mapping_collection, const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &q_collection, const UpdateFlags update_flags)
Definition: fe_values.cc:326
UpdateFlags update_flags
Definition: fe_values.h:3433
void reinit(const TriaIterator< DoFCellAccessor< DoFHandlerType, lda >> cell, const unsigned int face_no, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
Definition: fe_values.cc:241