From: Jean-Paul Pelteret Date: Thu, 10 Aug 2017 23:03:20 +0000 (-0600) Subject: Added Tensor::get_function_divergences_from_local_dof_values X-Git-Tag: v9.0.0-rc1~1291^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bec925eedf8aa43799d58cbcad718e9834c6c157;p=dealii.git Added Tensor::get_function_divergences_from_local_dof_values --- diff --git a/include/deal.II/fe/fe_values.h b/include/deal.II/fe/fe_values.h index 1e4fd53119..67d195f8ff 100644 --- a/include/deal.II/fe/fe_values.h +++ b/include/deal.II/fe/fe_values.h @@ -1592,6 +1592,14 @@ namespace FEValuesViews void get_function_divergences (const InputVector &fe_function, std::vector::type> &divergences) const; + /** + * @copydoc FEValuesViews::Tensor<2,dim,spacedim>::get_function_values_from_local_dof_values() + */ + template + void get_function_divergences_from_local_dof_values (const InputVector &dof_values, + std::vector::divergence_type> &values) const; + + private: /** * A pointer to the FEValuesBase object we operate on. diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 834f494f21..aa1cddfa30 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -1211,7 +1211,7 @@ namespace FEValuesViews do_function_divergences (const ArrayView &dof_values, const Table<2,dealii::Tensor<1,spacedim> > &shape_gradients, const std::vector::ShapeFunctionData> &shape_function_data, - std::vector >::type> &divergences) + std::vector::template OutputType::divergence_type> &divergences) { const unsigned int dofs_per_cell = dof_values.size(); const unsigned int n_quadrature_points = dofs_per_cell > 0 ? @@ -1219,7 +1219,7 @@ namespace FEValuesViews AssertDimension (divergences.size(), n_quadrature_points); std::fill (divergences.begin(), divergences.end(), - typename ProductType >::type()); + typename Tensor<2,dim,spacedim>::template OutputType::divergence_type()); for (unsigned int shape_function=0; shape_functionfinite_element_output.shape_gradients, shape_function_data, divergences); } + + + + template + template + void + Tensor<2, dim, spacedim>:: + get_function_divergences_from_local_dof_values (const InputVector &dof_values, + std::vector::divergence_type> &divergences) const + { + Assert(fe_values->update_flags & update_gradients, + (typename FEValuesBase::ExcAccessToUninitializedField("update_gradients"))); + Assert(fe_values->present_cell.get() != nullptr, + ExcMessage("FEValues object is not reinit'ed to any cell")); + AssertDimension (dof_values.size(), fe_values->dofs_per_cell); + + internal::do_function_divergences + (make_array_view(dof_values.begin(), dof_values.end()), + fe_values->finite_element_output.shape_gradients, shape_function_data, divergences); + } } diff --git a/source/fe/fe_values.inst.in b/source/fe/fe_values.inst.in index f9a258895b..aa9f5a8e08 100644 --- a/source/fe/fe_values.inst.in +++ b/source/fe/fe_values.inst.in @@ -230,6 +230,11 @@ for (VEC : GENERAL_CONTAINER_TYPES; void FEValuesViews::Tensor<2,deal_II_dimension, deal_II_space_dimension> ::get_function_values_from_local_dof_values> (const VEC&, std::vector::value_type>&) const; + + template + void FEValuesViews::Tensor<2,deal_II_dimension, deal_II_space_dimension> + ::get_function_divergences_from_local_dof_values> + (const VEC&, std::vector::divergence_type>&) const; #endif }