Reference documentation for deal.II version Git 3ff8e288d2 2020-02-20 17:30:00 +0100
\(\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>
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  AssertIndexRange(fe_index, fe_collection->size());
75  AssertIndexRange(mapping_index, mapping_collection->size());
76  AssertIndexRange(q_index, q_collection.size());
77 
78 
79  // set the triple of indices
80  // that we want to work with
81  present_fe_values_index = TableIndices<3>(fe_index, mapping_index, q_index);
82 
83  // first check whether we
84  // already have an object for
85  // this particular combination
86  // of indices
87  if (fe_values_table(present_fe_values_index).get() == nullptr)
88  fe_values_table(present_fe_values_index) =
89  std::make_shared<FEValuesType>((*mapping_collection)[mapping_index],
90  (*fe_collection)[fe_index],
91  q_collection[q_index],
92  update_flags);
93 
94  // now there definitely is one!
95  return *fe_values_table(present_fe_values_index);
96  }
97 } // namespace hp
98 
99 
100 namespace hp
101 {
102  // -------------------------- FEValues -------------------------
103 
104 
105  template <int dim, int spacedim>
108  const hp::FECollection<dim, spacedim> & fe_collection,
109  const hp::QCollection<dim> & q_collection,
112  fe_collection,
113  q_collection,
114  update_flags)
115  {}
116 
117 
118  template <int dim, int spacedim>
120  const hp::FECollection<dim, spacedim> &fe_collection,
121  const hp::QCollection<dim> & q_collection,
122  const UpdateFlags update_flags)
123  : hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>(fe_collection,
124  q_collection,
125  update_flags)
126  {}
127 
128 
129  template <int dim, int spacedim>
130  template <typename DoFHandlerType, bool lda>
131  void
134  const unsigned int q_index,
135  const unsigned int mapping_index,
136  const unsigned int fe_index)
137  {
138  // determine which indices we
139  // should actually use
140  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
141  real_fe_index = fe_index;
142 
143  if (real_q_index == numbers::invalid_unsigned_int)
144  {
145  if (this->q_collection.size() > 1)
146  real_q_index = cell->active_fe_index();
147  else
148  real_q_index = 0;
149  }
150 
151  if (real_mapping_index == numbers::invalid_unsigned_int)
152  {
153  if (this->mapping_collection->size() > 1)
154  real_mapping_index = cell->active_fe_index();
155  else
156  real_mapping_index = 0;
157  }
158 
159  if (real_fe_index == numbers::invalid_unsigned_int)
160  real_fe_index = cell->active_fe_index();
161 
162  // some checks
163  AssertIndexRange(real_q_index, this->q_collection.size());
164  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
165  AssertIndexRange(real_fe_index, this->fe_collection->size());
166 
167  // now finally actually get the
168  // corresponding object and
169  // initialize it
170  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
171  .reinit(cell);
172  }
173 
174 
175 
176  template <int dim, int spacedim>
177  void
179  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
180  const unsigned int q_index,
181  const unsigned int mapping_index,
182  const unsigned int fe_index)
183  {
184  // determine which indices we
185  // should actually use
186  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
187  real_fe_index = fe_index;
188 
189  if (real_q_index == numbers::invalid_unsigned_int)
190  real_q_index = 0;
191 
192  if (real_mapping_index == numbers::invalid_unsigned_int)
193  real_mapping_index = 0;
194 
195  if (real_fe_index == numbers::invalid_unsigned_int)
196  real_fe_index = 0;
197 
198  // some checks
199  AssertIndexRange(real_q_index, this->q_collection.size());
200  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
201  AssertIndexRange(real_fe_index, this->fe_collection->size());
202 
203  // now finally actually get the
204  // corresponding object and
205  // initialize it
206  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
207  .reinit(cell);
208  }
209 
210 
211  // -------------------------- FEFaceValues -------------------------
212 
213 
214  template <int dim, int spacedim>
217  const hp::FECollection<dim, spacedim> & fe_collection,
218  const hp::QCollection<dim - 1> & q_collection,
219  const UpdateFlags update_flags)
220  : hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>(
221  mapping,
222  fe_collection,
223  q_collection,
224  update_flags)
225  {}
226 
227 
228  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  fe_collection,
235  q_collection,
236  update_flags)
237  {}
238 
239 
240  template <int dim, int spacedim>
241  template <typename DoFHandlerType, bool lda>
242  void
245  const unsigned int face_no,
246  const unsigned int q_index,
247  const unsigned int mapping_index,
248  const unsigned int fe_index)
249  {
250  // determine which indices we
251  // should actually use
252  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
253  real_fe_index = fe_index;
254 
255  if (real_q_index == numbers::invalid_unsigned_int)
256  {
257  if (this->q_collection.size() > 1)
258  real_q_index = cell->active_fe_index();
259  else
260  real_q_index = 0;
261  }
262 
263  if (real_mapping_index == numbers::invalid_unsigned_int)
264  {
265  if (this->mapping_collection->size() > 1)
266  real_mapping_index = cell->active_fe_index();
267  else
268  real_mapping_index = 0;
269  }
270 
271  if (real_fe_index == numbers::invalid_unsigned_int)
272  real_fe_index = cell->active_fe_index();
273 
274  // some checks
275  AssertIndexRange(real_q_index, this->q_collection.size());
276  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
277  AssertIndexRange(real_fe_index, this->fe_collection->size());
278 
279  // now finally actually get the
280  // corresponding object and
281  // initialize it
282  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
283  .reinit(cell, face_no);
284  }
285 
286 
287 
288  template <int dim, int spacedim>
289  void
291  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
292  const unsigned int face_no,
293  const unsigned int q_index,
294  const unsigned int mapping_index,
295  const unsigned int fe_index)
296  {
297  // determine which indices we
298  // should actually use
299  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
300  real_fe_index = fe_index;
301 
302  if (real_q_index == numbers::invalid_unsigned_int)
303  real_q_index = 0;
304 
305  if (real_mapping_index == numbers::invalid_unsigned_int)
306  real_mapping_index = 0;
307 
308  if (real_fe_index == numbers::invalid_unsigned_int)
309  real_fe_index = 0;
310 
311  // some checks
312  AssertIndexRange(real_q_index, this->q_collection.size());
313  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
314  AssertIndexRange(real_fe_index, this->fe_collection->size());
315 
316  // now finally actually get the
317  // corresponding object and
318  // initialize it
319  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
320  .reinit(cell, face_no);
321  }
322 
323 
324  // -------------------------- FESubfaceValues -------------------------
325 
326 
327  template <int dim, int spacedim>
330  const hp::FECollection<dim, spacedim> & fe_collection,
331  const hp::QCollection<dim - 1> & q_collection,
332  const UpdateFlags update_flags)
333  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
334  mapping,
335  fe_collection,
336  q_collection,
337  update_flags)
338  {}
339 
340 
341  template <int dim, int spacedim>
343  const hp::FECollection<dim, spacedim> &fe_collection,
344  const hp::QCollection<dim - 1> & q_collection,
345  const UpdateFlags update_flags)
346  : hp::FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>(
347  fe_collection,
348  q_collection,
349  update_flags)
350  {}
351 
352 
353  template <int dim, int spacedim>
354  template <typename DoFHandlerType, bool lda>
355  void
358  const unsigned int face_no,
359  const unsigned int subface_no,
360  const unsigned int q_index,
361  const unsigned int mapping_index,
362  const unsigned int fe_index)
363  {
364  // determine which indices we
365  // should actually use
366  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
367  real_fe_index = fe_index;
368 
369  if (real_q_index == numbers::invalid_unsigned_int)
370  {
371  if (this->q_collection.size() > 1)
372  real_q_index = cell->active_fe_index();
373  else
374  real_q_index = 0;
375  }
376 
377  if (real_mapping_index == numbers::invalid_unsigned_int)
378  {
379  if (this->mapping_collection->size() > 1)
380  real_mapping_index = cell->active_fe_index();
381  else
382  real_mapping_index = 0;
383  }
384 
385  if (real_fe_index == numbers::invalid_unsigned_int)
386  real_fe_index = cell->active_fe_index();
387 
388  // some checks
389  AssertIndexRange(real_q_index, this->q_collection.size());
390  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
391  AssertIndexRange(real_fe_index, this->fe_collection->size());
392 
393  // now finally actually get the
394  // corresponding object and
395  // initialize it
396  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
397  .reinit(cell, face_no, subface_no);
398  }
399 
400 
401 
402  template <int dim, int spacedim>
403  void
405  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
406  const unsigned int face_no,
407  const unsigned int subface_no,
408  const unsigned int q_index,
409  const unsigned int mapping_index,
410  const unsigned int fe_index)
411  {
412  // determine which indices we
413  // should actually use
414  unsigned int real_q_index = q_index, real_mapping_index = mapping_index,
415  real_fe_index = fe_index;
416 
417  if (real_q_index == numbers::invalid_unsigned_int)
418  real_q_index = 0;
419 
420  if (real_mapping_index == numbers::invalid_unsigned_int)
421  real_mapping_index = 0;
422 
423  if (real_fe_index == numbers::invalid_unsigned_int)
424  real_fe_index = 0;
425 
426  // some checks
427  AssertIndexRange(real_q_index, this->q_collection.size());
428  AssertIndexRange(real_mapping_index, this->mapping_collection->size());
429  AssertIndexRange(real_fe_index, this->fe_collection->size());
430 
431  // now finally actually get the
432  // corresponding object and
433  // initialize it
434  this->select_fe_values(real_fe_index, real_mapping_index, real_q_index)
435  .reinit(cell, face_no, subface_no);
436  }
437 } // namespace hp
438 
439 
440 // explicit instantiations
441 #include "fe_values.inst"
442 
443 
444 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:356
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:1641
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:3401
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)
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:215
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:132
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:3880
FEValues(const Mapping< dim, spacedim > &mapping, const FiniteElement< dim, spacedim > &fe, const Quadrature< dim > &quadrature, const UpdateFlags update_flags)
Definition: fe_values.cc:4402
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:328
UpdateFlags update_flags
Definition: fe_values.h:3447
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:243