From 3c85a23805100915e3c5fb193145a95c6d30d64f Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Mon, 16 Nov 2020 18:28:38 -0700 Subject: [PATCH] Better document some FEValues functions. While there, mark some input arguments as 'const'. --- include/deal.II/fe/fe_values.h | 176 ++++++++++++++++++++------------- source/fe/fe_values.cc | 24 ++--- 2 files changed, 119 insertions(+), 81 deletions(-) diff --git a/include/deal.II/fe/fe_values.h b/include/deal.II/fe/fe_values.h index c0ad804b59..6bb2024f22 100644 --- a/include/deal.II/fe/fe_values.h +++ b/include/deal.II/fe/fe_values.h @@ -2387,20 +2387,58 @@ public: std::vector> &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 indices 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 value 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 indices 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 value 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 indices 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 value should allow for the same * multiple of the components of the finite element. * @@ -2461,10 +2481,6 @@ public: * values[p][i] if quadrature_points_fastest == false, and * values[i][p] 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 &indices, ArrayView< std::vector>> - 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>> - & 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 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 &indices, ArrayView< std::vector>> - 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> &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 & 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> &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 & indices, std::vector> &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>> - & 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 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 &indices, ArrayView< std::vector>> - third_derivatives, - bool quadrature_points_fastest = false) const; + third_derivatives, + const bool quadrature_points_fastest = false) const; //@} /// @name Cell degrees of freedom diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 34f1d8c6b8..ea659a9abf 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -3677,7 +3677,7 @@ FEValuesBase::get_function_values( const InputVector & fe_function, const ArrayView & indices, ArrayView> 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::get_function_gradients( const InputVector & fe_function, const ArrayView &indices, ArrayView>> - 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::get_function_hessians( const InputVector &fe_function, std::vector< std::vector>> - & 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::get_function_hessians( const InputVector & fe_function, const ArrayView &indices, ArrayView>> - 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::get_function_laplacians( const InputVector & fe_function, const ArrayView & indices, std::vector> &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::get_function_third_derivatives( const InputVector &fe_function, std::vector< std::vector>> - & 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::get_function_third_derivatives( const InputVector & fe_function, const ArrayView &indices, ArrayView>> - 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, -- 2.39.5