From: Wolfgang Bangerth Date: Mon, 16 Jan 2017 19:36:35 +0000 (-0700) Subject: Document FETools::compute_node_matrix. X-Git-Tag: v8.5.0-rc1~232^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F3809%2Fhead;p=dealii.git Document FETools::compute_node_matrix. --- diff --git a/include/deal.II/fe/fe_tools.h b/include/deal.II/fe/fe_tools.h index 203013761a..5eb1dcb976 100644 --- a/include/deal.II/fe/fe_tools.h +++ b/include/deal.II/fe/fe_tools.h @@ -232,31 +232,65 @@ namespace FETools FullMatrix &matrix); /** - * Compute the matrix of nodal values of a finite element applied to all its - * shape functions. - * - * This function is supposed to help building finite elements from - * polynomial spaces and should be called inside the constructor of an - * element. Applied to a completely initialized finite element, the result - * should be the unit matrix by definition of the node values. - * - * Using this matrix allows the construction of the basis of shape functions - * in two steps. - * - *
    - * - *
  1. Define the space of shape functions using an arbitrary basis - * wj and compute the matrix M of node functionals - * Ni applied to these basis functions. - * - *
  2. Compute the basis vj of the finite element shape - * function space by applying M-1 to the basis - * wj. - * - *
- * - * @note The FiniteElement must provide generalized support points and and - * interpolation functions. + * This is a rather specialized function used during the construction of + * finite element objects. It is used to build the basis of shape functions + * for an element given a set of polynomials and interpolation points. The + * function is only implemented for finite elements with exactly @p dim + * vector components. In particular, this applies to classes derived from + * the FE_PolyTensor class. + * + * Specifically, the purpose of this function is as follows: FE_PolyTensor + * receives, from its derived classes, an argument that describes a polynomial + * space. This space may be parameterized in terms of monomials, or in some + * other way, but is in general not in the form that we use for finite + * elements where we typically want to use a basis that is derived from + * some kind of node functional (e.g., the interpolation at specific points). + * Concretely, assume that the basis used by the polynomial space is + * $\{\tilde\varphi_j(\mathbf x)\}_{j=1}^N$, and that the node functionals + * of the finite element are $\{\Psi_i\}_{i=1}^N$. We then want to compute a + * basis $\{\varphi_j(\mathbf x)\}_{j=1}^N$ for the finite element space so + * that $\Psi_i[\varphi_j] = \delta_{ij}$. To do this, we can set + * $\varphi_j(\mathbf x) = \sum_{k=1}^N c_{jk} \tilde\varphi_k(\mathbf x)$ + * where we need to determine the expansion coefficients $c_{jk}$. We do this + * by applying $\Psi_i$ to both sides of the equation, to obtain + * @f{align*} + * \Psi_i [\varphi_j] = \sum_{k=1}^N c_{jk} \Psi_i[\tilde\varphi_k], + * @f} + * and we know that the left hand side equals $\delta_{ij}$. + * If you think of this as a system of $N\times N$ equations for the + * elements of a matrix on the left and on the right, then this can be + * written as + * @f{align*} + * I = C X^T + * @f} + * where $C$ is the matrix of coefficients $c_{jk}$ and + * $X_{ik} = \Psi_i[\tilde\varphi_k]$. Consequently, in order to compute + * the expansion coefficients $C=X^{-T}$, we need to apply the node + * functionals to all functions of the "raw" basis of the polynomial space. + * + * Until the finite element receives this matrix $X$ back, it describes its + * shape functions (e.g., in FiniteElement::shape_value()) in the form + * $\tilde\varphi_j$. After it calls this function, it has the expansion + * coefficients and can describe its shape functions as $\varphi_j$. + * + * This function therefore computes this matrix $X$, for the following + * specific circumstances: + * - That the node functionals $\Psi_i$ are point evaluations at points + * $\mathbf x_i$ that the finite element in question describes via its + * "generalized" support points (through + * FiniteElement::get_generalized_support_points()). These point + * evaluations need to necessarily evaluate the value of a shape + * function at that point (the shape function may be vector-valued, and + * so the functional may be a linear combination of the individual + * components of the values); but, in particular, the nodal functions may + * not be integrals over entire edges or faces, + * or other non-local functionals. In other words, we assume that + * $\Psi_i[\tilde\varphi_j] = f_j(\tilde\varphi_j(\mathbf x_i))$ + * where $f_j$ is a function of the (possibly vector-valued) argument + * that returns a scalar. + * - That the finite element has exactly @p dim vector components. + * - That the function $f_j$ is given by whatever the element implements + * through the FiniteElement::interpolate() function. */ template void compute_node_matrix(FullMatrix &M,