From 83caa02912f638e8b268e6328fd8660c3759f283 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Wed, 4 Oct 2017 08:41:32 -0600 Subject: [PATCH] Use only one FEValues object. > > Currently, there is one FEValues object that is reinit'd in the innermost > function, but that function may be called multiple times for the same > cell for different base elements. That is wasteful, so avoid it. > The solution is to reinit it at the outermost place. --- .../deal.II/numerics/vector_tools.templates.h | 34 +++++++------------ 1 file changed, 12 insertions(+), 22 deletions(-) diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 289cde3fe2..a4bf302661 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -109,14 +109,13 @@ namespace VectorTools // Input: // conformity: conformity of the finite element, used to select // appropriate type of transformation - // fe_values_jacobians, cell: used to reinitialize an fe_values object - // if values of jacobians (and inverses of - // jacobians) are needed + // fe_values_jacobians: used for jacobians (and inverses of + // jacobians). the object is supposed to be + // reinit()'d for the current cell // function_values, offset: function_values is manipulated in place // starting at position offset - template + template void transform(const typename FiniteElementData::Conformity conformity, - const T1 &cell, const unsigned int offset, T2 &fe_values_jacobians, T3 &function_values) @@ -129,7 +128,6 @@ namespace VectorTools // For given mapping F_K: \hat K \to K, we have to transform // \hat u = (dF_K)^T u\circ F_K - fe_values_jacobians.reinit(cell); for (unsigned int i = 0; i < function_values.size(); ++i) { const auto &jacobians = @@ -155,7 +153,6 @@ namespace VectorTools // For given mapping F_K: \hat K \to K, we have to transform // \hat w = det(dF_K) (dF_K)^{-1} w\circ F_K - fe_values_jacobians.reinit(cell); for (unsigned int i = 0; i < function_values.size(); ++i) { const auto &jacobians = @@ -209,10 +206,9 @@ namespace VectorTools // [ rest see above] // Output: the offset after we have handled the element at // a given offset - template + template unsigned int apply_transform(const FiniteElement &fe, - const T1 &cell, const unsigned int offset, T2 &fe_values_jacobians, T3 &function_values) @@ -233,7 +229,6 @@ namespace VectorTools // correctly handle nested fe systems. current_offset = apply_transform(base_fe, - cell, current_offset, fe_values_jacobians, function_values); @@ -244,7 +239,6 @@ namespace VectorTools else { transform(fe.conforming_space, - cell, offset, fe_values_jacobians, function_values); @@ -349,19 +343,15 @@ namespace VectorTools hp::MappingCollection mapping_collection(mapping); + // An FEValues object to evaluate (generalized) support point + // locations as well as Jacobians and their inverses. + // the latter are only needed for Hcurl or Hdiv conforming elements, + // but we'll just always include them. hp::FEValues fe_values( mapping_collection, fe, support_quadrature, - update_quadrature_points); - - // An extra FEValues object to compute jacobians. - // Only re-initialized in case of Hcurl or Hdiv conforming elements, - // i.e. if we really need the information. - hp::FEValues fe_values_jacobians( - mapping_collection, - fe, - support_quadrature, + update_quadrature_points | update_jacobians | update_inverse_jacobians); // @@ -415,8 +405,8 @@ namespace VectorTools const unsigned int offset = apply_transform( - fe[fe_index], cell, /* starting_offset = */ 0, - fe_values_jacobians, function_values); + fe[fe_index], /* starting_offset = */ 0, + fe_values, function_values); (void)offset; Assert(offset == n_components, ExcInternalError()); } -- 2.39.5