From 2ca5ba39e3a77afda898b04464aac21de5faaa3a Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Tue, 11 Sep 2007 15:08:43 +0000 Subject: [PATCH] start unifying names of gradients and hessians git-svn-id: https://svn.dealii.org/trunk@15191 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_update_flags.h | 12 +- deal.II/deal.II/include/fe/fe_values.h | 221 +++++++++++++++--- deal.II/deal.II/source/fe/fe.cc | 8 +- .../deal.II/source/fe/fe_dgp_nonparametric.cc | 16 +- deal.II/deal.II/source/fe/fe_values.cc | 44 ++-- .../deal.II/source/fe/fe_values.instance.h | 14 +- .../numerics/derivative_approximation.cc | 4 +- 7 files changed, 242 insertions(+), 77 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 df31c2bf09..51be7bd88f 100644 --- a/deal.II/deal.II/include/fe/fe_update_flags.h +++ b/deal.II/deal.II/include/fe/fe_update_flags.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 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 @@ -109,7 +109,7 @@ enum UpdateFlags * functions in coordinates of * the real cell. */ - update_second_derivatives = 0x0004, + update_hessians = 0x0004, //! Outter normal vector, not normalized /** * Vector product of tangential @@ -211,7 +211,13 @@ enum UpdateFlags * is used in conjunction with * H_div subspaces like RT and ABF. */ - update_cell_JxW_values = 0x2000 + update_cell_JxW_values = 0x2000, + /** + * @deprecated Update second + * derivatives. + */ + update_second_derivatives = update_hessians + }; diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index e081ead4aa..4162fdd2c0 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name: $ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 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 @@ -125,7 +125,7 @@ class FEValuesData * Likewise for second order * derivatives. */ - typedef std::vector > > GradGradVector; + typedef std::vector > > HessianVector; /** * Store the values of the shape @@ -154,7 +154,7 @@ class FEValuesData * for the layout of the data in * this field. */ - GradGradVector shape_2nd_derivatives; + HessianVector shape_hessians; /** * Store an array of weights @@ -238,7 +238,7 @@ class FEValuesData * given shape function occupies * in the #shape_values, * #shape_gradients and - * #shape_2nd_derivatives + * #shape_hessians * arrays. If all shape functions * are primitive, then this is * the identity mapping. If, on @@ -357,7 +357,7 @@ class FEValuesData * and you have to walk over all (or only the non-zero) components of * the shape function using this set of functions. * - *
  • get_function_values(), get_function_grads(), etc.: Compute a + *
  • get_function_values(), get_function_gradients(), etc.: Compute a * finite element function or its derivative in quadrature points. * *
  • reinit: initialize the FEValues object for a certain cell. @@ -591,6 +591,13 @@ class FEValuesBase : protected FEValuesData, * function. */ const Tensor<2,dim> & + shape_hessian (const unsigned int function_no, + const unsigned int point_no) const; + + /** + * @deprecated Wrapper for shape_hessian() + */ + const Tensor<2,dim> & shape_2nd_derivative (const unsigned int function_no, const unsigned int point_no) const; @@ -603,7 +610,7 @@ class FEValuesBase : protected FEValuesData, * is scalar, then only component * zero is allowed and the return * value equals that of the - * shape_2nd_derivative() + * shape_hessian() * function. If the finite * element is vector valued but * all shape functions are @@ -611,7 +618,7 @@ class FEValuesBase : protected FEValuesData, * non-zero in only one * component), then the value * returned by - * shape_2nd_derivative() + * shape_hessian() * equals that of this function * for exactly one * component. This function is @@ -622,6 +629,14 @@ class FEValuesBase : protected FEValuesData, * function cannot be used. */ Tensor<2,dim> + shape_hessian_component (const unsigned int function_no, + const unsigned int point_no, + const unsigned int component) const; + + /** + * @deprecated Wrapper for shape_hessian_component() + */ + Tensor<2,dim> shape_2nd_derivative_component (const unsigned int function_no, const unsigned int point_no, const unsigned int component) const; @@ -860,7 +875,7 @@ class FEValuesBase : protected FEValuesData, * finite element in use is a scalar one, * i.e. has only one vector component. If * it is a vector-valued one, then use - * the other get_function_grads() + * the other get_function_gradients() * function. * * The function assumes that the @@ -884,6 +899,14 @@ class FEValuesBase : protected FEValuesData, * on the unit cell). */ template + void get_function_gradients (const InputVector &fe_function, + std::vector > &gradients) const; + + /** + * @deprecated Use + * get_function_gradients() instead. + */ + template void get_function_grads (const InputVector &fe_function, std::vector > &gradients) const; @@ -923,6 +946,15 @@ class FEValuesBase : protected FEValuesData, * on the unit cell). */ template + void get_function_gradients (const InputVector &fe_function, + std::vector > > &gradients) const; + + + /** + * @deprecated Use + * get_function_gradients() instead. + */ + template void get_function_grads (const InputVector &fe_function, std::vector > > &gradients) const; @@ -933,6 +965,15 @@ class FEValuesBase : protected FEValuesData, * corresponding arguments. */ template + void get_function_gradients (const InputVector& fe_function, + const VectorSlice >& indices, + std::vector >& gradients) const; + + /** + * @deprecated Use + * get_function_gradients() instead. + */ + template void get_function_grads (const InputVector& fe_function, const VectorSlice >& indices, std::vector >& gradients) const; @@ -944,6 +985,16 @@ class FEValuesBase : protected FEValuesData, * corresponding arguments. */ template + void get_function_gradients (const InputVector& fe_function, + const VectorSlice >& indices, + std::vector > >& gradients, + bool quadrature_points_fastest = false) const; + + /** + * @deprecated Use + * get_function_gradients() instead. + */ + template void get_function_grads (const InputVector& fe_function, const VectorSlice >& indices, std::vector > >& gradients, @@ -958,7 +1009,7 @@ class FEValuesBase : protected FEValuesData, * points. * * The function assumes that the - * @p second_derivatives object + * @p hessians object * already has the correct size. * * This function may only be used if the @@ -966,7 +1017,7 @@ class FEValuesBase : protected FEValuesData, * i.e. has only one vector component. If * it is a vector-valued one, then use * the other - * get_function_2nd_derivatives() + * get_function_hessians() * function. * * The actual data type of the input @@ -987,8 +1038,8 @@ class FEValuesBase : protected FEValuesData, */ template void - get_function_2nd_derivatives (const InputVector& fe_function, - std::vector >& second_derivatives) const; + get_function_hessians (const InputVector& fe_function, + std::vector >& hessians) const; /** @@ -999,7 +1050,7 @@ class FEValuesBase : protected FEValuesData, * @p cell at the quadrature points. * * The function assumes that the - * @p second_derivatives object already has + * @p hessians object already has * the right size. * * This function does the same as @@ -1025,9 +1076,27 @@ class FEValuesBase : protected FEValuesData, */ template void - get_function_2nd_derivatives (const InputVector &fe_function, - std::vector > > &second_derivatives, - bool quadrature_points_fastest = false) const; + get_function_hessians (const InputVector &fe_function, + std::vector > > &hessians, + bool quadrature_points_fastest = false) const; + + /** + * @deprecated Wrapper for get_function_hessians() + */ + template + void + get_function_2nd_derivatives (const InputVector&, + std::vector >&) const; + + /** + * @deprecated Wrapper for get_function_hessians() + */ + template + void + get_function_2nd_derivatives (const InputVector&, + std::vector > >&, + bool = false) const; + //@} /** @@ -1768,7 +1837,7 @@ class FEValues : public FEValuesBase * about degrees of * freedom. These functions are, * above all, the - * get_function_value/grads/2nd_derivatives + * get_function_value/gradients/hessians * functions. If you want to call * these functions, you have to * call the @p reinit variants @@ -2071,7 +2140,7 @@ class FEFaceValues : public FEFaceValuesBase * information about degrees of * freedom. These functions are, * above all, the - * get_function_value/grads/2nd_derivatives + * get_function_value/gradients/hessians * functions. If you want to call * these functions, you have to * call the @p reinit variants @@ -2254,7 +2323,7 @@ class FESubfaceValues : public FEFaceValuesBase * information about degrees of * freedom. These functions are, * above all, the - * get_function_value/grads/2nd_derivatives + * get_function_value/gradients/hessians * functions. If you want to call * these functions, you have to * call the @p reinit variants @@ -2522,22 +2591,22 @@ FEValuesBase::shape_grad_component (const unsigned int i, template inline const Tensor<2,dim> & -FEValuesBase::shape_2nd_derivative (const unsigned int i, +FEValuesBase::shape_hessian (const unsigned int i, const unsigned int j) const { - Assert (this->update_flags & update_second_derivatives, + Assert (this->update_flags & update_hessians, ExcAccessToUninitializedField()); Assert (fe->is_primitive (i), ExcShapeFunctionNotPrimitive(i)); - Assert (ishape_2nd_derivatives.size(), - ExcIndexRange (i, 0, this->shape_2nd_derivatives.size())); - Assert (jshape_2nd_derivatives[0].size(), - ExcIndexRange (j, 0, this->shape_2nd_derivatives[0].size())); + Assert (ishape_hessians.size(), + ExcIndexRange (i, 0, this->shape_hessians.size())); + Assert (jshape_hessians[0].size(), + ExcIndexRange (j, 0, this->shape_hessians[0].size())); // if the entire FE is primitive, // then we can take a short-cut: if (fe->is_primitive()) - return this->shape_2nd_derivatives[i][j]; + return this->shape_hessians[i][j]; else // otherwise, use the mapping // between shape function numbers @@ -2548,19 +2617,29 @@ FEValuesBase::shape_2nd_derivative (const unsigned int i, // question to which vector // component the call of this // function refers - return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j]; + return this->shape_hessians[this->shape_function_to_row_table[i]][j]; } +template +inline +const Tensor<2,dim> & +FEValuesBase::shape_2nd_derivative (const unsigned int i, + const unsigned int j) const +{ + return shape_hessian(i,j); +} + + template inline Tensor<2,dim> -FEValuesBase::shape_2nd_derivative_component (const unsigned int i, +FEValuesBase::shape_hessian_component (const unsigned int i, const unsigned int j, const unsigned int component) const { - Assert (this->update_flags & update_second_derivatives, + Assert (this->update_flags & update_hessians, ExcAccessToUninitializedField()); Assert (component < fe->n_components(), ExcIndexRange(component, 0, fe->n_components())); @@ -2577,7 +2656,7 @@ FEValuesBase::shape_2nd_derivative_component (const unsigned int i, if (fe->is_primitive(i)) { if (component == fe->system_to_component_index(i).first) - return this->shape_2nd_derivatives[this->shape_function_to_row_table[i]][j]; + return this->shape_hessians[this->shape_function_to_row_table[i]][j]; else return Tensor<2,dim>(); } @@ -2609,11 +2688,21 @@ FEValuesBase::shape_2nd_derivative_component (const unsigned int i, std::count (fe->get_nonzero_components(i).begin(), fe->get_nonzero_components(i).begin()+component, true)); - return this->shape_2nd_derivatives[row][j]; + return this->shape_hessians[row][j]; }; } +template +inline +Tensor<2,dim> +FEValuesBase::shape_2nd_derivative_component (const unsigned int i, + const unsigned int j, + const unsigned int component) const +{ + return shape_hessian_component(i,j,component); +} + template inline @@ -2691,6 +2780,76 @@ FEValuesBase::JxW (const unsigned int i) const } +template +template +void +FEValuesBase::get_function_grads (const InputVector &fe_function, + std::vector > &gradients) const +{ + get_function_gradients(fe_function, gradients); +} + + +template +template +void FEValuesBase::get_function_grads ( + const InputVector& fe_function, + const VectorSlice >& indices, + std::vector > &values) const +{ + get_function_gradients(fe_function, indices, values); +} + + +template +template +void +FEValuesBase:: +get_function_grads (const InputVector &fe_function, + std::vector > > &gradients) const +{ + get_function_gradients(fe_function, gradients); +} + + +template +template +void FEValuesBase::get_function_grads ( + const InputVector& fe_function, + const VectorSlice >& indices, + std::vector > >& values, + bool q_points_fastest) const +{ + get_function_gradients(fe_function, indices, values, q_points_fastest); +} + + + +template +template +inline void +FEValuesBase:: +get_function_2nd_derivatives (const InputVector &fe_function, + std::vector > &hessians) const +{ + get_function_hessians(fe_function, hessians); +} + + + +template +template +void +FEValuesBase:: +get_function_2nd_derivatives (const InputVector &fe_function, + std::vector > > &hessians, + bool quadrature_points_fastest) const +{ + get_function_hessians(fe_function, hessians, quadrature_points_fastest); +} + + + /*------------------------ Inline functions: FEValues ----------------------------*/ diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 6def406a71..91abd88683 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -845,17 +845,17 @@ FiniteElement::compute_2nd ( FEValuesData &data) const { Assert ((fe_internal.update_each | fe_internal.update_once) - & update_second_derivatives, + & update_hessians, ExcInternalError()); // make sure we have as many entries as there are nonzero components -// Assert (data.shape_2nd_derivatives.size() == +// Assert (data.shape_hessians.size() == // std::accumulate (n_nonzero_components_table.begin(), // n_nonzero_components_table.end(), // 0U), // ExcInternalError()); // Number of quadrature points - const unsigned int n_q_points = data.shape_2nd_derivatives[0].size(); + const unsigned int n_q_points = data.shape_hessians[0].size(); // first reinit the fe_values // objects used for the finite @@ -971,7 +971,7 @@ FiniteElement::compute_2nd ( for (unsigned int q=0; q::update_each (const UpdateFlags flags) const { UpdateFlags out = flags; - if (flags & (update_values | update_gradients | update_second_derivatives)) + if (flags & (update_values | update_gradients | update_hessians)) out |= update_q_points ; return out; @@ -255,7 +255,7 @@ FE_DGPNonparametric::get_data ( data->grads.resize (this->dofs_per_cell); } - if (flags & update_second_derivatives) + if (flags & update_hessians) { data->grad_grads.resize (this->dofs_per_cell); } @@ -300,8 +300,8 @@ FE_DGPNonparametric::fill_fe_values ( data.shape_values[k][i] = fe_data.values[k]; if (flags & update_gradients) data.shape_gradients[k][i] = fe_data.grads[k]; - if (flags & update_second_derivatives) - data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k]; + if (flags & update_hessians) + data.shape_hessians[k][i] = fe_data.grad_grads[k]; } } } @@ -341,8 +341,8 @@ FE_DGPNonparametric::fill_fe_face_values ( data.shape_values[k][i] = fe_data.values[k]; if (flags & update_gradients) data.shape_gradients[k][i] = fe_data.grads[k]; - if (flags & update_second_derivatives) - data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k]; + if (flags & update_hessians) + data.shape_hessians[k][i] = fe_data.grad_grads[k]; } } } @@ -383,8 +383,8 @@ FE_DGPNonparametric::fill_fe_subface_values ( data.shape_values[k][i] = fe_data.values[k]; if (flags & update_gradients) data.shape_gradients[k][i] = fe_data.grads[k]; - if (flags & update_second_derivatives) - data.shape_2nd_derivatives[k][i] = fe_data.grad_grads[k]; + if (flags & update_hessians) + data.shape_hessians[k][i] = fe_data.grad_grads[k]; } } } diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 5bdcf6b55d..a87b8d0be8 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -153,7 +153,7 @@ FEValuesBase::TriaCellIterator::message_string = ("You have previously called the FEValues::reinit function with a\n" "cell iterator of type Triangulation::cell_iterator. However,\n" "when you do this, you cannot call some functions in the FEValues\n" - "class, such as the get_function_values/grads/2nd_derivatives\n" + "class, such as the get_function_values/gradients/hessians\n" "functions. If you need these functions, then you need to call\n" "FEValues::reinit with an iterator type that allows to extract\n" "degrees of freedom, such as DoFHandler::cell_iterator."); @@ -298,8 +298,8 @@ FEValuesData::initialize (const unsigned int n_quadrature_points, this->shape_gradients.resize (n_nonzero_shape_components, std::vector > (n_quadrature_points)); - if (flags & update_second_derivatives) - this->shape_2nd_derivatives.resize (n_nonzero_shape_components, + if (flags & update_hessians) + this->shape_hessians.resize (n_nonzero_shape_components, std::vector > (n_quadrature_points)); if (flags & update_q_points) @@ -663,7 +663,7 @@ template template void FEValuesBase:: -get_function_grads (const InputVector &fe_function, +get_function_gradients (const InputVector &fe_function, std::vector > &gradients) const { Assert (this->update_flags & update_gradients, ExcAccessToUninitializedField()); @@ -705,7 +705,7 @@ get_function_grads (const InputVector &fe_function, template template -void FEValuesBase::get_function_grads ( +void FEValuesBase::get_function_gradients ( const InputVector& fe_function, const VectorSlice >& indices, std::vector > &values) const @@ -745,7 +745,7 @@ template template void FEValuesBase:: -get_function_grads (const InputVector &fe_function, +get_function_gradients (const InputVector &fe_function, std::vector > > &gradients) const { Assert (gradients.size() == n_quadrature_points, @@ -800,7 +800,7 @@ get_function_grads (const InputVector &fe_function, template template -void FEValuesBase::get_function_grads ( +void FEValuesBase::get_function_gradients ( const InputVector& fe_function, const VectorSlice >& indices, std::vector > >& values, @@ -893,14 +893,14 @@ template template void FEValuesBase:: -get_function_2nd_derivatives (const InputVector &fe_function, - std::vector > &second_derivatives) const +get_function_hessians (const InputVector &fe_function, + std::vector > &hessians) const { Assert (fe->n_components() == 1, ExcDimensionMismatch(fe->n_components(), 1)); - Assert (second_derivatives.size() == n_quadrature_points, - ExcDimensionMismatch(second_derivatives.size(), n_quadrature_points)); - Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField()); + Assert (hessians.size() == n_quadrature_points, + ExcDimensionMismatch(hessians.size(), n_quadrature_points)); + Assert (this->update_flags & update_hessians, ExcAccessToUninitializedField()); Assert (present_cell.get() != 0, ExcMessage ("FEValues object is not reinit'ed to any cell")); Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(), @@ -913,7 +913,7 @@ get_function_2nd_derivatives (const InputVector &fe_function, present_cell->get_interpolated_dof_values(fe_function, dof_values); // initialize with zero - std::fill_n (second_derivatives.begin(), n_quadrature_points, Tensor<2,dim>()); + std::fill_n (hessians.begin(), n_quadrature_points, Tensor<2,dim>()); // add up contributions of trial // functions. note that here we @@ -924,9 +924,9 @@ get_function_2nd_derivatives (const InputVector &fe_function, for (unsigned int point=0; point tmp = this->shape_2nd_derivative(shape_func,point); + Tensor<2,dim> tmp = this->shape_hessian(shape_func,point); tmp *= dof_values(shape_func); - second_derivatives[point] += tmp; + hessians[point] += tmp; }; } @@ -936,7 +936,7 @@ template template void FEValuesBase:: -get_function_2nd_derivatives (const InputVector &fe_function, +get_function_hessians (const InputVector &fe_function, std::vector > > &second_derivs, bool quadrature_points_fastest) const { @@ -948,7 +948,7 @@ get_function_2nd_derivatives (const InputVector &fe_func Assert (second_derivs[i].size() == n_components, ExcDimensionMismatch(second_derivs[i].size(), n_components)); - Assert (this->update_flags & update_second_derivatives, ExcAccessToUninitializedField()); + Assert (this->update_flags & update_hessians, ExcAccessToUninitializedField()); Assert (present_cell.get() != 0, ExcMessage ("FEValues object is not reinit'ed to any cell")); Assert (fe_function.size() == present_cell->n_dofs_for_dof_handler(), @@ -971,7 +971,7 @@ get_function_2nd_derivatives (const InputVector &fe_func for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { - Tensor<2,dim> tmp(shape_2nd_derivative(shape_func,point)); + Tensor<2,dim> tmp(shape_hessian(shape_func,point)); tmp *= dof_values(shape_func); second_derivs[fe->system_to_component_index(shape_func).first][point] += tmp; @@ -979,7 +979,7 @@ get_function_2nd_derivatives (const InputVector &fe_func else for (unsigned int c=0; c tmp = this->shape_2nd_derivative_component(shape_func,point,c); + Tensor<2,dim> tmp = this->shape_hessian_component(shape_func,point,c); tmp *= dof_values(shape_func); second_derivs[c][point] += tmp; } @@ -988,7 +988,7 @@ get_function_2nd_derivatives (const InputVector &fe_func for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { - Tensor<2,dim> tmp(shape_2nd_derivative(shape_func,point)); + Tensor<2,dim> tmp(shape_hessian(shape_func,point)); tmp *= dof_values(shape_func); second_derivs[point][fe->system_to_component_index(shape_func).first] += tmp; @@ -996,7 +996,7 @@ get_function_2nd_derivatives (const InputVector &fe_func else for (unsigned int c=0; c tmp = this->shape_2nd_derivative_component(shape_func,point,c); + Tensor<2,dim> tmp = this->shape_hessian_component(shape_func,point,c); tmp *= dof_values(shape_func); second_derivs[point][c] += tmp; } @@ -1010,7 +1010,7 @@ FEValuesBase::memory_consumption () const { return (MemoryConsumption::memory_consumption (this->shape_values) + MemoryConsumption::memory_consumption (this->shape_gradients) + - MemoryConsumption::memory_consumption (this->shape_2nd_derivatives) + + MemoryConsumption::memory_consumption (this->shape_hessians) + MemoryConsumption::memory_consumption (this->JxW_values) + MemoryConsumption::memory_consumption (this->quadrature_points) + MemoryConsumption::memory_consumption (this->normal_vectors) + diff --git a/deal.II/deal.II/source/fe/fe_values.instance.h b/deal.II/deal.II/source/fe/fe_values.instance.h index 1cc5a87cf6..94f3f68484 100644 --- a/deal.II/deal.II/source/fe/fe_values.instance.h +++ b/deal.II/deal.II/source/fe/fe_values.instance.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 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 @@ -56,27 +56,27 @@ void FEValuesBase::get_function_values std::vector > &, bool) const; template -void FEValuesBase::get_function_grads +void FEValuesBase::get_function_gradients (const IN&, std::vector > &) const; template -void FEValuesBase::get_function_grads +void FEValuesBase::get_function_gradients (const IN&, const VectorSlice >&, std::vector > &) const; template -void FEValuesBase::get_function_grads +void FEValuesBase::get_function_gradients (const IN&, std::vector > > &) const; template -void FEValuesBase::get_function_grads +void FEValuesBase::get_function_gradients (const IN&, const VectorSlice >&, std::vector > > &, bool) const; template -void FEValuesBase::get_function_2nd_derivatives +void FEValuesBase::get_function_hessians (const IN&, std::vector > &) const; template -void FEValuesBase::get_function_2nd_derivatives +void FEValuesBase::get_function_hessians (const IN&, std::vector > > &, bool) const; DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/numerics/derivative_approximation.cc b/deal.II/deal.II/source/numerics/derivative_approximation.cc index df45dce7f4..bc170f5aa3 100644 --- a/deal.II/deal.II/source/numerics/derivative_approximation.cc +++ b/deal.II/deal.II/source/numerics/derivative_approximation.cc @@ -406,14 +406,14 @@ get_projected_derivative (const FEValues &fe_values, if (fe_values.get_fe().n_components() == 1) { std::vector values (1); - fe_values.get_function_2nd_derivatives (solution, values); + fe_values.get_function_hessians (solution, values); return values[0]; } else { std::vector > values (1, std::vector(fe_values.get_fe().n_components())); - fe_values.get_function_2nd_derivatives (solution, values); + fe_values.get_function_hessians (solution, values); return values[0][component]; }; } -- 2.39.5