From: Maximilian Bergbauer Date: Thu, 19 May 2022 11:58:06 +0000 (+0200) Subject: Implement class NonMatching::MappingInfo X-Git-Tag: v9.4.0-rc1~78^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F13835%2Fhead;p=dealii.git Implement class NonMatching::MappingInfo --- diff --git a/include/deal.II/matrix_free/fe_point_evaluation.h b/include/deal.II/matrix_free/fe_point_evaluation.h index e1f7e5e474..c722e0e158 100644 --- a/include/deal.II/matrix_free/fe_point_evaluation.h +++ b/include/deal.II/matrix_free/fe_point_evaluation.h @@ -33,6 +33,8 @@ #include #include +#include + DEAL_II_NAMESPACE_OPEN namespace internal @@ -364,48 +366,6 @@ namespace internal template std::vector> get_polynomial_space(const FiniteElement &fe); - - /** - * Compute the mapping related data for the given @p mapping, - * @p cell and @p unit_points that is required by the FEPointEvaluation - * class. - */ - template - void - compute_mapping_data_for_generic_points( - const Mapping & mapping, - const typename Triangulation::cell_iterator &cell, - const ArrayView> & unit_points, - const UpdateFlags update_flags, - internal::FEValuesImplementation::MappingRelatedData - &mapping_data) - { - UpdateFlags update_flags_mapping = update_default; - // translate update flags - if (update_flags & update_jacobians) - update_flags_mapping |= update_jacobians; - if (update_flags & update_gradients || - update_flags & update_inverse_jacobians) - update_flags_mapping |= update_inverse_jacobians; - if (update_flags & update_quadrature_points) - update_flags_mapping |= update_quadrature_points; - - if (const MappingQ *mapping_q = - dynamic_cast *>(&mapping)) - { - mapping_q->fill_mapping_data_for_generic_points(cell, - unit_points, - update_flags_mapping, - mapping_data); - } - else if (const MappingCartesian *mapping_cartesian = - dynamic_cast *>( - &mapping)) - { - mapping_cartesian->fill_mapping_data_for_generic_points( - cell, unit_points, update_flags_mapping, mapping_data); - } - } } // namespace FEPointEvaluation } // namespace internal @@ -476,55 +436,39 @@ public: const unsigned int first_selected_component = 0); /** - * Set up the mapping information for the given cell, e.g., by computing the - * Jacobian of the mapping for the given points if gradients of the functions - * are requested. + * Constructor to make the present class able to re-use the geometry + * data also used by other `FEPointEvaluation` objects. * - * @param[in] cell An iterator to the current cell + * @param mapping_info The MappingInfo class describes the geometry-related + * data for evaluating finite-element solutions. This object enables to + * construct such an object on the outside, possibly re-using it between + * several objects or between several calls to the same cell and unit points. * - * @param[in] unit_points List of points in the reference locations of the - * current cell where the FiniteElement object should be - * evaluated/integrated in the evaluate() and integrate() functions. + * @param fe The FiniteElement object that is used for the evaluation, which + * is typically the same on all cells to be evaluated. + * + * @param first_selected_component For multi-component FiniteElement + * objects, this parameter allows to select a range of `n_components` + * components starting from this parameter. */ - void - reinit(const typename Triangulation::cell_iterator &cell, - const ArrayView> &unit_points); + FEPointEvaluation(NonMatching::MappingInfo &mapping_info, + const FiniteElement & fe, + const unsigned int first_selected_component = 0); /** - * Set up the mapping information for the given cell. This function is - * an alternative to the function with the same name above and uses - * 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 get the mapping data from an initialized - * FEPointEvaluation object by calling the function get_mapping_data(). + * Set up the mapping information for the given cell, e.g., by computing the + * Jacobian of the mapping for the given points if gradients of the functions + * are requested. * * @param[in] cell An iterator to the current cell * * @param[in] unit_points List of points in the reference locations of the * current cell where the FiniteElement object should be * evaluated/integrated in the evaluate() and integrate() functions. - * - * @param[in] mapping_data A mapping data object that is precomputed - * for the given cell and positions. */ void reinit(const typename Triangulation::cell_iterator &cell, - const ArrayView> &unit_points, - const dealii::internal::FEValuesImplementation:: - MappingRelatedData &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 & - get_mapping_data() const; + const ArrayView> &unit_points); /** * This function interpolates the finite element solution, represented by @@ -656,6 +600,17 @@ public: unit_point(const unsigned int point_index) const; private: + /** + * Common setup function for both constructors. Does the setup for both fast + * and slow path. + * + * @param first_selected_component For multi-component FiniteElement + * objects, this parameter allows to select a range of `n_components` + * components starting from this parameter. + */ + void + setup(const unsigned int first_selected_component); + /** * Pointer to the Mapping object passed to the constructor. */ @@ -739,29 +694,34 @@ private: /** * The desired update flags for the evaluation. */ - UpdateFlags update_flags; + const UpdateFlags update_flags; /** - * The update flags specific for the mapping in the fast evaluation path. + * The FEValues object underlying the slow evaluation path. */ - UpdateFlags update_flags_mapping; + std::shared_ptr> fe_values; /** - * The FEValues object underlying the slow evaluation path. + * Pointer to mapping info on the fly computed during reinit. */ - std::shared_ptr> fe_values; + std::unique_ptr> + mapping_info_on_the_fly; /** - * Array to store temporary data computed by the mapping for the fast - * evaluation path. + * Pointer to currently used mapping info (either on the fly or external + * precomputed). */ - internal::FEValuesImplementation::MappingRelatedData - mapping_data; + SmartPointer> mapping_info; /** * The reference points specified at reinit(). */ std::vector> unit_points; + + /** + * Bool indicating if fast path is chosen. + */ + bool fast_path; }; // ----------------------- template and inline function ---------------------- @@ -776,73 +736,97 @@ FEPointEvaluation::FEPointEvaluation( : mapping(&mapping) , fe(&fe) , update_flags(update_flags) - , update_flags_mapping(update_default) + , mapping_info_on_the_fly( + std::make_unique>(mapping, + update_flags)) + , mapping_info(mapping_info_on_the_fly.get()) +{ + setup(first_selected_component); +} + + + +template +FEPointEvaluation::FEPointEvaluation( + NonMatching::MappingInfo &mapping_info, + const FiniteElement & fe, + const unsigned int first_selected_component) + : mapping(&mapping_info.get_mapping()) + , fe(&fe) + , update_flags(mapping_info.get_update_flags()) + , mapping_info(&mapping_info) +{ + setup(first_selected_component); +} + + + +template +void +FEPointEvaluation::setup( + const unsigned int first_selected_component) { AssertIndexRange(first_selected_component + n_components, - fe.n_components() + 1); + fe->n_components() + 1); bool same_base_element = true; unsigned int base_element_number = 0; component_in_base_element = 0; unsigned int component = 0; - for (; base_element_number < fe.n_base_elements(); ++base_element_number) - if (component + fe.element_multiplicity(base_element_number) > + for (; base_element_number < fe->n_base_elements(); ++base_element_number) + if (component + fe->element_multiplicity(base_element_number) > first_selected_component) { if (first_selected_component + n_components > - component + fe.element_multiplicity(base_element_number)) + component + fe->element_multiplicity(base_element_number)) same_base_element = false; component_in_base_element = first_selected_component - component; break; } else - component += fe.element_multiplicity(base_element_number); + component += fe->element_multiplicity(base_element_number); - if (internal::FEPointEvaluation::is_fast_path_supported(mapping) && + if (internal::FEPointEvaluation::is_fast_path_supported(*mapping) && internal::FEPointEvaluation::is_fast_path_supported( - fe, base_element_number) && + *fe, base_element_number) && same_base_element) { internal::MatrixFreeFunctions::ShapeInfo shape_info; - shape_info.reinit(QMidpoint<1>(), fe, base_element_number); + shape_info.reinit(QMidpoint<1>(), *fe, base_element_number); renumber = shape_info.lexicographic_numbering; dofs_per_component = shape_info.dofs_per_component_on_cell; poly = internal::FEPointEvaluation::get_polynomial_space( - fe.base_element(base_element_number)); + fe->base_element(base_element_number)); polynomials_are_hat_functions = (poly.size() == 2 && poly[0].value(0.) == 1. && poly[0].value(1.) == 0. && poly[1].value(0.) == 0. && poly[1].value(1.) == 1.); + + fast_path = true; } else { - nonzero_shape_function_component.resize(fe.n_dofs_per_cell()); + nonzero_shape_function_component.resize(fe->n_dofs_per_cell()); for (unsigned int d = 0; d < n_components; ++d) { const unsigned int component = first_selected_component + d; - for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i) + for (unsigned int i = 0; i < fe->n_dofs_per_cell(); ++i) { - const bool is_primitive = fe.is_primitive() || fe.is_primitive(i); + const bool is_primitive = + fe->is_primitive() || fe->is_primitive(i); if (is_primitive) nonzero_shape_function_component[i][d] = - (component == fe.system_to_component_index(i).first); + (component == fe->system_to_component_index(i).first); else nonzero_shape_function_component[i][d] = - (fe.get_nonzero_components(i)[component] == true); + (fe->get_nonzero_components(i)[component] == true); } } - } - // translate update flags - if (update_flags & update_jacobians) - update_flags_mapping |= update_jacobians; - if (update_flags & update_gradients || - update_flags & update_inverse_jacobians) - update_flags_mapping |= update_inverse_jacobians; - if (update_flags & update_quadrature_points) - update_flags_mapping |= update_quadrature_points; + fast_path = false; + } } @@ -853,73 +837,25 @@ FEPointEvaluation::reinit( const typename Triangulation::cell_iterator &cell, const ArrayView> & unit_points) { - // If using the fast path, we need to precompute the mapping data. - if (!poly.empty()) - internal::FEPointEvaluation::compute_mapping_data_for_generic_points( - *mapping, cell, unit_points, update_flags_mapping, mapping_data); - - // then call the other version of this function with the precomputed data - reinit(cell, unit_points, mapping_data); -} - + // reinit is only allowed for mapping computation on the fly + AssertThrow(mapping_info_on_the_fly.get() != nullptr, ExcNotImplemented()); + mapping_info->reinit(cell, unit_points); -template -void -FEPointEvaluation::reinit( - const typename Triangulation::cell_iterator &cell, - const ArrayView> & unit_points, - const dealii::internal::FEValuesImplementation::MappingRelatedData - &precomputed_mapping_data) -{ - this->unit_points.resize(unit_points.size()); - std::copy(unit_points.begin(), unit_points.end(), this->unit_points.begin()); - - if (!poly.empty()) - { - // Check the mapping data for consistency. - if ((update_flags_mapping & update_jacobians) != 0u) - Assert(precomputed_mapping_data.jacobians.size() == unit_points.size(), - ExcDimensionMismatch(precomputed_mapping_data.jacobians.size(), - unit_points.size())); - if ((update_flags_mapping & update_inverse_jacobians) != 0u) - Assert(precomputed_mapping_data.inverse_jacobians.size() == - unit_points.size(), - ExcDimensionMismatch( - precomputed_mapping_data.inverse_jacobians.size(), - unit_points.size())); - - if ((update_flags_mapping & update_quadrature_points) != 0u) - Assert(precomputed_mapping_data.inverse_jacobians.size() == - unit_points.size(), - ExcDimensionMismatch( - precomputed_mapping_data.quadrature_points.size(), - unit_points.size())); - - mapping_data = precomputed_mapping_data; - } - else + if (!fast_path) { fe_values = std::make_shared>( *mapping, *fe, Quadrature( std::vector>(unit_points.begin(), unit_points.end())), - update_flags | update_flags_mapping); + update_flags); fe_values->reinit(cell); - mapping_data.initialize(unit_points.size(), update_flags_mapping); - if (update_flags_mapping & update_jacobians) - for (unsigned int q = 0; q < unit_points.size(); ++q) - mapping_data.jacobians[q] = fe_values->jacobian(q); - if (update_flags_mapping & update_inverse_jacobians) - for (unsigned int q = 0; q < unit_points.size(); ++q) - mapping_data.inverse_jacobians[q] = fe_values->inverse_jacobian(q); - if (update_flags_mapping & update_quadrature_points) - for (unsigned int q = 0; q < unit_points.size(); ++q) - mapping_data.quadrature_points[q] = fe_values->quadrature_point(q); } + this->unit_points = + std::vector>(unit_points.begin(), unit_points.end()); + if (update_flags & update_values) values.resize(unit_points.size(), numbers::signaling_nan()); if (update_flags & update_gradients) @@ -929,29 +865,31 @@ FEPointEvaluation::reinit( -template -const dealii::internal::FEValuesImplementation::MappingRelatedData & -FEPointEvaluation::get_mapping_data() const -{ - return mapping_data; -} - - - template void FEPointEvaluation::evaluate( const ArrayView & solution_values, const EvaluationFlags::EvaluationFlags &evaluation_flag) { + const bool precomputed_mapping = mapping_info_on_the_fly.get() == nullptr; + if (precomputed_mapping) + { + unit_points = mapping_info->get_unit_points(); + + if (update_flags & update_values) + values.resize(unit_points.size(), numbers::signaling_nan()); + if (update_flags & update_gradients) + gradients.resize(unit_points.size(), + numbers::signaling_nan()); + } + if (unit_points.empty()) return; AssertDimension(solution_values.size(), fe->dofs_per_cell); if (((evaluation_flag & EvaluationFlags::values) || (evaluation_flag & EvaluationFlags::gradients)) && - !poly.empty()) + fast_path) { // fast path with tensor product evaluation if (solution_renumbered.size() != dofs_per_component) @@ -1002,20 +940,21 @@ FEPointEvaluation::evaluate( ExcNotInitialized()); for (unsigned int j = 0; j < n_lanes && i + j < n_points; ++j) { - Assert(update_flags_mapping & update_inverse_jacobians, - ExcNotInitialized()); internal::FEPointEvaluation::EvaluatorTypeTraits< dim, n_components, Number>::set_gradient(val_and_grad.second, j, unit_gradients[i + j]); - gradients[i + j] = static_cast< - typename internal::FEPointEvaluation:: - EvaluatorTypeTraits:: - gradient_type>(apply_transformation( - mapping_data.inverse_jacobians[i + j].transpose(), - unit_gradients[i + j])); + gradients[i + j] = + static_cast::gradient_type>( + apply_transformation(mapping_info->get_mapping_data() + .inverse_jacobians[i + j] + .transpose(), + unit_gradients[i + j])); } } } @@ -1084,6 +1023,18 @@ FEPointEvaluation::integrate( const ArrayView & solution_values, const EvaluationFlags::EvaluationFlags &integration_flags) { + const bool precomputed_mapping = mapping_info_on_the_fly.get() == nullptr; + if (precomputed_mapping) + { + unit_points = mapping_info->get_unit_points(); + + if (update_flags & update_values) + values.resize(unit_points.size(), numbers::signaling_nan()); + if (update_flags & update_gradients) + gradients.resize(unit_points.size(), + numbers::signaling_nan()); + } + if (unit_points.size() == 0) // no evaluation points provided { std::fill(solution_values.begin(), solution_values.end(), 0.0); @@ -1093,7 +1044,7 @@ FEPointEvaluation::integrate( AssertDimension(solution_values.size(), fe->dofs_per_cell); if (((integration_flags & EvaluationFlags::values) || (integration_flags & EvaluationFlags::gradients)) && - !poly.empty()) + fast_path) { // fast path with tensor product integration @@ -1139,14 +1090,12 @@ FEPointEvaluation::integrate( if (integration_flags & EvaluationFlags::gradients) for (unsigned int j = 0; j < n_lanes && i + j < n_points; ++j) { - Assert(update_flags_mapping & update_inverse_jacobians, - ExcNotInitialized()); gradients[i + j] = static_cast:: - gradient_type>( - apply_transformation(mapping_data.inverse_jacobians[i + j], - gradients[i + j])); + gradient_type>(apply_transformation( + mapping_info->get_mapping_data().inverse_jacobians[i + j], + gradients[i + j])); internal::FEPointEvaluation:: EvaluatorTypeTraits::get_gradient( gradient, j, gradients[i + j]); @@ -1269,7 +1218,7 @@ inline const typename FEPointEvaluation:: FEPointEvaluation::get_unit_gradient( const unsigned int point_index) const { - Assert(!poly.empty(), + Assert(fast_path, ExcMessage("Unit gradients are currently only implemented for tensor " "product finite elements combined with MappingQ " "mappings")); @@ -1308,9 +1257,9 @@ inline DerivativeForm<1, dim, spacedim> FEPointEvaluation::jacobian( const unsigned int point_index) const { - Assert(update_flags_mapping & update_jacobians, ExcNotInitialized()); - AssertIndexRange(point_index, mapping_data.jacobians.size()); - return mapping_data.jacobians[point_index]; + AssertIndexRange(point_index, + mapping_info->get_mapping_data().jacobians.size()); + return mapping_info->get_mapping_data().jacobians[point_index]; } @@ -1320,11 +1269,9 @@ inline DerivativeForm<1, spacedim, dim> FEPointEvaluation::inverse_jacobian( const unsigned int point_index) const { - Assert(update_flags_mapping & update_inverse_jacobians || - update_flags_mapping & update_gradients, - ExcNotInitialized()); - AssertIndexRange(point_index, mapping_data.inverse_jacobians.size()); - return mapping_data.inverse_jacobians[point_index]; + AssertIndexRange(point_index, + mapping_info->get_mapping_data().inverse_jacobians.size()); + return mapping_info->get_mapping_data().inverse_jacobians[point_index]; } @@ -1334,9 +1281,9 @@ inline Point FEPointEvaluation::real_point( const unsigned int point_index) const { - Assert(update_flags_mapping & update_quadrature_points, ExcNotInitialized()); - AssertIndexRange(point_index, mapping_data.quadrature_points.size()); - return mapping_data.quadrature_points[point_index]; + AssertIndexRange(point_index, + mapping_info->get_mapping_data().quadrature_points.size()); + return mapping_info->get_mapping_data().quadrature_points[point_index]; } diff --git a/include/deal.II/non_matching/mapping_info.h b/include/deal.II/non_matching/mapping_info.h new file mode 100644 index 0000000000..1a70680da9 --- /dev/null +++ b/include/deal.II/non_matching/mapping_info.h @@ -0,0 +1,252 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2011 - 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +#ifndef dealii_non_matching_mapping_info_h +#define dealii_non_matching_mapping_info_h + + +#include + +#include +#include +#include + +#include +#include +#include +#include + +#include + + +DEAL_II_NAMESPACE_OPEN + +namespace NonMatching +{ + /** + * This class provides the mapping information computation and mapping data + * storage to be used together with FEPointEvaluation. + */ + template + class MappingInfo : public Subscriptor + { + public: + /** + * Constructor. + * + * @param mapping The Mapping class describing the geometry of a cell. + * + * @param update_flags Specify the quantities to be computed by the mapping + * during the call of reinit(). These update flags are also handed to a + * FEEvaluation object if you construct it with this MappingInfo object. + */ + MappingInfo(const Mapping &mapping, const UpdateFlags update_flags); + + /** + * Reinitialize the mapping information for the incoming cell and unit + * points. + */ + void + reinit(const typename Triangulation::cell_iterator &cell, + const ArrayView> &unit_points); + + /** + * Getter function for current unit points. + */ + const std::vector> & + get_unit_points() const; + + /** + * Getter function for computed mapping data. This function accesses + * internal data and is therefore not a stable interface. + */ + const dealii::internal::FEValuesImplementation::MappingRelatedData + & + get_mapping_data() const; + + /** + * Getter function for underlying mapping. + */ + const Mapping & + get_mapping() const; + + /** + * Getter function for the update flags. + */ + UpdateFlags + get_update_flags() const; + + private: + /** + * Compute the mapping related data for the given @p mapping, + * @p cell and @p unit_points that is required by the FEPointEvaluation + * class. + */ + void + compute_mapping_data_for_generic_points( + const typename Triangulation::cell_iterator &cell, + const ArrayView> & unit_points); + + /** + * The reference points specified at reinit(). + */ + std::vector> unit_points; + + /** + * A pointer to the underlying mapping. + */ + const SmartPointer> mapping; + + /** + * The desired update flags for the evaluation. + */ + const UpdateFlags update_flags; + + /** + * The update flags for the desired mapping information. + */ + UpdateFlags update_flags_mapping; + + /** + * The internal data container for mapping information. The implementation + * is subject to future changes. + */ + dealii::internal::FEValuesImplementation::MappingRelatedData + mapping_data; + }; + + // ----------------------- template functions ---------------------- + + + template + MappingInfo::MappingInfo(const Mapping &mapping, + const UpdateFlags update_flags) + : mapping(&mapping) + , update_flags(update_flags) + { + update_flags_mapping = update_default; + // translate update flags + if (update_flags & update_jacobians) + update_flags_mapping |= update_jacobians; + if (update_flags & update_gradients || + update_flags & update_inverse_jacobians) + update_flags_mapping |= update_inverse_jacobians; + + // always save quadrature points for now + update_flags_mapping |= update_quadrature_points; + } + + + + template + void + MappingInfo::reinit( + const typename Triangulation::cell_iterator &cell, + const ArrayView> & unit_points) + { + this->unit_points = + std::vector>(unit_points.begin(), unit_points.end()); + compute_mapping_data_for_generic_points(cell, unit_points); + } + + + + template + const std::vector> & + MappingInfo::get_unit_points() const + { + return unit_points; + } + + + + template + const dealii::internal::FEValuesImplementation::MappingRelatedData & + MappingInfo::get_mapping_data() const + { + return mapping_data; + } + + + + template + const Mapping & + MappingInfo::get_mapping() const + { + return *mapping; + } + + + + template + UpdateFlags + MappingInfo::get_update_flags() const + { + return update_flags; + } + + + + template + void + MappingInfo::compute_mapping_data_for_generic_points( + const typename Triangulation::cell_iterator &cell, + const ArrayView> & unit_points) + { + if (const MappingQ *mapping_q = + dynamic_cast *>(&(*mapping))) + { + mapping_q->fill_mapping_data_for_generic_points(cell, + unit_points, + update_flags_mapping, + mapping_data); + } + else if (const MappingCartesian *mapping_cartesian = + dynamic_cast *>( + &(*mapping))) + { + mapping_cartesian->fill_mapping_data_for_generic_points( + cell, unit_points, update_flags_mapping, mapping_data); + } + else + { + FE_DGQ dummy_fe(1); + dealii::FEValues fe_values( + *mapping, + dummy_fe, + Quadrature( + std::vector>(unit_points.begin(), unit_points.end())), + update_flags_mapping); + fe_values.reinit(cell); + mapping_data.initialize(unit_points.size(), update_flags_mapping); + if (update_flags_mapping & update_jacobians) + for (unsigned int q = 0; q < unit_points.size(); ++q) + mapping_data.jacobians[q] = fe_values.jacobian(q); + if (update_flags_mapping & update_inverse_jacobians) + for (unsigned int q = 0; q < unit_points.size(); ++q) + mapping_data.inverse_jacobians[q] = fe_values.inverse_jacobian(q); + if (update_flags_mapping & update_quadrature_points) + for (unsigned int q = 0; q < unit_points.size(); ++q) + mapping_data.quadrature_points[q] = fe_values.quadrature_point(q); + } + } +} // namespace NonMatching + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/tests/matrix_free/point_evaluation_15.cc b/tests/matrix_free/point_evaluation_15.cc index 18c4476e26..8dc2df4be1 100644 --- a/tests/matrix_free/point_evaluation_15.cc +++ b/tests/matrix_free/point_evaluation_15.cc @@ -84,16 +84,15 @@ test() dof_handler.distribute_dofs(fe); Vector vector(dof_handler.n_dofs()); - // Use the first evaluator to compute the values of the function - FEPointEvaluation evaluator1(mapping, - fe, - update_values | update_gradients); + // mapping info object for precomputed mapping + NonMatching::MappingInfo mapping_info(mapping, + update_values | + update_gradients); + // First evaluator to use precomputed mapping + FEPointEvaluation evaluator1(mapping_info, fe); - // The second evaluator will use the precomputed mapping information - // from evaluator1 to compute the values of the function - FEPointEvaluation evaluator2(mapping, - fe, - update_values | update_gradients); + // Second evaluator to use precomputed mapping + FEPointEvaluation evaluator2(mapping_info, fe); VectorTools::interpolate(mapping, dof_handler, MyFunction(), vector); @@ -108,11 +107,10 @@ test() solution_values.begin(), solution_values.end()); - evaluator1.reinit(cell, unit_points); + mapping_info.reinit(cell, unit_points); + evaluator1.evaluate(solution_values, EvaluationFlags::values | EvaluationFlags::gradients); - - evaluator2.reinit(cell, unit_points, evaluator1.get_mapping_data()); evaluator2.evaluate(solution_values, EvaluationFlags::values | EvaluationFlags::gradients); diff --git a/tests/matrix_free/point_evaluation_16.cc b/tests/matrix_free/point_evaluation_16.cc index 9f9fbaaa2e..e2c12e0766 100644 --- a/tests/matrix_free/point_evaluation_16.cc +++ b/tests/matrix_free/point_evaluation_16.cc @@ -88,12 +88,12 @@ test(const unsigned int degree) dof_handler.distribute_dofs(fe); Vector vector(dof_handler.n_dofs()); - FEPointEvaluation evaluator(mapping, - fe, - update_values | update_gradients); - FEPointEvaluation evaluator2(mapping, - fe, - update_values | update_gradients); + NonMatching::MappingInfo mapping_info(mapping, + update_values | + update_gradients); + + FEPointEvaluation evaluator(mapping_info, fe); + FEPointEvaluation evaluator2(mapping_info, fe); VectorTools::interpolate(mapping, dof_handler, MyFunction(), vector); @@ -108,11 +108,11 @@ test(const unsigned int degree) solution_values.begin(), solution_values.end()); - evaluator.reinit(cell, unit_points); + mapping_info.reinit(cell, unit_points); + evaluator.evaluate(solution_values, EvaluationFlags::values | EvaluationFlags::gradients); - evaluator2.reinit(cell, unit_points, evaluator.get_mapping_data()); evaluator2.evaluate(solution_values, EvaluationFlags::values | EvaluationFlags::gradients); diff --git a/tests/matrix_free/point_evaluation_17.cc b/tests/matrix_free/point_evaluation_17.cc new file mode 100644 index 0000000000..1a2c5ad76c --- /dev/null +++ b/tests/matrix_free/point_evaluation_17.cc @@ -0,0 +1,163 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Check precomputed mapping information against mapping +// information on the fly for a vector valued FE and MappingCartesian. +// This is a modified version of point_evaluation_13.cc. + +#include + +#include + +#include +#include +#include + +#include +#include + +#include + +#include + +#include + +#include + +#include "../tests.h" + + + +template +class MyFunction : public Function +{ +public: + MyFunction() + : Function(dim) + {} + + double + value(const Point &p, const unsigned int component) const override + { + return component + p[component]; + } +}; + + + +template +void +test() +{ + using namespace dealii; + Triangulation tria; + + GridGenerator::subdivided_hyper_cube(tria, 2, 0, 1); + + MappingCartesian mapping; + deallog << "Cartesian linear mapping" << std::endl; + + std::vector> unit_points; + for (unsigned int i = 0; i < 13; ++i) + { + Point p; + for (unsigned int d = 0; d < dim; ++d) + p[d] = static_cast(i) / 17. + 0.015625 * d; + unit_points.push_back(p); + } + + FESystem fe(FE_Q(2), dim); + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + Vector vector(dof_handler.n_dofs()); + + // mapping info object for precomputed mapping + NonMatching::MappingInfo mapping_info(mapping, + update_values | + update_gradients); + // The first evaluator will use the precomputed mapping + FEPointEvaluation evaluator1(mapping_info, fe); + + // The second evaluator will compute the mapping on the fly + FEPointEvaluation evaluator2(mapping, + fe, + update_values | update_gradients); + + VectorTools::interpolate(mapping, dof_handler, MyFunction(), vector); + + std::vector solution_values(fe.dofs_per_cell); + std::vector solution_values2(fe.dofs_per_cell); + + FEValuesExtractors::Vector extractor(0); + + for (const auto &cell : dof_handler.active_cell_iterators()) + { + cell->get_dof_values(vector, + solution_values.begin(), + solution_values.end()); + + mapping_info.reinit(cell, unit_points); + evaluator1.evaluate(solution_values, + EvaluationFlags::values | EvaluationFlags::gradients); + + evaluator2.reinit(cell, unit_points); + evaluator2.evaluate(solution_values, + EvaluationFlags::values | EvaluationFlags::gradients); + + deallog << "Cell with center " << cell->center(true) << std::endl; + for (unsigned int i = 0; i < unit_points.size(); ++i) + deallog + << mapping.transform_unit_to_real_cell(cell, unit_points[i]) << ": " + << evaluator1.get_value(i) << " error value " + << (evaluator1.get_value(i) - evaluator2.get_value(i)).norm() + << " error grad " + << (evaluator1.get_gradient(i) - evaluator2.get_gradient(i)).norm() + << std::endl; + deallog << std::endl; + + for (unsigned int i = 0; i < unit_points.size(); ++i) + { + evaluator1.submit_value(evaluator1.get_value(i), i); + evaluator1.submit_gradient(evaluator1.get_gradient(i), i); + + evaluator2.submit_value(evaluator2.get_value(i), i); + evaluator2.submit_gradient(evaluator2.get_gradient(i), i); + } + + evaluator1.integrate(solution_values, + EvaluationFlags::values | + EvaluationFlags::gradients); + evaluator2.integrate(solution_values2, + EvaluationFlags::values | + EvaluationFlags::gradients); + + for (unsigned int i = 0; i < solution_values.size(); ++i) + deallog << solution_values[i] - solution_values2[i] << ' '; + deallog << std::endl; + } +} + + + +int +main() +{ + initlog(); + deallog << std::setprecision(10); + + test<2>(); + test<3>(); +} diff --git a/tests/matrix_free/point_evaluation_17.output b/tests/matrix_free/point_evaluation_17.output new file mode 100644 index 0000000000..0a2facb0c2 --- /dev/null +++ b/tests/matrix_free/point_evaluation_17.output @@ -0,0 +1,195 @@ + +DEAL::Cartesian linear mapping +DEAL::Cell with center 0.2500000000 0.2500000000 +DEAL::0.000000000 0.007812500000: 0.000000000 1.007812500 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.03722426471: 0.02941176471 1.037224265 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.06663602941: 0.05882352941 1.066636029 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.09604779412: 0.08823529412 1.096047794 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.1254595588: 0.1176470588 1.125459559 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.1548713235: 0.1470588235 1.154871324 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.1842830882: 0.1764705882 1.184283088 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.2136948529: 0.2058823529 1.213694853 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.2431066176: 0.2352941176 1.243106618 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.2725183824: 0.2647058824 1.272518382 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.3019301471: 0.2941176471 1.301930147 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.3313419118: 0.3235294118 1.331341912 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.3607536765: 0.3529411765 1.360753676 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.2500000000 +DEAL::0.5000000000 0.007812500000: 0.5000000000 1.007812500 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.03722426471: 0.5294117647 1.037224265 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.06663602941: 0.5588235294 1.066636029 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.09604779412: 0.5882352941 1.096047794 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.1254595588: 0.6176470588 1.125459559 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.1548713235: 0.6470588235 1.154871324 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.1842830882: 0.6764705882 1.184283088 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.2136948529: 0.7058823529 1.213694853 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.2431066176: 0.7352941176 1.243106618 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.2725183824: 0.7647058824 1.272518382 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.3019301471: 0.7941176471 1.301930147 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.3313419118: 0.8235294118 1.331341912 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.3607536765: 0.8529411765 1.360753676 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.2500000000 0.7500000000 +DEAL::0.000000000 0.5078125000: 0.000000000 1.507812500 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.5372242647: 0.02941176471 1.537224265 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.5666360294: 0.05882352941 1.566636029 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.5960477941: 0.08823529412 1.596047794 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.6254595588: 0.1176470588 1.625459559 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.6548713235: 0.1470588235 1.654871324 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.6842830882: 0.1764705882 1.684283088 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.7136948529: 0.2058823529 1.713694853 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.7431066176: 0.2352941176 1.743106618 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.7725183824: 0.2647058824 1.772518382 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.8019301471: 0.2941176471 1.801930147 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.8313419118: 0.3235294118 1.831341912 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.8607536765: 0.3529411765 1.860753676 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.7500000000 +DEAL::0.5000000000 0.5078125000: 0.5000000000 1.507812500 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.5372242647: 0.5294117647 1.537224265 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.5666360294: 0.5588235294 1.566636029 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.5960477941: 0.5882352941 1.596047794 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.6254595588: 0.6176470588 1.625459559 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.6548713235: 0.6470588235 1.654871324 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.6842830882: 0.6764705882 1.684283088 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.7136948529: 0.7058823529 1.713694853 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.7431066176: 0.7352941176 1.743106618 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.7725183824: 0.7647058824 1.772518382 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.8019301471: 0.7941176471 1.801930147 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.8313419118: 0.8235294118 1.831341912 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.8607536765: 0.8529411765 1.860753676 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cartesian linear mapping +DEAL::Cell with center 0.2500000000 0.2500000000 0.2500000000 +DEAL::0.000000000 0.007812500000 0.01562500000: 0.000000000 1.007812500 2.015625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.03722426471 0.04503676471: 0.02941176471 1.037224265 2.045036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.06663602941 0.07444852941: 0.05882352941 1.066636029 2.074448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.09604779412 0.1038602941: 0.08823529412 1.096047794 2.103860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.1254595588 0.1332720588: 0.1176470588 1.125459559 2.133272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.1548713235 0.1626838235: 0.1470588235 1.154871324 2.162683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.1842830882 0.1920955882: 0.1764705882 1.184283088 2.192095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.2136948529 0.2215073529: 0.2058823529 1.213694853 2.221507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.2431066176 0.2509191176: 0.2352941176 1.243106618 2.250919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.2725183824 0.2803308824: 0.2647058824 1.272518382 2.280330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.3019301471 0.3097426471: 0.2941176471 1.301930147 2.309742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.3313419118 0.3391544118: 0.3235294118 1.331341912 2.339154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.3607536765 0.3685661765: 0.3529411765 1.360753676 2.368566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.2500000000 0.2500000000 +DEAL::0.5000000000 0.007812500000 0.01562500000: 0.5000000000 1.007812500 2.015625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.03722426471 0.04503676471: 0.5294117647 1.037224265 2.045036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.06663602941 0.07444852941: 0.5588235294 1.066636029 2.074448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.09604779412 0.1038602941: 0.5882352941 1.096047794 2.103860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.1254595588 0.1332720588: 0.6176470588 1.125459559 2.133272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.1548713235 0.1626838235: 0.6470588235 1.154871324 2.162683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.1842830882 0.1920955882: 0.6764705882 1.184283088 2.192095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.2136948529 0.2215073529: 0.7058823529 1.213694853 2.221507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.2431066176 0.2509191176: 0.7352941176 1.243106618 2.250919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.2725183824 0.2803308824: 0.7647058824 1.272518382 2.280330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.3019301471 0.3097426471: 0.7941176471 1.301930147 2.309742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.3313419118 0.3391544118: 0.8235294118 1.331341912 2.339154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.3607536765 0.3685661765: 0.8529411765 1.360753676 2.368566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.2500000000 0.7500000000 0.2500000000 +DEAL::0.000000000 0.5078125000 0.01562500000: 0.000000000 1.507812500 2.015625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.5372242647 0.04503676471: 0.02941176471 1.537224265 2.045036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.5666360294 0.07444852941: 0.05882352941 1.566636029 2.074448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.5960477941 0.1038602941: 0.08823529412 1.596047794 2.103860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.6254595588 0.1332720588: 0.1176470588 1.625459559 2.133272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.6548713235 0.1626838235: 0.1470588235 1.654871324 2.162683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.6842830882 0.1920955882: 0.1764705882 1.684283088 2.192095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.7136948529 0.2215073529: 0.2058823529 1.713694853 2.221507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.7431066176 0.2509191176: 0.2352941176 1.743106618 2.250919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.7725183824 0.2803308824: 0.2647058824 1.772518382 2.280330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.8019301471 0.3097426471: 0.2941176471 1.801930147 2.309742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.8313419118 0.3391544118: 0.3235294118 1.831341912 2.339154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.8607536765 0.3685661765: 0.3529411765 1.860753676 2.368566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.7500000000 0.2500000000 +DEAL::0.5000000000 0.5078125000 0.01562500000: 0.5000000000 1.507812500 2.015625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.5372242647 0.04503676471: 0.5294117647 1.537224265 2.045036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.5666360294 0.07444852941: 0.5588235294 1.566636029 2.074448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.5960477941 0.1038602941: 0.5882352941 1.596047794 2.103860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.6254595588 0.1332720588: 0.6176470588 1.625459559 2.133272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.6548713235 0.1626838235: 0.6470588235 1.654871324 2.162683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.6842830882 0.1920955882: 0.6764705882 1.684283088 2.192095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.7136948529 0.2215073529: 0.7058823529 1.713694853 2.221507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.7431066176 0.2509191176: 0.7352941176 1.743106618 2.250919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.7725183824 0.2803308824: 0.7647058824 1.772518382 2.280330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.8019301471 0.3097426471: 0.7941176471 1.801930147 2.309742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.8313419118 0.3391544118: 0.8235294118 1.831341912 2.339154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.8607536765 0.3685661765: 0.8529411765 1.860753676 2.368566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.2500000000 0.2500000000 0.7500000000 +DEAL::0.000000000 0.007812500000 0.5156250000: 0.000000000 1.007812500 2.515625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.03722426471 0.5450367647: 0.02941176471 1.037224265 2.545036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.06663602941 0.5744485294: 0.05882352941 1.066636029 2.574448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.09604779412 0.6038602941: 0.08823529412 1.096047794 2.603860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.1254595588 0.6332720588: 0.1176470588 1.125459559 2.633272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.1548713235 0.6626838235: 0.1470588235 1.154871324 2.662683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.1842830882 0.6920955882: 0.1764705882 1.184283088 2.692095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.2136948529 0.7215073529: 0.2058823529 1.213694853 2.721507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.2431066176 0.7509191176: 0.2352941176 1.243106618 2.750919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.2725183824 0.7803308824: 0.2647058824 1.272518382 2.780330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.3019301471 0.8097426471: 0.2941176471 1.301930147 2.809742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.3313419118 0.8391544118: 0.3235294118 1.331341912 2.839154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.3607536765 0.8685661765: 0.3529411765 1.360753676 2.868566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.2500000000 0.7500000000 +DEAL::0.5000000000 0.007812500000 0.5156250000: 0.5000000000 1.007812500 2.515625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.03722426471 0.5450367647: 0.5294117647 1.037224265 2.545036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.06663602941 0.5744485294: 0.5588235294 1.066636029 2.574448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.09604779412 0.6038602941: 0.5882352941 1.096047794 2.603860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.1254595588 0.6332720588: 0.6176470588 1.125459559 2.633272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.1548713235 0.6626838235: 0.6470588235 1.154871324 2.662683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.1842830882 0.6920955882: 0.6764705882 1.184283088 2.692095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.2136948529 0.7215073529: 0.7058823529 1.213694853 2.721507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.2431066176 0.7509191176: 0.7352941176 1.243106618 2.750919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.2725183824 0.7803308824: 0.7647058824 1.272518382 2.780330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.3019301471 0.8097426471: 0.7941176471 1.301930147 2.809742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.3313419118 0.8391544118: 0.8235294118 1.331341912 2.839154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.3607536765 0.8685661765: 0.8529411765 1.360753676 2.868566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.2500000000 0.7500000000 0.7500000000 +DEAL::0.000000000 0.5078125000 0.5156250000: 0.000000000 1.507812500 2.515625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.02941176471 0.5372242647 0.5450367647: 0.02941176471 1.537224265 2.545036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.05882352941 0.5666360294 0.5744485294: 0.05882352941 1.566636029 2.574448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.08823529412 0.5960477941 0.6038602941: 0.08823529412 1.596047794 2.603860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.1176470588 0.6254595588 0.6332720588: 0.1176470588 1.625459559 2.633272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.1470588235 0.6548713235 0.6626838235: 0.1470588235 1.654871324 2.662683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.1764705882 0.6842830882 0.6920955882: 0.1764705882 1.684283088 2.692095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.2058823529 0.7136948529 0.7215073529: 0.2058823529 1.713694853 2.721507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.2352941176 0.7431066176 0.7509191176: 0.2352941176 1.743106618 2.750919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.2647058824 0.7725183824 0.7803308824: 0.2647058824 1.772518382 2.780330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.2941176471 0.8019301471 0.8097426471: 0.2941176471 1.801930147 2.809742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.3235294118 0.8313419118 0.8391544118: 0.3235294118 1.831341912 2.839154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.3529411765 0.8607536765 0.8685661765: 0.3529411765 1.860753676 2.868566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 +DEAL::Cell with center 0.7500000000 0.7500000000 0.7500000000 +DEAL::0.5000000000 0.5078125000 0.5156250000: 0.5000000000 1.507812500 2.515625000 error value 0.000000000 error grad 0.000000000 +DEAL::0.5294117647 0.5372242647 0.5450367647: 0.5294117647 1.537224265 2.545036765 error value 0.000000000 error grad 0.000000000 +DEAL::0.5588235294 0.5666360294 0.5744485294: 0.5588235294 1.566636029 2.574448529 error value 0.000000000 error grad 0.000000000 +DEAL::0.5882352941 0.5960477941 0.6038602941: 0.5882352941 1.596047794 2.603860294 error value 0.000000000 error grad 0.000000000 +DEAL::0.6176470588 0.6254595588 0.6332720588: 0.6176470588 1.625459559 2.633272059 error value 0.000000000 error grad 0.000000000 +DEAL::0.6470588235 0.6548713235 0.6626838235: 0.6470588235 1.654871324 2.662683824 error value 0.000000000 error grad 0.000000000 +DEAL::0.6764705882 0.6842830882 0.6920955882: 0.6764705882 1.684283088 2.692095588 error value 0.000000000 error grad 0.000000000 +DEAL::0.7058823529 0.7136948529 0.7215073529: 0.7058823529 1.713694853 2.721507353 error value 0.000000000 error grad 0.000000000 +DEAL::0.7352941176 0.7431066176 0.7509191176: 0.7352941176 1.743106618 2.750919118 error value 0.000000000 error grad 0.000000000 +DEAL::0.7647058824 0.7725183824 0.7803308824: 0.7647058824 1.772518382 2.780330882 error value 0.000000000 error grad 0.000000000 +DEAL::0.7941176471 0.8019301471 0.8097426471: 0.7941176471 1.801930147 2.809742647 error value 0.000000000 error grad 0.000000000 +DEAL::0.8235294118 0.8313419118 0.8391544118: 0.8235294118 1.831341912 2.839154412 error value 0.000000000 error grad 0.000000000 +DEAL::0.8529411765 0.8607536765 0.8685661765: 0.8529411765 1.860753676 2.868566176 error value 0.000000000 error grad 0.000000000 +DEAL:: +DEAL::0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000