Reference documentation for deal.II version Git ac854ef673 2021-03-06 10:09:45 +0100
\(\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\}}\)
mapping_q.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 2019 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 
22 #include <deal.II/base/utilities.h>
23 
25 
26 #include <deal.II/fe/fe_q.h>
27 #include <deal.II/fe/fe_values.h>
28 #include <deal.II/fe/mapping_q.h>
29 
31 
33 
34 #include <memory>
35 #include <numeric>
36 
38 
39 
40 template <int dim, int spacedim>
42  : use_mapping_q1_on_current_cell(false)
43 {}
44 
45 
46 
47 template <int dim, int spacedim>
48 std::size_t
50 {
51  return (
56 }
57 
58 
59 
60 template <int dim, int spacedim>
61 MappingQ<dim, spacedim>::MappingQ(const unsigned int degree,
62  const bool use_mapping_q_on_all_cells)
63  : polynomial_degree(degree)
64  ,
65 
66  // see whether we want to use *this* mapping objects on *all* cells,
67  // or defer to an explicit Q1 mapping on interior cells. if
68  // degree==1, then we are already that Q1 mapping, so we don't need
69  // it; if dim!=spacedim, there is also no need for anything because
70  // we're most likely on a curved manifold
71  use_mapping_q_on_all_cells(degree == 1 || use_mapping_q_on_all_cells ||
72  (dim != spacedim))
73  ,
74  // create a Q1 mapping for use on interior cells (if necessary)
75  // or to create a good initial guess in transform_real_to_unit_cell()
76  q1_mapping(std::make_shared<MappingQGeneric<dim, spacedim>>(1))
77  ,
78 
79  // create a Q_p mapping; if p=1, simply share the Q_1 mapping already
80  // created via the shared_ptr objects
81  qp_mapping(this->polynomial_degree > 1 ?
82  std::make_shared<const MappingQGeneric<dim, spacedim>>(degree) :
83  q1_mapping)
84 {}
85 
86 
87 
88 template <int dim, int spacedim>
92 {
93  // Note that we really do have to use clone() here, since mapping.q1_mapping
94  // may be MappingQ1Eulerian and mapping.qp_mapping may be MappingQEulerian.
95  std::shared_ptr<const Mapping<dim, spacedim>> other_q1_map =
96  mapping.q1_mapping->clone();
97  q1_mapping = std::dynamic_pointer_cast<const MappingQGeneric<dim, spacedim>>(
98  other_q1_map);
99  Assert(q1_mapping != nullptr, ExcInternalError());
100  Assert(q1_mapping->get_degree() == 1, ExcInternalError());
101 
102  // Same as the other constructor: if possible reuse the Q1 mapping
103  if (this->polynomial_degree == 1)
104  {
106  }
107  else
108  {
109  std::shared_ptr<const Mapping<dim, spacedim>> other_qp_map =
110  mapping.qp_mapping->clone();
111  qp_mapping =
112  std::dynamic_pointer_cast<const MappingQGeneric<dim, spacedim>>(
113  other_qp_map);
114  Assert(qp_mapping != nullptr, ExcInternalError());
115  }
116 }
117 
118 
119 
120 template <int dim, int spacedim>
121 unsigned int
123 {
124  return polynomial_degree;
125 }
126 
127 
128 
129 template <int dim, int spacedim>
130 inline bool
132 {
133  return true;
134 }
135 
136 
137 
138 template <int dim, int spacedim>
141 {
142  return (q1_mapping->requires_update_flags(in) |
143  qp_mapping->requires_update_flags(in));
144 }
145 
146 
147 
148 template <int dim, int spacedim>
149 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
151  const Quadrature<dim> &quadrature) const
152 {
153  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
154  std::make_unique<InternalData>();
155  auto &data = dynamic_cast<InternalData &>(*data_ptr);
156 
157  // build the Q1 and Qp internal data objects in parallel
159  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
161  *qp_mapping,
162  update_flags,
163  quadrature);
164 
166  data.mapping_q1_data = Utilities::dynamic_unique_cast<
168  std::move(q1_mapping->get_data(update_flags, quadrature)));
169 
170  // wait for the task above to finish and use returned value
171  data.mapping_qp_data = Utilities::dynamic_unique_cast<
172  typename MappingQGeneric<dim, spacedim>::InternalData>(
173  std::move(do_get_data.return_value()));
174  return data_ptr;
175 }
176 
177 
178 
179 template <int dim, int spacedim>
180 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
182  const UpdateFlags update_flags,
183  const hp::QCollection<dim - 1> &quadrature) const
184 {
185  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
186  std::make_unique<InternalData>();
187  auto &data = dynamic_cast<InternalData &>(*data_ptr);
188 
189  std::unique_ptr<typename MappingQGeneric<dim, spacedim>::InternalDataBase> (
190  MappingQGeneric<dim, spacedim>::*mapping_get_face_data)(
191  const UpdateFlags, const hp::QCollection<dim - 1> &) const =
193 
194  // build the Q1 and Qp internal data objects in parallel
196  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
197  do_get_data = Threads::new_task(mapping_get_face_data,
198  *qp_mapping,
199  update_flags,
200  quadrature);
201 
203  data.mapping_q1_data = Utilities::dynamic_unique_cast<
205  std::move(q1_mapping->get_face_data(update_flags, quadrature)));
206 
207  // wait for the task above to finish and use returned value
208  data.mapping_qp_data = Utilities::dynamic_unique_cast<
209  typename MappingQGeneric<dim, spacedim>::InternalData>(
210  std::move(do_get_data.return_value()));
211  return data_ptr;
212 }
213 
214 
215 
216 template <int dim, int spacedim>
217 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
219  const UpdateFlags update_flags,
220  const Quadrature<dim - 1> &quadrature) const
221 {
222  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
223  std::make_unique<InternalData>();
224  auto &data = dynamic_cast<InternalData &>(*data_ptr);
225 
226  // build the Q1 and Qp internal data objects in parallel
228  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>>
229  do_get_data =
231  *qp_mapping,
232  update_flags,
233  quadrature);
234 
236  data.mapping_q1_data = Utilities::dynamic_unique_cast<
238  std::move(q1_mapping->get_subface_data(update_flags, quadrature)));
239 
240  // wait for the task above to finish and use returned value
241  data.mapping_qp_data = Utilities::dynamic_unique_cast<
242  typename MappingQGeneric<dim, spacedim>::InternalData>(
243  std::move(do_get_data.return_value()));
244  return data_ptr;
245 }
246 
247 
248 // Note that the CellSimilarity flag is modifiable, since MappingQ can need to
249 // recalculate data even when cells are similar.
250 template <int dim, int spacedim>
253  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
254  const CellSimilarity::Similarity cell_similarity,
255  const Quadrature<dim> & quadrature,
256  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
258  &output_data) const
259 {
260  // convert data object to internal data for this class. fails with an
261  // exception if that is not possible
262  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
263  ExcInternalError());
264  const InternalData &data = static_cast<const InternalData &>(internal_data);
265 
266  // check whether this cell needs the full mapping or can be treated by a
267  // reduced Q1 mapping, e.g. if the cell is in the interior of the domain
269  !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
270 
271 
272  // call the base class. we need to ensure that the flag indicating whether
273  // we can use some similarity has to be modified - for a general MappingQ,
274  // the data needs to be recomputed anyway since then the mapping changes the
275  // data. this needs to be known also for later operations, so modify the
276  // variable here. this also affects the calculation of the next cell -- if
277  // we use Q1 data on the next cell, the data will still be invalid.
278  const CellSimilarity::Similarity updated_cell_similarity =
279  ((data.use_mapping_q1_on_current_cell == false) &&
280  (this->polynomial_degree > 1) ?
282  cell_similarity);
283 
284  // depending on the results above, decide whether the Q1 mapping or
285  // the Qp mapping needs to handle this cell
287  q1_mapping->fill_fe_values(cell,
288  updated_cell_similarity,
289  quadrature,
290  *data.mapping_q1_data,
291  output_data);
292  else
293  qp_mapping->fill_fe_values(cell,
294  updated_cell_similarity,
295  quadrature,
296  *data.mapping_qp_data,
297  output_data);
298 
299  return updated_cell_similarity;
300 }
301 
302 
303 
304 template <int dim, int spacedim>
305 void
307  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
308  const unsigned int face_no,
309  const hp::QCollection<dim - 1> & quadrature,
310  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
312  &output_data) const
313 {
314  // convert data object to internal data for this class. fails with an
315  // exception if that is not possible
316  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
317  ExcInternalError());
318  const InternalData &data = static_cast<const InternalData &>(internal_data);
319 
320  // check whether this cell needs the full mapping or can be treated by a
321  // reduced Q1 mapping, e.g. if the cell is entirely in the interior of the
322  // domain. note that it is not sufficient to ask whether the present _face_
323  // is in the interior, as the mapping on the face depends on the mapping of
324  // the cell, which in turn depends on the fact whether _any_ of the faces of
325  // this cell is at the boundary, not only the present face
327  !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
328 
329  // depending on the results above, decide whether the Q1 mapping or
330  // the Qp mapping needs to handle this cell
332  q1_mapping->fill_fe_face_values(
333  cell, face_no, quadrature, *data.mapping_q1_data, output_data);
334  else
335  qp_mapping->fill_fe_face_values(
336  cell, face_no, quadrature, *data.mapping_qp_data, output_data);
337 }
338 
339 
340 template <int dim, int spacedim>
341 void
343  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
344  const unsigned int face_no,
345  const unsigned int subface_no,
346  const Quadrature<dim - 1> & quadrature,
347  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
349  &output_data) const
350 {
351  // convert data object to internal data for this class. fails with an
352  // exception if that is not possible
353  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
354  ExcInternalError());
355  const InternalData &data = static_cast<const InternalData &>(internal_data);
356 
357  // check whether this cell needs the full mapping or can be treated by a
358  // reduced Q1 mapping, e.g. if the cell is entirely in the interior of the
359  // domain. note that it is not sufficient to ask whether the present _face_
360  // is in the interior, as the mapping on the face depends on the mapping of
361  // the cell, which in turn depends on the fact whether _any_ of the faces of
362  // this cell is at the boundary, not only the present face
364  !(use_mapping_q_on_all_cells || cell->has_boundary_lines());
365 
366  // depending on the results above, decide whether the Q1 mapping or
367  // the Qp mapping needs to handle this cell
369  q1_mapping->fill_fe_subface_values(cell,
370  face_no,
371  subface_no,
372  quadrature,
373  *data.mapping_q1_data,
374  output_data);
375  else
376  qp_mapping->fill_fe_subface_values(cell,
377  face_no,
378  subface_no,
379  quadrature,
380  *data.mapping_qp_data,
381  output_data);
382 }
383 
384 
385 
386 template <int dim, int spacedim>
387 void
389  const ArrayView<const Tensor<1, dim>> & input,
390  const MappingKind mapping_kind,
391  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
392  const ArrayView<Tensor<1, spacedim>> & output) const
393 {
394  AssertDimension(input.size(), output.size());
395 
396  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
397  Assert(data != nullptr, ExcInternalError());
398 
399  // check whether we should in fact work on the Q1 portion of it
400  if (data->use_mapping_q1_on_current_cell)
401  q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
402  else
403  // otherwise use the full mapping
404  qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
405 }
406 
407 
408 
409 template <int dim, int spacedim>
410 void
412  const ArrayView<const DerivativeForm<1, dim, spacedim>> &input,
413  const MappingKind mapping_kind,
414  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
415  const ArrayView<Tensor<2, spacedim>> & output) const
416 {
417  AssertDimension(input.size(), output.size());
418  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
419  &mapping_data) != nullptr),
420  ExcInternalError());
421  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
422 
423  // check whether we should in fact work on the Q1 portion of it
425  q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
426  else
427  // otherwise use the full mapping
428  qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
429 }
430 
431 
432 
433 template <int dim, int spacedim>
434 void
436  const ArrayView<const Tensor<2, dim>> & input,
437  const MappingKind mapping_kind,
438  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
439  const ArrayView<Tensor<2, spacedim>> & output) const
440 {
441  AssertDimension(input.size(), output.size());
442  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
443  &mapping_data) != nullptr),
444  ExcInternalError());
445  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
446 
447  // check whether we should in fact work on the Q1 portion of it
449  q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
450  else
451  // otherwise use the full mapping
452  qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
453 }
454 
455 
456 
457 template <int dim, int spacedim>
458 void
460  const ArrayView<const DerivativeForm<2, dim, spacedim>> &input,
461  const MappingKind mapping_kind,
462  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
463  const ArrayView<Tensor<3, spacedim>> & output) const
464 {
465  AssertDimension(input.size(), output.size());
466  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
467  &mapping_data) != nullptr),
468  ExcInternalError());
469  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
470 
471  // check whether we should in fact work on the Q1 portion of it
473  q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
474  else
475  // otherwise use the full mapping
476  qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
477 }
478 
479 
480 
481 template <int dim, int spacedim>
482 void
484  const ArrayView<const Tensor<3, dim>> & input,
485  const MappingKind mapping_kind,
486  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
487  const ArrayView<Tensor<3, spacedim>> & output) const
488 {
489  AssertDimension(input.size(), output.size());
490  Assert((dynamic_cast<const typename MappingQ<dim, spacedim>::InternalData *>(
491  &mapping_data) != nullptr),
492  ExcInternalError());
493  const InternalData *data = dynamic_cast<const InternalData *>(&mapping_data);
494 
495  // check whether we should in fact work on the Q1 portion of it
497  q1_mapping->transform(input, mapping_kind, *data->mapping_q1_data, output);
498  else
499  // otherwise use the full mapping
500  qp_mapping->transform(input, mapping_kind, *data->mapping_qp_data, output);
501 }
502 
503 
504 
505 template <int dim, int spacedim>
508  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
509  const Point<dim> & p) const
510 {
511  // first see, whether we want to use a linear or a higher order
512  // mapping, then either use our own facilities or that of the Q1
513  // mapping we store
514  if (use_mapping_q_on_all_cells || cell->has_boundary_lines())
515  return qp_mapping->transform_unit_to_real_cell(cell, p);
516  else
517  return q1_mapping->transform_unit_to_real_cell(cell, p);
518 }
519 
520 
521 
522 template <int dim, int spacedim>
525  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
526  const Point<spacedim> & p) const
527 {
528  if (cell->has_boundary_lines() || use_mapping_q_on_all_cells ||
529  (dim != spacedim))
530  return qp_mapping->transform_real_to_unit_cell(cell, p);
531  else
532  return q1_mapping->transform_real_to_unit_cell(cell, p);
533 }
534 
535 
536 
537 template <int dim, int spacedim>
538 std::unique_ptr<Mapping<dim, spacedim>>
540 {
541  return std::make_unique<MappingQ<dim, spacedim>>(
543 }
544 
545 
546 
547 template <int dim, int spacedim>
550  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
551 {
552  if (cell->has_boundary_lines() || use_mapping_q_on_all_cells ||
553  (dim != spacedim))
554  return BoundingBox<spacedim>(
555  qp_mapping->compute_mapping_support_points(cell));
556  else
557  return BoundingBox<spacedim>(q1_mapping->get_vertices(cell));
558 }
559 
560 
561 
562 template <int dim, int spacedim>
563 bool
565  const ReferenceCell &reference_cell) const
566 {
567  Assert(dim == reference_cell.get_dimension(),
568  ExcMessage("The dimension of your mapping (" +
569  Utilities::to_string(dim) +
570  ") and the reference cell cell_type (" +
571  Utilities::to_string(reference_cell.get_dimension()) +
572  " ) do not agree."));
573 
574  return reference_cell.is_hyper_cube();
575 }
576 
577 
578 
579 // explicit instantiations
580 #include "mapping_q.inst"
581 
582 
unsigned int get_degree() const
Definition: mapping_q.cc:122
virtual std::size_t memory_consumption() const override
Definition: mapping_q.cc:49
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
std::shared_ptr< const MappingQGeneric< dim, spacedim > > qp_mapping
Definition: mapping_q.h:377
bool is_hyper_cube() const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_q.cc:150
Task< RT > new_task(const std::function< RT()> &function)
bool use_mapping_q1_on_current_cell
Definition: mapping_q.h:256
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const override
Definition: mapping_q.cc:388
STL namespace.
const bool use_mapping_q_on_all_cells
Definition: mapping_q.h:342
MappingKind
Definition: mapping.h:64
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:482
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:181
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:342
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1466
UpdateFlags
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
Abstract base class for mapping classes.
Definition: mapping.h:303
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
const unsigned int polynomial_degree
Definition: mapping_q.h:336
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_q.cc:549
std::unique_ptr< typename MappingQGeneric< dim, spacedim >::InternalData > mapping_qp_data
Definition: mapping_q.h:270
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:218
unsigned int get_dimension() const
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_q.cc:564
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
Definition: mapping_q.cc:507
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:306
std::unique_ptr< To > dynamic_unique_cast(std::unique_ptr< From > &&p)
Definition: utilities.h:1403
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
virtual bool preserves_vertex_locations() const override
Definition: mapping_q.cc:131
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_q.cc:539
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:252
std::unique_ptr< typename MappingQGeneric< dim, spacedim >::InternalData > mapping_q1_data
Definition: mapping_q.h:263
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_q.cc:140
MappingQ(const unsigned int polynomial_degree, const bool use_mapping_q_on_all_cells=false)
Definition: mapping_q.cc:61
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
Definition: mapping_q.cc:524
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static ::ExceptionBase & ExcInternalError()
std::shared_ptr< const MappingQGeneric< dim, spacedim > > q1_mapping
Definition: mapping_q.h:361