From: Martin Kronbichler Date: Wed, 7 Oct 2020 11:25:59 +0000 (+0200) Subject: New test case X-Git-Tag: v9.3.0-rc1~1021^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a0609430c1335b76defaf53d3693340298671ea1;p=dealii.git New test case --- diff --git a/tests/mappings/mapping_points_real_to_unit.cc b/tests/mappings/mapping_points_real_to_unit.cc new file mode 100644 index 0000000000..43b81a56de --- /dev/null +++ b/tests/mappings/mapping_points_real_to_unit.cc @@ -0,0 +1,148 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + + + +// on a test case similar to mapping_real_to_unit_q4_curved, check the +// implementation of the many-point interface +// Mapping::transform_points_real_to_unit_cell for both a MappingQGeneric and +// MappingFEField + +#include + +#include + +#include +#include +#include +#include + +#include +#include +#include + +#include + +#include "../tests.h" + + +template +void +test_real_to_unit_cell(const Mapping &mapping) +{ + // define a boundary that fits the the vertices of the hyper cube we're + // going to create below + SphericalManifold boundary; + + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -1, 1); + + // set the boundary indicator for one face of the single cell + triangulation.set_manifold(1, boundary); + triangulation.begin_active()->face(0)->set_boundary_id(1); + + const unsigned int n_points = 5; + std::vector> unit_points(Utilities::pow(n_points, dim)); + + switch (dim) + { + case 1: + for (unsigned int x = 0; x < n_points; ++x) + unit_points[x][0] = static_cast(x) / n_points; + break; + + case 2: + for (unsigned int x = 0; x < n_points; ++x) + for (unsigned int y = 0; y < n_points; ++y) + { + unit_points[y * n_points + x][0] = + static_cast(x) / n_points; + unit_points[y * n_points + x][1] = + static_cast(y) / n_points; + } + break; + + case 3: + for (unsigned int x = 0; x < n_points; ++x) + for (unsigned int y = 0; y < n_points; ++y) + for (unsigned int z = 0; z < n_points; ++z) + { + unit_points[z * n_points * n_points + y * n_points + x][0] = + static_cast(x) / n_points; + unit_points[z * n_points * n_points + y * n_points + x][1] = + static_cast(y) / n_points; + unit_points[z * n_points * n_points + y * n_points + x][2] = + static_cast(z) / n_points; + } + break; + } + + typename Triangulation::active_cell_iterator cell = + triangulation.begin_active(); + std::vector> real_points(unit_points.size()); + for (unsigned int i = 0; i < unit_points.size(); ++i) + real_points[i] = mapping.transform_unit_to_real_cell(cell, unit_points[i]); + std::vector> new_points(unit_points.size()); + ArrayView> points_view(new_points); + mapping.transform_points_real_to_unit_cell(cell, real_points, points_view); + for (unsigned int i = 0; i < unit_points.size(); ++i) + { + // for each of the points, verify that applying the forward map and + // then pull back get the same point again + AssertThrow(unit_points[i].distance(new_points[i]) < 1e-10, + ExcInternalError()); + } + deallog << "OK" << std::endl; +} + + + +template +void +test() +{ + deallog << "dim=" << dim << ", spacedim=" << spacedim << std::endl; + deallog << "MappingQ(1): "; + test_real_to_unit_cell(MappingQGeneric(1)); + deallog << "MappingQ(4): "; + test_real_to_unit_cell(MappingQGeneric(4)); + + deallog << "MappingFEField(FESystem(FE_Q(4))): "; + Triangulation triangulation; + GridGenerator::hyper_cube(triangulation, -1, 1); + + FESystem fe(FE_Q(4), spacedim); + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs(fe); + const ComponentMask mask(spacedim, true); + Vector location_vector(dof_handler.n_dofs()); + VectorTools::get_position_vector(dof_handler, location_vector, mask); + MappingFEField mapping(dof_handler, location_vector, mask); + test_real_to_unit_cell(mapping); +} + + + +int +main() +{ + initlog(); + + test<1, 1>(); + test<2, 2>(); + test<3, 3>(); + + test<1, 2>(); +} diff --git a/tests/mappings/mapping_points_real_to_unit.output b/tests/mappings/mapping_points_real_to_unit.output new file mode 100644 index 0000000000..8ddfa4c656 --- /dev/null +++ b/tests/mappings/mapping_points_real_to_unit.output @@ -0,0 +1,17 @@ + +DEAL::dim=1, spacedim=1 +DEAL::MappingQ(1): OK +DEAL::MappingQ(4): OK +DEAL::MappingFEField(FESystem(FE_Q(4))): OK +DEAL::dim=2, spacedim=2 +DEAL::MappingQ(1): OK +DEAL::MappingQ(4): OK +DEAL::MappingFEField(FESystem(FE_Q(4))): OK +DEAL::dim=3, spacedim=3 +DEAL::MappingQ(1): OK +DEAL::MappingQ(4): OK +DEAL::MappingFEField(FESystem(FE_Q(4))): OK +DEAL::dim=1, spacedim=2 +DEAL::MappingQ(1): OK +DEAL::MappingQ(4): OK +DEAL::MappingFEField(FESystem(FE_Q(4))): OK