From 1a5ee42636ca7f673f9480cf450dd4d3914047f7 Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 13 Nov 2000 13:33:56 +0000 Subject: [PATCH] Implement the second FEValues::get_function_2nd_derivatives function, for vector-valued finite elements. git-svn-id: https://svn.dealii.org/trunk@3483 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe_values.h | 27 +++++++++++++ deal.II/deal.II/source/fe/fe_values.cc | 53 ++++++++++++++++++++++++++ 2 files changed, 80 insertions(+) diff --git a/deal.II/deal.II/include/fe/fe_values.h b/deal.II/deal.II/include/fe/fe_values.h index d85ab1d674..e32fdc8ba6 100644 --- a/deal.II/deal.II/include/fe/fe_values.h +++ b/deal.II/deal.II/include/fe/fe_values.h @@ -468,6 +468,33 @@ class FEValuesBase void get_function_2nd_derivatives (const InputVector &fe_function, vector > &second_derivatives) const; + + /** + * Return the tensor of second + * derivatives of the finite + * element function characterized + * by @p{fe_function} restricted to + * @p{cell} at the quadrature points. + * + * The function assumes that the + * @p{second_derivatives} object already has + * the right size. + * + * This function does the same as + * the other one with the same + * name, but applies to + * vector-valued finite elements. + * + * The actual data type of the + * input vector may be either a + * @p{Vector}, + * @p{Vector}, or + * @p{BlockVector}. + */ + template + void get_function_2nd_derivatives (const InputVector &fe_function, + vector > > &second_derivatives) const; + /** * Return the position of the @p{i}th * quadrature point in real space. diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 29ed57a3cb..3e02dbb2f1 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -286,6 +286,48 @@ void FEValuesBase::get_function_2nd_derivatives (const InputVector &fe +template +template +void +FEValuesBase:: +get_function_2nd_derivatives (const InputVector &fe_function, + vector > > &second_derivs) const +{ + Assert (n_quadrature_points == second_derivs.size(), + ExcWrongNoOfComponents()); + Assert (selected_datasetn_components(), + ExcWrongVectorSize(second_derivs[i].size(), fe->n_components())); + Assert (update_flags & update_second_derivatives, ExcAccessToUninitializedField()); + + // get function values of dofs + // on this cell + Vector dof_values (dofs_per_cell); + if (present_cell->active()) + present_cell->get_dof_values (fe_function, dof_values); + else + present_cell->get_interpolated_dof_values(fe_function, dof_values); + + // initialize with zero + for (unsigned i=0;i()); + + // add up contributions of trial + // functions + for (unsigned int point=0; point tmp(shape_2nd_derivatives[shape_func][point]); + tmp *= dof_values(shape_func); + second_derivs[point][fe->system_to_component_index(shape_func).first] + += tmp; + }; +}; + + + template const Point & FEValuesBase::quadrature_point (const unsigned int i) const @@ -929,3 +971,14 @@ void FEValuesBase::get_function_2nd_derivatives (const Vector template void FEValuesBase::get_function_2nd_derivatives (const BlockVector &, vector > &) const; +//----------------------------------------------------------------------------- + +template +void FEValuesBase::get_function_2nd_derivatives (const Vector &, + vector > > &) const; +template +void FEValuesBase::get_function_2nd_derivatives (const Vector &, + vector > > &) const; +template +void FEValuesBase::get_function_2nd_derivatives (const BlockVector &, + vector > > &) const; -- 2.39.5