]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Better document some FEValues functions. 11194/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 17 Nov 2020 01:28:38 +0000 (18:28 -0700)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 17 Nov 2020 01:28:38 +0000 (18:28 -0700)
While there, mark some input arguments as 'const'.

include/deal.II/fe/fe_values.h
source/fe/fe_values.cc

index c0ad804b591a57030f12ed78c6009026e3cf29da..6bb2024f22137ad6747f58fc0edf2b0667171c7e 100644 (file)
@@ -2387,20 +2387,58 @@ public:
     std::vector<Vector<typename InputVector::value_type>> &values) const;
 
   /**
-   * Generate function values from an arbitrary vector.
-   *
-   * This function offers the possibility to extract function values in
-   * quadrature points from vectors not corresponding to a whole
-   * discretization.
-   *
-   * The vector <tt>indices</tt> corresponds to the degrees of freedom on a
-   * single cell. Its length may even be a multiple of the number of dofs per
-   * cell. Then, the vectors in <tt>value</tt> should allow for the same
-   * multiple of the components of the finite element.
+   * Generate function values from an arbitrary vector. This function
+   * does in essence the same as the first function of this name above,
+   * except that it does not make the assumption that the input vector
+   * corresponds to a DoFHandler that describes the unknowns of a finite
+   * element field (and for which we would then assume that
+   * `fe_function.size() == dof_handler.n_dofs()`). Rather, the nodal
+   * values corresponding to the current cell are elements of an otherwise
+   * arbitrary vector, and these elements are indexed by the second
+   * argument to this function. What the rest of the `fe_function` input
+   * argument corresponds to is of no consequence to this function.
+   *
+   * Given this, the function above corresponds to passing `fe_function`
+   * as first argument to the current function, and using the
+   * `local_dof_indices` array that results from the following call as
+   * second argument to the current function:
+   * @code
+   *   cell->get_dof_indices (local_dof_indices);
+   * @endcode
+   * (See DoFCellAccessor::get_dof_indices() for more information.)
    *
-   * You may want to use this function, if you want to access just a single
-   * block from a BlockVector, if you have a multi-level vector or if you
-   * already have a local representation of your finite element data.
+   * Likewise, the function above is equivalent to calling
+   * @code
+   *   cell->get_dof_values (fe_function, local_dof_values);
+   * @endcode
+   * and then calling the current function with `local_dof_values` as
+   * first argument, and an array with indices `{0,...,fe.dofs_per_cell-1}`
+   * as second argument.
+   *
+   * The point of the current function is that one sometimes wants to
+   * evaluate finite element functions at quadrature points with nodal
+   * values that are not stored in a global vector -- for example, one could
+   * modify these local values first, such as by applying a limiter
+   * or by ensuring that all nodal values are positive, before evaluating
+   * the finite element field that corresponds to these local values on the
+   * current cell. Another application is where one wants to postprocess
+   * the solution on a cell into a different finite element space on every
+   * cell, without actually creating a corresponding DoFHandler -- in that
+   * case, all one would compute is a local representation of that
+   * postprocessed function, characterized by its nodal values; this function
+   * then allows the evaluation of that representation at quadrature points.
+   *
+   * @param[in] fe_function A vector of nodal values. This vector can have
+   *   an arbitrary size, as long as all elements index by `indices` can
+   *   actually be accessed.
+   *
+   * @param[in] indices A vector of indices into `fe_function`. This vector
+   *   must have length equal to the number of degrees of freedom on the
+   *   current cell, and must identify elements in `fe_function` in the
+   *   order in which degrees of freedom are indexed on the reference cell.
+   *
+   * @param[out] values A vector of values of the given finite element field,
+   *   at the quadrature points on the current object.
    *
    * @dealiiRequiresUpdateFlags{update_values}
    */
@@ -2414,21 +2452,8 @@ public:
   /**
    * Generate vector function values from an arbitrary vector.
    *
-   * This function offers the possibility to extract function values in
-   * quadrature points from vectors not corresponding to a whole
-   * discretization.
-   *
-   * The vector <tt>indices</tt> corresponds to the degrees of freedom on a
-   * single cell. Its length may even be a multiple of the number of dofs per
-   * cell. Then, the vectors in <tt>value</tt> should allow for the same
-   * multiple of the components of the finite element.
-   *
-   * You may want to use this function, if you want to access just a single
-   * block from a BlockVector, if you have a multi-level vector or if you
-   * already have a local representation of your finite element data.
-   *
-   * Since this function allows for fairly general combinations of argument
-   * sizes, be aware that the checks on the arguments may not detect errors.
+   * This function corresponds to the previous one, just for the vector-valued
+   * case.
    *
    * @dealiiRequiresUpdateFlags{update_values}
    */
