From b579135b5f776974d8fed3db6ba84077ae8a39a3 Mon Sep 17 00:00:00 2001 From: bangerth Date: Wed, 7 Nov 2007 22:59:17 +0000 Subject: [PATCH] Add a bit to Guido's essay on update flags. Link to from other places. git-svn-id: https://svn.dealii.org/trunk@15458 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_update_flags.h | 3 + deal.II/deal.II/include/fe/fe_values.h | 6 + deal.II/doc/doxygen/headers/fe.h | 7 +- deal.II/doc/doxygen/headers/update_flags.h | 177 ++++++++++++++----- 4 files changed, 143 insertions(+), 50 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe_update_flags.h b/deal.II/deal.II/include/fe/fe_update_flags.h index 0d60c43819..5d6711e033 100644 --- a/deal.II/deal.II/include/fe/fe_update_flags.h +++ b/deal.II/deal.II/include/fe/fe_update_flags.h @@ -72,6 +72,9 @@ DEAL_II_NAMESPACE_OPEN * computation once or on each cell. Subsequent calls to the functions * @p update_once and @p update_each should just select among these * flags, but should not add new flags. + * + * The mechanism by which all this is accomplished is also discussed + * on the page on @ref UpdateFlagsEssay. */ enum UpdateFlags { diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index 45d814b666..42840193c7 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -59,6 +59,9 @@ template class Quadrature; * @note All data fields are public, but this is not * critical, because access to this object is private in FEValues. * + * The purpose of this class is discussed + * on the page on @ref UpdateFlagsEssay. + * * @author Guido Kanschat, 2000 */ template @@ -400,6 +403,9 @@ class FEValuesData * auxiliary data fields for updating. Still, it is recommended to * give all needed update flags to FEValues. * + * The mechanisms by which this class works is also discussed + * on the page on @ref UpdateFlagsEssay. + * * @author Wolfgang Bangerth, 1998, 2003, Guido Kanschat, 2001 */ template diff --git a/deal.II/doc/doxygen/headers/fe.h b/deal.II/doc/doxygen/headers/fe.h index 6db1fdeec8..5fbb1c2024 100644 --- a/deal.II/doc/doxygen/headers/fe.h +++ b/deal.II/doc/doxygen/headers/fe.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2003, 2005, 2006 by the deal.II authors +// Copyright (C) 2003, 2005, 2006, 2007 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -92,6 +92,11 @@ * All these classes are used in all tutorial programs from step-3 onward, and * are described there in significant detail. * + * The actual workings of the FEValues class and friends is + * complicated because it has to be general yet efficient. The page on + * @ref UpdateFlagsEssay attempts to give an overview of how this + * works. + * * @ingroup feall */ diff --git a/deal.II/doc/doxygen/headers/update_flags.h b/deal.II/doc/doxygen/headers/update_flags.h index 60d9d8ab61..5e2d97ae76 100644 --- a/deal.II/doc/doxygen/headers/update_flags.h +++ b/deal.II/doc/doxygen/headers/update_flags.h @@ -15,33 +15,99 @@ @page UpdateFlagsEssay The interplay of UpdateFlags, Mapping and FiniteElement in FEValues -@section update_flags_Intro Introduction - -The class FEValues was introduced because most of the time values of -shape functions etc. are used in Quadrature points only and their -computation might be sped up by doing it all at once. Furthermore, -some quantities, like the shape function values of standard scalar -finite elements do not depend on the actual mesh cell at all and do -not have to be reevaluated on each cell. The interface between -FEValues and FiniteElement should be aware of such things and do the -best thing possible. - -Now, in order not to compute too much, FEValues has to be told what to -update. This is where UpdateFlags comes into play. And each FEValues -uses two sets of UpdateFlags internally: one that describes the values -that can be computed on the reference cell only and one requiring -recomputation on each cell. Furthermore, since the actual computations -are performed by Mapping and FiniteElement objects, both of them have -UpdateFlags. - -@section update_flags_Once Update once or each - -Your desired set of values for computation given, the FiniteElement -and Mapping involved have to perform the following tasks: +

