} // namespace FEPointEvaluation
} // namespace internal
+
+
/**
* This class provides an interface to the evaluation of interpolated solution
* values and gradients on cells on arbitrary reference point positions. These
* a precomputed @p mapping_data object. This function can be used
* to avoid duplicated evaluation of the mapping if multiple
* FEPointEvaluation objects for the different components of the same FESystem
- * are used. You can precompute the mapping data
- * using the function internal::FEPointEvaluation::compute_mapping_data().
+ * are used. You can get the mapping data from an initialized FEPointEvaluation
+ * object by calling the function get_mapping_data().
*
* @param[in] cell An iterator to the current cell
*
const dealii::internal::FEValuesImplementation::
MappingRelatedData<dim, spacedim> &mapping_data);
+ /**
+ * Returns the mapping data that was computed during the last call to
+ * the reinit() function. This can be useful if multiple FEPointEvaluation
+ * objects are used for multiple components of a FESystem. The mapping data
+ * can be retrieved from the first FEPointEvaluation object and subsequently
+ * passed to the other objects, avoiding the need to recompute the mapping
+ * data.
+ */
+ const dealii::internal::FEValuesImplementation::MappingRelatedData<dim,
+ spacedim> &
+ get_mapping_data() const;
+
/**
* This function interpolates the finite element solution, represented by
* `solution_values`, on the cell and `unit_points` passed to reinit().
+template <int n_components, int dim, int spacedim, typename Number>
+const dealii::internal::FEValuesImplementation::MappingRelatedData<dim,
+ spacedim> &
+FEPointEvaluation<n_components, dim, spacedim, Number>::get_mapping_data() const
+{
+ return mapping_data;
+}
+
+
+
template <int n_components, int dim, int spacedim, typename Number>
void
FEPointEvaluation<n_components, dim, spacedim, Number>::evaluate(
evaluator1.evaluate(solution_values,
EvaluationFlags::values | EvaluationFlags::gradients);
- internal::FEValuesImplementation::MappingRelatedData<dim> mapping_data;
- internal::FEPointEvaluation::compute_mapping_data_for_generic_points<dim>(
- mapping,
- cell,
- unit_points,
- update_values | update_gradients,
- mapping_data);
-
- evaluator2.reinit(cell, unit_points, mapping_data);
+ evaluator2.reinit(cell, unit_points, evaluator1.get_mapping_data());
evaluator2.evaluate(solution_values,
EvaluationFlags::values | EvaluationFlags::gradients);
evaluator.evaluate(solution_values,
EvaluationFlags::values | EvaluationFlags::gradients);
- internal::FEValuesImplementation::MappingRelatedData<dim> mapping_data;
- internal::FEPointEvaluation::compute_mapping_data_for_generic_points<dim>(
- mapping,
- cell,
- unit_points,
- update_values | update_gradients,
- mapping_data);
-
- evaluator2.reinit(cell, unit_points, mapping_data);
+ evaluator2.reinit(cell, unit_points, evaluator.get_mapping_data());
evaluator2.evaluate(solution_values,
EvaluationFlags::values | EvaluationFlags::gradients);