@@ -2441,14 +2466,9 @@ public:
 
 
   /**
-   * Generate vector function values from an arbitrary vector.
-   *
-   * This function offers the possibility to extract function values in
-   * quadrature points from vectors not corresponding to a whole
-   * discretization.
-   *
-   * The vector <tt>indices</tt> corresponds to the degrees of freedom on a
-   * single cell. Its length may even be a multiple of the number of dofs per
+   * Generate vector function values from an arbitrary vector. This
+   * function is similar to the previous one, but the `indices`
+   * vector may also be a multiple of the number of dofs per
    * cell. Then, the vectors in <tt>value</tt> should allow for the same
    * multiple of the components of the finite element.
    *
@@ -2461,10 +2481,6 @@ public:
    * <tt>values[p][i]</tt> if <tt>quadrature_points_fastest == false</tt>, and
    * <tt>values[i][p]</tt> otherwise.
    *
-   * You may want to use this function, if you want to access just a single
-   * block from a BlockVector, if you have a multi-level vector or if you
-   * already have a local representation of your finite element data.
-   *
    * Since this function allows for fairly general combinations of argument
    * sizes, be aware that the checks on the arguments may not detect errors.
    *
@@ -2550,8 +2566,10 @@ public:
       &gradients) const;
 
   /**
-   * Function gradient access with more flexibility. See get_function_values()
-   * with corresponding arguments.
+   * This function relates to the first of the get_function_gradients() function
+   * above in the same way as the get_function_values() with similar arguments
+   * relates to the first of the get_function_values() functions. See there for
+   * more information.
    *
    * @dealiiRequiresUpdateFlags{update_gradients}
    */
@@ -2564,8 +2582,10 @@ public:
       &gradients) const;
 
   /**
-   * Function gradient access with more flexibility. See get_function_values()
-   * with corresponding arguments.
+   * This function relates to the first of the get_function_gradients() function
+   * above in the same way as the get_function_values() with similar arguments
+   * relates to the first of the get_function_values() functions. See there for
+   * more information.
    *
    * @dealiiRequiresUpdateFlags{update_gradients}
    */
@@ -2576,8 +2596,8 @@ public:
     const ArrayView<const types::global_dof_index> &indices,
     ArrayView<
       std::vector<Tensor<1, spacedim, typename InputVector::value_type>>>
-         gradients,
-    bool quadrature_points_fastest = false) const;
+               gradients,
+    const bool quadrature_points_fastest = false) const;
 
   //@}
   /// @name Access to second derivatives
@@ -2652,12 +2672,16 @@ public:
     const InputVector &fe_function,
     std::vector<
       std::vector<Tensor<2, spacedim, typename InputVector::value_type>>>
-      &  hessians,
-    bool quadrature_points_fastest = false) const;
+      &        hessians,
+    const bool quadrature_points_fastest = false) const;
 
   /**
-   * Access to the second derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_hessians() function
+   * above in the same way as the get_function_values() with similar arguments
+   * relates to the first of the get_function_values() functions. See there for
+   * more information.
+   *
+   * @dealiiRequiresUpdateFlags{update_hessians}
    */
   template <class InputVector>
   void
@@ -2668,8 +2692,10 @@ public:
       &hessians) const;
 
   /**
-   * Access to the second derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_hessians() function
+   * above in the same way as the get_function_values() with similar arguments
+   * relates to the first of the get_function_values() functions. See there for
+   * more information.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2680,8 +2706,8 @@ public:
     const ArrayView<const types::global_dof_index> &indices,
     ArrayView<
       std::vector<Tensor<2, spacedim, typename InputVector::value_type>>>
-         hessians,
-    bool quadrature_points_fastest = false) const;
+               hessians,
+    const bool quadrature_points_fastest = false) const;
 
   /**
    * Compute the (scalar) Laplacian (i.e. the trace of the tensor of second
@@ -2755,8 +2781,10 @@ public:
     std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
 
   /**
-   * Access to the second derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_laplacians()
+   * function above in the same way as the get_function_values() with similar
+   * arguments relates to the first of the get_function_values() functions. See
+   * there for more information.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2768,8 +2796,10 @@ public:
     std::vector<typename InputVector::value_type> & laplacians) const;
 
   /**
-   * Access to the second derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_laplacians()
+   * function above in the same way as the get_function_values() with similar
+   * arguments relates to the first of the get_function_values() functions. See
+   * there for more information.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2781,8 +2811,10 @@ public:
     std::vector<Vector<typename InputVector::value_type>> &laplacians) const;
 
   /**
-   * Access to the second derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_laplacians()
+   * function above in the same way as the get_function_values() with similar
+   * arguments relates to the first of the get_function_values() functions. See
+   * there for more information.
    *
    * @dealiiRequiresUpdateFlags{update_hessians}
    */
@@ -2792,7 +2824,7 @@ public:
     const InputVector &                                         fe_function,
     const ArrayView<const types::global_dof_index> &            indices,
     std::vector<std::vector<typename InputVector::value_type>> &laplacians,
