From 2ef5db7cb3927958e4dbc474f092d63b43be0413 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Sat, 9 Jan 2021 15:32:43 -0700 Subject: [PATCH] Add some comments to VectorTools::interpolate_boundary_values(). --- doc/doxygen/references.bib | 14 +++ .../deal.II/numerics/vector_tools_boundary.h | 86 +++++++++++-------- 2 files changed, 63 insertions(+), 37 deletions(-) diff --git a/doc/doxygen/references.bib b/doc/doxygen/references.bib index 4aa9c31852..739844d763 100644 --- a/doc/doxygen/references.bib +++ b/doc/doxygen/references.bib @@ -1049,3 +1049,17 @@ year = {2008}, volume = {40}, doi = {10.1002/(SICI)1097-0207(19970228)40:4<727::AID-NME86>3.0.CO;2-N} } + +@article{Bartels2004, + doi = {10.1007/s00211-004-0548-3}, + url = {https://doi.org/10.1007/s00211-004-0548-3}, + year = {2004}, + month = sep, + publisher = {Springer Science and Business Media {LLC}}, + volume = {99}, + number = {1}, + pages = {1--24}, + author = {S. Bartels and C. Carstensen and G. Dolzmann}, + title = {Inhomogeneous Dirichlet conditions in a priori and a posteriori finite element error analysis}, + journal = {Numerische Mathematik} +} \ No newline at end of file diff --git a/include/deal.II/numerics/vector_tools_boundary.h b/include/deal.II/numerics/vector_tools_boundary.h index 69e7818119..332cdc0aad 100644 --- a/include/deal.II/numerics/vector_tools_boundary.h +++ b/include/deal.II/numerics/vector_tools_boundary.h @@ -47,7 +47,8 @@ namespace VectorTools //@{ /** - * Compute Dirichlet boundary conditions. This function makes up a map of + * Compute constraints on the solution that corresponds to the imposition + * of Dirichlet boundary conditions. This function creates a map of * degrees of freedom subject to Dirichlet boundary conditions and the * corresponding values to be assigned to them, by interpolation around the * boundary. For each degree of freedom at the boundary, if its index @@ -95,6 +96,42 @@ namespace VectorTools * of these non-primitive shape functions must be @p false. * * See the general documentation of this namespace for more information. + * + * @note When solving a partial differential equation with boundary + * conditions $u|_{\partial\Omega}=g$ (or on *parts* of the boundary), + * then this boundary condition is in general not satisfiable exactly + * using finite elements in the form $u_h|_{\partial\Omega}=g$. That is + * because the function $g$ is generally not a polynomial, whereas + * $u_h|_{\partial\Omega}$ *is* a polynomial on each face of the + * mesh that is located at the boundary. In other words, it is in + * general not possible to *impose* such boundary condition; what one + * *can* do, however, is to impose + * @f[ u_h|_{\partial\Omega}=I_h^{\partial\Omega} g, @f] + * where $I_h^{\partial\Omega} g$ is a function that equals $g$ at each node + * of the finite element space located on the boundary, and is piecewise + * polynomial in between. In other words, $I_h^{\partial\Omega}$ is an + * *interpolation operator* and $I_h^{\partial\Omega} g$ are the + * interpolated boundary values -- thus the name. The use of + * $I_h^{\partial\Omega} g$ instead of $g$ as boundary values imposes + * an additional error (in the same spirit as using quadrature introduces + * an additional error compared to being able to compute the integrals of + * the weak form exactly). In most cases, this additional error is of the + * same order as the other error terms in the finite element method, + * though there are some subtle differences when measuring the error in + * the $L^2$ norm. For some details, see @cite Bartels2004 . + * + * @note An alternative to using the interpolant, + * @f[ u_h|_{\partial\Omega}=I_h^{\partial\Omega} g @f] + * is to use the *projection* of the boundary values $g$ onto the + * finite element space on the boundary: + * @f[ u_h|_{\partial\Omega}=\Pi_h^{\partial\Omega} g. @f] + * The projection is available using the project_boundary_values() + * function. Using the projection may have some theoretical advantages + * (see again @cite Bartels2004) but has the practical disadvantage + * that computing the projection is far more expensive than computing + * the interpolation because the latter can be done one face at a time + * whereas the projection requires the solution of a problem on the entire + * boundary. */ template void @@ -210,42 +247,9 @@ namespace VectorTools * the module on * @ref constraints. * - * The parameter @p boundary_component corresponds to the number @p - * boundary_id of the face. - * - * The flags in the last parameter, @p component_mask denote which - * components of the finite element space shall be interpolated. If it is - * left as specified by the default value (i.e. an empty array), all - * components are interpolated. If it is different from the default value, - * it is assumed that the number of entries equals the number of components - * in the boundary functions and the finite element, and those components in - * the given boundary function will be used for which the respective flag - * was set in the component mask. See also - * @ref GlossComponentMask. - * As an example, assume that you are solving the Stokes equations in 2d, - * with variables $(u,v,p)$ and that you only want to interpolate boundary - * values for the pressure, then the component mask should correspond to - * (true,true,false). - * - * @note Whether a component mask has been specified or not, the number of - * components of the functions in @p function_map must match that of the - * finite element used by @p dof. In other words, for the example above, you - * need to provide a Function object that has 3 components (the two - * velocities and the pressure), even though you are only interested in the - * first two of them. interpolate_boundary_values() will then call this - * function to obtain a vector of 3 values at each interpolation point but - * only take the first two and discard the third. In other words, you are - * free to return whatever you like in the third component of the vector - * returned by Function::vector_value, but the Function object must state - * that it has 3 components. - * - * If the finite element used has shape functions that are non-zero in more - * than one component (in deal.II speak: they are non-primitive), then these - * components can presently not be used for interpolating boundary values. - * Thus, the elements in the component mask corresponding to the components - * of these non-primitive shape functions must be @p false. - * - * See the general documentation of this namespace for more information. + * This function is fundamentally equivalent to the ones above except that it + * puts its results into an AffineConstraint object rather than a `std::map`. + * See the functions above for more comments. * * @ingroup constraints */ @@ -410,6 +414,14 @@ namespace VectorTools * vector component of the finite element used in @p dof. This entry is the * component number in @p boundary_functions that should be used for this * component in @p dof. By default, no remapping is applied. + * + * @note Using the *projection* rather than the *interpolation* of boundary + * values makes relatively little difference in practice. That said, + * it is far more computationally expensive to compute projections because + * the require the solution of a problem that couples all unknowns on the + * boundary, whereas interpolation works on one face at a time. For + * some more theoretical considerations, see the documentation of the first + * interpolate_boundary_values() function above. */ template void -- 2.39.5