Introduction

+ +In order to compute local contributions of an individual to the global +matrix and right hand side, integrals are usually transformed to the +reference cell. For example, for the Laplace matrix, on each cell we +have to compute +@f[ + A^K_{ij} = \sum_{q}J^{-1}(\hat{\bf x}_q) \hat \nabla \varphi_i(\hat{\bf x}_q) \cdot + J^{-1}(\hat{\bf x}_q) \hat \nabla \varphi_j(\hat{\bf x}_q)\ |\textrm{det}\ J(\hat{\bf x}_q)| + w_q, +@f] +where a hat indicates reference coordinates, and $J(\hat{\bf +x}_q)$ is the Jacobian +$\frac{\partial F_K(\hat{\bf x}_q)}{\partial\bf \hat x}$ of the mapping, +evaluated at a quadrature point $\hat{\bf x}_q$ on the reference cell. + +In order to evaluate such an expression in an application code, we +have to access three different kinds of objects: a quadrature +object that describes locations $\hat{\bf x}_q$ and weights $w_q$ of +quadrature points on the reference cell; a finite element object that +describes the gradients $\hat\nabla \varphi_i(\hat{\bf x}_q)$ of shape +functions on the unit cell; and a mapping object that provides the +Jacobian as well as its determinant. Dealing with all these +objects would be cumbersome and error prone. + +On the other hand, these three kinds of objects almost always appear together, +and it is in fact very rare for deal.II application codes to do anything with +quadrature, finite element, or mapping objects besides using them together. +For this reason, we have introduced the FEValues abstraction +combining information on the shape functions, the geometry of the actual mesh +cell and a quadrature rule on a reference cell. Upon construction it takes one +object of each of the three mentioned categories. Later, it can be +"re-initialized" for a concrete grid cell and then provides mapped quadrature +points and weights, mapped shape function values and derivatives as well as +some properties of the transformation from the reference cell to the actual +mesh cell. + +Since computation of any of these values is potentially expensive (for +example when using high order mappings with high order elements), the +FEValues class only computes what is explicitly asked for. To this +end, it takes a list of flags of type UpdateFlags at construction time +specifying which quantities should be updated each time a cell is +visited. In addition, allowing further optimizations, the functions +filling the data fields of FEValues are able to distinguish between +values that have to be recomputed on each cell (for example mapped +gradients) and quantities that do not change from cell to cell (for +example the values of shape functions of the usual $Q_k$ +finite elements at the same quadrature points on different cells; this +property does not hold for the shape functions of Raviart-Thomas +elements, however, which must be rotated with the local cell). + +At construction time, each FEValues object splits the given +UpdateFlags into two sets: those flags that describe quantities that +can be computed on the reference cell once at the beginning +("update_once"), and those that require recomputation on each cell +("update_each"). + +Furthermore, each FEValues object is associated with a Mapping and a +FiniteElement object that do the actual computations. It therefore +keeps a set up UpdateFlags for both update_once and update_each +operations for both the Mapping and the FiniteElement object it is +associated with. + + +

Update once or each