-    bool quadrature_points_fastest = false) const;
+    const bool quadrature_points_fastest = false) const;
 
   //@}
   /// @name Access to third derivatives of global finite element fields
@@ -2867,12 +2899,16 @@ public:
     const InputVector &fe_function,
     std::vector<
       std::vector<Tensor<3, spacedim, typename InputVector::value_type>>>
-      &  third_derivatives,
-    bool quadrature_points_fastest = false) const;
+      &        third_derivatives,
+    const bool quadrature_points_fastest = false) const;
 
   /**
-   * Access to the third derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_third_derivatives()
+   * function above in the same way as the get_function_values() with similar
+   * arguments relates to the first of the get_function_values() functions. See
+   * there for more information.
+   *
+   * @dealiiRequiresUpdateFlags{update_3rd_derivatives}
    */
   template <class InputVector>
   void
@@ -2883,8 +2919,10 @@ public:
       &third_derivatives) const;
 
   /**
-   * Access to the third derivatives of a function with more flexibility. See
-   * get_function_values() with corresponding arguments.
+   * This function relates to the first of the get_function_third_derivatives()
+   * function above in the same way as the get_function_values() with similar
+   * arguments relates to the first of the get_function_values() functions. See
+   * there for more information.
    *
    * @dealiiRequiresUpdateFlags{update_3rd_derivatives}
    */
@@ -2895,8 +2933,8 @@ public:
     const ArrayView<const types::global_dof_index> &indices,
     ArrayView<
       std::vector<Tensor<3, spacedim, typename InputVector::value_type>>>
-         third_derivatives,
-    bool quadrature_points_fastest = false) const;
+               third_derivatives,
+    const bool quadrature_points_fastest = false) const;
   //@}
 
   /// @name Cell degrees of freedom
index 34f1d8c6b842bc8a8edd783f924140c79622bb75..ea659a9abf79a7f73fb078bdb4627e050649ca35 100644 (file)
@@ -3677,7 +3677,7 @@ FEValuesBase<dim, spacedim>::get_function_values(
   const InputVector &                                      fe_function,
   const ArrayView<const types::global_dof_index> &         indices,
   ArrayView<std::vector<typename InputVector::value_type>> values,
-  bool quadrature_points_fastest) const
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_values,
@@ -3790,8 +3790,8 @@ FEValuesBase<dim, spacedim>::get_function_gradients(
   const InputVector &                             fe_function,
   const ArrayView<const types::global_dof_index> &indices,
   ArrayView<std::vector<Tensor<1, spacedim, typename InputVector::value_type>>>
-       gradients,
-  bool quadrature_points_fastest) const
+             gradients,
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   // Size of indices must be a multiple of dofs_per_cell such that an integer
@@ -3874,8 +3874,8 @@ FEValuesBase<dim, spacedim>::get_function_hessians(
   const InputVector &fe_function,
   std::vector<
     std::vector<Tensor<2, spacedim, typename InputVector::value_type>>>
-    &  hessians,
-  bool quadrature_points_fastest) const
+    &        hessians,
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_hessians,
@@ -3905,8 +3905,8 @@ FEValuesBase<dim, spacedim>::get_function_hessians(
   const InputVector &                             fe_function,
   const ArrayView<const types::global_dof_index> &indices,
   ArrayView<std::vector<Tensor<2, spacedim, typename InputVector::value_type>>>
-       hessians,
-  bool quadrature_points_fastest) const
+             hessians,
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_hessians,
@@ -4043,7 +4043,7 @@ FEValuesBase<dim, spacedim>::get_function_laplacians(
   const InputVector &                                         fe_function,
   const ArrayView<const types::global_dof_index> &            indices,
   std::vector<std::vector<typename InputVector::value_type>> &laplacians,
-  bool quadrature_points_fastest) const
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(indices.size() % dofs_per_cell == 0,
@@ -4126,8 +4126,8 @@ FEValuesBase<dim, spacedim>::get_function_third_derivatives(
   const InputVector &fe_function,
   std::vector<
     std::vector<Tensor<3, spacedim, typename InputVector::value_type>>>
-    &  third_derivatives,
-  bool quadrature_points_fastest) const
+    &        third_derivatives,
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_3rd_derivatives,
@@ -4157,8 +4157,8 @@ FEValuesBase<dim, spacedim>::get_function_third_derivatives(
   const InputVector &                             fe_function,
   const ArrayView<const types::global_dof_index> &indices,
   ArrayView<std::vector<Tensor<3, spacedim, typename InputVector::value_type>>>
-       third_derivatives,
-  bool quadrature_points_fastest) const
+             third_derivatives,
+  const bool quadrature_points_fastest) const
 {
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_3rd_derivatives,

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.