From: Marco Feder Date: Wed, 3 May 2023 08:52:55 +0000 (+0200) Subject: fix empty quad in MappingInfo::reinit_cells() X-Git-Tag: v9.5.0-rc1~269^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7381a46466bf164d42b0bdda80792d9186fc4019;p=dealii.git fix empty quad in MappingInfo::reinit_cells() --- diff --git a/include/deal.II/non_matching/mapping_info.h b/include/deal.II/non_matching/mapping_info.h index 7a106c18ac..6312d22c24 100644 --- a/include/deal.II/non_matching/mapping_info.h +++ b/include/deal.II/non_matching/mapping_info.h @@ -402,7 +402,7 @@ namespace NonMatching std::vector> quadrature_vector(n_cells); for (unsigned int cell_index = 0; cell_index < n_cells; ++cell_index) quadrature_vector[cell_index] = - Quadrature(quadrature_vector[cell_index].get_points()); + Quadrature(unit_points_vector[cell_index]); reinit_cells(cell_iterator_range, quadrature_vector, n_unfiltered_cells); } diff --git a/tests/non_matching/mapping_info_02.cc b/tests/non_matching/mapping_info_02.cc new file mode 100644 index 0000000000..fb0ca11143 --- /dev/null +++ b/tests/non_matching/mapping_info_02.cc @@ -0,0 +1,120 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 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 the NonMatching::MappingInfo class together with FEPointEvaluation when + * `reinit_cells` is called with a vector of vector of qpoints, rather than a + * vector of quadratures. + */ + +#include + +#include + +#include + +#include +#include + +#include + +#include +#include +#include +#include + +#include + +#include "../tests.h" + +void +test_vec_of_qpoints() +{ + constexpr unsigned int dim = 2; + constexpr unsigned int degree = 1; + + FE_Q fe_q(degree); + + Triangulation tria; + + MappingQ mapping(degree); + + GridGenerator::subdivided_hyper_cube(tria, 4); + + DoFHandler dof_handler(tria); + + dof_handler.distribute_dofs(fe_q); + + Functions::SignedDistance::Sphere level_set; + + Vector level_set_vec(dof_handler.n_dofs()); + + VectorTools::interpolate(dof_handler, level_set, level_set_vec); + + NonMatching::MeshClassifier mesh_classifier(dof_handler, level_set_vec); + mesh_classifier.reclassify(); + + hp::QCollection<1> q_collection((QGauss<1>(degree))); + + NonMatching::DiscreteQuadratureGenerator quadrature_generator( + q_collection, dof_handler, level_set_vec); + + // FEPointEvaluation + NonMatching::MappingInfo mapping_info_cell( + mapping, update_values | update_gradients | update_JxW_values); + + std::vector>> unit_points_vector; + for (const auto &cell : tria.active_cell_iterators()) + { + quadrature_generator.generate(cell); + unit_points_vector.emplace_back( + quadrature_generator.get_inside_quadrature().get_points()); + } + + mapping_info_cell.reinit_cells(tria.active_cell_iterators(), + unit_points_vector); + + Vector src(dof_handler.n_dofs()); + std::vector solution_values_in(fe_q.dofs_per_cell); + + for (auto &v : src) + v = random_value(); + FEPointEvaluation<1, dim, dim, double> fe_point_cell(mapping_info_cell, fe_q); + for (const auto &cell : dof_handler.active_cell_iterators()) + { + cell->get_dof_values(src, + solution_values_in.begin(), + solution_values_in.end()); + + fe_point_cell.reinit(cell->active_cell_index()); + fe_point_cell.evaluate(solution_values_in, + EvaluationFlags::values | + EvaluationFlags::gradients); + for (const auto q : fe_point_cell.quadrature_point_indices()) + { + deallog << fe_point_cell.get_value(q) << std::endl; + deallog << fe_point_cell.get_gradient(q) << std::endl; + } + } +} + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); + initlog(); + test_vec_of_qpoints(); +} diff --git a/tests/non_matching/mapping_info_02.output b/tests/non_matching/mapping_info_02.output new file mode 100644 index 0000000000..28c8091576 --- /dev/null +++ b/tests/non_matching/mapping_info_02.output @@ -0,0 +1,33 @@ + +DEAL::0.704027 +DEAL::-0.860928 0.693937 +DEAL::0.575505 +DEAL::-0.167248 -0.620078 +DEAL::0.553163 +DEAL::-0.0114928 -0.562178 +DEAL::0.488064 +DEAL::-0.543415 1.43809 +DEAL::0.671952 +DEAL::0.333629 -0.950543 +DEAL::0.497412 +DEAL::-1.72995 -0.00467201 +DEAL::0.460992 +DEAL::1.43859 -0.175191 +DEAL::0.678670 +DEAL::0.449138 -0.138598 +DEAL::0.664544 +DEAL::-0.258019 0.891278 +DEAL::0.586666 +DEAL::-0.365003 0.718707 +DEAL::0.436108 +DEAL::-0.302437 -0.0385718 +DEAL::0.471293 +DEAL::-1.75106 0.226205 +DEAL::0.451239 +DEAL::1.78818 -1.47462 +DEAL::0.473064 +DEAL::-0.171464 -2.58544 +DEAL::0.512506 +DEAL::0.0738393 -1.94578 +DEAL::0.515943 +DEAL::-1.71419 -0.781521