1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2000 - 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 16 #ifndef dealii_mapping_h 17 #define dealii_mapping_h 18 19 20 #include <deal.II/base/config.h> 21 22 #include <deal.II/base/array_view.h> 23 #include <deal.II/base/derivative_form.h> 24 25 #include <deal.II/fe/fe_update_flags.h> 26 27 #include <deal.II/grid/tria.h> 28 29 #include <array> 30 #include <cmath> 31 #include <memory> 32 33 DEAL_II_NAMESPACE_OPEN 34 35 template <typename ElementType, typename MemorySpaceType> 36 class ArrayView; 37 template <int dim> 38 class Quadrature; 39 template <int dim, int spacedim> 40 class FEValues; 41 template <int dim, int spacedim> 42 class FEValuesBase; 43 template <int dim, int spacedim> 44 class FEValues; 45 template <int dim, int spacedim> 46 class FEFaceValues; 47 template <int dim, int spacedim> 48 class FESubfaceValues; 49 50 51 /** 52 * The transformation kind used for the Mapping::transform() functions. 53 * 54 * Special finite elements may need special Mapping from the reference cell to 55 * the actual mesh cell. In order to be most flexible, this enum provides an 56 * extensible interface for arbitrary transformations. Nevertheless, these 57 * must be implemented in the transform() functions of inheriting classes in 58 * order to work. 59 * 60 * @ingroup mapping 61 */ 62 enum MappingKind 63 { 64 /** 65 * No mapping, i.e., shape functions are not mapped from a reference cell 66 * but instead are defined right on the real-space cell. 67 */ 68 mapping_none = 0x0000, 69 70 /** 71 * Covariant mapping (see Mapping::transform() for details). 72 */ 73 mapping_covariant = 0x0001, 74 75 /** 76 * Contravariant mapping (see Mapping::transform() for details). 77 */ 78 mapping_contravariant = 0x0002, 79 80 /** 81 * Mapping of the gradient of a covariant vector field (see 82 * Mapping::transform() for details). 83 */ 84 mapping_covariant_gradient = 0x0003, 85 86 /** 87 * Mapping of the gradient of a contravariant vector field (see 88 * Mapping::transform() for details). 89 */ 90 mapping_contravariant_gradient = 0x0004, 91 92 /** 93 * The Piola transform usually used for Hdiv elements. Piola transform is 94 * the standard transformation of vector valued elements in H<sup>div</sup>. 95 * It amounts to a contravariant transformation scaled by the inverse of the 96 * volume element. 97 */ 98 mapping_piola = 0x0100, 99 100 /** 101 * Transformation for the gradient of a vector field corresponding to a 102 * mapping_piola transformation (see Mapping::transform() for details). 103 */ 104 mapping_piola_gradient = 0x0101, 105 106 /** 107 * The mapping used for Nedelec elements. 108 * 109 * Curl-conforming elements are mapped as covariant vectors. Nevertheless, 110 * we introduce a separate mapping kind, such that we can use the same flag 111 * for the vector and its gradient (see Mapping::transform() for details). 112 */ 113 mapping_nedelec = 0x0200, 114 115 /** 116 * The mapping used for Raviart-Thomas elements. 117 */ 118 mapping_raviart_thomas = 0x0300, 119 120 /** 121 * The mapping used for BDM elements. 122 */ 123 mapping_bdm = mapping_raviart_thomas, 124 125 /** 126 * The mappings for 2-forms and third order tensors. 127 * 128 * These are mappings typpically applied to hessians transformed to the 129 * reference cell. 130 * 131 * Mapping of the hessian of a covariant vector field (see 132 * Mapping::transform() for details). 133 */ 134 mapping_covariant_hessian, 135 136 /** 137 * Mapping of the hessian of a contravariant vector field (see 138 * Mapping::transform() for details). 139 */ 140 mapping_contravariant_hessian, 141 142 /** 143 * Mapping of the hessian of a piola vector field (see Mapping::transform() 144 * for details). 145 */ 146 mapping_piola_hessian 147 }; 148 149 150 /** 151 * @short Abstract base class for mapping classes. 152 * 153 * This class declares the interface for the functionality to describe 154 * mappings from the reference (unit) cell to a cell in real space, as well as 155 * for filling the information necessary to use the FEValues, FEFaceValues, 156 * and FESubfaceValues classes. Concrete implementations of these interfaces 157 * are provided in derived classes. 158 * 159 * <h3>Mathematics of the mapping</h3> 160 * 161 * The mapping is a transformation $\mathbf x = \mathbf F_K(\hat{\mathbf x})$ 162 * which maps points $\hat{\mathbf x}$ in the reference cell 163 * $[0,1]^\text{dim}$ to points $\mathbf x$ in the actual grid cell 164 * $K\subset{\mathbb R}^\text{spacedim}$. Many of the applications of such 165 * mappings require the Jacobian of this mapping, $J(\hat{\mathbf x}) = 166 * \hat\nabla {\mathbf F}_K(\hat{\mathbf x})$. For instance, if 167 * dim=spacedim=2, we have 168 * @f[ 169 * J(\hat{\mathbf x}) = \left(\begin{matrix} 170 * \frac{\partial x}{\partial \hat x} & \frac{\partial x}{\partial \hat y} 171 * \\ 172 * \frac{\partial y}{\partial \hat x} & \frac{\partial y}{\partial \hat y} 173 * \end{matrix}\right) 174 * @f] 175 * 176 * <h4>%Mapping of scalar functions</h4> 177 * 178 * The shape functions of scalar finite elements are typically defined on a 179 * reference cell and are then simply mapped according to the rule 180 * @f[ 181 * \varphi(\mathbf x) = \varphi\bigl(\mathbf F_K(\hat{\mathbf x})\bigr) 182 * = \hat \varphi(\hat{\mathbf x}). 183 * @f] 184 * 185 * 186 * <h4>%Mapping of integrals</h4> 187 * 188 * Using simply a change of variables, integrals of scalar functions over a 189 * cell $K$ can be expressed as an integral over the reference cell $\hat K$. 190 * Specifically, The volume form $d\hat x$ is transformed so that 191 * @f[ 192 * \int_K u(\mathbf x)\,dx = \int_{\hat K} \hat 193 * u(\hat{\mathbf x}) \left|\text{det}J(\hat{\mathbf x})\right| 194 * \,d\hat x. 195 * @f] 196 * 197 * In expressions where such integrals are approximated by quadrature, this 198 * then leads to terms of the form 199 * @f[ 200 * \int_K u(\mathbf x)\,dx 201 * \approx 202 * \sum_{q} 203 * \hat u(\hat{\mathbf x}_q) 204 * \underbrace{\left|\text{det}J(\hat{\mathbf x}_q)\right| w_q}_{=: 205 * \text{JxW}_q}. 206 * @f] 207 * Here, the weights $\text{JxW}_q$ of each quadrature point (where <i>JxW</i> 208 * mnemonically stands for <i>Jacobian times Quadrature Weights</i>) take the 209 * role of the $dx$ in the original integral. Consequently, they appear in all 210 * code that computes integrals approximated by quadrature, and are accessed 211 * by FEValues::JxW(). 212 * 213 * @todo Document what happens in the codimension-1 case. 214 * 215 * 216 * <h4>%Mapping of vector fields, differential forms and gradients of vector 217 * fields</h4> 218 * 219 * The transformation of vector fields or differential forms (gradients of 220 * scalar functions) $\mathbf v$, and gradients of vector fields $\mathbf T$ 221 * follows the general form 222 * 223 * @f[ 224 * \mathbf v(\mathbf x) = \mathbf A(\hat{\mathbf x}) 225 * \hat{\mathbf v}(\hat{\mathbf x}), 226 * \qquad 227 * \mathbf T(\mathbf x) = \mathbf A(\hat{\mathbf x}) 228 * \hat{\mathbf T}(\hat{\mathbf x}) \mathbf B(\hat{\mathbf x}). 229 * @f] 230 * The differential forms <b>A</b> and <b>B</b> are determined by the kind of 231 * object being transformed. These transformations are performed through the 232 * transform() functions, and the type of object being transformed is 233 * specified by their MappingKind argument. See the documentation there for 234 * possible choices. 235 * 236 * <h4>Derivatives of the mapping</h4> 237 * 238 * Some applications require the derivatives of the mapping, of which the 239 * first order derivative is the mapping Jacobian, $J_{iJ}(\hat{\mathbf 240 * x})=\frac{\partial x_i}{\partial \hat x_J}$, described above. Higher order 241 * derivatives of the mapping are similarly defined, for example the Jacobian 242 * derivative, $\hat H_{iJK}(\hat{\mathbf x}) = \frac{\partial^2 243 * x_i}{\partial \hat x_J \partial \hat x_K}$, and the Jacobian second 244 * derivative, $\hat K_{iJKL}(\hat{\mathbf x}) = \frac{\partial^3 245 * x_i}{\partial \hat x_J \partial \hat x_K \partial \hat x_L}$. It is also 246 * useful to define the "pushed-forward" versions of the higher order 247 * derivatives: the Jacobian pushed-forward derivative, $H_{ijk}(\hat{\mathbf 248 * x}) = \frac{\partial^2 x_i}{\partial \hat x_J \partial \hat 249 * x_K}(J_{jJ})^{-1}(J_{kK})^{-1}$, and the Jacobian pushed-forward second 250 * derivative, $K_{ijkl}(\hat{\mathbf x}) = \frac{\partial^3 x_i}{\partial 251 * \hat x_J \partial \hat x_K \partial \hat 252 * x_L}(J_{jJ})^{-1}(J_{kK})^{-1}(J_{lL})^{-1}$. These pushed-forward versions 253 * can be used to compute the higher order derivatives of functions defined on 254 * the reference cell with respect to the real cell coordinates. For instance, 255 * the Jacobian derivative with respect to the real cell coordinates is given 256 * by: 257 * 258 * @f[ 259 * \frac{\partial}{\partial x_j}\left[J_{iJ}(\hat{\mathbf x})\right] = 260 * H_{ikn}(\hat{\mathbf x})J_{nJ}(\hat{\mathbf x}), 261 * @f] 262 * and the derivative of the Jacobian inverse with respect to the real cell 263 * coordinates is similarly given by: 264 * @f[ 265 * \frac{\partial}{\partial x_j}\left[\left(J_{iJ}(\hat{\mathbf 266 * x})\right)^{-1}\right] = -H_{nik}(\hat{\mathbf x})\left(J_{nJ}(\hat{\mathbf 267 * x})\right)^{-1}. 268 * @f] 269 * 270 * In a similar fashion, higher order derivatives, with respect to the real 271 * cell coordinates, of functions defined on the reference cell can be defined 272 * using the Jacobian pushed-forward higher-order derivatives. For example, 273 * the derivative, with respect to the real cell coordinates, of the Jacobian 274 * pushed-forward derivative is given by: 275 * 276 * @f[ 277 * \frac{\partial}{\partial x_l}\left[H_{ijk}(\hat{\mathbf x})\right] = 278 * K_{ijkl}(\hat{\mathbf x}) -H_{mjl}(\hat{\mathbf x})H_{imk}(\hat{\mathbf 279 * x})-H_{mkl}(\hat{\mathbf x})H_{imj}(\hat{\mathbf x}). 280 * @f] 281 * 282 * <h3>References</h3> 283 * 284 * A general publication on differential geometry and finite elements is the 285 * survey 286 * <ul> 287 * <li>Douglas N. Arnold, Richard S. Falk, and Ragnar Winther. <i>Finite 288 * element exterior calculus: from Hodge theory to numerical stability.</i> 289 * Bull. Amer. Math. Soc. (N.S.), 47:281-354, 2010. <a 290 * href="http://dx.doi.org/10.1090/S0273-0979-10-01278-4">DOI: 291 * 10.1090/S0273-0979-10-01278-4</a>. 292 * </ul> 293 * 294 * The description of the Piola transform has been taken from the <a 295 * href="http://www.math.uh.edu/~rohop/spring_05/downloads/">lecture notes</a> 296 * by Ronald H. W. Hoppe, University of Houston, Chapter 7. 297 * 298 * @ingroup mapping 299 */ 300 template <int dim, int spacedim = dim> 301 class Mapping : public Subscriptor 302 { 303 public: 304 /** 305 * Virtual destructor. 306 */ 307 virtual ~Mapping() override = default; 308 309 /** 310 * Return a pointer to a copy of the present object. The caller of this copy 311 * then assumes ownership of it. 312 * 313 * The function is declared abstract virtual in this base class, and derived 314 * classes will have to implement it. 315 * 316 * This function is mainly used by the hp::MappingCollection class. 317 */ 318 virtual std::unique_ptr<Mapping<dim, spacedim>> 319 clone() const = 0; 320 321 /** 322 * Return the mapped vertices of a cell. 323 * 324 * Most of the time, these values will simply be the coordinates of the 325 * vertices of a cell as returned by <code>cell-@>vertex(v)</code> for 326 * vertex <code>v</code>, i.e., information stored by the triangulation. 327 * However, there are also mappings that add displacements or choose 328 * completely different locations, e.g., MappingQEulerian, 329 * MappingQ1Eulerian, or MappingFEField. 330 * 331 * The default implementation of this function simply returns the 332 * information stored by the triangulation, i.e., 333 * <code>cell-@>vertex(v)</code>. 334 */ 335 virtual std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> 336 get_vertices( 337 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const; 338 339 /** 340 * Return the mapped center of a cell. 341 * 342 * If you are using a (bi-,tri-)linear mapping that preserves vertex 343 * locations, this function simply returns the value also produced by 344 * `cell->center()`. However, there are also mappings that add displacements 345 * or choose completely different locations, e.g., MappingQEulerian, 346 * MappingQ1Eulerian, or MappingFEField, and mappings based on high order 347 * polynomials, for which the center may not coincide with the average of 348 * the vertex locations. 349 * 350 * By default, this function returns the push forward of the center of the 351 * reference cell. If the parameter 352 * @p map_center_of_reference_cell is set to false, than the return value 353 * will be the average of the vertex locations, as returned by the 354 * get_vertices() method. 355 * 356 * @param[in] cell The cell for which you want to compute the center 357 * @param[in] map_center_of_reference_cell A flag that switches the algorithm 358 * for the computation of the cell center from 359 * transform_unit_to_real_cell() applied to the center of the reference cell 360 * to computing the vertex averages. 361 */ 362 virtual Point<spacedim> 363 get_center(const typename Triangulation<dim, spacedim>::cell_iterator &cell, 364 const bool map_center_of_reference_cell = true) const; 365 366 /** 367 * Return the bounding box of a mapped cell. 368 * 369 * If you are using a (bi-,tri-)linear mapping that preserves vertex 370 * locations, this function simply returns the value also produced by 371 * `cell->bounding_box()`. However, there are also mappings that add 372 * displacements or choose completely different locations, e.g., 373 * MappingQEulerian, MappingQ1Eulerian, or MappingFEField. 374 * 375 * This function returns the bounding box containing all the vertices of the 376 * cell, as returned by the get_vertices() method. Beware of the fact that 377 * for higher order mappings this bounding box is only an approximation of the 378 * true bounding box, since it does not take into account curved faces, and it 379 * may be smaller than the true bounding box. 380 * 381 * @param[in] cell The cell for which you want to compute the bounding box 382 */ 383 virtual BoundingBox<spacedim> 384 get_bounding_box( 385 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const; 386 387 /** 388 * Return whether the mapping preserves vertex locations. In other words, 389 * this function returns whether the mapped location of the reference cell 390 * vertices (given by GeometryInfo::unit_cell_vertex()) equals the result of 391 * <code>cell-@>vertex()</code> (i.e., information stored by the 392 * triangulation). 393 * 394 * For example, implementations in derived classes return @p true for 395 * MappingQ, MappingQGeneric, MappingCartesian, but @p false for 396 * MappingQEulerian, MappingQ1Eulerian, and MappingFEField. 397 */ 398 virtual bool 399 preserves_vertex_locations() const = 0; 400 401 /** 402 * @name Mapping points between reference and real cells 403 * @{ 404 */ 405 406 /** 407 * Map the point @p p on the unit cell to the corresponding point on the 408 * real cell @p cell. 409 * 410 * @param cell Iterator to the cell that will be used to define the mapping. 411 * @param p Location of a point on the reference cell. 412 * @return The location of the reference point mapped to real space using 413 * the mapping defined by the class derived from the current one that 414 * implements the mapping, and the coordinates of the cell identified by the 415 * first argument. 416 */ 417 virtual Point<spacedim> 418 transform_unit_to_real_cell( 419 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 420 const Point<dim> & p) const = 0; 421 422 /** 423 * Map the point @p p on the real @p cell to the corresponding point on the 424 * unit cell, and return its coordinates. This function provides the inverse 425 * of the mapping provided by transform_unit_to_real_cell(). 426 * 427 * In the codimension one case, this function returns the normal projection 428 * of the real point @p p on the curve or surface identified by the @p cell. 429 * 430 * @note Polynomial mappings from the reference (unit) cell coordinates to 431 * the coordinate system of a real cell are not always invertible if the 432 * point for which the inverse mapping is to be computed lies outside the 433 * cell's boundaries. In such cases, the current function may fail to 434 * compute a point on the reference cell whose image under the mapping 435 * equals the given point @p p. If this is the case then this function 436 * throws an exception of type Mapping::ExcTransformationFailed . Whether 437 * the given point @p p lies outside the cell can therefore be determined by 438 * checking whether the returned reference coordinates lie inside or outside 439 * the reference cell (e.g., using GeometryInfo::is_inside_unit_cell()) or 440 * whether the exception mentioned above has been thrown. 441 * 442 * @param cell Iterator to the cell that will be used to define the mapping. 443 * @param p Location of a point on the given cell. 444 * @return The reference cell location of the point that when mapped to real 445 * space equals the coordinates given by the second argument. This mapping 446 * uses the mapping defined by the class derived from the current one that 447 * implements the mapping, and the coordinates of the cell identified by the 448 * first argument. 449 */ 450 virtual Point<dim> 451 transform_real_to_unit_cell( 452 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 453 const Point<spacedim> & p) const = 0; 454 455 /** 456 * Transform the point @p p on the real @p cell to the corresponding point 457 * on the reference cell, and then project this point to a (dim-1)-dimensional 458 * point in the coordinate system of the face with 459 * the given face number @p face_no. Ideally the point @p p is near the face 460 * @p face_no, but any point in the cell can technically be projected. 461 * 462 * This function does not make physical sense when dim=1, so it throws an 463 * exception in this case. 464 */ 465 Point<dim - 1> 466 project_real_point_to_unit_point_on_face( 467 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 468 const unsigned int face_no, 469 const Point<spacedim> & p) const; 470 471 /** 472 * @} 473 */ 474 475 476 /** 477 * @name Exceptions 478 * @{ 479 */ 480 481 /** 482 * Exception 483 */ 484 DeclException0(ExcInvalidData); 485 486 487 /** 488 * Computing the mapping between a real space point and a point in reference 489 * space failed, typically because the given point lies outside the cell 490 * where the inverse mapping is not unique. 491 * 492 * @ingroup Exceptions 493 */ 494 DeclExceptionMsg( 495 ExcTransformationFailed, 496 "Computing the mapping between a real space point and a point in reference " 497 "space failed, typically because the given point lies outside the cell " 498 "where the inverse mapping is not unique."); 499 500 /** 501 * deal.II assumes the Jacobian determinant to be positive. When the cell 502 * geometry is distorted under the image of the mapping, the mapping becomes 503 * invalid and this exception is thrown. 504 * 505 * @ingroup Exceptions 506 */ 507 DeclException3(ExcDistortedMappedCell, 508 Point<spacedim>, 509 double, 510 int, 511 << "The image of the mapping applied to cell with center [" 512 << arg1 << "] is distorted. The cell geometry or the " 513 << "mapping are invalid, giving a non-positive volume " 514 << "fraction of " << arg2 << " in quadrature point " << arg3 515 << "."); 516 517 /** 518 * @} 519 */ 520 521 /** 522 * @name Interface with FEValues 523 * @{ 524 */ 525 526 public: 527 /** 528 * Base class for internal data of mapping objects. The internal mechanism 529 * is that upon construction of a FEValues object, it asks the mapping and 530 * finite element classes that are to be used to allocate memory for their 531 * own purpose in which they may store data that only needs to be computed 532 * once. For example, most finite elements will store the values of the 533 * shape functions at the quadrature points in this object, since they do 534 * not change from cell to cell and only need to be computed once. The same 535 * may be true for Mapping classes that want to only evaluate the shape 536 * functions used for mapping once at the quadrature points. 537 * 538 * Since different FEValues objects using different quadrature rules might 539 * access the same mapping object at the same time, it is necessary to 540 * create one such object per FEValues object. FEValues does this by calling 541 * Mapping::get_data(), or in reality the implementation of the 542 * corresponding function in derived classes. Ownership of the object 543 * created by Mapping::get_data() is then transferred to the FEValues 544 * object, but a reference to this object is passed to the mapping object 545 * every time it is asked to compute information on a concrete cell. This 546 * happens when FEValues::reinit() (or the corresponding classes in 547 * FEFaceValues and FESubfaceValues) call Mapping::fill_fe_values() (and 548 * similarly via Mapping::fill_fe_face_values() and 549 * Mapping::fill_fe_subface_values()). 550 * 551 * The purpose of this class is for mapping objects to store information 552 * that can be computed once at the beginning, on the reference cell, and to 553 * access it later when computing information on a concrete cell. As such, 554 * the object handed to Mapping::fill_fe_values() is marked as 555 * <code>const</code>, because the assumption is that at the time this 556 * information is used, it will not need to modified again. However, classes 557 * derived from Mapping can also use such objects for two other purposes: 558 * 559 * - To provide scratch space for computations that are done in 560 * Mapping::fill_fe_values() and similar functions. Some of the derived 561 * classes would like to use scratch arrays and it would be a waste of time 562 * to allocate these arrays every time this function is called, just to de- 563 * allocate it again at the end of the function. Rather, one could allocate 564 * this memory once as a member variable of the current class, and simply 565 * use it in Mapping::fill_fe_values(). 566 * - After calling Mapping::fill_fe_values(), FEValues::reinit() 567 * calls FiniteElement::fill_fe_values() where the finite element computes 568 * values, gradients, etc of the shape functions using both information 569 * computed once at the beginning using a mechanism similar to the one 570 * described here (see FiniteElement::InternalDataBase) as well as the data 571 * already computed by Mapping::fill_fe_values(). As part of its work, some 572 * implementations of FiniteElement::fill_fe_values() need to transform 573 * shape function data, and they do so by calling Mapping::transform(). The 574 * call to the latter function also receives a reference to the 575 * Mapping::InternalDataBase object. Since Mapping::transform() may be 576 * called many times on each cell, it is sometimes worth for derived classes 577 * to compute some information only once in Mapping::fill_fe_values() and 578 * reuse it in Mapping::transform(). This information can also be stored in 579 * the classes that derived mapping classes derive from InternalDataBase. 580 * 581 * In both of these cases, the InternalDataBase object being passed around 582 * is "morally const", i.e., no external observer can tell whether a scratch 583 * array or some intermediate data for Mapping::transform() is being 584 * modified by Mapping::fill_fe_values() or not. Consequently, the 585 * InternalDataBase objects are always passed around as <code>const</code> 586 * objects. Derived classes that would like to make use of the two 587 * additional uses outlined above therefore need to mark the member 588 * variables they want to use for these purposes as <code>mutable</code> to 589 * allow for their modification despite the fact that the surrounding object 590 * is marked as <code>const</code>. 591 */ 592 class InternalDataBase 593 { 594 public: 595 /** 596 * Constructor. Sets update_flags to @p update_default and @p first_cell 597 * to @p true. 598 */ 599 InternalDataBase(); 600 601 /** 602 * Copy construction is forbidden. 603 */ 604 InternalDataBase(const InternalDataBase &) = delete; 605 606 /** 607 * Virtual destructor for derived classes 608 */ 609 virtual ~InternalDataBase() = default; 610 611 /** 612 * A set of update flags specifying the kind of information that an 613 * implementation of the Mapping interface needs to compute on each cell 614 * or face, i.e., in Mapping::fill_fe_values() and friends. 615 * 616 * This set of flags is stored here by implementations of 617 * Mapping::get_data(), Mapping::get_face_data(), or 618 * Mapping::get_subface_data(), and is that subset of the update flags 619 * passed to those functions that require re-computation on every cell. 620 * (The subset of the flags corresponding to information that can be 621 * computed once and for all already at the time of the call to 622 * Mapping::get_data() -- or an implementation of that interface -- need 623 * not be stored here because it has already been taken care of.) 624 */ 625 UpdateFlags update_each; 626 627 /** 628 * Return an estimate (in bytes) for the memory consumption of this object. 629 */ 630 virtual std::size_t 631 memory_consumption() const; 632 }; 633 634 635 protected: 636 /** 637 * Given a set of update flags, compute which other quantities <i>also</i> 638 * need to be computed in order to satisfy the request by the given flags. 639 * Then return the combination of the original set of flags and those just 640 * computed. 641 * 642 * As an example, if @p update_flags contains update_JxW_values (i.e., the 643 * product of the determinant of the Jacobian and the weights provided by 644 * the quadrature formula), a mapping may require the computation of the 645 * full Jacobian matrix in order to compute its determinant. They would then 646 * return not just update_JxW_values, but also update_jacobians. (This is 647 * not how it is actually done internally in the derived classes that 648 * compute the JxW values -- they set update_contravariant_transformation 649 * instead, from which the determinant can also be computed -- but this does 650 * not take away from the instructiveness of the example.) 651 * 652 * An extensive discussion of the interaction between this function and 653 * FEValues can be found in the 654 * @ref FE_vs_Mapping_vs_FEValues 655 * documentation module. 656 * 657 * @see UpdateFlags 658 */ 659 virtual UpdateFlags 660 requires_update_flags(const UpdateFlags update_flags) const = 0; 661 662 /** 663 * Create and return a pointer to an object into which mappings can store 664 * data that only needs to be computed once but that can then be used 665 * whenever the mapping is applied to a concrete cell (e.g., in the various 666 * transform() functions, as well as in the fill_fe_values(), 667 * fill_fe_face_values() and fill_fe_subface_values() that form the 668 * interface of mappings with the FEValues class). 669 * 670 * Derived classes will return pointers to objects of a type derived from 671 * Mapping::InternalDataBase (see there for more information) and may pre- 672 * compute some information already (in accordance with what will be asked 673 * of the mapping in the future, as specified by the update flags) and for 674 * the given quadrature object. Subsequent calls to transform() or 675 * fill_fe_values() and friends will then receive back the object created 676 * here (with the same set of update flags and for the same quadrature 677 * object). Derived classes can therefore pre-compute some information in 678 * their get_data() function and store it in the internal data object. 679 * 680 * The mapping classes do not keep track of the objects created by this 681 * function. Ownership will therefore rest with the caller. 682 * 683 * An extensive discussion of the interaction between this function and 684 * FEValues can be found in the 685 * @ref FE_vs_Mapping_vs_FEValues 686 * documentation module. 687 * 688 * @param update_flags A set of flags that define what is expected of the 689 * mapping class in future calls to transform() or the fill_fe_values() 690 * group of functions. This set of flags may contain flags that mappings do 691 * not know how to deal with (e.g., for information that is in fact computed 692 * by the finite element classes, such as UpdateFlags::update_values). 693 * Derived classes will need to store these flags, or at least that subset 694 * of flags that will require the mapping to perform any actions in 695 * fill_fe_values(), in InternalDataBase::update_each. 696 * @param quadrature The quadrature object for which mapping information 697 * will have to be computed. This includes the locations and weights of 698 * quadrature points. 699 * @return A pointer to a newly created object of type InternalDataBase (or 700 * a derived class). Ownership of this object passes to the calling 701 * function. 702 * 703 * @note C++ allows that virtual functions in derived classes may return 704 * pointers to objects not of type InternalDataBase but in fact pointers to 705 * objects of classes <i>derived</i> from InternalDataBase. (This feature is 706 * called "covariant return types".) This is useful in some contexts where 707 * the calling is within the derived class and will immediately make use of 708 * the returned object, knowing its real (derived) type. 709 */ 710 virtual std::unique_ptr<InternalDataBase> 711 get_data(const UpdateFlags update_flags, 712 const Quadrature<dim> &quadrature) const = 0; 713 714 /** 715 * Like get_data(), but in preparation for later calls to transform() or 716 * fill_fe_face_values() that will need information about mappings from the 717 * reference face to a face of a concrete cell. 718 * 719 * @param update_flags A set of flags that define what is expected of the 720 * mapping class in future calls to transform() or the fill_fe_values() 721 * group of functions. This set of flags may contain flags that mappings do 722 * not know how to deal with (e.g., for information that is in fact computed 723 * by the finite element classes, such as UpdateFlags::update_values). 724 * Derived classes will need to store these flags, or at least that subset 725 * of flags that will require the mapping to perform any actions in 726 * fill_fe_values(), in InternalDataBase::update_each. 727 * @param quadrature The quadrature object for which mapping information 728 * will have to be computed. This includes the locations and weights of 729 * quadrature points. 730 * @return A pointer to a newly created object of type InternalDataBase (or 731 * a derived class). Ownership of this object passes to the calling 732 * function. 733 * 734 * @note C++ allows that virtual functions in derived classes may return 735 * pointers to objects not of type InternalDataBase but in fact pointers to 736 * objects of classes <i>derived</i> from InternalDataBase. (This feature is 737 * called "covariant return types".) This is useful in some contexts where 738 * the calling is within the derived class and will immediately make use of 739 * the returned object, knowing its real (derived) type. 740 */ 741 virtual std::unique_ptr<InternalDataBase> 742 get_face_data(const UpdateFlags update_flags, 743 const Quadrature<dim - 1> &quadrature) const = 0; 744 745 /** 746 * Like get_data() and get_face_data(), but in preparation for later calls 747 * to transform() or fill_fe_subface_values() that will need information 748 * about mappings from the reference face to a child of a face (i.e., 749 * subface) of a concrete cell. 750 * 751 * @param update_flags A set of flags that define what is expected of the 752 * mapping class in future calls to transform() or the fill_fe_values() 753 * group of functions. This set of flags may contain flags that mappings do 754 * not know how to deal with (e.g., for information that is in fact computed 755 * by the finite element classes, such as UpdateFlags::update_values). 756 * Derived classes will need to store these flags, or at least that subset 757 * of flags that will require the mapping to perform any actions in 758 * fill_fe_values(), in InternalDataBase::update_each. 759 * @param quadrature The quadrature object for which mapping information 760 * will have to be computed. This includes the locations and weights of 761 * quadrature points. 762 * @return A pointer to a newly created object of type InternalDataBase (or 763 * a derived class). Ownership of this object passes to the calling 764 * function. 765 * 766 * @note C++ allows that virtual functions in derived classes may return 767 * pointers to objects not of type InternalDataBase but in fact pointers to 768 * objects of classes <i>derived</i> from InternalDataBase. (This feature is 769 * called "covariant return types".) This is useful in some contexts where 770 * the calling is within the derived class and will immediately make use of 771 * the returned object, knowing its real (derived) type. 772 */ 773 virtual std::unique_ptr<InternalDataBase> 774 get_subface_data(const UpdateFlags update_flags, 775 const Quadrature<dim - 1> &quadrature) const = 0; 776 777 /** 778 * Compute information about the mapping from the reference cell to the real 779 * cell indicated by the first argument to this function. Derived classes 780 * will have to implement this function based on the kind of mapping they 781 * represent. It is called by FEValues::reinit(). 782 * 783 * Conceptually, this function's represents the application of the mapping 784 * $\mathbf x=\mathbf F_K(\hat {\mathbf x})$ from reference coordinates 785 * $\mathbf\in [0,1]^d$ to real space coordinates $\mathbf x$ for a given 786 * cell $K$. Its purpose is to compute the following kinds of data: 787 * 788 * - Data that results from the application of the mapping itself, e.g., 789 * computing the location $\mathbf x_q = \mathbf F_K(\hat{\mathbf x}_q)$ of 790 * quadrature points on the real cell, and that is directly useful to users 791 * of FEValues, for example during assembly. 792 * - Data that is necessary for finite element implementations to compute 793 * their shape functions on the real cell. To this end, the 794 * FEValues::reinit() function calls FiniteElement::fill_fe_values() after 795 * the current function, and the output of this function serves as input to 796 * FiniteElement::fill_fe_values(). Examples of information that needs to be 797 * computed here for use by the finite element classes is the Jacobian of 798 * the mapping, $\hat\nabla \mathbf F_K(\hat{\mathbf x})$ or its inverse, 799 * for example to transform the gradients of shape functions on the 800 * reference cell to the gradients of shape functions on the real cell. 801 * 802 * The information computed by this function is used to fill the various 803 * member variables of the output argument of this function. Which of the 804 * member variables of that structure should be filled is determined by the 805 * update flags stored in the Mapping::InternalDataBase object passed to 806 * this function. 807 * 808 * An extensive discussion of the interaction between this function and 809 * FEValues can be found in the 810 * @ref FE_vs_Mapping_vs_FEValues 811 * documentation module. 812 * 813 * @param[in] cell The cell of the triangulation for which this function is 814 * to compute a mapping from the reference cell to. 815 * @param[in] cell_similarity Whether or not the cell given as first 816 * argument is simply a translation, rotation, etc of the cell for which 817 * this function was called the most recent time. This information is 818 * computed simply by matching the vertices (as stored by the Triangulation) 819 * between the previous and the current cell. The value passed here may be 820 * modified by implementations of this function and should then be returned 821 * (see the discussion of the return value of this function). 822 * @param[in] quadrature A reference to the quadrature formula in use for 823 * the current evaluation. This quadrature object is the same as the one 824 * used when creating the @p internal_data object. The object is used both 825 * to map the location of quadrature points, as well as to compute the JxW 826 * values for each quadrature point (which involves the quadrature weights). 827 * @param[in] internal_data A reference to an object previously created by 828 * get_data() and that may be used to store information the mapping can 829 * compute once on the reference cell. See the documentation of the 830 * Mapping::InternalDataBase class for an extensive description of the 831 * purpose of these objects. 832 * @param[out] output_data A reference to an object whose member variables 833 * should be computed. Not all of the members of this argument need to be 834 * filled; which ones need to be filled is determined by the update flags 835 * stored inside the @p internal_data object. 836 * @return An updated value of the @p cell_similarity argument to this 837 * function. The returned value will be used for the corresponding argument 838 * when FEValues::reinit() calls FiniteElement::fill_fe_values(). In most 839 * cases, derived classes will simply want to return the value passed for @p 840 * cell_similarity. However, implementations of this function may downgrade 841 * the level of cell similarity. This is, for example, the case for classes 842 * that take not only into account the locations of the vertices of a cell 843 * (as reported by the Triangulation), but also other information specific 844 * to the mapping. The purpose is that FEValues::reinit() can compute 845 * whether a cell is similar to the previous one only based on the cell's 846 * vertices, whereas the mapping may also consider displacement fields 847 * (e.g., in the MappingQ1Eulerian and MappingFEField classes). In such 848 * cases, the mapping may conclude that the previously computed cell 849 * similarity is too optimistic, and invalidate it for subsequent use in 850 * FiniteElement::fill_fe_values() by returning a less optimistic cell 851 * similarity value. 852 * 853 * @note FEValues ensures that this function is always called with the same 854 * pair of @p internal_data and @p output_data objects. In other words, if 855 * an implementation of this function knows that it has written a piece of 856 * data into the output argument in a previous call, then there is no need 857 * to copy it there again in a later call if the implementation knows that 858 * this is the same value. 859 */ 860 virtual CellSimilarity::Similarity 861 fill_fe_values( 862 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 863 const CellSimilarity::Similarity cell_similarity, 864 const Quadrature<dim> & quadrature, 865 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data, 866 dealii::internal::FEValuesImplementation::MappingRelatedData<dim, spacedim> 867 &output_data) const = 0; 868 869 /** 870 * This function is the equivalent to Mapping::fill_fe_values(), but for 871 * faces of cells. See there for an extensive discussion of its purpose. It 872 * is called by FEFaceValues::reinit(). 873 * 874 * @param[in] cell The cell of the triangulation for which this function is 875 * to compute a mapping from the reference cell to. 876 * @param[in] face_no The number of the face of the given cell for which 877 * information is requested. 878 * @param[in] quadrature A reference to the quadrature formula in use for 879 * the current evaluation. This quadrature object is the same as the one 880 * used when creating the @p internal_data object. The object is used both 881 * to map the location of quadrature points, as well as to compute the JxW 882 * values for each quadrature point (which involves the quadrature weights). 883 * @param[in] internal_data A reference to an object previously created by 884 * get_data() and that may be used to store information the mapping can 885 * compute once on the reference cell. See the documentation of the 886 * Mapping::InternalDataBase class for an extensive description of the 887 * purpose of these objects. 888 * @param[out] output_data A reference to an object whose member variables 889 * should be computed. Not all of the members of this argument need to be 890 * filled; which ones need to be filled is determined by the update flags 891 * stored inside the @p internal_data object. 892 */ 893 virtual void 894 fill_fe_face_values( 895 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 896 const unsigned int face_no, 897 const Quadrature<dim - 1> & quadrature, 898 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data, 899 dealii::internal::FEValuesImplementation::MappingRelatedData<dim, spacedim> 900 &output_data) const = 0; 901 902 /** 903 * This function is the equivalent to Mapping::fill_fe_values(), but for 904 * subfaces (i.e., children of faces) of cells. See there for an extensive 905 * discussion of its purpose. It is called by FESubfaceValues::reinit(). 906 * 907 * @param[in] cell The cell of the triangulation for which this function is 908 * to compute a mapping from the reference cell to. 909 * @param[in] face_no The number of the face of the given cell for which 910 * information is requested. 911 * @param[in] subface_no The number of the child of a face of the given cell 912 * for which information is requested. 913 * @param[in] quadrature A reference to the quadrature formula in use for 914 * the current evaluation. This quadrature object is the same as the one 915 * used when creating the @p internal_data object. The object is used both 916 * to map the location of quadrature points, as well as to compute the JxW 917 * values for each quadrature point (which involves the quadrature weights). 918 * @param[in] internal_data A reference to an object previously created by 919 * get_data() and that may be used to store information the mapping can 920 * compute once on the reference cell. See the documentation of the 921 * Mapping::InternalDataBase class for an extensive description of the 922 * purpose of these objects. 923 * @param[out] output_data A reference to an object whose member variables 924 * should be computed. Not all of the members of this argument need to be 925 * filled; which ones need to be filled is determined by the update flags 926 * stored inside the @p internal_data object. 927 */ 928 virtual void 929 fill_fe_subface_values( 930 const typename Triangulation<dim, spacedim>::cell_iterator &cell, 931 const unsigned int face_no, 932 const unsigned int subface_no, 933 const Quadrature<dim - 1> & quadrature, 934 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data, 935 dealii::internal::FEValuesImplementation::MappingRelatedData<dim, spacedim> 936 &output_data) const = 0; 937 938 /** 939 * @} 940 */ 941 942 public: 943 /** 944 * @name Functions to transform tensors from reference to real coordinates 945 * @{ 946 */ 947 948 /** 949 * Transform a field of vectors or 1-differential forms according to the 950 * selected MappingKind. 951 * 952 * @note Normally, this function is called by a finite element, filling 953 * FEValues objects. For this finite element, there should be an alias 954 * MappingKind like @p mapping_bdm, @p mapping_nedelec, etc. This alias 955 * should be preferred to using the kinds below. 956 * 957 * The mapping kinds currently implemented by derived classes are: 958 * <ul> 959 * <li> @p mapping_contravariant: maps a vector field on the reference cell 960 * to the physical cell through the Jacobian: 961 * @f[ 962 * \mathbf u(\mathbf x) = J(\hat{\mathbf x})\hat{\mathbf u}(\hat{\mathbf 963 * x}). 964 * @f] 965 * In physics, this is usually referred to as the contravariant 966 * transformation. Mathematically, it is the push forward of a vector field. 967 * 968 * <li> @p mapping_covariant: maps a field of one-forms on the reference 969 * cell to a field of one-forms on the physical cell. (Theoretically this 970 * would refer to a DerivativeForm<1,dim,1> but we canonically identify this 971 * type with a Tensor<1,dim>). Mathematically, it is the pull back of the 972 * differential form 973 * @f[ 974 * \mathbf u(\mathbf x) = J(\hat{\mathbf x})(J(\hat{\mathbf x})^{T} 975 * J(\hat{\mathbf x}))^{-1}\hat{\mathbf u}(\hat{\mathbf x}). 976 * @f] 977 * Gradients of scalar differentiable functions are transformed this way. 978 * 979 * In the case when dim=spacedim the previous formula reduces to 980 * @f[ 981 * \mathbf u(\mathbf x) = J(\hat{\mathbf x})^{-T}\hat{\mathbf 982 * u}(\hat{\mathbf x}) 983 * @f] 984 * because we assume that the mapping $\mathbf F_K$ is always invertible, 985 * and consequently its Jacobian $J$ is an invertible matrix. 986 * 987 * <li> @p mapping_piola: A field of <i>dim-1</i>-forms on the reference 988 * cell is also represented by a vector field, but again transforms 989 * differently, namely by the Piola transform 990 * @f[ 991 * \mathbf u(\mathbf x) = \frac{1}{\text{det}\;J(\hat{\mathbf x})} 992 * J(\hat{\mathbf x}) \hat{\mathbf u}(\hat{\mathbf x}). 993 * @f] 994 * </ul> 995 * 996 * @param[in] input An array (or part of an array) of input objects that 997 * should be mapped. 998 * @param[in] kind The kind of mapping to be applied. 999 * @param[in] internal A pointer to an object of type 1000 * Mapping::InternalDataBase that contains information previously stored by 1001 * the mapping. The object pointed to was created by the get_data(), 1002 * get_face_data(), or get_subface_data() function, and will have been 1003 * updated as part of a call to fill_fe_values(), fill_fe_face_values(), or 1004 * fill_fe_subface_values() for the current cell, before calling the current 1005 * function. In other words, this object also represents with respect to 1006 * which cell the transformation should be applied to. 1007 * @param[out] output An array (or part of an array) into which the 1008 * transformed objects should be placed. (Note that the array view is @p 1009 * const, but the tensors it points to are not.) 1010 */ 1011 virtual void 1012 transform(const ArrayView<const Tensor<1, dim>> & input, 1013 const MappingKind kind, 1014 const typename Mapping<dim, spacedim>::InternalDataBase &internal, 1015 const ArrayView<Tensor<1, spacedim>> &output) const = 0; 1016 1017 /** 1018 * Transform a field of differential forms from the reference cell to the 1019 * physical cell. It is useful to think of $\mathbf{T} = \nabla \mathbf u$ 1020 * and $\hat{\mathbf T} = \hat \nabla \hat{\mathbf u}$, with $\mathbf u$ a 1021 * vector field. The mapping kinds currently implemented by derived classes 1022 * are: 1023 * <ul> 1024 * <li> @p mapping_covariant: maps a field of forms on the reference cell to 1025 * a field of forms on the physical cell. Mathematically, it is the pull 1026 * back of the differential form 1027 * @f[ 1028 * \mathbf T(\mathbf x) = \hat{\mathbf T}(\hat{\mathbf x}) 1029 * J(\hat{\mathbf x})(J(\hat{\mathbf x})^{T} 1030 * J(\hat{\mathbf x}))^{-1}. 1031 * @f] 1032 * Jacobians of spacedim-vector valued differentiable functions are 1033 * transformed this way. 1034 * 1035 * In the case when dim=spacedim the previous formula reduces to 1036 * @f[ 1037 * \mathbf T(\mathbf x) = \hat{\mathbf u}(\hat{\mathbf x}) 1038 * J(\hat{\mathbf x})^{-1}. 1039 * @f] 1040 * </ul> 1041 * 1042 * @note It would have been more reasonable to make this transform a 1043 * template function with the rank in <code>DerivativeForm@<1, dim, 1044 * rank@></code>. Unfortunately C++ does not allow templatized virtual 1045 * functions. This is why we identify <code>DerivativeForm@<1, dim, 1046 * 1@></code> with a <code>Tensor@<1,dim@></code> when using 1047 * mapping_covariant() in the function transform() above this one. 1048 * 1049 * @param[in] input An array (or part of an array) of input objects that 1050 * should be mapped. 1051 * @param[in] kind The kind of mapping to be applied. 1052 * @param[in] internal A pointer to an object of type 1053 * Mapping::InternalDataBase that contains information previously stored by 1054 * the mapping. The object pointed to was created by the get_data(), 1055 * get_face_data(), or get_subface_data() function, and will have been 1056 * updated as part of a call to fill_fe_values(), fill_fe_face_values(), or 1057 * fill_fe_subface_values() for the current cell, before calling the current 1058 * function. In other words, this object also represents with respect to 1059 * which cell the transformation should be applied to. 1060 * @param[out] output An array (or part of an array) into which the 1061 * transformed objects should be placed. (Note that the array view is @p 1062 * const, but the tensors it points to are not.) 1063 */ 1064 virtual void 1065 transform(const ArrayView<const DerivativeForm<1, dim, spacedim>> &input, 1066 const MappingKind kind, 1067 const typename Mapping<dim, spacedim>::InternalDataBase &internal, 1068 const ArrayView<Tensor<2, spacedim>> &output) const = 0; 1069 1070 /** 1071 * Transform a tensor field from the reference cell to the physical cell. 1072 * These tensors are usually the Jacobians in the reference cell of vector 1073 * fields that have been pulled back from the physical cell. The mapping 1074 * kinds currently implemented by derived classes are: 1075 * <ul> 1076 * <li> @p mapping_contravariant_gradient: it assumes $\mathbf u(\mathbf x) 1077 * = J \hat{\mathbf u}$ so that 1078 * @f[ 1079 * \mathbf T(\mathbf x) = 1080 * J(\hat{\mathbf x}) \hat{\mathbf T}(\hat{\mathbf x}) 1081 * J(\hat{\mathbf x})^{-1}. 1082 * @f] 1083 * <li> @p mapping_covariant_gradient: it assumes $\mathbf u(\mathbf x) = 1084 * J^{-T} \hat{\mathbf u}$ so that 1085 * @f[ 1086 * \mathbf T(\mathbf x) = 1087 * J(\hat{\mathbf x})^{-T} \hat{\mathbf T}(\hat{\mathbf x}) 1088 * J(\hat{\mathbf x})^{-1}. 1089 * @f] 1090 * <li> @p mapping_piola_gradient: it assumes $\mathbf u(\mathbf x) = 1091 * \frac{1}{\text{det}\;J(\hat{\mathbf x})} J(\hat{\mathbf x}) \hat{\mathbf 1092 * u}(\hat{\mathbf x})$ so that 1093 * @f[ 1094 * \mathbf T(\mathbf x) = 1095 * \frac{1}{\text{det}\;J(\hat{\mathbf x})} 1096 * J(\hat{\mathbf x}) \hat{\mathbf T}(\hat{\mathbf x}) 1097 * J(\hat{\mathbf x})^{-1}. 1098 * @f] 1099 * </ul> 1100 * 1101 * @todo The formulas for mapping_covariant_gradient, 1102 * mapping_contravariant_gradient and mapping_piola_gradient are only true 1103 * as stated for linear mappings. If, for example, the mapping is bilinear 1104 * (or has a higher order polynomial degree) then there is a missing term 1105 * associated with the derivative of $J$. 1106 * 1107 * @param[in] input An array (or part of an array) of input objects that 1108 * should be mapped. 1109 * @param[in] kind The kind of mapping to be applied. 1110 * @param[in] internal A pointer to an object of type 1111 * Mapping::InternalDataBase that contains information previously stored by 1112 * the mapping. The object pointed to was created by the get_data(), 1113 * get_face_data(), or get_subface_data() function, and will have been 1114 * updated as part of a call to fill_fe_values(), fill_fe_face_values(), or 1115 * fill_fe_subface_values() for the current cell, before calling the current 1116 * function. In other words, this object also represents with respect to 1117 * which cell the transformation should be applied to. 1118 * @param[out] output An array (or part of an array) into which the 1119 * transformed objects should be placed. (Note that the array view is @p 1120 * const, but the tensors it points to are not.) 1121 */ 1122 virtual void 1123 transform(const ArrayView<const Tensor<2, dim>> & input, 1124 const MappingKind kind, 1125 const typename Mapping<dim, spacedim>::InternalDataBase &internal, 1126 const ArrayView<Tensor<2, spacedim>> &output) const = 0; 1127 1128 /** 1129 * Transform a tensor field from the reference cell to the physical cell. 1130 * This tensors are most of times the hessians in the reference cell of 1131 * vector fields that have been pulled back from the physical cell. 1132 * 1133 * The mapping kinds currently implemented by derived classes are: 1134 * <ul> 1135 * <li> @p mapping_covariant_gradient: maps a field of forms on the 1136 * reference cell to a field of forms on the physical cell. Mathematically, 1137 * it is the pull back of the differential form 1138 * @f[ 1139 * \mathbf T_{ijk}(\mathbf x) = \hat{\mathbf T}_{iJK}(\hat{\mathbf x}) 1140 * J_{jJ}^{\dagger} J_{kK}^{\dagger}@f], 1141 * 1142 * where @f[ J^{\dagger} = J(\hat{\mathbf x})(J(\hat{\mathbf x})^{T} 1143 * J(\hat{\mathbf x}))^{-1}. 1144 * @f] 1145 * </ul> 1146 * 1147 * Hessians of spacedim-vector valued differentiable functions are 1148 * transformed this way (After subtraction of the product of the derivative 1149 * with the Jacobian gradient). 1150 * 1151 * In the case when dim=spacedim the previous formula reduces to 1152 * @f[J^{\dagger} = J^{-1}@f] 1153 * 1154 * @param[in] input An array (or part of an array) of input objects that 1155 * should be mapped. 1156 * @param[in] kind The kind of mapping to be applied. 1157 * @param[in] internal A pointer to an object of type 1158 * Mapping::InternalDataBase that contains information previously stored by 1159 * the mapping. The object pointed to was created by the get_data(), 1160 * get_face_data(), or get_subface_data() function, and will have been 1161 * updated as part of a call to fill_fe_values(), fill_fe_face_values(), or 1162 * fill_fe_subface_values() for the current cell, before calling the current 1163 * function. In other words, this object also represents with respect to 1164 * which cell the transformation should be applied to. 1165 * @param[out] output An array (or part of an array) into which the 1166 * transformed objects should be placed. (Note that the array view is @p 1167 * const, but the tensors it points to are not.) 1168 */ 1169 virtual void 1170 transform(const ArrayView<const DerivativeForm<2, dim, spacedim>> &input, 1171 const MappingKind kind, 1172 const typename Mapping<dim, spacedim>::InternalDataBase &internal, 1173 const ArrayView<Tensor<3, spacedim>> &output) const = 0; 1174 1175 /** 1176 * Transform a field of 3-differential forms from the reference cell to the 1177 * physical cell. It is useful to think of $\mathbf{T}_{ijk} = D^2_{jk} 1178 * \mathbf u_i$ and $\mathbf{\hat T}_{IJK} = \hat D^2_{JK} \mathbf{\hat 1179 * u}_I$, with $\mathbf u_i$ a vector field. 1180 * 1181 * The mapping kinds currently implemented by derived classes are: 1182 * <ul> 1183 * <li> @p mapping_contravariant_hessian: it assumes $\mathbf u_i(\mathbf x) 1184 * = J_{iI} \hat{\mathbf u}_I$ so that 1185 * @f[ 1186 * \mathbf T_{ijk}(\mathbf x) = 1187 * J_{iI}(\hat{\mathbf x}) \hat{\mathbf T}_{IJK}(\hat{\mathbf x}) 1188 * J_{jJ}(\hat{\mathbf x})^{-1} J_{kK}(\hat{\mathbf x})^{-1}. 1189 * @f] 1190 * <li> @p mapping_covariant_hessian: it assumes $\mathbf u_i(\mathbf x) = 1191 * J_{iI}^{-T} \hat{\mathbf u}_I$ so that 1192 * @f[ 1193 * \mathbf T_{ijk}(\mathbf x) = 1194 * J_iI(\hat{\mathbf x})^{-1} \hat{\mathbf T}_{IJK}(\hat{\mathbf x}) 1195 * J_{jJ}(\hat{\mathbf x})^{-1} J_{kK}(\hat{\mathbf x})^{-1}. 1196 * @f] 1197 * <li> @p mapping_piola_hessian: it assumes $\mathbf u_i(\mathbf x) = 1198 * \frac{1}{\text{det}\;J(\hat{\mathbf x})} J_{iI}(\hat{\mathbf x}) 1199 * \hat{\mathbf u}(\hat{\mathbf x})$ so that 1200 * @f[ 1201 * \mathbf T_{ijk}(\mathbf x) = 1202 * \frac{1}{\text{det}\;J(\hat{\mathbf x})} 1203 * J_{iI}(\hat{\mathbf x}) \hat{\mathbf T}_{IJK}(\hat{\mathbf x}) 1204 * J_{jJ}(\hat{\mathbf x})^{-1} J_{kK}(\hat{\mathbf x})^{-1}. 1205 * @f] 1206 * </ul> 1207 * 1208 * @param[in] input An array (or part of an array) of input objects that 1209 * should be mapped. 1210 * @param[in] kind The kind of mapping to be applied. 1211 * @param[in] internal A pointer to an object of type 1212 * Mapping::InternalDataBase that contains information previously stored by 1213 * the mapping. The object pointed to was created by the get_data(), 1214 * get_face_data(), or get_subface_data() function, and will have been 1215 * updated as part of a call to fill_fe_values(), fill_fe_face_values(), or 1216 * fill_fe_subface_values() for the current cell, before calling the current 1217 * function. In other words, this object also represents with respect to 1218 * which cell the transformation should be applied to. 1219 * @param[out] output An array (or part of an array) into which the 1220 * transformed objects should be placed. 1221 */ 1222 virtual void 1223 transform(const ArrayView<const Tensor<3, dim>> & input, 1224 const MappingKind kind, 1225 const typename Mapping<dim, spacedim>::InternalDataBase &internal, 1226 const ArrayView<Tensor<3, spacedim>> &output) const = 0; 1227 1228 /** 1229 * @} 1230 */ 1231 1232 1233 // Give class @p FEValues access to the private <tt>get_...data</tt> and 1234 // <tt>fill_fe_...values</tt> functions. 1235 friend class FEValuesBase<dim, spacedim>; 1236 friend class FEValues<dim, spacedim>; 1237 friend class FEFaceValues<dim, spacedim>; 1238 friend class FESubfaceValues<dim, spacedim>; 1239 }; 1240 1241 1242 DEAL_II_NAMESPACE_CLOSE 1243 1244 #endif 1245