From 5bcb2d64858f7de55ce7541a100035587ac602bb Mon Sep 17 00:00:00 2001 From: David Wells Date: Sat, 25 Feb 2023 15:28:13 -0500 Subject: [PATCH] Make Mapping::get_center() work with simplices. --- include/deal.II/fe/mapping.h | 14 +- source/fe/mapping.cc | 12 +- tests/simplex/mapping_01.cc | 105 +++++++++++ tests/simplex/mapping_01.output | 324 ++++++++++++++++++++++++++++++++ 4 files changed, 441 insertions(+), 14 deletions(-) create mode 100644 tests/simplex/mapping_01.cc create mode 100644 tests/simplex/mapping_01.output diff --git a/include/deal.II/fe/mapping.h b/include/deal.II/fe/mapping.h index 99c3fc2109..051406b2c4 100644 --- a/include/deal.II/fe/mapping.h +++ b/include/deal.II/fe/mapping.h @@ -361,7 +361,7 @@ public: const unsigned int face_no) const; /** - * Return the mapped center of a cell. + * Return one of two possible mapped centers of a cell. * * If you are using a (bi-,tri-)linear mapping that preserves vertex * locations, this function simply returns the value also produced by @@ -371,21 +371,21 @@ public: * polynomials, for which the center may not coincide with the average of * the vertex locations. * - * By default, this function returns the push forward of the center of the + * By default, this function returns the push forward of the barycenter of the * reference cell. If the parameter - * @p map_center_of_reference_cell is set to false, than the return value - * will be the average of the vertex locations, as returned by the + * @p map_barycenter_of_reference_cell is set to false, then the returned + * value will be the average of the vertex locations, as returned by the * get_vertices() method. * * @param[in] cell The cell for which you want to compute the center - * @param[in] map_center_of_reference_cell A flag that switches the algorithm - * for the computation of the cell center from + * @param[in] map_barycenter_of_reference_cell A flag that switches the + * algorithm for the computation of the cell center from * transform_unit_to_real_cell() applied to the center of the reference cell * to computing the vertex averages. */ virtual Point get_center(const typename Triangulation::cell_iterator &cell, - const bool map_center_of_reference_cell = true) const; + const bool map_barycenter_of_reference_cell = true) const; /** * Return the bounding box of a mapped cell. diff --git a/source/fe/mapping.cc b/source/fe/mapping.cc index 1da11f2088..7d95cc3caf 100644 --- a/source/fe/mapping.cc +++ b/source/fe/mapping.cc @@ -82,14 +82,12 @@ template Point Mapping::get_center( const typename Triangulation::cell_iterator &cell, - const bool map_center_of_reference_cell) const + const bool map_barycenter_of_reference_cell) const { - if (map_center_of_reference_cell) + if (map_barycenter_of_reference_cell) { - Point reference_center; - for (unsigned int d = 0; d < dim; ++d) - reference_center[d] = .5; - return transform_unit_to_real_cell(cell, reference_center); + return transform_unit_to_real_cell( + cell, cell->reference_cell().template barycenter()); } else { @@ -97,7 +95,7 @@ Mapping::get_center( Point center; for (const auto &v : vertices) center += v; - return center / GeometryInfo::vertices_per_cell; + return center / cell->n_vertices(); } } diff --git a/tests/simplex/mapping_01.cc b/tests/simplex/mapping_01.cc new file mode 100644 index 0000000000..0a9f053278 --- /dev/null +++ b/tests/simplex/mapping_01.cc @@ -0,0 +1,105 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2023 - 2023 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. +// +// --------------------------------------------------------------------- + +// Test mapping centers with simplices + +#include + +#include + +#include +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +template +class Position : public Function +{ +public: + Position() + : Function(dim) + {} + + double + value(const Point &point, const unsigned int component) const + { + if (component == 0) + return point[component] + std::sin(point[dim - 1] * numbers::PI); + return point[component]; + } +}; + +void +test(const unsigned int mapping_degree) +{ + const int dim = 3; + + Triangulation tria; + GridGenerator::subdivided_hyper_cube_with_simplices(tria, 2); + + FE_SimplexP fe(mapping_degree); + FESystem position_fe(fe, dim); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + DoFHandler position_dof_handler(tria); + position_dof_handler.distribute_dofs(position_fe); + + Vector position_vector(position_dof_handler.n_dofs()); + MappingFE mapping_interpolation(FE_SimplexP(1)); + VectorTools::interpolate(mapping_interpolation, + position_dof_handler, + Position(), + position_vector); + + MappingFEField mapping(position_dof_handler, position_vector); + + for (const auto &cell : tria.active_cell_iterators()) + { + deallog << "cell = " << cell << std::endl; + deallog << "mapped barycenter = " << mapping.get_center(cell, true) + << std::endl; + Point manually_mapped_center; + for (unsigned int d = 0; d < dim; ++d) + manually_mapped_center[d] = + Position().value(cell->barycenter(), d); + deallog << "exact mapped barycenter = " << manually_mapped_center + << std::endl; + deallog << "mapped vertex mean = " << mapping.get_center(cell, false) + << std::endl; + } +} + +int +main() +{ + initlog(); + + test(1); // linear mapping + + deallog << std::endl << "quadratic mapping" << std::endl << std::endl; + + test(2); // quadratic mapping +} diff --git a/tests/simplex/mapping_01.output b/tests/simplex/mapping_01.output new file mode 100644 index 0000000000..88e2683c8c --- /dev/null +++ b/tests/simplex/mapping_01.output @@ -0,0 +1,324 @@ + +DEAL::cell = 0.0 +DEAL::mapped barycenter = 0.375000 0.125000 0.125000 +DEAL::exact mapped barycenter = 0.507683 0.125000 0.125000 +DEAL::mapped vertex mean = 0.375000 0.125000 0.125000 +DEAL::cell = 0.1 +DEAL::mapped barycenter = 0.625000 0.375000 0.125000 +DEAL::exact mapped barycenter = 0.757683 0.375000 0.125000 +DEAL::mapped vertex mean = 0.625000 0.375000 0.125000 +DEAL::cell = 0.2 +DEAL::mapped barycenter = 1.12500 0.125000 0.375000 +DEAL::exact mapped barycenter = 1.29888 0.125000 0.375000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.375000 +DEAL::cell = 0.3 +DEAL::mapped barycenter = 0.875000 0.375000 0.375000 +DEAL::exact mapped barycenter = 1.04888 0.375000 0.375000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.375000 +DEAL::cell = 0.4 +DEAL::mapped barycenter = 0.750000 0.250000 0.250000 +DEAL::exact mapped barycenter = 0.957107 0.250000 0.250000 +DEAL::mapped vertex mean = 0.750000 0.250000 0.250000 +DEAL::cell = 0.5 +DEAL::mapped barycenter = 1.12500 0.125000 0.125000 +DEAL::exact mapped barycenter = 1.25768 0.125000 0.125000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.125000 +DEAL::cell = 0.6 +DEAL::mapped barycenter = 0.875000 0.375000 0.125000 +DEAL::exact mapped barycenter = 1.00768 0.375000 0.125000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.125000 +DEAL::cell = 0.7 +DEAL::mapped barycenter = 1.37500 0.125000 0.375000 +DEAL::exact mapped barycenter = 1.54888 0.125000 0.375000 +DEAL::mapped vertex mean = 1.37500 0.125000 0.375000 +DEAL::cell = 0.8 +DEAL::mapped barycenter = 1.62500 0.375000 0.375000 +DEAL::exact mapped barycenter = 1.79888 0.375000 0.375000 +DEAL::mapped vertex mean = 1.62500 0.375000 0.375000 +DEAL::cell = 0.9 +DEAL::mapped barycenter = 1.25000 0.250000 0.250000 +DEAL::exact mapped barycenter = 1.45711 0.250000 0.250000 +DEAL::mapped vertex mean = 1.25000 0.250000 0.250000 +DEAL::cell = 0.10 +DEAL::mapped barycenter = 0.625000 0.625000 0.125000 +DEAL::exact mapped barycenter = 0.757683 0.625000 0.125000 +DEAL::mapped vertex mean = 0.625000 0.625000 0.125000 +DEAL::cell = 0.11 +DEAL::mapped barycenter = 0.375000 0.875000 0.125000 +DEAL::exact mapped barycenter = 0.507683 0.875000 0.125000 +DEAL::mapped vertex mean = 0.375000 0.875000 0.125000 +DEAL::cell = 0.12 +DEAL::mapped barycenter = 0.875000 0.625000 0.375000 +DEAL::exact mapped barycenter = 1.04888 0.625000 0.375000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.375000 +DEAL::cell = 0.13 +DEAL::mapped barycenter = 1.12500 0.875000 0.375000 +DEAL::exact mapped barycenter = 1.29888 0.875000 0.375000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.375000 +DEAL::cell = 0.14 +DEAL::mapped barycenter = 0.750000 0.750000 0.250000 +DEAL::exact mapped barycenter = 0.957107 0.750000 0.250000 +DEAL::mapped vertex mean = 0.750000 0.750000 0.250000 +DEAL::cell = 0.15 +DEAL::mapped barycenter = 0.875000 0.625000 0.125000 +DEAL::exact mapped barycenter = 1.00768 0.625000 0.125000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.125000 +DEAL::cell = 0.16 +DEAL::mapped barycenter = 1.12500 0.875000 0.125000 +DEAL::exact mapped barycenter = 1.25768 0.875000 0.125000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.125000 +DEAL::cell = 0.17 +DEAL::mapped barycenter = 1.62500 0.625000 0.375000 +DEAL::exact mapped barycenter = 1.79888 0.625000 0.375000 +DEAL::mapped vertex mean = 1.62500 0.625000 0.375000 +DEAL::cell = 0.18 +DEAL::mapped barycenter = 1.37500 0.875000 0.375000 +DEAL::exact mapped barycenter = 1.54888 0.875000 0.375000 +DEAL::mapped vertex mean = 1.37500 0.875000 0.375000 +DEAL::cell = 0.19 +DEAL::mapped barycenter = 1.25000 0.750000 0.250000 +DEAL::exact mapped barycenter = 1.45711 0.750000 0.250000 +DEAL::mapped vertex mean = 1.25000 0.750000 0.250000 +DEAL::cell = 0.20 +DEAL::mapped barycenter = 1.12500 0.125000 0.625000 +DEAL::exact mapped barycenter = 1.29888 0.125000 0.625000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.625000 +DEAL::cell = 0.21 +DEAL::mapped barycenter = 0.875000 0.375000 0.625000 +DEAL::exact mapped barycenter = 1.04888 0.375000 0.625000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.625000 +DEAL::cell = 0.22 +DEAL::mapped barycenter = 0.375000 0.125000 0.875000 +DEAL::exact mapped barycenter = 0.507683 0.125000 0.875000 +DEAL::mapped vertex mean = 0.375000 0.125000 0.875000 +DEAL::cell = 0.23 +DEAL::mapped barycenter = 0.625000 0.375000 0.875000 +DEAL::exact mapped barycenter = 0.757683 0.375000 0.875000 +DEAL::mapped vertex mean = 0.625000 0.375000 0.875000 +DEAL::cell = 0.24 +DEAL::mapped barycenter = 0.750000 0.250000 0.750000 +DEAL::exact mapped barycenter = 0.957107 0.250000 0.750000 +DEAL::mapped vertex mean = 0.750000 0.250000 0.750000 +DEAL::cell = 0.25 +DEAL::mapped barycenter = 1.37500 0.125000 0.625000 +DEAL::exact mapped barycenter = 1.54888 0.125000 0.625000 +DEAL::mapped vertex mean = 1.37500 0.125000 0.625000 +DEAL::cell = 0.26 +DEAL::mapped barycenter = 1.62500 0.375000 0.625000 +DEAL::exact mapped barycenter = 1.79888 0.375000 0.625000 +DEAL::mapped vertex mean = 1.62500 0.375000 0.625000 +DEAL::cell = 0.27 +DEAL::mapped barycenter = 1.12500 0.125000 0.875000 +DEAL::exact mapped barycenter = 1.25768 0.125000 0.875000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.875000 +DEAL::cell = 0.28 +DEAL::mapped barycenter = 0.875000 0.375000 0.875000 +DEAL::exact mapped barycenter = 1.00768 0.375000 0.875000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.875000 +DEAL::cell = 0.29 +DEAL::mapped barycenter = 1.25000 0.250000 0.750000 +DEAL::exact mapped barycenter = 1.45711 0.250000 0.750000 +DEAL::mapped vertex mean = 1.25000 0.250000 0.750000 +DEAL::cell = 0.30 +DEAL::mapped barycenter = 0.875000 0.625000 0.625000 +DEAL::exact mapped barycenter = 1.04888 0.625000 0.625000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.625000 +DEAL::cell = 0.31 +DEAL::mapped barycenter = 1.12500 0.875000 0.625000 +DEAL::exact mapped barycenter = 1.29888 0.875000 0.625000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.625000 +DEAL::cell = 0.32 +DEAL::mapped barycenter = 0.625000 0.625000 0.875000 +DEAL::exact mapped barycenter = 0.757683 0.625000 0.875000 +DEAL::mapped vertex mean = 0.625000 0.625000 0.875000 +DEAL::cell = 0.33 +DEAL::mapped barycenter = 0.375000 0.875000 0.875000 +DEAL::exact mapped barycenter = 0.507683 0.875000 0.875000 +DEAL::mapped vertex mean = 0.375000 0.875000 0.875000 +DEAL::cell = 0.34 +DEAL::mapped barycenter = 0.750000 0.750000 0.750000 +DEAL::exact mapped barycenter = 0.957107 0.750000 0.750000 +DEAL::mapped vertex mean = 0.750000 0.750000 0.750000 +DEAL::cell = 0.35 +DEAL::mapped barycenter = 1.62500 0.625000 0.625000 +DEAL::exact mapped barycenter = 1.79888 0.625000 0.625000 +DEAL::mapped vertex mean = 1.62500 0.625000 0.625000 +DEAL::cell = 0.36 +DEAL::mapped barycenter = 1.37500 0.875000 0.625000 +DEAL::exact mapped barycenter = 1.54888 0.875000 0.625000 +DEAL::mapped vertex mean = 1.37500 0.875000 0.625000 +DEAL::cell = 0.37 +DEAL::mapped barycenter = 0.875000 0.625000 0.875000 +DEAL::exact mapped barycenter = 1.00768 0.625000 0.875000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.875000 +DEAL::cell = 0.38 +DEAL::mapped barycenter = 1.12500 0.875000 0.875000 +DEAL::exact mapped barycenter = 1.25768 0.875000 0.875000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.875000 +DEAL::cell = 0.39 +DEAL::mapped barycenter = 1.25000 0.750000 0.750000 +DEAL::exact mapped barycenter = 1.45711 0.750000 0.750000 +DEAL::mapped vertex mean = 1.25000 0.750000 0.750000 +DEAL:: +DEAL::quadratic mapping +DEAL:: +DEAL::cell = 0.0 +DEAL::mapped barycenter = 0.530330 0.125000 0.125000 +DEAL::exact mapped barycenter = 0.507683 0.125000 0.125000 +DEAL::mapped vertex mean = 0.375000 0.125000 0.125000 +DEAL::cell = 0.1 +DEAL::mapped barycenter = 0.780330 0.375000 0.125000 +DEAL::exact mapped barycenter = 0.757683 0.375000 0.125000 +DEAL::mapped vertex mean = 0.625000 0.375000 0.125000 +DEAL::cell = 0.2 +DEAL::mapped barycenter = 1.28033 0.125000 0.375000 +DEAL::exact mapped barycenter = 1.29888 0.125000 0.375000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.375000 +DEAL::cell = 0.3 +DEAL::mapped barycenter = 1.03033 0.375000 0.375000 +DEAL::exact mapped barycenter = 1.04888 0.375000 0.375000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.375000 +DEAL::cell = 0.4 +DEAL::mapped barycenter = 0.957107 0.250000 0.250000 +DEAL::exact mapped barycenter = 0.957107 0.250000 0.250000 +DEAL::mapped vertex mean = 0.750000 0.250000 0.250000 +DEAL::cell = 0.5 +DEAL::mapped barycenter = 1.28033 0.125000 0.125000 +DEAL::exact mapped barycenter = 1.25768 0.125000 0.125000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.125000 +DEAL::cell = 0.6 +DEAL::mapped barycenter = 1.03033 0.375000 0.125000 +DEAL::exact mapped barycenter = 1.00768 0.375000 0.125000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.125000 +DEAL::cell = 0.7 +DEAL::mapped barycenter = 1.53033 0.125000 0.375000 +DEAL::exact mapped barycenter = 1.54888 0.125000 0.375000 +DEAL::mapped vertex mean = 1.37500 0.125000 0.375000 +DEAL::cell = 0.8 +DEAL::mapped barycenter = 1.78033 0.375000 0.375000 +DEAL::exact mapped barycenter = 1.79888 0.375000 0.375000 +DEAL::mapped vertex mean = 1.62500 0.375000 0.375000 +DEAL::cell = 0.9 +DEAL::mapped barycenter = 1.45711 0.250000 0.250000 +DEAL::exact mapped barycenter = 1.45711 0.250000 0.250000 +DEAL::mapped vertex mean = 1.25000 0.250000 0.250000 +DEAL::cell = 0.10 +DEAL::mapped barycenter = 0.780330 0.625000 0.125000 +DEAL::exact mapped barycenter = 0.757683 0.625000 0.125000 +DEAL::mapped vertex mean = 0.625000 0.625000 0.125000 +DEAL::cell = 0.11 +DEAL::mapped barycenter = 0.530330 0.875000 0.125000 +DEAL::exact mapped barycenter = 0.507683 0.875000 0.125000 +DEAL::mapped vertex mean = 0.375000 0.875000 0.125000 +DEAL::cell = 0.12 +DEAL::mapped barycenter = 1.03033 0.625000 0.375000 +DEAL::exact mapped barycenter = 1.04888 0.625000 0.375000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.375000 +DEAL::cell = 0.13 +DEAL::mapped barycenter = 1.28033 0.875000 0.375000 +DEAL::exact mapped barycenter = 1.29888 0.875000 0.375000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.375000 +DEAL::cell = 0.14 +DEAL::mapped barycenter = 0.957107 0.750000 0.250000 +DEAL::exact mapped barycenter = 0.957107 0.750000 0.250000 +DEAL::mapped vertex mean = 0.750000 0.750000 0.250000 +DEAL::cell = 0.15 +DEAL::mapped barycenter = 1.03033 0.625000 0.125000 +DEAL::exact mapped barycenter = 1.00768 0.625000 0.125000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.125000 +DEAL::cell = 0.16 +DEAL::mapped barycenter = 1.28033 0.875000 0.125000 +DEAL::exact mapped barycenter = 1.25768 0.875000 0.125000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.125000 +DEAL::cell = 0.17 +DEAL::mapped barycenter = 1.78033 0.625000 0.375000 +DEAL::exact mapped barycenter = 1.79888 0.625000 0.375000 +DEAL::mapped vertex mean = 1.62500 0.625000 0.375000 +DEAL::cell = 0.18 +DEAL::mapped barycenter = 1.53033 0.875000 0.375000 +DEAL::exact mapped barycenter = 1.54888 0.875000 0.375000 +DEAL::mapped vertex mean = 1.37500 0.875000 0.375000 +DEAL::cell = 0.19 +DEAL::mapped barycenter = 1.45711 0.750000 0.250000 +DEAL::exact mapped barycenter = 1.45711 0.750000 0.250000 +DEAL::mapped vertex mean = 1.25000 0.750000 0.250000 +DEAL::cell = 0.20 +DEAL::mapped barycenter = 1.28033 0.125000 0.625000 +DEAL::exact mapped barycenter = 1.29888 0.125000 0.625000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.625000 +DEAL::cell = 0.21 +DEAL::mapped barycenter = 1.03033 0.375000 0.625000 +DEAL::exact mapped barycenter = 1.04888 0.375000 0.625000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.625000 +DEAL::cell = 0.22 +DEAL::mapped barycenter = 0.530330 0.125000 0.875000 +DEAL::exact mapped barycenter = 0.507683 0.125000 0.875000 +DEAL::mapped vertex mean = 0.375000 0.125000 0.875000 +DEAL::cell = 0.23 +DEAL::mapped barycenter = 0.780330 0.375000 0.875000 +DEAL::exact mapped barycenter = 0.757683 0.375000 0.875000 +DEAL::mapped vertex mean = 0.625000 0.375000 0.875000 +DEAL::cell = 0.24 +DEAL::mapped barycenter = 0.957107 0.250000 0.750000 +DEAL::exact mapped barycenter = 0.957107 0.250000 0.750000 +DEAL::mapped vertex mean = 0.750000 0.250000 0.750000 +DEAL::cell = 0.25 +DEAL::mapped barycenter = 1.53033 0.125000 0.625000 +DEAL::exact mapped barycenter = 1.54888 0.125000 0.625000 +DEAL::mapped vertex mean = 1.37500 0.125000 0.625000 +DEAL::cell = 0.26 +DEAL::mapped barycenter = 1.78033 0.375000 0.625000 +DEAL::exact mapped barycenter = 1.79888 0.375000 0.625000 +DEAL::mapped vertex mean = 1.62500 0.375000 0.625000 +DEAL::cell = 0.27 +DEAL::mapped barycenter = 1.28033 0.125000 0.875000 +DEAL::exact mapped barycenter = 1.25768 0.125000 0.875000 +DEAL::mapped vertex mean = 1.12500 0.125000 0.875000 +DEAL::cell = 0.28 +DEAL::mapped barycenter = 1.03033 0.375000 0.875000 +DEAL::exact mapped barycenter = 1.00768 0.375000 0.875000 +DEAL::mapped vertex mean = 0.875000 0.375000 0.875000 +DEAL::cell = 0.29 +DEAL::mapped barycenter = 1.45711 0.250000 0.750000 +DEAL::exact mapped barycenter = 1.45711 0.250000 0.750000 +DEAL::mapped vertex mean = 1.25000 0.250000 0.750000 +DEAL::cell = 0.30 +DEAL::mapped barycenter = 1.03033 0.625000 0.625000 +DEAL::exact mapped barycenter = 1.04888 0.625000 0.625000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.625000 +DEAL::cell = 0.31 +DEAL::mapped barycenter = 1.28033 0.875000 0.625000 +DEAL::exact mapped barycenter = 1.29888 0.875000 0.625000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.625000 +DEAL::cell = 0.32 +DEAL::mapped barycenter = 0.780330 0.625000 0.875000 +DEAL::exact mapped barycenter = 0.757683 0.625000 0.875000 +DEAL::mapped vertex mean = 0.625000 0.625000 0.875000 +DEAL::cell = 0.33 +DEAL::mapped barycenter = 0.530330 0.875000 0.875000 +DEAL::exact mapped barycenter = 0.507683 0.875000 0.875000 +DEAL::mapped vertex mean = 0.375000 0.875000 0.875000 +DEAL::cell = 0.34 +DEAL::mapped barycenter = 0.957107 0.750000 0.750000 +DEAL::exact mapped barycenter = 0.957107 0.750000 0.750000 +DEAL::mapped vertex mean = 0.750000 0.750000 0.750000 +DEAL::cell = 0.35 +DEAL::mapped barycenter = 1.78033 0.625000 0.625000 +DEAL::exact mapped barycenter = 1.79888 0.625000 0.625000 +DEAL::mapped vertex mean = 1.62500 0.625000 0.625000 +DEAL::cell = 0.36 +DEAL::mapped barycenter = 1.53033 0.875000 0.625000 +DEAL::exact mapped barycenter = 1.54888 0.875000 0.625000 +DEAL::mapped vertex mean = 1.37500 0.875000 0.625000 +DEAL::cell = 0.37 +DEAL::mapped barycenter = 1.03033 0.625000 0.875000 +DEAL::exact mapped barycenter = 1.00768 0.625000 0.875000 +DEAL::mapped vertex mean = 0.875000 0.625000 0.875000 +DEAL::cell = 0.38 +DEAL::mapped barycenter = 1.28033 0.875000 0.875000 +DEAL::exact mapped barycenter = 1.25768 0.875000 0.875000 +DEAL::mapped vertex mean = 1.12500 0.875000 0.875000 +DEAL::cell = 0.39 +DEAL::mapped barycenter = 1.45711 0.750000 0.750000 +DEAL::exact mapped barycenter = 1.45711 0.750000 0.750000 +DEAL::mapped vertex mean = 1.25000 0.750000 0.750000 -- 2.39.5