Reference documentation for deal.II version Git 840703ecf3 2019-10-14 17:15:20 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
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>
28  const ::hp::MappingCollection<dim, FEValuesType::space_dimension>
29  &mapping_collection,
30  const ::hp::FECollection<dim, FEValuesType::space_dimension>
31  & fe_collection,
32  const ::hp::QCollection<q_dim> &q_collection,
33  const UpdateFlags update_flags)
34  : fe_collection(&fe_collection)
35  , mapping_collection(&mapping_collection)
36  , q_collection(q_collection)
37  , fe_values_table(fe_collection.size(),
38  mapping_collection.size(),
39  q_collection.size())
40  , present_fe_values_index(numbers::invalid_unsigned_int,
41  numbers::invalid_unsigned_int,
42  numbers::invalid_unsigned_int)
43  , update_flags(update_flags)
44  {}
45 
46 
47  template <int dim, int q_dim, class FEValuesType>
49  const ::hp::FECollection<dim, FEValuesType::space_dimension>
50  & fe_collection,
51  const ::hp::QCollection<q_dim> &q_collection,
53  : fe_collection(&fe_collection)
54  , mapping_collection(
55  &::hp::StaticMappingQ1<dim, FEValuesType::space_dimension>::
56  mapping_collection)
57  , q_collection(q_collection)
58  , fe_values_table(fe_collection.size(), 1, q_collection.size())
59  , present_fe_values_index(numbers::invalid_unsigned_int,
60  numbers::invalid_unsigned_int,
61  numbers::invalid_unsigned_int)
62  , update_flags(update_flags)
63  {}
64 
65 
66 
67  template <int dim, int q_dim, class FEValuesType>
68  FEValuesType &
70  const unsigned int fe_index,
71  const unsigned int mapping_index,
72  const unsigned int q_index)
73  {
74  Assert(fe_index < fe_collection->size(),
75  ExcIndexRange(fe_index, 0, fe_collection->size()));
76  Assert(mapping_index < mapping_collection->size(),
77  ExcIndexRange(mapping_index, 0, mapping_collection->size()));
78  Assert(q_index < q_collection.size(),
79  ExcIndexRange(q_index, 0, q_collection.size()));
80 
81 
82  // set the triple of indices
83  // that we want to work with
84  present_fe_values_index = TableIndices<3>(fe_index, mapping_index, q_index);
85 
86  // first check whether we
87  // already have an object for
88  // this particular combination
89  // of indices
90  if (fe_values_table(present_fe_values_index).get() == nullptr)
91  fe_values_table(present_fe_values_index) =
92  std::make_shared<FEValuesType>((*mapping_collection)[mapping_index],
93  (*fe_collection)[fe_index],
94  q_collection[q_index],
95  update_flags);
96 
97  // now there definitely is one!
98  return *fe_values_table(present_fe_values_index);
99  }
100 } // namespace hp
101 
102 
103 namespace hp
104 {
105  // -------------------------- FEValues -------------------------
106 
107 
108  template <int dim, int spacedim>
111  const hp::FECollection<dim, spacedim> & fe_collection,
112  const hp::QCollection<dim> & q_collection,
115  fe_collection,
116  q_collection,
117  update_flags)
118  {}
119 
120 
121  template <int dim, int spacedim>
123  const hp::FECollection<dim, spacedim> &fe_collection,
124  const hp::QCollection<dim> & q_collection,
125  const UpdateFlags update_flags)
126  : hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>(fe_collection,
127  q_collection,
128  update_flags)
129  {}
130 
131 
132  template <int dim, int spacedim>
133  template <typename DoFHandlerType, bool lda>
134  void
137  const unsigned int q_index,
138  const unsigned int mapping_index,
139  const unsigned int fe_index)
140  {
141  // determine which indices we
142  // should actually use
143  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
144  real_fe_index = fe_index;
145 
146  if (real_q_index == numbers::invalid_unsigned_int)
147  {
148  if (this->q_collection.size() > 1)
149  real_q_index = cell->active_fe_index();
150  else
151  real_q_index = 0;
152  }
153 
154  if (real_mapping_index == numbers::invalid_unsigned_int)
155  {
156  if (this->mapping_collection->size() > 1)
157  real_mapping_index = cell->active_fe_index();
158  else
159  real_mapping_index = 0;
160  }
161 
162  if (real_fe_index == numbers::invalid_unsigned_int)
163  real_fe_index = cell->active_fe_index();
164 
165  // some checks
166  Assert(real_q_index < this->q_collection.size(),
167  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
168  Assert(real_mapping_index < this->mapping_collection->size(),
169  ExcIndexRange(real_mapping_index,
170  0,
171  this->mapping_collection->size()));
172  Assert(real_fe_index < this->fe_collection->size(),
173  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
174 
175  // now finally actually get the
176  // corresponding object and
177  // initialize it
178  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
179  .reinit(cell);
180  }
181 
182 
183 
184  template <int dim, int spacedim>
185  void
187  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
188  const unsigned int q_index,
189  const unsigned int mapping_index,
190  const unsigned int fe_index)
191  {
192  // determine which indices we
193  // should actually use
194  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
195  real_fe_index = fe_index;
196 
197  if (real_q_index == numbers::invalid_unsigned_int)
198  real_q_index = 0;
199 
200  if (real_mapping_index == numbers::invalid_unsigned_int)
201  real_mapping_index = 0;
202 
203  if (real_fe_index == numbers::invalid_unsigned_int)
204  real_fe_index = 0;
205 
206  // some checks
207  Assert(real_q_index < this->q_collection.size(),
208  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
209  Assert(real_mapping_index < this->mapping_collection->size(),
210  ExcIndexRange(real_mapping_index,
211  0,
212  this->mapping_collection->size()));
213  Assert(real_fe_index < this->fe_collection->size(),
214  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
215 
216  // now finally actually get the
217  // corresponding object and
218  // initialize it
219  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
220  .reinit(cell);
221  }
222 
223 
224  // -------------------------- FEFaceValues -------------------------
225 
226 
227  template <int dim, int spacedim>
230  const hp::FECollection<dim, spacedim> & fe_collection,
231  const hp::QCollection<dim - 1> & q_collection,
232  const UpdateFlags update_flags)
233  : hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>(
234  mapping,
235  fe_collection,
236  q_collection,
237  update_flags)
238  {}
239 
240 
241  template <int dim, int spacedim>
243  const hp::FECollection<dim, spacedim> &fe_collection,
244  const hp::QCollection<dim - 1> & q_collection,
245  const UpdateFlags update_flags)
246  : hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>(
247  fe_collection,
248  q_collection,
249  update_flags)
250  {}
251 
252 
253  template <int dim, int spacedim>
254  template <typename DoFHandlerType, bool lda>
255  void
258  const unsigned int face_no,
259  const unsigned int q_index,
260  const unsigned int mapping_index,
261  const unsigned int fe_index)
262  {
263  // determine which indices we
264  // should actually use
265  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
266  real_fe_index = fe_index;
267 
268  if (real_q_index == numbers::invalid_unsigned_int)
269  {
270  if (this->q_collection.size() > 1)
271  real_q_index = cell->active_fe_index();
272  else
273  real_q_index = 0;
274  }
275 
276  if (real_mapping_index == numbers::invalid_unsigned_int)
277  {
278  if (this->mapping_collection->size() > 1)
279  real_mapping_index = cell->active_fe_index();
280  else
281  real_mapping_index = 0;
282  }
283 
284  if (real_fe_index == numbers::invalid_unsigned_int)
285  real_fe_index = cell->active_fe_index();
286 
287  // some checks
288  Assert(real_q_index < this->q_collection.size(),
289  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
290  Assert(real_mapping_index < this->mapping_collection->size(),
291  ExcIndexRange(real_mapping_index,
292  0,
293  this->mapping_collection->size()));
294  Assert(real_fe_index < this->fe_collection->size(),
295  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
296 
297  // now finally actually get the
298  // corresponding object and
299  // initialize it
300  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
301  .reinit(cell, face_no);
302  }
303 
304 
305 
306  template <int dim, int spacedim>
307  void
309  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
310  const unsigned int face_no,
311  const unsigned int q_index,
312  const unsigned int mapping_index,
313  const unsigned int fe_index)
314  {
315  // determine which indices we
316  // should actually use
317  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
318  real_fe_index = fe_index;
319 
320  if (real_q_index == numbers::invalid_unsigned_int)
321  real_q_index = 0;
322 
323  if (real_mapping_index == numbers::invalid_unsigned_int)
324  real_mapping_index = 0;
325 
326  if (real_fe_index == numbers::invalid_unsigned_int)
327  real_fe_index = 0;
328 
329  // some checks
330  Assert(real_q_index < this->q_collection.size(),
331  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
332  Assert(real_mapping_index < this->mapping_collection->size(),
333  ExcIndexRange(real_mapping_index,
334  0,
335  this->mapping_collection->size()));
336  Assert(real_fe_index < this->fe_collection->size(),
337  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
338 
339  // now finally actually get the
340  // corresponding object and
341  // initialize it
342  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
343  .reinit(cell, face_no);
344  }
345 
346 
347  // -------------------------- FESubfaceValues -------------------------
348 
349 
350  template <int dim, int spacedim>
353  const hp::FECollection<dim, spacedim> & fe_collection,
354  const hp::QCollection<dim - 1> & q_collection,
355  const UpdateFlags update_flags)
356  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
357  mapping,
358  fe_collection,
359  q_collection,
360  update_flags)
361  {}
362 
363 
364  template <int dim, int spacedim>
366  const hp::FECollection<dim, spacedim> &fe_collection,
367  const hp::QCollection<dim - 1> & q_collection,
368  const UpdateFlags update_flags)
369  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
370  fe_collection,
371  q_collection,
372  update_flags)
373  {}
374 
375 
376  template <int dim, int spacedim>
377  template <typename DoFHandlerType, bool lda>
378  void
381  const unsigned int face_no,
382  const unsigned int subface_no,
383  const unsigned int q_index,
384  const unsigned int mapping_index,
385  const unsigned int fe_index)
386  {
387  // determine which indices we
388  // should actually use
389  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
390  real_fe_index = fe_index;
391 
392  if (real_q_index == numbers::invalid_unsigned_int)
393  {
394  if (this->q_collection.size() > 1)
395  real_q_index = cell->active_fe_index();
396  else
397  real_q_index = 0;
398  }
399 
400  if (real_mapping_index == numbers::invalid_unsigned_int)
401  {
402  if (this->mapping_collection->size() > 1)
403  real_mapping_index = cell->active_fe_index();
404  else
405  real_mapping_index = 0;
406  }
407 
408  if (real_fe_index == numbers::invalid_unsigned_int)
409  real_fe_index = cell->active_fe_index();
410 
411  // some checks
412  Assert(real_q_index < this->q_collection.size(),
413  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
414  Assert(real_mapping_index < this->mapping_collection->size(),
415  ExcIndexRange(real_mapping_index,
416  0,
417  this->mapping_collection->size()));
418  Assert(real_fe_index < this->fe_collection->size(),
419  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
420 
421  // now finally actually get the
422  // corresponding object and
423  // initialize it
424  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
425  .reinit(cell, face_no, subface_no);
426  }
427 
428 
429 
430  template <int dim, int spacedim>
431  void
433  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
434  const unsigned int face_no,
435  const unsigned int subface_no,
436  const unsigned int q_index,
437  const unsigned int mapping_index,
438  const unsigned int fe_index)
439  {
440  // determine which indices we
441  // should actually use
442  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
443  real_fe_index = fe_index;
444 
445  if (real_q_index == numbers::invalid_unsigned_int)
446  real_q_index = 0;
447 
448  if (real_mapping_index == numbers::invalid_unsigned_int)
449  real_mapping_index = 0;
450 
451  if (real_fe_index == numbers::invalid_unsigned_int)
452  real_fe_index = 0;
453 
454  // some checks
455  Assert(real_q_index < this->q_collection.size(),
456  ExcIndexRange(real_q_index, 0, this->q_collection.size()));
457  Assert(real_mapping_index < this->mapping_collection->size(),
458  ExcIndexRange(real_mapping_index,
459  0,
460  this->mapping_collection->size()));
461  Assert(real_fe_index < this->fe_collection->size(),
462  ExcIndexRange(real_fe_index, 0, this->fe_collection->size()));
463 
464  // now finally actually get the
465  // corresponding object and
466  // initialize it
467  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
468  .reinit(cell, face_no, subface_no);
469  }
470 } // namespace hp
471 
472 
473 // explicit instantiations
474 #include "fe_values.inst"
475 
476 
477 DEAL_II_NAMESPACE_CLOSE
static const unsigned int invalid_unsigned_int
Definition: types.h:187
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:379
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:3083
const SmartPointer< const MappingCollection< dim, ::FEValues< dim, spacedim > ::space_dimension >, FEValuesBase< dim, q_dim, ::FEValues< dim, spacedim > > > mapping_collection
Definition: fe_values.h:147
const SmartPointer< const Mapping< dim, spacedim >, FEValuesBase< dim, spacedim > > mapping
Definition: fe_values.h:3393
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
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)
#define Assert(cond, exc)
Definition: exceptions.h:1407
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:228
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:135
Definition: hp.h:117
::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:69
static const unsigned int space_dimension
Definition: fe_values.h:3841
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4409
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:351
UpdateFlags update_flags
Definition: fe_values.h:3439
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:256