Reference documentation for deal.II version Git 73c87d96ef 2021-11-30 22:54:44 +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_manifold.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2021 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 
24 
26 
27 #include <deal.II/fe/fe.h>
28 #include <deal.II/fe/fe_tools.h>
29 #include <deal.II/fe/fe_values.h>
31 #include <deal.II/fe/mapping_q1.h>
32 
33 #include <deal.II/grid/manifold.h>
34 #include <deal.II/grid/tria.h>
36 
38 
39 #include <algorithm>
40 #include <array>
41 #include <cmath>
42 #include <memory>
43 #include <numeric>
44 
45 
47 
48 template <int dim, int spacedim>
49 std::size_t
51 {
52  return (
65 }
66 
67 
68 template <int dim, int spacedim>
69 void
71  const UpdateFlags update_flags,
72  const Quadrature<dim> &q,
73  const unsigned int n_original_q_points)
74 {
75  // store the flags in the internal data object so we can access them
76  // in fill_fe_*_values()
77  this->update_each = update_flags;
78 
79  // Store the quadrature
80  this->quad = q;
81 
82  // Resize the weights
84 
85  // see if we need the (transformation) shape function values
86  // and/or gradients and resize the necessary arrays
87  if (this->update_each &
90 
92  covariant.resize(n_original_q_points);
93 
95  contravariant.resize(n_original_q_points);
96 
98  volume_elements.resize(n_original_q_points);
99 }
100 
101 
102 template <int dim, int spacedim>
103 void
105  const UpdateFlags update_flags,
106  const Quadrature<dim> &q,
107  const unsigned int n_original_q_points)
108 {
109  initialize(update_flags, q, n_original_q_points);
110 
111  if (dim > 1)
112  {
114  {
115  aux.resize(dim - 1,
116  std::vector<Tensor<1, spacedim>>(n_original_q_points));
117 
118  // Compute tangentials to the unit cell.
119  for (const unsigned int i : GeometryInfo<dim>::face_indices())
120  {
121  unit_tangentials[i].resize(n_original_q_points);
122  std::fill(unit_tangentials[i].begin(),
123  unit_tangentials[i].end(),
125  if (dim > 2)
126  {
128  .resize(n_original_q_points);
129  std::fill(
131  .begin(),
133  .end(),
135  }
136  }
137  }
138  }
139 }
140 
141 
142 
143 template <int dim, int spacedim>
146 {}
147 
148 
149 
150 template <int dim, int spacedim>
151 std::unique_ptr<Mapping<dim, spacedim>>
153 {
154  return std::make_unique<MappingManifold<dim, spacedim>>(*this);
155 }
156 
157 
158 
159 template <int dim, int spacedim>
163  const Point<spacedim> &) const
164 {
165  Assert(false, ExcNotImplemented());
166  return {};
167 }
168 
169 
170 
171 template <int dim, int spacedim>
175  const Point<dim> & p) const
176 {
177  std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> vertices;
178  std::array<double, GeometryInfo<dim>::vertices_per_cell> weights;
179 
180  for (const unsigned int v : GeometryInfo<dim>::vertex_indices())
181  {
182  vertices[v] = cell->vertex(v);
184  }
185  return cell->get_manifold().get_new_point(
186  make_array_view(vertices.begin(), vertices.end()),
187  make_array_view(weights.begin(), weights.end()));
188 }
189 
190 
191 
192 // In the code below, GCC tries to instantiate MappingManifold<3,4> when
193 // seeing which of the overloaded versions of
194 // do_transform_real_to_unit_cell_internal() to call. This leads to bad
195 // error messages and, generally, nothing very good. Avoid this by ensuring
196 // that this class exists, but does not have an inner InternalData
197 // type, thereby ruling out the codim-1 version of the function
198 // below when doing overload resolution.
199 template <>
200 class MappingManifold<3, 4>
201 {};
202 
203 
204 
205 template <int dim, int spacedim>
208  const UpdateFlags in) const
209 {
210  // add flags if the respective quantities are necessary to compute
211  // what we need. note that some flags appear in both the conditions
212  // and in subsequent set operations. this leads to some circular
213  // logic. the only way to treat this is to iterate. since there are
214  // 5 if-clauses in the loop, it will take at most 5 iterations to
215  // converge. do them:
216  UpdateFlags out = in;
217  for (unsigned int i = 0; i < 5; ++i)
218  {
219  // The following is a little incorrect:
220  // If not applied on a face,
221  // update_boundary_forms does not
222  // make sense. On the other hand,
223  // it is necessary on a
224  // face. Currently,
225  // update_boundary_forms is simply
226  // ignored for the interior of a
227  // cell.
229  out |= update_boundary_forms;
230 
235 
236  if (out &
241 
242  // The contravariant transformation used in the Piola
243  // transformation, which requires the determinant of the Jacobi
244  // matrix of the transformation. Because we have no way of
245  // knowing here whether the finite elements wants to use the
246  // contravariant of the Piola transforms, we add the JxW values
247  // to the list of flags to be updated for each cell.
249  out |= update_JxW_values;
250 
251  if (out & update_normal_vectors)
252  out |= update_JxW_values;
253  }
254 
255  // Now throw an exception if we stumble upon something that was not
256  // implemented yet
261 
262  return out;
263 }
264 
265 
266 
267 template <int dim, int spacedim>
268 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
270  const Quadrature<dim> &q) const
271 {
272  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
273  std::make_unique<InternalData>();
274  auto &data = dynamic_cast<InternalData &>(*data_ptr);
275  data.initialize(this->requires_update_flags(update_flags), q, q.size());
276 
277  return data_ptr;
278 }
279 
280 
281 
282 template <int dim, int spacedim>
283 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
285  const UpdateFlags update_flags,
286  const hp::QCollection<dim - 1> &quadrature) const
287 {
288  AssertDimension(quadrature.size(), 1);
289 
290  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
291  std::make_unique<InternalData>();
292  auto &data = dynamic_cast<InternalData &>(*data_ptr);
293  data.initialize_face(this->requires_update_flags(update_flags),
295  ReferenceCells::get_hypercube<dim>(), quadrature[0]),
296  quadrature[0].size());
297 
298  return data_ptr;
299 }
300 
301 
302 
303 template <int dim, int spacedim>
304 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
306  const UpdateFlags update_flags,
307  const Quadrature<dim - 1> &quadrature) const
308 {
309  std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
310  std::make_unique<InternalData>();
311  auto &data = dynamic_cast<InternalData &>(*data_ptr);
312  data.initialize_face(this->requires_update_flags(update_flags),
314  ReferenceCells::get_hypercube<dim>(), quadrature),
315  quadrature.size());
316 
317  return data_ptr;
318 }
319 
320 
321 
322 namespace internal
323 {
324  namespace MappingManifoldImplementation
325  {
326  namespace
327  {
334  template <int dim, int spacedim>
335  void
337  const typename QProjector<dim>::DataSetDescriptor data_set,
338  const typename ::MappingManifold<dim, spacedim>::InternalData
339  & data,
340  std::vector<Point<spacedim>> &quadrature_points)
341  {
342  const UpdateFlags update_flags = data.update_each;
343 
344  AssertDimension(data.vertices.size(),
346 
347  if (update_flags & update_quadrature_points)
348  {
349  for (unsigned int point = 0; point < quadrature_points.size();
350  ++point)
351  {
352  quadrature_points[point] = data.manifold->get_new_point(
353  make_array_view(data.vertices),
355  data.cell_manifold_quadrature_weights[point + data_set]));
356  }
357  }
358  }
359 
360 
361 
368  template <int dim, int spacedim>
369  void
371  const typename ::QProjector<dim>::DataSetDescriptor data_set,
372  const typename ::MappingManifold<dim, spacedim>::InternalData
373  &data)
374  {
375  const UpdateFlags update_flags = data.update_each;
376 
377  if (update_flags & update_contravariant_transformation)
378  {
379  const unsigned int n_q_points = data.contravariant.size();
380 
381  std::fill(data.contravariant.begin(),
382  data.contravariant.end(),
384 
386  data.vertices.size());
387  for (unsigned int point = 0; point < n_q_points; ++point)
388  {
389  // Start by figuring out how to compute the direction in
390  // the reference space:
391  const Point<dim> &p = data.quad.point(point + data_set);
392 
393  // And get its image on the manifold:
394  const Point<spacedim> P = data.manifold->get_new_point(
395  make_array_view(data.vertices),
397  data.cell_manifold_quadrature_weights[point + data_set]));
398 
399  // To compute the Jacobian, we choose dim points aligned
400  // with the dim reference axes, which are still in the
401  // given cell, and ask for the tangent vector in these
402  // directions. Choosing the points is somewhat arbitrary,
403  // so we try to be smart and we pick points which are
404  // on the opposite quadrant w.r.t. the evaluation
405  // point.
406  for (unsigned int i = 0; i < dim; ++i)
407  {
408  const Point<dim> ei = Point<dim>::unit_vector(i);
409  const double pi = p[i];
410  Assert(pi >= 0 && pi <= 1.0,
412  "Was expecting a quadrature point "
413  "inside the unit reference element."));
414 
415  // In the length L, we store also the direction sign,
416  // which is positive, if the coordinate is < .5,
417  const double L = pi > .5 ? -pi : 1 - pi;
418 
419  const Point<dim> np(p + L * ei);
420 
421  // Get the weights to compute the np point in real space
422  for (const unsigned int j :
424  data.vertex_weights[j] =
426 
427  const Point<spacedim> NP = data.manifold->get_new_point(
428  make_array_view(data.vertices),
429  make_array_view(data.vertex_weights));
430 
431  const Tensor<1, spacedim> T =
432  data.manifold->get_tangent_vector(P, NP);
433 
434  for (unsigned int d = 0; d < spacedim; ++d)
435  data.contravariant[point][d][i] = T[d] / L;
436  }
437  }
438 
439  if (update_flags & update_covariant_transformation)
440  {
441  const unsigned int n_q_points = data.contravariant.size();
442  for (unsigned int point = 0; point < n_q_points; ++point)
443  {
444  data.covariant[point] =
445  (data.contravariant[point]).covariant_form();
446  }
447  }
448 
449  if (update_flags & update_volume_elements)
450  {
451  const unsigned int n_q_points = data.contravariant.size();
452  for (unsigned int point = 0; point < n_q_points; ++point)
453  data.volume_elements[point] =
454  data.contravariant[point].determinant();
455  }
456  }
457  }
458  } // namespace
459  } // namespace MappingManifoldImplementation
460 } // namespace internal
461 
462 
463 
464 template <int dim, int spacedim>
469  const Quadrature<dim> & quadrature,
470  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
472  &output_data) const
473 {
474  // ensure that the following static_cast is really correct:
475  Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
476  ExcInternalError());
477  const InternalData &data = static_cast<const InternalData &>(internal_data);
478 
479  const unsigned int n_q_points = quadrature.size();
480 
481  data.store_vertices(cell);
482  data.manifold = &(cell->get_manifold());
483 
485  spacedim>(
487  data,
488  output_data.quadrature_points);
489 
491  spacedim>(
493 
494  const UpdateFlags update_flags = data.update_each;
495  const std::vector<double> &weights = quadrature.get_weights();
496 
497  // Multiply quadrature weights by absolute value of Jacobian determinants or
498  // the area element g=sqrt(DX^t DX) in case of codim > 0
499 
500  if (update_flags & (update_normal_vectors | update_JxW_values))
501  {
502  AssertDimension(output_data.JxW_values.size(), n_q_points);
503 
504  Assert(!(update_flags & update_normal_vectors) ||
505  (output_data.normal_vectors.size() == n_q_points),
506  ExcDimensionMismatch(output_data.normal_vectors.size(),
507  n_q_points));
508 
509 
510  for (unsigned int point = 0; point < n_q_points; ++point)
511  {
512  if (dim == spacedim)
513  {
514  const double det = data.contravariant[point].determinant();
515 
516  // check for distorted cells.
517 
518  // TODO: this allows for anisotropies of up to 1e6 in 3D and
519  // 1e12 in 2D. might want to find a finer
520  // (dimension-independent) criterion
521  Assert(det > 1e-12 * Utilities::fixed_power<dim>(
522  cell->diameter() / std::sqrt(double(dim))),
524  cell->center(), det, point)));
525 
526  output_data.JxW_values[point] = weights[point] * det;
527  }
528  // if dim==spacedim, then there is no cell normal to
529  // compute. since this is for FEValues (and not FEFaceValues),
530  // there are also no face normals to compute
531  else // codim>0 case
532  {
533  Tensor<1, spacedim> DX_t[dim];
534  for (unsigned int i = 0; i < spacedim; ++i)
535  for (unsigned int j = 0; j < dim; ++j)
536  DX_t[j][i] = data.contravariant[point][i][j];
537 
538  Tensor<2, dim> G; // First fundamental form
539  for (unsigned int i = 0; i < dim; ++i)
540  for (unsigned int j = 0; j < dim; ++j)
541  G[i][j] = DX_t[i] * DX_t[j];
542 
543  output_data.JxW_values[point] =
544  std::sqrt(determinant(G)) * weights[point];
545 
546  if (update_flags & update_normal_vectors)
547  {
548  Assert(spacedim == dim + 1,
549  ExcMessage(
550  "There is no (unique) cell normal for " +
552  "-dimensional cells in " +
553  Utilities::int_to_string(spacedim) +
554  "-dimensional space. This only works if the "
555  "space dimension is one greater than the "
556  "dimensionality of the mesh cells."));
557 
558  if (dim == 1)
559  output_data.normal_vectors[point] =
560  cross_product_2d(-DX_t[0]);
561  else // dim == 2
562  output_data.normal_vectors[point] =
563  cross_product_3d(DX_t[0], DX_t[1]);
564 
565  output_data.normal_vectors[point] /=
566  output_data.normal_vectors[point].norm();
567 
568  if (cell->direction_flag() == false)
569  output_data.normal_vectors[point] *= -1.;
570  }
571  } // codim>0 case
572  }
573  }
574 
575 
576 
577  // copy values from InternalData to vector given by reference
578  if (update_flags & update_jacobians)
579  {
580  AssertDimension(output_data.jacobians.size(), n_q_points);
581  for (unsigned int point = 0; point < n_q_points; ++point)
582  output_data.jacobians[point] = data.contravariant[point];
583  }
584 
585  // copy values from InternalData to vector given by reference
586  if (update_flags & update_inverse_jacobians)
587  {
588  AssertDimension(output_data.inverse_jacobians.size(), n_q_points);
589  for (unsigned int point = 0; point < n_q_points; ++point)
590  output_data.inverse_jacobians[point] =
591  data.covariant[point].transpose();
592  }
593 
595 }
596 
597 
598 
599 namespace internal
600 {
601  namespace MappingManifoldImplementation
602  {
603  namespace
604  {
614  template <int dim, int spacedim>
615  void
617  const ::MappingManifold<dim, spacedim> &mapping,
618  const typename ::Triangulation<dim, spacedim>::cell_iterator
619  & cell,
620  const unsigned int face_no,
621  const unsigned int subface_no,
622  const unsigned int n_q_points,
623  const std::vector<double> &weights,
624  const typename ::MappingManifold<dim, spacedim>::InternalData
625  &data,
627  &output_data)
628  {
629  const UpdateFlags update_flags = data.update_each;
630 
631  if (update_flags & update_boundary_forms)
632  {
633  AssertDimension(output_data.boundary_forms.size(), n_q_points);
634  if (update_flags & update_normal_vectors)
635  AssertDimension(output_data.normal_vectors.size(), n_q_points);
636  if (update_flags & update_JxW_values)
637  AssertDimension(output_data.JxW_values.size(), n_q_points);
638 
639  // map the unit tangentials to the real cell. checking for d!=dim-1
640  // eliminates compiler warnings regarding unsigned int expressions <
641  // 0.
642  for (unsigned int d = 0; d != dim - 1; ++d)
643  {
645  data.unit_tangentials.size(),
646  ExcInternalError());
647  Assert(
648  data.aux[d].size() <=
649  data
650  .unit_tangentials[face_no +
652  .size(),
653  ExcInternalError());
654 
655  mapping.transform(
657  data
658  .unit_tangentials[face_no +
661  data,
662  make_array_view(data.aux[d]));
663  }
664 
665  // if dim==spacedim, we can use the unit tangentials to compute the
666  // boundary form by simply taking the cross product
667  if (dim == spacedim)
668  {
669  for (unsigned int i = 0; i < n_q_points; ++i)
670  switch (dim)
671  {
672  case 1:
673  // in 1d, we don't have access to any of the data.aux
674  // fields (because it has only dim-1 components), but we
675  // can still compute the boundary form by simply
676  // looking at the number of the face
677  output_data.boundary_forms[i][0] =
678  (face_no == 0 ? -1 : +1);
679  break;
680  case 2:
681  output_data.boundary_forms[i] =
682  cross_product_2d(data.aux[0][i]);
683  break;
684  case 3:
685  output_data.boundary_forms[i] =
686  cross_product_3d(data.aux[0][i], data.aux[1][i]);
687  break;
688  default:
689  Assert(false, ExcNotImplemented());
690  }
691  }
692  else //(dim < spacedim)
693  {
694  // in the codim-one case, the boundary form results from the
695  // cross product of all the face tangential vectors and the cell
696  // normal vector
697  //
698  // to compute the cell normal, use the same method used in
699  // fill_fe_values for cells above
700  AssertDimension(data.contravariant.size(), n_q_points);
701 
702  for (unsigned int point = 0; point < n_q_points; ++point)
703  {
704  switch (dim)
705  {
706  case 1:
707  {
708  // J is a tangent vector
709  output_data.boundary_forms[point] =
710  data.contravariant[point].transpose()[0];
711  output_data.boundary_forms[point] /=
712  (face_no == 0 ? -1. : +1.) *
713  output_data.boundary_forms[point].norm();
714 
715  break;
716  }
717 
718  case 2:
719  {
721  data.contravariant[point].transpose();
722 
723  Tensor<1, spacedim> cell_normal =
724  cross_product_3d(DX_t[0], DX_t[1]);
725  cell_normal /= cell_normal.norm();
726 
727  // then compute the face normal from the face
728  // tangent and the cell normal:
729  output_data.boundary_forms[point] =
730  cross_product_3d(data.aux[0][point], cell_normal);
731 
732  break;
733  }
734 
735  default:
736  Assert(false, ExcNotImplemented());
737  }
738  }
739  }
740 
741  if (update_flags & (update_normal_vectors | update_JxW_values))
742  for (unsigned int i = 0; i < output_data.boundary_forms.size();
743  ++i)
744  {
745  if (update_flags & update_JxW_values)
746  {
747  output_data.JxW_values[i] =
748  output_data.boundary_forms[i].norm() * weights[i];
749 
750  if (subface_no != numbers::invalid_unsigned_int)
751  {
752  const double area_ratio =
754  cell->subface_case(face_no), subface_no);
755  output_data.JxW_values[i] *= area_ratio;
756  }
757  }
758 
759  if (update_flags & update_normal_vectors)
760  output_data.normal_vectors[i] =
761  Point<spacedim>(output_data.boundary_forms[i] /
762  output_data.boundary_forms[i].norm());
763  }
764 
765  if (update_flags & update_jacobians)
766  for (unsigned int point = 0; point < n_q_points; ++point)
767  output_data.jacobians[point] = data.contravariant[point];
768 
769  if (update_flags & update_inverse_jacobians)
770  for (unsigned int point = 0; point < n_q_points; ++point)
771  output_data.inverse_jacobians[point] =
772  data.covariant[point].transpose();
773  }
774  }
775 
776 
783  template <int dim, int spacedim>
784  void
786  const ::MappingManifold<dim, spacedim> &mapping,
787  const typename ::Triangulation<dim, spacedim>::cell_iterator
788  & cell,
789  const unsigned int face_no,
790  const unsigned int subface_no,
791  const typename QProjector<dim>::DataSetDescriptor data_set,
792  const Quadrature<dim - 1> & quadrature,
793  const typename ::MappingManifold<dim, spacedim>::InternalData
794  &data,
796  &output_data)
797  {
798  data.store_vertices(cell);
799 
800  data.manifold = &cell->face(face_no)->get_manifold();
801 
802  maybe_compute_q_points<dim, spacedim>(data_set,
803  data,
804  output_data.quadrature_points);
805  maybe_update_Jacobians<dim, spacedim>(data_set, data);
806 
807  maybe_compute_face_data(mapping,
808  cell,
809  face_no,
810  subface_no,
811  quadrature.size(),
812  quadrature.get_weights(),
813  data,
814  output_data);
815  }
816 
817  template <int dim, int spacedim, int rank>
818  void
820  const ArrayView<const Tensor<rank, dim>> & input,
821  const MappingKind mapping_kind,
822  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
823  const ArrayView<Tensor<rank, spacedim>> & output)
824  {
825  AssertDimension(input.size(), output.size());
826  Assert((dynamic_cast<const typename ::
827  MappingManifold<dim, spacedim>::InternalData *>(
828  &mapping_data) != nullptr),
829  ExcInternalError());
830  const typename ::MappingManifold<dim, spacedim>::InternalData
831  &data =
832  static_cast<const typename ::MappingManifold<dim, spacedim>::
833  InternalData &>(mapping_data);
834 
835  switch (mapping_kind)
836  {
838  {
839  Assert(
840  data.update_each & update_contravariant_transformation,
842  "update_contravariant_transformation"));
843 
844  for (unsigned int i = 0; i < output.size(); ++i)
845  output[i] =
846  apply_transformation(data.contravariant[i], input[i]);
847 
848  return;
849  }
850 
851  case mapping_piola:
852  {
853  Assert(
854  data.update_each & update_contravariant_transformation,
856  "update_contravariant_transformation"));
857  Assert(
858  data.update_each & update_volume_elements,
860  "update_volume_elements"));
861  Assert(rank == 1, ExcMessage("Only for rank 1"));
862  if (rank != 1)
863  return;
864 
865  for (unsigned int i = 0; i < output.size(); ++i)
866  {
867  output[i] =
868  apply_transformation(data.contravariant[i], input[i]);
869  output[i] /= data.volume_elements[i];
870  }
871  return;
872  }
873  // We still allow this operation as in the
874  // reference cell Derivatives are Tensor
875  // rather than DerivativeForm
876  case mapping_covariant:
877  {
878  Assert(
879  data.update_each & update_contravariant_transformation,
881  "update_covariant_transformation"));
882 
883  for (unsigned int i = 0; i < output.size(); ++i)
884  output[i] = apply_transformation(data.covariant[i], input[i]);
885 
886  return;
887  }
888 
889  default:
890  Assert(false, ExcNotImplemented());
891  }
892  }
893 
894 
895  template <int dim, int spacedim, int rank>
896  void
898  const ArrayView<const Tensor<rank, dim>> & input,
899  const MappingKind mapping_kind,
900  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
901  const ArrayView<Tensor<rank, spacedim>> & output)
902  {
903  AssertDimension(input.size(), output.size());
904  Assert((dynamic_cast<const typename ::
905  MappingManifold<dim, spacedim>::InternalData *>(
906  &mapping_data) != nullptr),
907  ExcInternalError());
908  const typename ::MappingManifold<dim, spacedim>::InternalData
909  &data =
910  static_cast<const typename ::MappingManifold<dim, spacedim>::
911  InternalData &>(mapping_data);
912 
913  switch (mapping_kind)
914  {
916  {
917  Assert(
918  data.update_each & update_covariant_transformation,
920  "update_covariant_transformation"));
921  Assert(
922  data.update_each & update_contravariant_transformation,
924  "update_contravariant_transformation"));
925  Assert(rank == 2, ExcMessage("Only for rank 2"));
926 
927  for (unsigned int i = 0; i < output.size(); ++i)
928  {
930  apply_transformation(data.contravariant[i],
931  transpose(input[i]));
932  output[i] =
933  apply_transformation(data.covariant[i], A.transpose());
934  }
935 
936  return;
937  }
938 
940  {
941  Assert(
942  data.update_each & update_covariant_transformation,
944  "update_covariant_transformation"));
945  Assert(rank == 2, ExcMessage("Only for rank 2"));
946 
947  for (unsigned int i = 0; i < output.size(); ++i)
948  {
950  apply_transformation(data.covariant[i],
951  transpose(input[i]));
952  output[i] =
953  apply_transformation(data.covariant[i], A.transpose());
954  }
955 
956  return;
957  }
958 
960  {
961  Assert(
962  data.update_each & update_covariant_transformation,
964  "update_covariant_transformation"));
965  Assert(
966  data.update_each & update_contravariant_transformation,
968  "update_contravariant_transformation"));
969  Assert(
970  data.update_each & update_volume_elements,
972  "update_volume_elements"));
973  Assert(rank == 2, ExcMessage("Only for rank 2"));
974 
975  for (unsigned int i = 0; i < output.size(); ++i)
976  {
978  apply_transformation(data.covariant[i], input[i]);
980  apply_transformation(data.contravariant[i],
981  A.transpose());
982 
983  output[i] = transpose(T);
984  output[i] /= data.volume_elements[i];
985  }
986 
987  return;
988  }
989 
990  default:
991  Assert(false, ExcNotImplemented());
992  }
993  }
994 
995 
996 
997  template <int dim, int spacedim>
998  void
1000  const ArrayView<const Tensor<3, dim>> & input,
1001  const MappingKind mapping_kind,
1002  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1003  const ArrayView<Tensor<3, spacedim>> & output)
1004  {
1005  AssertDimension(input.size(), output.size());
1006  Assert((dynamic_cast<const typename ::
1007  MappingManifold<dim, spacedim>::InternalData *>(
1008  &mapping_data) != nullptr),
1009  ExcInternalError());
1010  const typename ::MappingManifold<dim, spacedim>::InternalData
1011  &data =
1012  static_cast<const typename ::MappingManifold<dim, spacedim>::
1013  InternalData &>(mapping_data);
1014 
1015  switch (mapping_kind)
1016  {
1018  {
1019  Assert(
1020  data.update_each & update_covariant_transformation,
1022  "update_covariant_transformation"));
1023  Assert(
1024  data.update_each & update_contravariant_transformation,
1026  "update_contravariant_transformation"));
1027 
1028  for (unsigned int q = 0; q < output.size(); ++q)
1029  for (unsigned int i = 0; i < spacedim; ++i)
1030  {
1031  double tmp1[dim][dim];
1032  for (unsigned int J = 0; J < dim; ++J)
1033  for (unsigned int K = 0; K < dim; ++K)
1034  {
1035  tmp1[J][K] =
1036  data.contravariant[q][i][0] * input[q][0][J][K];
1037  for (unsigned int I = 1; I < dim; ++I)
1038  tmp1[J][K] +=
1039  data.contravariant[q][i][I] * input[q][I][J][K];
1040  }
1041  for (unsigned int j = 0; j < spacedim; ++j)
1042  {
1043  double tmp2[dim];
1044  for (unsigned int K = 0; K < dim; ++K)
1045  {
1046  tmp2[K] = data.covariant[q][j][0] * tmp1[0][K];
1047  for (unsigned int J = 1; J < dim; ++J)
1048  tmp2[K] += data.covariant[q][j][J] * tmp1[J][K];
1049  }
1050  for (unsigned int k = 0; k < spacedim; ++k)
1051  {
1052  output[q][i][j][k] =
1053  data.covariant[q][k][0] * tmp2[0];
1054  for (unsigned int K = 1; K < dim; ++K)
1055  output[q][i][j][k] +=
1056  data.covariant[q][k][K] * tmp2[K];
1057  }
1058  }
1059  }
1060  return;
1061  }
1062 
1064  {
1065  Assert(
1066  data.update_each & update_covariant_transformation,
1068  "update_covariant_transformation"));
1069 
1070  for (unsigned int q = 0; q < output.size(); ++q)
1071  for (unsigned int i = 0; i < spacedim; ++i)
1072  {
1073  double tmp1[dim][dim];
1074  for (unsigned int J = 0; J < dim; ++J)
1075  for (unsigned int K = 0; K < dim; ++K)
1076  {
1077  tmp1[J][K] =
1078  data.covariant[q][i][0] * input[q][0][J][K];
1079  for (unsigned int I = 1; I < dim; ++I)
1080  tmp1[J][K] +=
1081  data.covariant[q][i][I] * input[q][I][J][K];
1082  }
1083  for (unsigned int j = 0; j < spacedim; ++j)
1084  {
1085  double tmp2[dim];
1086  for (unsigned int K = 0; K < dim; ++K)
1087  {
1088  tmp2[K] = data.covariant[q][j][0] * tmp1[0][K];
1089  for (unsigned int J = 1; J < dim; ++J)
1090  tmp2[K] += data.covariant[q][j][J] * tmp1[J][K];
1091  }
1092  for (unsigned int k = 0; k < spacedim; ++k)
1093  {
1094  output[q][i][j][k] =
1095  data.covariant[q][k][0] * tmp2[0];
1096  for (unsigned int K = 1; K < dim; ++K)
1097  output[q][i][j][k] +=
1098  data.covariant[q][k][K] * tmp2[K];
1099  }
1100  }
1101  }
1102 
1103  return;
1104  }
1105 
1106  case mapping_piola_hessian:
1107  {
1108  Assert(
1109  data.update_each & update_covariant_transformation,
1111  "update_covariant_transformation"));
1112  Assert(
1113  data.update_each & update_contravariant_transformation,
1115  "update_contravariant_transformation"));
1116  Assert(
1117  data.update_each & update_volume_elements,
1119  "update_volume_elements"));
1120 
1121  for (unsigned int q = 0; q < output.size(); ++q)
1122  for (unsigned int i = 0; i < spacedim; ++i)
1123  {
1124  double factor[dim];
1125  for (unsigned int I = 0; I < dim; ++I)
1126  factor[I] =
1127  data.contravariant[q][i][I] / data.volume_elements[q];
1128  double tmp1[dim][dim];
1129  for (unsigned int J = 0; J < dim; ++J)
1130  for (unsigned int K = 0; K < dim; ++K)
1131  {
1132  tmp1[J][K] = factor[0] * input[q][0][J][K];
1133  for (unsigned int I = 1; I < dim; ++I)
1134  tmp1[J][K] += factor[I] * input[q][I][J][K];
1135  }
1136  for (unsigned int j = 0; j < spacedim; ++j)
1137  {
1138  double tmp2[dim];
1139  for (unsigned int K = 0; K < dim; ++K)
1140  {
1141  tmp2[K] = data.covariant[q][j][0] * tmp1[0][K];
1142  for (unsigned int J = 1; J < dim; ++J)
1143  tmp2[K] += data.covariant[q][j][J] * tmp1[J][K];
1144  }
1145  for (unsigned int k = 0; k < spacedim; ++k)
1146  {
1147  output[q][i][j][k] =
1148  data.covariant[q][k][0] * tmp2[0];
1149  for (unsigned int K = 1; K < dim; ++K)
1150  output[q][i][j][k] +=
1151  data.covariant[q][k][K] * tmp2[K];
1152  }
1153  }
1154  }
1155 
1156  return;
1157  }
1158 
1159  default:
1160  Assert(false, ExcNotImplemented());
1161  }
1162  }
1163 
1164 
1165 
1166  template <int dim, int spacedim, int rank>
1167  void
1169  const ArrayView<const DerivativeForm<rank, dim, spacedim>> &input,
1170  const MappingKind mapping_kind,
1171  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1172  const ArrayView<Tensor<rank + 1, spacedim>> & output)
1173  {
1174  AssertDimension(input.size(), output.size());
1175  Assert((dynamic_cast<const typename ::
1176  MappingManifold<dim, spacedim>::InternalData *>(
1177  &mapping_data) != nullptr),
1178  ExcInternalError());
1179  const typename ::MappingManifold<dim, spacedim>::InternalData
1180  &data =
1181  static_cast<const typename ::MappingManifold<dim, spacedim>::
1182  InternalData &>(mapping_data);
1183 
1184  switch (mapping_kind)
1185  {
1186  case mapping_covariant:
1187  {
1188  Assert(
1189  data.update_each & update_contravariant_transformation,
1191  "update_covariant_transformation"));
1192 
1193  for (unsigned int i = 0; i < output.size(); ++i)
1194  output[i] = apply_transformation(data.covariant[i], input[i]);
1195 
1196  return;
1197  }
1198  default:
1199  Assert(false, ExcNotImplemented());
1200  }
1201  }
1202  } // namespace
1203  } // namespace MappingManifoldImplementation
1204 } // namespace internal
1205 
1206 
1207 
1208 template <int dim, int spacedim>
1209 void
1212  const unsigned int face_no,
1213  const hp::QCollection<dim - 1> & quadrature,
1214  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
1216  &output_data) const
1217 {
1218  AssertDimension(quadrature.size(), 1);
1219 
1220  // ensure that the following cast is really correct:
1221  Assert((dynamic_cast<const InternalData *>(&internal_data) != nullptr),
1222  ExcInternalError());
1223  const InternalData &data = static_cast<const InternalData &>(internal_data);
1224 
1226  *this,
1227  cell,
1228  face_no,
1231  ReferenceCells::get_hypercube<dim>(),
1232  face_no,
1233  cell->face_orientation(face_no),
1234  cell->face_flip(face_no),
1235  cell->face_rotation(face_no),
1236  quadrature[0].size()),
1237  quadrature[0],
1238  data,
1239  output_data);
1240 }
1241 
1242 
1243 
1244 template <int dim, int spacedim>
1245 void
1248  const unsigned int face_no,
1249  const unsigned int subface_no,
1250  const Quadrature<dim - 1> & quadrature,
1251  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
1253  &output_data) const
1254 {
1255  // ensure that the following cast is really correct:
1256  Assert((dynamic_cast<const InternalData *>(&internal_data) != nullptr),
1257  ExcInternalError());
1258  const InternalData &data = static_cast<const InternalData &>(internal_data);
1259 
1261  *this,
1262  cell,
1263  face_no,
1264  subface_no,
1266  ReferenceCells::get_hypercube<dim>(),
1267  face_no,
1268  subface_no,
1269  cell->face_orientation(face_no),
1270  cell->face_flip(face_no),
1271  cell->face_rotation(face_no),
1272  quadrature.size(),
1273  cell->subface_case(face_no)),
1274  quadrature,
1275  data,
1276  output_data);
1277 }
1278 
1279 
1280 
1281 template <int dim, int spacedim>
1282 void
1284  const ArrayView<const Tensor<1, dim>> & input,
1285  const MappingKind mapping_kind,
1286  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1287  const ArrayView<Tensor<1, spacedim>> & output) const
1288 {
1290  mapping_kind,
1291  mapping_data,
1292  output);
1293 }
1294 
1295 
1296 
1297 template <int dim, int spacedim>
1298 void
1300  const ArrayView<const DerivativeForm<1, dim, spacedim>> &input,
1301  const MappingKind mapping_kind,
1302  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1303  const ArrayView<Tensor<2, spacedim>> & output) const
1304 {
1306  input, mapping_kind, mapping_data, output);
1307 }
1308 
1309 
1310 
1311 template <int dim, int spacedim>
1312 void
1314  const ArrayView<const Tensor<2, dim>> & input,
1315  const MappingKind mapping_kind,
1316  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1317  const ArrayView<Tensor<2, spacedim>> & output) const
1318 {
1319  switch (mapping_kind)
1320  {
1321  case mapping_contravariant:
1323  mapping_kind,
1324  mapping_data,
1325  output);
1326  return;
1327 
1332  input, mapping_kind, mapping_data, output);
1333  return;
1334  default:
1335  Assert(false, ExcNotImplemented());
1336  }
1337 }
1338 
1339 
1340 
1341 template <int dim, int spacedim>
1342 void
1344  const ArrayView<const DerivativeForm<2, dim, spacedim>> &input,
1345  const MappingKind mapping_kind,
1346  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1347  const ArrayView<Tensor<3, spacedim>> & output) const
1348 {
1349  AssertDimension(input.size(), output.size());
1350  Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
1351  ExcInternalError());
1352  const InternalData &data = static_cast<const InternalData &>(mapping_data);
1353 
1354  switch (mapping_kind)
1355  {
1357  {
1360  "update_covariant_transformation"));
1361 
1362  for (unsigned int q = 0; q < output.size(); ++q)
1363  for (unsigned int i = 0; i < spacedim; ++i)
1364  for (unsigned int j = 0; j < spacedim; ++j)
1365  {
1366  double tmp[dim];
1367  for (unsigned int K = 0; K < dim; ++K)
1368  {
1369  tmp[K] = data.covariant[q][j][0] * input[q][i][0][K];
1370  for (unsigned int J = 1; J < dim; ++J)
1371  tmp[K] += data.covariant[q][j][J] * input[q][i][J][K];
1372  }
1373  for (unsigned int k = 0; k < spacedim; ++k)
1374  {
1375  output[q][i][j][k] = data.covariant[q][k][0] * tmp[0];
1376  for (unsigned int K = 1; K < dim; ++K)
1377  output[q][i][j][k] += data.covariant[q][k][K] * tmp[K];
1378  }
1379  }
1380  return;
1381  }
1382 
1383  default:
1384  Assert(false, ExcNotImplemented());
1385  }
1386 }
1387 
1388 
1389 
1390 template <int dim, int spacedim>
1391 void
1393  const ArrayView<const Tensor<3, dim>> & input,
1394  const MappingKind mapping_kind,
1395  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
1396  const ArrayView<Tensor<3, spacedim>> & output) const
1397 {
1398  switch (mapping_kind)
1399  {
1400  case mapping_piola_hessian:
1404  input, mapping_kind, mapping_data, output);
1405  return;
1406  default:
1407  Assert(false, ExcNotImplemented());
1408  }
1409 }
1410 
1411 //--------------------------- Explicit instantiations -----------------------
1412 #include "mapping_manifold.inst"
1413 
1414 
Transformed quadrature weights.
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
constexpr Number determinant(const SymmetricTensor< 2, dim, Number > &)
static const unsigned int invalid_unsigned_int
Definition: types.h:196
static const char L
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
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
SmartPointer< const Manifold< dim, spacedim > > manifold
Contravariant transformation.
unsigned int size() const
Definition: collection.h:263
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
const std::vector< double > & get_weights() const
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
Definition: generators.cc:485
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
std::vector< Tensor< 1, spacedim > > boundary_forms
Volume element.
static Point< dim, Number > unit_vector(const unsigned int i)
Outer normal vector, not normalized.
Tensor< 1, spacedim, typename ProductType< Number1, Number2 >::type > apply_transformation(const DerivativeForm< 1, dim, spacedim, Number1 > &grad_F, const Tensor< 1, dim, Number2 > &d_x)
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
Determinant of the Jacobian.
Transformed quadrature points.
void store_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
virtual std::size_t memory_consumption() const override
void transform_differential_forms(const ArrayView< const DerivativeForm< rank, dim, spacedim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank+1, spacedim >> &output)
MappingKind
Definition: mapping.h:71
static DataSetDescriptor cell()
Definition: qprojector.h:563
void transform_gradients(const ArrayView< const Tensor< rank, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank, spacedim >> &output)
void maybe_update_Jacobians(const CellSimilarity::Similarity cell_similarity, const typename ::QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data)
Triangulation< dim, spacedim >::cell_iterator cell
constexpr Tensor< 1, dim, typename ProductType< Number1, Number2 >::type > cross_product_3d(const Tensor< 1, dim, Number1 > &src1, const Tensor< 1, dim, Number2 > &src2)
Definition: tensor.h:2650
std::vector< Point< spacedim > > vertices
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
std::vector< std::vector< Tensor< 1, spacedim > > > aux
static ::ExceptionBase & ExcMessage(std::string arg1)
static const char T
#define Assert(cond, exc)
Definition: exceptions.h:1461
UpdateFlags
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
Definition: mapping.h:310
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition: array_view.h:698
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
VectorType::value_type * end(VectorType &V)
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
DerivativeForm< 1, spacedim, dim, Number > transpose() const
Gradient of volume element.
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:473
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:185
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::vector< double > volume_elements
unsigned int size() const
std::vector< double > vertex_weights
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
void compute_manifold_quadrature_weights(const Quadrature< dim > &quadrature)
std::vector< Point< spacedim > > quadrature_points
static const char A
void do_fill_fe_face_values(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const typename QProjector< dim >::DataSetDescriptor data_set, const Quadrature< dim - 1 > &quadrature, const typename ::MappingQ< dim, spacedim >::InternalData &data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
VectorType::value_type * begin(VectorType &V)
Normal vectors.
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
void maybe_compute_face_data(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int n_q_points, const std::vector< double > &weights, const typename ::MappingQ< dim, spacedim >::InternalData &data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
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
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
static ::ExceptionBase & ExcNotImplemented()
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
void maybe_compute_q_points(const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Point< spacedim >> &quadrature_points)
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
numbers::NumberTraits< Number >::real_type norm() const
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
void transform_hessians(const ArrayView< const Tensor< 3, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< 3, spacedim >> &output)
constexpr Tensor< 1, dim, Number > cross_product_2d(const Tensor< 1, dim, Number > &src)
Definition: tensor.h:2625
void transform_fields(const ArrayView< const Tensor< rank, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank, spacedim >> &output)
static double subface_ratio(const internal::SubfaceCase< dim > &subface_case, const unsigned int subface_no)
Covariant transformation.
std::vector< std::vector< double > > cell_manifold_quadrature_weights
std::vector< Tensor< 1, spacedim > > normal_vectors
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
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
MappingManifold()=default
UpdateFlags update_each
Definition: mapping.h:659
static ::ExceptionBase & ExcInternalError()