+ +Sometimes in order to compute one quantity, something else also has to be +computed. For example, if only update_values is requested by the user of +a FEValues object and if the associated finite element is of type FE_Q, +then update_once=update_values and update_each=0 +for the FiniteElement computations and update_once=update_each=0 +for the Mapping object: we can compute the values of the shape function +at the quadrature points of each cell once on the reference cell, and they +will have the same value at the quadrature point of a real cell. There is +nothing the finite element object has to do when we move to the next cell +(these would be update_each calculations) and there is nothing the mapping +has to do, either up front (update_once) or on each cell (update_each). + +However, this is not the case if we used a FE_RaviartThomas element: there, +computing the values of the shape functions on a cell involves knowing the +Jacobian of the mapping, and so if a user requests update_values +of a FEValues object that uses a FE_RaviartThomas element, then we can set +update_once=update_values and update_each=0 +for the FiniteElement, but need to set update_once=0 +update_each=update_jacobians for the Mapping object. + +To accomodate this structure, at the time a FEValues object is constructed, +it asks both the FiniteElement and the Mapping object it uses the following:
    -
  1. Are any additional values required in order to compute these -values? If so, add these flags to the current set. For instance, the -derivative od a standard scalar element requires the inverse of the +
  2. Are any additional values required in order to compute the +currently required values? If so, add these flags to the current set. +Another example to the one above would be that the +derivative of a standard scalar element requires the inverse of the Jacobian of the Mapping.
  3. Given the enhanced set, determine the subsets of values that is performed on the reference cell only and on each cell, respectively. @@ -64,32 +130,29 @@ actual set of flags from the desired one looks like this: That is, a FiniteElement can set additional flags which are honored by the Mapping. -@section update_flags_Functions Generation of the actual data -The computation of the fields accessible through FEValues proceeds in -two different steps: -@subsection update_flags_Init Initialization +

    Generation of the actual data

    -Functions of FEValues classes involved: -
      -
    • FEValues::initialize() -
    • FEFaceValues::initialize() -
    • FESubfaceValues::initialize() -
    +As outlined above, data is computed at two different times: once at +the beginning on the reference cell, and once whenever we move to an +actual cell. The functions involved in each of these steps are +discussed next: + + +

    Initialization

    -These function is called by the constructor of FEValues, FEFaceValues -and FESubfaceValues, respectively, and allows the +Computing data on the reference cell before we even visit the first +real cell is a two-step process. First, the constructor of FEValues, +FEFaceValues and FESubfaceValues, respectively, need to allow the Mapping and FiniteElement objects to set up internal data -structures. These structures are internal in the sense that FEValues -only forwards them to Mapping and FiniteElement if necessary, but does -not access their contents. FEValues just stores a pointer to their -base class Mapping::InternalDataBase. Concrete Mapping and -FiniteElement classes will derive their own sutructures containing all -the data needed to speed up updating the fields of FEValues on a cell -or face. - -The functions actually providing these internal data structures are: +structures. These structures are internal in the following sense: the +FEValues object asks the finite element and mapping objects to create +an object of type FiniteElement::InternalDataBase and +Mapping::InternalDataBase each; the actual finite element and mapping +class may in fact create objects of a derived type if they wish to +store some data beyond what these base classes already provide. The +functions involved in this are
    • Mapping::get_data()
    • Mapping::get_face_data() @@ -99,9 +162,25 @@ The functions actually providing these internal data structures are:
    • FiniteElement::get_subface_data()
    -@subsection update_flags_Reinit Reinitialization for a mesh cell +The FEValues object then takes over ownership of these objects and will +destroy them at the end of the FEValues object's lifetime. After this, +the FEValues object asks the FiniteElement and Mapping objects to fill +these InternalDataBase objects with the data that pertains to what +can and needs to be computed on the reference cell. This is done in these +functions: +
      +
    • FEValues::initialize() +
    • FEFaceValues::initialize() +
    • FESubfaceValues::initialize() +
    + + +

    Reinitialization for a mesh cell

    -Functions of FEValues classes involved: +Once initialization is over and we call FEValues::reinit, FEFaceValues::reinit +or FESubfaceValues::reinit to move to a concrete cell or face, we need +to calculate the update_each kinds of data. This done in the following +functions:
    • FEValues::reinit() calls Mapping::fill_fe_values(), then FiniteElement::fill_fe_values()
    • FEFaceValues::reinit() calls Mapping::fill_fe_face_values(), then FiniteElement::fill_fe_face_values() -- 2.39.5