From: Jean-Paul Pelteret Date: Wed, 5 Jan 2022 18:25:04 +0000 (+0100) Subject: Add test for filtered iterators and FEInterfaceValues X-Git-Tag: v9.4.0-rc1~667^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=04ef4844c5eceaf41e07c52242faec2f9fbe8bf4;p=dealii.git Add test for filtered iterators and FEInterfaceValues --- diff --git a/tests/grid/filtered_iterator_08.cc b/tests/grid/filtered_iterator_08.cc new file mode 100644 index 0000000000..b9e3ea00c3 --- /dev/null +++ b/tests/grid/filtered_iterator_08.cc @@ -0,0 +1,106 @@ +// --------------------------------------------------------------------- +// +// 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 that filtered iterators can be used with using FEInterfaceValues. + +#include + +#include + +#include +#include + +#include +#include +#include +#include + +#include + +#include "../tests.h" + + +template +void +test() +{ + Triangulation tria; + GridGenerator::hyper_cube(tria, 0, 1, true); + tria.refine_global(1); + + const FE_Q fe(1); + const QGauss quadrature(1); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + FEInterfaceValues fe_interface_values(fe, + quadrature, + update_default); + const unsigned int invalid_subface = numbers::invalid_unsigned_int; + const unsigned int face_index = 3; + + for (const auto &cell : + filter_iterators(dof_handler.active_cell_iterators(), + IteratorFilters::ActiveFEIndexEqualTo(0))) + { + if (cell == dof_handler.begin_active()) + { + const auto cell_neighbor = cell->neighbor(face_index); + + // This commented out code block identifies that the filtered iterator + // is of a different type to that returned by cell->neighbor() + // + // using C1 = typename std::decay::type; + // using C2 = typename std::decay::type; + // static_assert(std::is_same::value, + // "Cell iterator types are not identical."); + + fe_interface_values.reinit(cell, + face_index, + invalid_subface, + cell->neighbor(face_index), + cell->neighbor_of_neighbor(face_index), + invalid_subface); + } + } + + for (const auto &cell : dof_handler.active_cell_iterators() | + IteratorFilters::ActiveFEIndexEqualTo(0)) + { + if (cell == dof_handler.begin_active()) + { + fe_interface_values.reinit(cell, + face_index, + invalid_subface, + cell->neighbor(face_index), + cell->neighbor_of_neighbor(face_index), + invalid_subface); + } + } +} + + +int +main() +{ + initlog(); + + test<2, 2>(); + test<2, 3>(); + test<3, 3>(); + + deallog << "OK" << std::endl; +} diff --git a/tests/grid/filtered_iterator_08.output b/tests/grid/filtered_iterator_08.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/grid/filtered_iterator_08.output @@ -0,0 +1,2 @@ + +DEAL::OK